# differential equations – Problem with stopping NDSolve after a condition is met

I’m trying to write a piece of Mathematica code that is essentially a differential equation solver that needs to take a specified function $$V(t,q)$$, and then numerically solves the differential equation $$ddot{phi}(t) = -3sqrt{frac{1}{3}(frac{1}{2}dot{phi}^2 + V)}dot{phi} – partial_phi V$$.
The issue I’m having is that I need NDSolve to stop trying to solve when a certain condition is met, specifically when either of the functions $$varepsilon_i = frac{1}{2}(frac{partial_phi V}{V})^2$$ or $$|eta_i| = |frac{partial_phi^2 V}{V}|$$ gets greater than 1. Currently I’m running the code with a test function $$V = Lambda(1-lambda^{-q}phi(t)^q)$$

``````Clear parameters used in the potential

Clear((Lambda), q)

q = 4;

Define potential and derivatives

V(t_, q_) := (CapitalLambda) (1 - (Lambda)^-q (Phi)(t)^q)^2

Subscript(V, (Phi))(t_, q_) := !(
*SubscriptBox(((PartialD)), ((Phi)(t)))(V(t, q)))

Subscript(V, (Phi)(Phi))(t_, q_) := !(
*SubscriptBox(((PartialD)), ((Phi)(t), (Phi)(t)))(V(t, q)))

Define Hubble parameter

H(t_, q_) := Sqrt(1/3 (1/2 (Phi)dot(t)^2 + V(t, q)))

Slow-roll parameters

Subscript((Epsilon), i)(t_, q_) := 1/2 (Subscript(V, (Phi))(t, q)/V(t, q))^2

Subscript((Eta), i)(t_, q_) := Subscript(V, (Phi)(Phi))(t, q)/V(t, q)

Setup the Klein-Gordon equation and the e-folds

(Phi)ddot(t_) := -3 H(t, q) (Phi)dot(t) - Subscript(V, (Phi))(t, q)

(Lambda)range = 2;
(Phi)dotrange = 2;

For((Lambda)set = 1, (Lambda)set < (Lambda)range, (Lambda)set++,
(Lambda) = (Lambda)set*1;
For((Phi)dotset = 1, (Phi)dotset < (Phi)dotrange, (Phi)dotset++,
(Phi)dotinit = (Phi)dotset/100;
sol = NDSolve({(Phi)'(t) == (Phi)dot(t), (Phi)dot'(t) == (Phi)ddot(t),
lna'(t) == lnadot(t), (Phi)(0) == 0, (Phi)dot(0) == (Phi)dotinit,
lna(0) == 0, (CapitalLambda) == 1}, {(Phi), (Phi)dot,
lna, (CapitalLambda)}, {t, 0, 10^16},
Method -> {"EventLocator",
"Event" :> {Subscript((Epsilon), i)(t, q) > 1,
Abs(Subscript((Eta), i)(t, q)) > 1},
"EventAction" :> Throw(tend = t, "StopIntegration")});
))
``````

The issue is that Mathematica does not seem happy with me using the conditions \$varepsilon_i > 1\$ and \$|eta_i|>1\$ and returns the error

``````NDSolve::nbnum1: The function value Abs((32 (CapitalLambda) (Phi)(0.000307756)^6-24 (CapitalLambda) (Phi)(0.000307756)^2 (1-Power(<<2>>)))/((CapitalLambda) (1-(Phi)(<<1>>)^4)^2))>1 is not True or False when the arguments are {0.000307756,0.000177687,3.07592*10^-6,0.00999467,0.577365,0.00999467,-0.0173117}.
``````

but works perfectly fine if I explicitly evaluate \$varepsilon_i \$ and \$eta_i\$ and instead use

``````sol = NDSolve({(Phi)'(t) == (Phi)dot(t), (Phi)dot'(
t) == (Phi)ddot(t),
0, (Phi)dot(0) == (Phi)dotinit, lna(0) == 0,
(CapitalLambda) == 1}, {(Phi), (Phi)dot, lna, (CapitalLambda)}, {t, 0, 10^16},
Method -> {"EventLocator",
"Event" :> {(
2 q^2 (Lambda)^(-2 q) (Phi)(t)^(-2 +
2 q))/(1 - (Lambda)^-q (Phi)(t)^q)^2 > 1,
Abs((2 q^2 (Lambda)^(-2 q) (CapitalLambda) (Phi)(t)^(-2 + 2 q) -
2 (-1 + q) q (Lambda)^-q (CapitalLambda) (Phi)(t)^(-2 + q) (1 - (Lambda)^-q (Phi)(t)^q))/(
(CapitalLambda) (1 - (Lambda)^-q (Phi)(t)^q)^2)) > 1},
"EventAction" :> Throw(tend = t, "StopIntegration")});
``````

I’d like to know if there’s a workaround so that I can just specify $$V$$ and any parameters and not have to explicitly evaluate $$varepsilon_i$$ and $$eta_i$$ . Apologies if I haven’t included all the relevant information to the problem or if I’ve made some basic errors, my experience with Mathematica is quite limited.

Posted on Categories Articles