I want to solve a differential equation in Mathematica using NDsolve, the differential equation is given by:
$y(t)''+3H(t)y(t)'+\frac{\partial V(u(y))}{\partial y}=0$
where $H(t)$ is given by:
$H=\frac{1}{\sqrt{3}a}\sqrt{(y(t)'^2+V(u(y)))}$
The problem is that the potential $V(u(y))=\frac{u^4}{(1+bu^2/a^2)^2}$ doesn't depend explicitly on $y$, it depends on $u$, an the relation between $y(u)$ is given by (and the code that I have, where $yr=y(u)$):
yr[y_] := a*Sqrt[(1 + 6 b)/b]*ArcSinh[(Sqrt[b (1 + 6 b)] u[y])/a] - Sqrt[6]*a*ArcTanh[(Sqrt[6]*b*u[y])/Sqrt[a^2 + b (1 + 6 b) u[y]^2]]
V[y_] := u[y]^4/(1 + (b*u[y]^2)/a^2)^2;
H[t_] = Sqrt[1/(3 a^2) (1/2 y'[t]^2 + V[y])];
eom[t_] =
Simplify[y''[t] + 3 H[t]*y'[t] + D[V[y], u[y]]*(D[yr[y], u[y]])^-1, a^2 + b*y[t]^2 > 0 && a > 0 && b > 0] // Numerator
a = 10^18;
b = 10^4;
solution =
NDSolve[{eom[t] == 0, y[0] == 100*a, y'[0] == 0}, y, {t, Subscript[t, initial], Subscript[t, final]}, Method -> "StiffnessSwitching", MaxSteps -> Automatic]
So the problem comes precisely here, I can not inverse the previous relation $u(y)$ (if it would be the case, then I can introduce this in the potential function $V$ to get $V(y)$ and the differential equation is easily solved), so how can be done to solve the differential equation in terms of the variable $y(t)$?
Also I don't know exactly how define a variable in terms of others to be solved in NDSolve, since actually I have $u(y(t))$.
I put the following numerical conditions: $a=10^{18}$, $b=10^4$, $y(t=0)=100*a$ and $y(t=0)'=0$
Comment:
Maybe can be useful to take into account the chain rule for the potential function $V(u(y))$:
$\frac{\partial V(u(y))}{\partial y}=\frac{\partial V(u(y))}{\partial u}\frac{\partial u}{\partial y}=\frac{\partial V(u(y))}{\partial u} * f(h)$
But still this functions depends on $u$ not in $y$
y == yr[y]
, maybe? $\endgroup$ – Michael E2 May 30 '18 at 17:14a
andb
, if you ploty
vsu
, you will see that over a very large range that essentiallyu[y] = y
which vastly simplifies the problem. If those values are not typical, you can still convert all your equations to u[t] and try to NDSolve them that way. $\endgroup$ – Bill Watts May 31 '18 at 1:45