在 python 中使用有限差分
Using finite difference in python
我正在尝试将 Python 与 Numpy 结合使用,使用有限差分法求解基本方程。该代码为我提供了 a 的正确第一个值,即它给了我 a[1];但是,之后的所有其他值都为零?
我不知道我做错了什么,因为它显然适用于第一个值,所以我该如何解决这个问题?
任何想法都会很有帮助。
from numpy import *
import numpy as np
import matplotlib.pyplot as plt
import scipy as sp
from scipy.integrate import odeint
def solver(omega_m, dt):
#t_0, H_0, a_0, dt, n, T; always the same; omega's change
t_0 = 0.0004
a_0 = 0.001
H_0 = 1./13.7
T = 13.7
dt = float(dt)
n = int(round((T - t_0)/dt))
x = zeros(n+1)
t = linspace(t_0, T, n+1)
x[0] = a_0
for i in range (0, n):
x[i+1] = x[i] + (H_0 * ((omega_m)**(1./2.)) * ((x[i])**(-1./2.)) * dt)
return x, t
a, t = solver(omega_m =1, dt=0.001)
print a, t
您的函数在第一次迭代后 returns,因为您的 return 语句在 for 循环内。您应该缩进 return 语句,这样您的循环就不会过早终止:
for i in range (0, n):
x[i+1] = x[i] + (H_0 * ((omega_m)**(1./2.)) * ((x[i])**(-1./2.)) * dt)
return x, t
我正在尝试将 Python 与 Numpy 结合使用,使用有限差分法求解基本方程。该代码为我提供了 a 的正确第一个值,即它给了我 a[1];但是,之后的所有其他值都为零? 我不知道我做错了什么,因为它显然适用于第一个值,所以我该如何解决这个问题? 任何想法都会很有帮助。
from numpy import *
import numpy as np
import matplotlib.pyplot as plt
import scipy as sp
from scipy.integrate import odeint
def solver(omega_m, dt):
#t_0, H_0, a_0, dt, n, T; always the same; omega's change
t_0 = 0.0004
a_0 = 0.001
H_0 = 1./13.7
T = 13.7
dt = float(dt)
n = int(round((T - t_0)/dt))
x = zeros(n+1)
t = linspace(t_0, T, n+1)
x[0] = a_0
for i in range (0, n):
x[i+1] = x[i] + (H_0 * ((omega_m)**(1./2.)) * ((x[i])**(-1./2.)) * dt)
return x, t
a, t = solver(omega_m =1, dt=0.001)
print a, t
您的函数在第一次迭代后 returns,因为您的 return 语句在 for 循环内。您应该缩进 return 语句,这样您的循环就不会过早终止:
for i in range (0, n):
x[i+1] = x[i] + (H_0 * ((omega_m)**(1./2.)) * ((x[i])**(-1./2.)) * dt)
return x, t