Python:以数字方式找到积分的主值 [英] Python: Find principal value of an integral numerically
问题描述
我正在使用python数值求解积分:
其中 a(x)可以取任何值；在[-1; 1]的内部或外部的正数，负数和eta是无穷小正数.第二个外部积分会改变a(x)的值
我正在尝试使用 Sokhotski–Plemelj解决此问题定理:
但是，这涉及确定原理值，我在python中找不到任何方法.我知道它是在Matlab中实现的，但是有人知道图书馆还是通过其他方法确定python中的主要值(如果存在主要值)?
您可以使用sympy直接评估积分.它的eta-> 0的实部是主值:
from sympy import *
x, y, eta = symbols('x y eta', real=True)
re(integrate(1/(x - y + I*eta), (x, -1, 1))).simplify().subs({eta: 0})
# -> log(Abs(-y + 1)/Abs(y + 1))
Matlab的符号工具箱int当然会为您提供相同的结果(我不知道Matlab中与此相关的其他工具---请指定您是否知道特定的工具.)
您询问了有关主值的数值计算的问题.答案是，如果您只有一个函数f(y)，而该函数的解析形式或行为未知，则通常无法进行数值计算.您需要了解诸如被积物的极点在哪里以及它们的阶数之类的东西.
另一方面，如果您知道积分的形式为f(y) / (y - y_0)，则scipy.integrate.quad可以为您计算本金，例如:
import numpy as np
from scipy import integrate, special
# P \int_{-1}^1 dx 1/(x - wvar) * (1 + sin(x))
print(integrate.quad(lambda x: 1 + np.sin(x), -1, 1, weight='cauchy', wvar=0))
# -> (1.8921661407343657, 2.426947531830592e-13)
# Check against known result
print(2*special.sici(1)[0])
# -> 1.89216614073
有关详情.请参见此处. /p>
I'm solving the integral numerically using python:
where a(x) can take on any value; positive, negative, inside or outside the the [-1;1] and eta is an infinitesimal positive quantity. There is a second outer integral of which changes the value of a(x)
I'm trying to solve this using the Sokhotski–Plemelj theorem:
However this involves determining the principle value, which I can't find any method to in python. I know it's implemented in Matlab, but does anyone know of either a library or some other way of the determining the principal value in python (if a principle value exists)?
You can use sympy to evaluate the integral directly. Its real part with eta->0 is the principal value:
from sympy import *
x, y, eta = symbols('x y eta', real=True)
re(integrate(1/(x - y + I*eta), (x, -1, 1))).simplify().subs({eta: 0})
# -> log(Abs(-y + 1)/Abs(y + 1))
Matlab's symbolic toolbox int gives you the same result, of course (I'm not aware of other relevant tools in Matlab for this --- please specify if you know a specific one).
You asked about numerical computation of a principal value. The answer there is that if you only have a function f(y) whose analytical form or behavior you don't know, it's in general impossible to compute them numerically. You need to know things such as where the poles of the integrand are and what order they are.
If you on the other hand know your integral is of the form f(y) / (y - y_0), scipy.integrate.quad can compute the principal value for you, for example:
import numpy as np
from scipy import integrate, special
# P \int_{-1}^1 dx 1/(x - wvar) * (1 + sin(x))
print(integrate.quad(lambda x: 1 + np.sin(x), -1, 1, weight='cauchy', wvar=0))
# -> (1.8921661407343657, 2.426947531830592e-13)
# Check against known result
print(2*special.sici(1)[0])
# -> 1.89216614073
See here for details.
这篇关于Python:以数字方式找到积分的主值的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！
