定义用于从数据文件求解方程的函数

Question

我对 python 完全陌生，事实上任何基础编程语言，我都使用 Mathematica 进行我所有的符号和数字计算。我正在学习使用 python 并发现它真的很棒！这是一个我试图解决但毫无头绪的问题！例如我有一个数据文件

0.      1.    
0.01    0.9998000066665778    
0.02    0.9992001066609779    
...     ..

这就是{t, Cos[2t]}。我想根据这些数据定义一个函数，并用它来求解 python 中的方程。我的 Mathematica 直觉告诉我应该像这样定义函数：

    iFunc[x_] = Interpolation[iData, x]

剩下的工作就很简单了。例如

NDSolve[{y''[x] + iFunc[x] y[x] == 0, y[0] == 1, y[1] == 0}, y, {x, 0, 1}]

轻松解方程。（虽然我没有尝试过更复杂的案例）。现在如何在 python 中完成工作以及准确性对我来说是一个重要的问题。那么，现在我想问两个问题。

1.这是 Mathematica 中最准确的方法吗？

2。在 python 中更准确的解决问题的方法是什么？

这是我解决问题的拙劣尝试（使用来自 Whosebug 的大量输入），其中 cos(2t) 的定义有效：

from scipy.integrate import odeint
import numpy as np
import matplotlib.pyplot as plt
from math import cos
from scipy import interpolate

data = np.genfromtxt('cos2t.dat')

T = data[:,0]  #first column
phi = data[:,1]  #second column

f = interpolate.interp1d(T, phi)

tmin = 0.0# There should be a better way to define from the data
dt = 0.01
tmax = 2*np.pi
t = np.arange(tmin, tmax, dt) 

phinew = f(t)   # use interpolation function returned by `interp1d`

"""
def fun(z, t):
    x, y = z
    return np.array([y, -(cos(2*t))*x ])
"""    
def fun(z, t):
    x, y = z
    return np.array([y, -(phinew(t))*x ])    

sol1 = odeint(fun, [1, 0], t)[..., 0]
    
# for checking the plots        
plt.plot(t, sol1, label='sol')
plt.show()

*当我运行具有来自 cos(2t) 数据的插值函数的代码不起作用时...错误消息告诉

Traceback (most recent call last): File "testde.py", line 30, 
     in <module> sol1 = odeint(fun, [1, 0], t)[..., 0] 
File "/home/archimedes/anaconda3/lib/python3.6/site-packages/scip‌y/integrate/odepack.‌py", 
    line 215, in odeint ixpr, mxstep, mxhnil, mxordn, mxords) 
File "testde.py", 
    line 28, in fun return np.array([y, -(phinew(t))*x ])
TypeError: 'numpy.ndarray' object is not callable.

我真的看不懂。请帮忙...

Answer 1

转分题2

有

t = np.arange(tmin, tmax, dt) 

phinew = f(t)   # use interpolation function returned by `interp1d`

相当于

phinew = np.array([ f(s) for s in t])

你构造 phinew 不是作为可调用函数而是作为值数组，将圆数组关闭到数组的插值函数。在导数函数中直接使用标量函数f，

def fun(z, t):
    x, y = z
    return np.array([y, -f(t)*x ])

Answer 2

在 Mathematica 中，通常的方法很简单

iFunc = Interpolation[iData]

Interpolation[iData] 已经 returns 一个函数。

定义用于从数据文件求解方程的函数

defining a fuction to be used in solving equation from a data file

python

wolfram-mathematica

numerical-methods

differential-equations