如何将大型三角方程与python集成？答案

【问题标题】：How to integrate large trigonometric equations with python?如何将大型三角方程与python集成？
【发布时间】：2020-04-07 16:23:05
【问题描述】：

在我的程序中，我必须重新创建galerkin-method，我必须使用它进行一些集成，我现在使用sympy.integrate 进行。

我的问题是，如果方程包含多个三角函数或附加的 e 函数（它应该能够解决所有问题），sympy 将永远计算。我也在使用sympy.symbols，因为通过这些集成，我正在创建一个必须求解的方程组。我需要一种反衍生物，而不仅仅是一个值的解决方案。

是否有数字积分方法或其他可以返回准确值的方法？我尝试使用sympy.Integral(equation).evalf(1)，但错误太高或返回的十进制数字太长，这让我回到了太高的运行时间。

函数可能如下所示：

(-4*x**3*sin(2*x) + 2*cos(x)*cos(2*x))*(cos(x) + 1)

我必须整合其中的 20 个。

【问题讨论】：

对我来说，您实际上在寻找哪种类型的“解决方案功能”有点不清楚。您需要符号集成吗？那么数值积分有什么帮助呢？您想使用 galerkin 方法最终以数值方式求解系统吗？
Plus- 你确定存在符号积分吗？
我不知道我要找的词，它是一种反衍生物，我在编辑中添加了它。是的，猜猜是的

标签： python integration sympy

【解决方案1】：

如果您不需要符号积分，并且您知道积分界限，那么您可以进行数值积分。一个务实的起点将是梯形规则。 https://en.wikipedia.org/wiki/Trapezoidal_rule。精度可以通过更精细的步骤（在界限内）任意提高。

通过使用更高阶来提高准确性在编程上会更复杂一些，但在数值上更有效。实现通常只有在计算时间高于编程时间时才有意义）

import math
import matplotlib

sin=math.sin  
cos=math.cos

def f(x):
    return (-4*x**3*sin(2*x) + 2*cos(x)*cos(2*x))*(cos(x) + 1)

x0 = 0              # lower bound
x1 = 2.8            # upper bound (not exactly resolved in this example)
dx = 0.01           # integration step size

F_int = 0
x=x0

while x < x1:
    f_n = f(x)
    f_n_plusdx = f(x+dx)
    dF = 0.5*(f_n+f_n_plusdx)*dx     # calc the trapezoid area
    F_int += dF                      # sum of the integration
    x+=dx                            # next dx.

print(F_int)

【讨论】：

【解决方案2】：

SymPy 的集成函数尝试了多种不同的集成方法，其中一种是 Risch 算法，在某些情况下可能会非常慢。还有一种“手动”集成方法，它不如 Risch 完整，但偶尔会出现极度缓慢的情况。这里有一些描述： https://docs.sympy.org/latest/modules/integrals/integrals.html#internals

您给出的示例中的问题是它陷入了启发式。所以让我们试试“手动”吧：

In [1]: expr = (-4*x**3*sin(2*x) + 2*cos(x)*cos(2*x))*(cos(x) + 1)                                                                

In [2]: expr                                                                                                                      
Out[2]: 
⎛     3                             ⎞             
⎝- 4⋅x ⋅sin(2⋅x) + 2⋅cos(x)⋅cos(2⋅x)⎠⋅(cos(x) + 1)

In [3]: %time anti1 = expr.integrate(x, manual=True)                                                                              
CPU times: user 39.7 s, sys: 232 ms, total: 39.9 s
Wall time: 43.1 s

In [4]: anti1                                                                                                                     
Out[4]: 
   3    3                                                              3                                           ⌠            
8⋅x ⋅cos (x)      3               2                           x   4⋅sin (x)                           sin(4⋅x)     ⎮  2    3    
──────────── + 2⋅x ⋅cos(2⋅x) - 3⋅x ⋅sin(2⋅x) - 3⋅x⋅cos(2⋅x) + ─ - ───────── + 2⋅sin(x) + 2⋅sin(2⋅x) + ──────── - 8⋅⎮ x ⋅cos (x) dx
     3                                                        2       3                                  8         ⌡

这需要 40 秒，但结果并未完全积分：manualintegrate 在其中留下了一个未计算的积分。我们可以通过调用doit 使用普通的integrate 来完成它：

In [5]: %time anti1.doit()                                                                                                        
CPU times: user 4.46 s, sys: 142 ms, total: 4.61 s
Wall time: 4.81 s
Out[5]: 
   3    3                          2    3                                                    2                      3           
8⋅x ⋅cos (x)      3            16⋅x ⋅sin (x)      2           2         2            32⋅x⋅sin (x)⋅cos(x)   112⋅x⋅cos (x)        
──────────── + 2⋅x ⋅cos(2⋅x) - ───────────── - 8⋅x ⋅sin(x)⋅cos (x) - 3⋅x ⋅sin(2⋅x) - ─────────────────── - ───────────── - 3⋅x⋅c
     3                               3                                                        3                  9              

                     3                    2                                      
          x   284⋅sin (x)   112⋅sin(x)⋅cos (x)                           sin(4⋅x)
os(2⋅x) + ─ + ─────────── + ────────────────── + 2⋅sin(x) + 2⋅sin(2⋅x) + ────────
          2        27               9                                       8

所以又花了几秒钟才得到这个结果。正如我们可以验证的那样，这现在是一个完整的反导数：

In [6]: simplify(expr - _.diff(x))                                                                                                
Out[6]: 0

这意味着我们可以使用expr.integrate(x, manual=True).doit() 在大约 50 秒内完成这个特殊的积分。

实际上，如果将这个特定示例从 sin/cos 重写为 exp，则可以在 5 秒内完成：

In [1]: expr = (-4*x**3*sin(2*x) + 2*cos(x)*cos(2*x))*(cos(x) + 1)                                                                

In [2]: %time expr.rewrite(exp).expand().integrate(x).expand().rewrite(sin).simplify()                                            
CPU times: user 5.3 s, sys: 21.2 ms, total: 5.32 s
Wall time: 5.33 s
Out[2]: 
                                 3                                             2                                                
   3             3            2⋅x ⋅cos(3⋅x)      2             2            2⋅x ⋅sin(3⋅x)                                4⋅x⋅cos
2⋅x ⋅cos(x) + 2⋅x ⋅cos(2⋅x) + ───────────── - 6⋅x ⋅sin(x) - 3⋅x ⋅sin(2⋅x) - ───────────── - 12⋅x⋅cos(x) - 3⋅x⋅cos(2⋅x) - ───────
                                    3                                             3                                           9 


(3⋅x)   x                            13⋅sin(3⋅x)   sin(4⋅x)
───── + ─ + 13⋅sin(x) + 2⋅sin(2⋅x) + ─────────── + ────────
        2                                 27          8    

In [3]: simplify(expr - _.diff(x))                                                                                                
Out[3]: 0

虽然这个答案看起来与前一个答案不同，但使用三角恒等式重写三角表达式的方法有无数种，但它们应该等价于一个加法常数（根据反导数的要求）。

【讨论】：