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)
lnadot(t_) := H(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),
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" :> {(
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.