Nquad Python невероятно медленный для этого интеграла гладкой функции? - PullRequest
Новая
       56

Nquad Python невероятно медленный для этого интеграла гладкой функции?

0 голосов
/ 26 февраля 2019

Я надеялся максимально использовать функциональность Python для нового проекта, но столкнулся с серьезной проблемой эффективности вычислений, надеюсь, что она может быть решена, поэтому мне не нужно начинать ссылаться на C ++ на этом этапе.Ранняя стадия.

У меня есть следующий код, который связан с интегрированием по двум треугольникам, которые имеют общую вершину и, таким образом, подынтегральное выражение, которое может стать единичным.

Это, однако, не имеет значения,приведенный ниже код «регуляризован», так что нет особого поведения, и теперь интегрирование происходит по области (0,1) ^ 4.Поскольку интеграция гладкая, я ожидал, что nquad интегрирует это очень эффективно, однако на моем компьютере в ноутбуке Jupyter потребуется около 7 секунд (i7-8700K, 3,70 ГГц × 6, 32 ГБ ОЗУ).

Это кажетсяневероятно медленный для интеграции гладкой функции.Я использую Python только 2 месяца, поэтому, возможно, мой код плохо написан с точки зрения производительности?

Есть ли способ ускорить эту интеграцию или мне нужно прибегнуть к C ++ для быстрой интеграции? 

import numpy as np
import numpy.linalg as LA
from scipy.integrate import nquad
import time

##################################################
# Triangle vertices
##################################################
v_a_1 = np.array([[4,0,0]]).T
v_a_2 = np.array([[5,1,0]]).T
v_a_3 = np.array([[4,2,0]]).T

v_b_1 = np.array([[4,0,0]]).T
v_b_2 = np.array([[5,-1,0]]).T
v_b_3 = np.array([[4,-2,0]]).T

##################################################
# g_tau
##################################################
def g_tau():
    J_tau = v_a_2-v_a_1
    J_tau = np.append(J_tau, v_a_3-v_a_2,axis=1)
    G = np.dot(J_tau.T,J_tau)
    return np.sqrt(LA.det(G))

##################################################
# g_t
##################################################
def g_t():
    J_t = v_b_2-v_b_1
    J_t = np.append(J_t, v_b_3-v_b_2,axis=1)
    G = np.dot(J_t.T,J_t)
    return np.sqrt(LA.det(G))

##################################################
# chi_tau
##################################################
def chi_tau(x):
    return v_a_1 + (v_a_2-v_a_1)*x[0] + (v_a_3-v_a_2)*x[1]

##################################################
# chi_t
##################################################
def chi_t(y):
    return v_b_1 + (v_b_2-v_b_1)*y[0] + (v_b_3-v_b_2)*y[1]

##################################################
# k_
##################################################
def k_(x,y):
    return LA.norm(x+y)

##################################################
# k
##################################################
def k(x,y):
    return k_(chi_tau(x),chi_t(y))*g_tau()*g_t()

##################################################
# func
##################################################
def func(z, y, x, w):
    a_1_1 = z;
    a_1_2 = z * w;
    a_2_1 = z * x;
    a_2_2 = z * x * y;

    a_1 = np.array([a_1_1,a_1_2]).T
    a_2 = np.array([a_2_1,a_2_2]).T
    res = k(a_1,a_2)

    return res

start=time.time()
[I_app_b_3_b_3, err] = nquad(func, [[0,1], [0,1], [0,1], [0,1]])
end=time.time()

print(end-start)
...