SymPy 中更难的方程（带有导数和积分）和 ConditionSet

Question

我打算计算 b（它也是 Xo 轴上的 x），从 0 到 x 的曲线（函数）长度等于 1。

通过了解：https://www.mathsisfun.com/calculus/arc-length.html

（从 0 到 b 的积分） ∫ (1 + ((f'(x))^2)^(1/2) dx = 1

还有那个：

（从 a 到 b 的积分） ∫ f(x)dx = F(b) - F(a)

我们可以通过

计算

1 - F(0) + F(b) = 0 ，现在这是一个关于 x 的方程式，因为正如我所说，b 是 Xo 轴上的 x。

所以现在我尝试了 f(x) = x**3（完整代码将在下面）

F(b)等于这个怪物：https://www.wolframalpha.com/input/?i=integral&assumption=%7B%22C%22%2C+%22integral%22%7D+-%3E+%7B%22Calculator%22%7D&assumption=%7B%22F%22%2C+%22Integral%22%2C+%22integrand%22%7D+-%3E%22%281+%2B+9x%5E4%29%5E%281%2F2%29%22

我从 SymPy 得到的只是 ConditionSet 而不是数字。当然ConditionSet不能被evalf()

求值

所以这是我的问题：

我数学有错吗？
我的代码有错吗？如何改进？
SymPy 足以计算这个吗？
我是不是误解了文档？

代码：

from __future__ import division
import matplotlib.pyplot as plt
from sympy import *

x, y, z = symbols('x y z', real=True)

function1 = x**3

Antiderivative1 = integrate((1+(diff(function1))**2)**(1/2), x)

b = solveset(Eq(1 + Antiderivative1.subs(x, 0).evalf() - Antiderivative1, 0), x)

print(b)

这就是输出：

ConditionSet(x, Eq(x*hyper((-0.5, 1/4), (5/4,), 9*x**4*exp_polar(I*pi)) - 4.0*gamma(5/4)/gamma(1/4), 0), Complexes)

提前致谢，对语法错误深表歉意。

Answer 1

请注意，您应该使用 S(1)/2 或 Rational(1, 2)（或 sqrt）而不是 1/2，这将使您在 [=28 中得到 float =].这样我们就有了

In [16]: integrand = sqrt(1 + ((x**3).diff(x))**2)                                                                                

In [17]: integrand                                                                                                                
Out[17]: 
   __________
  ╱    4     
╲╱  9⋅x  + 1 

In [18]: antiderivative = integrand.integrate(x)                                                                                  

In [19]: antiderivative                                                                                                           
Out[19]: 
          ┌─  ⎛-1/2, 1/4 │    4  ⅈ⋅π⎞
x⋅Γ(1/4)⋅ ├─  ⎜          │ 9⋅x ⋅ℯ   ⎟
         2╵ 1 ⎝   5/4    │          ⎠
─────────────────────────────────────
               4⋅Γ(5/4)

虽然这与 Wolfram Alpha 的结果形式不同，但它很可能是相同的函数（直到加法常量）。从这个结果或 Wolfram Alpha 上的结果，我非常怀疑你会找到一个分析解决方案（使用 SymPy 或其他任何东西）。

但是您可以找到数值解。不幸的是，SymPy 的 lambdify 函数中存在一个错误，这意味着 nsolve 不适用于此函数：

In [22]: nsolve(equation, x, 1)                                                                                                   
...
NameError: name 'exp_polar' is not defined

不过我们可以自己用牛顿步来做：

In [76]: f = equation.lhs                                                                                                         

In [77]: fd = f.diff(x)                                                                                                           

In [78]: newton = lambda xi: (xi - f.subs(x, xi)/fd.subs(x, xi)).evalf()                                                          

In [79]: xj = 1.0                                                                                                                 

In [80]: xj = newton(xj); print(xj)                                                                                               
0.826749667942050

In [81]: xj = newton(xj); print(xj)                                                                                               
0.791950624620750

In [82]: xj = newton(xj); print(xj)                                                                                               
0.790708415511451

In [83]: xj = newton(xj); print(xj)                                                                                               
0.790706893629886

In [84]: xj = newton(xj); print(xj)                                                                                               
0.790706893627605

In [85]: xj = newton(xj); print(xj)                                                                                               
0.790706893627605

SymPy 中更难的方程（带有导数和积分）和 ConditionSet

Harder equations (with Derivatives and Integrals) and ConditionSet in SymPy

python

math

symbolic-math

sympy

integral