Skip to main content
 首页 » 编程设计

python之Mathematica 和 Python 中超越方程的解不匹配

2025年05月04日13itcoder

我正在求解超越方程组:

  1. cos(x)/x = 0.48283 + a*3.46891
  2. cos(y)/y = 0.47814 + b*28.6418
  3. a + b = 1
  4. 1.02 * sinc(x) = 1.03 * sinc(y)

碰巧我试图用两种不同的编程语言(Mathematica 和 Python)来解决上述系统

数学

运行代码

FindRoot[{Cos[x]/x == 0.482828 + a*3.46891,  
  Cos[y]/y == 0.47814 + b*28.6418, a + b == 1,  
  1.02*Sinc[x] == 1.03*Sinc[y]}, {{x, .2}, {y, .2}, {a, 0.3}, {b,  
   0.3}}, PrecisionGoal -> 6] 

返回

{x -> 0.261727, y -> 0.355888, a -> 0.924737, b -> 0.0752628} 

python

运行代码:

import numpy as np 
from scipy.optimize import root 
 
def fuuu(X, en,dv,tri,sti): 
        x, y, a, b = X 
 
        F = [np.cos(x) / x - en-a*dv, 
             np.cos(y) / y - tri-b*sti, 
             a + b - 1, 
             1.02 * np.sinc(x) - 1.03 * np.sinc(y)] 
 
        return F 
    root(fuuu, [0.2, 0.2, 0.3, 0.3], args=(0.482828,3.46891,0.47814,28.6418)).x 

返回

array([ 0.26843418,  0.27872813,  0.89626625,  0.10373375]) 

比较

假设“x”值相同。让我们忽略微小的差异。但是 y 值相差英里!物理意义完全改变。出于某种原因,我相信 Mathematica 的值(value)胜过 Python 的值(value)。

问题:

  1. 为什么计算结果不同?
  2. 现在哪一个是正确的?我必须在 python 中更改什么(假设 python 是有问题的)?

请您参考如下方法:

计算因 sinc 函数不同而不同。

(* Mathematica *) 
In[1] := Sinc[0.26843418] 
Out[1] = 0.988034 
 
# Python 
>>> np.sinc(0.26843418) 
0.88561519683835599 
>>> np.sin(0.26843418) / 0.26843418 
0.98803370932709034 

嗯?嗯let's RTFM

numpy.sinc(x)

Return the sinc function.

The sinc function is sin(πx)/(πx).

糟糕。 NumPy 的 sinc 有一个 different definition比 Mathematica 的 Sinc。

  • Mathematica 的 Sinc 使用非规范化 定义 sin(x)/x。这个定义通常用于数学和物理。
  • NumPy 的sinc 使用规范化 版本sin(πx)/(πx)。这个定义通常用于数字信号处理和信息论。它被称为归一化是因为
    ∫-∞∞ sin(πx)/(πx) dx = 1。

因此,如果您希望 NumPy 产生与 Mathematica 相同的结果,您需要将 x 和 y 除以 np.pi。

def fuuu(X, en,dv,tri,sti): 
    x, y, a, b = X 
    F = [np.cos(x) / x - en-a*dv, 
         np.cos(y) / y - tri-b*sti, 
         a + b - 1, 
         1.02 * np.sinc(x/np.pi) - 1.03 * np.sinc(y/np.pi)]    # <--- 
    return F 
>>> root(fuuu, [0.2, 0.2, 0.3, 0.3], args=(0.482828,3.46891,0.47814,28.6418)).x 
array([ 0.26172691,  0.3558877 ,  0.92473722,  0.07526278])