问题描述
我正在尝试在 python 中解决不同电位的耦合 ODE 系统。它适用于特定类型的势(指数),但是一旦势由幂律描述,python 生成的图就完全不连贯,它经常只为所有参数分配零值。我的编码适用于:
div(class='sm:grid sm:grid-cols-3 sm:gap-4 sm:items-start sm:border-t sm:border-gray-200 sm:pt-5')
label.block.text-sm.font-medium.text-gray-700(for='userPhotos' class='sm:mt-px sm:pt-2')
| Review Photos
.mt-1(class='sm:mt-0 sm:col-span-2')
.max-w-lg.flex.justify-center.px-6.pt-5.pb-6.border-2.border-gray-300.border-dashed.rounded-md
.space-y-1.text-center
svg.mx-auto.h-12.w-12.text-gray-400(stroke='currentColor' fill='none' viewBox='0 0 48 48' aria-hidden='true')
path(d='M28 8H12a4 4 0 00-4 4v20m32-12v8m0 0v8a4 4 0 01-4 4H12a4 4 0 01-4-4v-4m32-4l-3.172-3.172a4 4 0 00-5.656 0L28 28M8 32l9.172-9.172a4 4 0 015.656 0L28 28m0 0l4 4m4-24h8m-4-4v8m-12 4h.02' stroke-width='2' stroke-linecap='round' stroke-linejoin='round')
.flex.text-sm.text-gray-600
label.relative.cursor-pointer.bg-white.rounded-md.font-medium.text-indigo-600(for='userPhotos' class='hover:text-indigo-500 focus-within:outline-none focus-within:ring-2 focus-within:ring-offset-2 focus-within:ring-indigo-500')
span Upload files
input#userPhotos.sr-only(type='file' multiple,accept='image/*',name='userPhotos')
p.pl-1 or drag and drop
p.text-xs.text-gray-500
| PNG,JPG,GIF up to 10MB
它不适用于(并打印 RuntimeWarning 消息):
kr1 = 8*np.pi
#rho_m = a**(-3)
#V = np.e**(-kr1*x_1)
#dVdx = -kr1*np.e**(-kr1*x_1)
def quintessence (x,t):
a = x[0]
x_1 = x[1]
x_2 = x[2]
dadt = (kr1*a/np.sqrt(3))*np.sqrt(1/2 * x_2**2 + np.e**(-kr1*x_1) + a**(-3))
dx_1dt = x_2
dx_2dt = -np.sqrt(3)*kr1*np.sqrt(1/2 * x_2**2 + np.e**(-kr1*x_1) + a**(-3))*x_2 + kr1*np.e**(-kr1*x_1)
return[dadt,dx_1dt,dx_2dt]
x0 = [0.0001,0]
t = np.linspace(0,80,1000)
x = odeint(quintessence,x0,t)
a = x[:,0]
x_1 = x[:,1]
x_2 = x[:,2]
plt.plot(t,a)
plt.show()
知道第二段代码有什么问题吗?我对 python 比较陌生,不知道它的局限性。
解决方法
精华暗能量模型(见physics.SE为任意参考)是
x'' + 3*H*x' + V'(x) = 0
H = a'/a
这里,根据所写,
H = kr1/sqrt(3) * sqrt(E) with
E = 0.5*x'^2 + V(x) + a^(-3)
(( which leads to
E' = x'*(x''+V'(x)) - 3*a^(-4)*a' = - 3*H*x'^2 -3*a^(-3)*H
= -3*H*(x'^2+a^(-3))
))
我会在这些部分实施这个系统,因为
def V(x): return ...
def dV(x): return ...
def quintessence(u,t):
a,x,dx = u
E = 0.5*dx**2 + V(x) + a**-3
H = kr1/3**0.5 * E**0.5
da = a*H
ddx = -3*H*dx -dV(x)
return [da,dx,ddx]