求解短时间内包含整数表达式的方程式

Question

我正在尝试使用SymPy求解以下包含积分的方程：我试图使用下面的代码来计算整数部分，但是在r

中生成表达式需要花费很长时间。

from sympy import *
mean,std =0,1
Q=250
#defining Cumulative distribution function
def cdf(mean,std):
  t,x = symbols('t,x')
  cdf_eqn = (1/(std*sqrt(2*pi)))*exp(-(((t-mean)**2)/(2*std**2)))
  cdf = Integral(cdf_eqn, (t,-oo,x)).doit()
  return cdf
#defining Probability density function
def pdf(mean,std):
  x = symbols('x')
  pdf = (1/(std*sqrt(2*pi)))*exp(-((( (x - mean)**2)/(2*std**2)))).doit()
  return pdf
#multiplying cdf and pdf
r,x = symbols('r,x')
equation = cdf(mean=0,std=1).subs(x,x)*pdf(mean=0,std=1).subs(x,(r + Q -x))
#getting interating equation over the limits [0,r]
final_equation = Integral(equation, (x,0,r))
#solving the equation
final_equation.doit()

求解方程需要花费大量时间。如何使用SymPy或任何其他包/库（scipy？）在短时间内解决整个方程式]

代表我的朋友发帖。

Answer 1

SymPy似乎很难做到这一点。我在机器上等待了大约2分钟，但没有弹出任何消息。也许在分析上不能解决。

因此，我使用SciPy's root finding algorithm寻求一种数值方法。

import sympy as sp
from scipy.optimize import root_scalar
import time

start_time = time.time()

mean, std = 0, 1
Q = 250
p = 5
w = 2

x = sp.symbols("x")
r_symbol = sp.symbols("r")

pdf = (1 / (std * sp.sqrt(2 * sp.pi))) * sp.exp(-(((x - mean) ** 2) / (2 * std ** 2)))
cdf = sp.erf(sp.sqrt(2) * x / 2) / 2 + 1 / 2  # pre-calculated with integrate(pdf, (x, -oo, t))


def f(r: float) -> float:
    result = sp.N(-p + (p + w * cdf.subs(x, Q)) * cdf.subs(x, r) + \
                  w * sp.Integral(cdf * pdf.subs(x, (r + Q - x)), (x, 0, r)))
    return result


r0 = 1  # initial estimate for the root
bracket = (-10, 10)  # the upper and lower bounds of where the root is
solution = root_scalar(f, x0=r0, bracket=bracket)
print(solution)  # info about the convergence
print(solution.root)  # the actual number

end_time = time.time()
print("Time taken:", end_time - start_time)

哪个为我产生以下输出：
      converged: True
           flag: 'converged'
 function_calls: 14
     iterations: 13
           root: 0.5659488219328516
0.5659488219328516
Time taken: 26.701611518859863

也可以使用MatPlotLib或Desmos上的图在视觉上看到根：
我认为所花费的时间是合理的，因为它必须评估14个非常困难的积分。但是，Desmos几乎没有时间这样做，所以也许根本上有错误。