I have nonlinear differential equations that require outputs of other functions and these functions require time and some states as inputs.

$ p = sqrt {1- left (1- frac {h} {25} right) ^ 2} $

$ q = p left ( text {x1} ^ 2 sin (t) + text {x2} ^ 2 cos (t) + sin (2t) right) $

$ r = p left ( text {x1} ^ 2 cos (t) + text {x2} ^ 2 sin (t) + cos (2t) right) $

$ alp = tan ^ {- 1} (q, r) $ if q not zero and r not zero

$ f1 = alp * p * q * sqrt {q ^ 2 + r ^ 2} $

$ f2 = alp * p * r * sqrt {q ^ 2 + r ^ 2} $

$ g1 = int_0 ^ { text {hmax}} text {f1} , dh $

$ g2 = int_0 ^ { text {hmax}} text {f2} , dh $

equation1 =$ frac { text {y1} (10 sin ( text {x1}) + cos ( text {x2}))} {1000000} + frac {1} {10} text {y2} ( 0,1 sin ( text {x1}) + 6 cos ( text {x2})) + frac {1} {100} text {y3} ( sin ( text {x1}) + 4 cos ( text {x2})) = int_0 ^ { text {hmax}} text {f1} , dh + text {x1} ^ 2 sin (t) + cos (t) +2.5 $

equation2 =$ frac { text {y1} (5 cos ( text {x1}) + sin ( text {x2}))} {1000} + frac {1} {100} text {y3} ( cos ( text {x1}) + 4 sin ( text {x2})) + frac { text {y2} (2 sin ( text {x2}) + 40 cos ( text {x2} }))}} {1000000} = text {x2} ^ 2 left ( int_0 ^ { text {hmax}} text {f2} , dh right) + sin (t) + 2 $

equation3 =$ frac {1} {10} text {y3} (20 cos ( text {x1} (t)) + sin ( text {x2} (t))) + frac { text {y1 } (4 sin ( text {x1}) + 0.5 sin ( text {x2}))}} {1000000} + frac { text {y2} (4 sin ( text {x1}) + 0.1 cos ( text {x2}))} {1000} = text {x1} sin (t) + cos (t) + text {x2} +1.5 $

$ text {x1} = text {x1} * text {y1} + text {x2} * text {y3} $

$ text {x2} = text {x2} * text {y3} – text {x2} * text {y1} $

$ text {x3} = text {x1} * text {x2} * cos (t) + sin (t) + text {x4} $

$ text {x4} = text {x3} ^ 2 * sin (t) + text {x4} $

$ x1 (0) = 0.01, x2 (0) = 0.1, x3 (0) = 0, x4 (0) = 0 $

$ x $ implies a derivative with time & # 39; t & # 39;

Here is the code I wrote in Mathematica.

```
p(h_) := Sqrt(1 - (1 - 2*(h/50))^2)
q(t_, x2_, x1_, h_) :=(x2^2*Cos(t) + x1^2*Sin(t) + Sin(2*t))*p(h)
r(t_, x2_, x1_, h_) :=(x1^2*Cos(t) + x2^2*Sin(t) + Cos(2*t))*p(h)
alp(t_, x2_, x1_, h_) :=If(q(t, h, x1, x2) == 0 &&r(t, h, x1, x2) == 0,0,ArcTan(q(t, h, x1, x2),r(t, h, x1, x2)))
g(t1_, x11_, x21_) :=
Module({t = t1, x1 = x11,x2 = x21},
h = Range(0, 10, 0.1);
cc = ConstantArray(0,Length(h) - 1);
f1 = cc;
f2 = cc;
For(i = 1,i < Length(h), i++,
h1 = (h((i)) + h((i + 1)))*0.5;
p = Sqrt(1 - (1 - 2*(h1/50))^2);
f = p*Sin(alp(t, x1, x2, h1))*
Sqrt(q(t, x1, x2, h1)^2 + r(t, x1, x2, h1)^2)*
{q(t, x1, x2, h1),r(t, x1, x2, h1), 0}*(-h((i)) + h((i + 1)));
f1((i)) = f((1));
f2((i)) = f((2));
);
g1 = Total(f1,{1});
g2 = Total(f2, {1});
{g1, g2}
)
l1(t_, x1_, x2_) :=(y1*(10*Sin(x1) + Cos(x2)))/
10^6 + (y2*(Sin(x1)*0.1 +6*Cos(x2)))/10 +
(y3*(Sin(x1) + 4*Cos(x2)))/10^2;
l2(t_, x1_, x2_) :=
(y1*(Sin(x2) + 5*Cos(x1)))/10^3 +
(y2*(40*Cos(x2) + 2*Sin(x2)))/
10^6 + (y3*(4*Sin(x2) +
Cos(x1)))/10^2;
l3(t_, x1_, x2_) :=
(y1*(0.5*Sin(x2) + 4*Sin(x1)))/
10^6 + (y2*(0.1*Cos(x2) +
4*Sin(x1)))/10^3 +
(y3*(Sin(x2) + 20*Cos(x1)))/10;
eq(t_, x1_, x2_) := {l1(t_, x1_, x2_) == l4(t_, x1_, x2_),
l2(t_, x1_, x2_) == l5(t_, x1_, x2_), l3(t_, x1_, x2_) == l6(t_, x1_, x2_)};
y11(t_, x1_, x2_) = NSolve(eq(t, x1, x2), {y1, y2,
y3})((1)) /. {Rule -> Set};
ynew1(t_, x1_, x2_) := y11(t, x1, x2)((1));
ynew3(t_, x1_, x2_) := y11(t, x1, x2)((3));
diffeq = {
Derivative(1)(x1)(t) ==x1(t)*ynew1(t, x1(t), x2(t)) -
Cross({ynew1(t, x1(t), x2(t)), 0, ynew3(t, x1(t), x2(t))},
{x1(t), x2(t), 0})((1)),
Derivative(1)(x2)(t) ==
x2(t)*ynew3(t, x1(t), x2(t)) -
Cross({ynew1(t, x1(t), x2(t)), 0, ynew3(t, x1(t), x2(t))},
{x1(t), x2(t), 0})((2)),
Derivative(1)(x3)(t) ==
Sin(t) + x1(t)*x2(t)*Cos(t),
Derivative(1)(x4)(t) == x3(t)^2*Sin(t),
x1(0) == 0.01,
x2(0) == 0.01,
x3(0) == 0,
x4(0) == 0
};
xsol = NDSolve(diffeq, {x1, x2},
{t, 0, 10});
Plot(xsol(t), {t, 0, 10},
PlotRange -> All)
```

I have tried a lot but this code does not work. Please help me with this code.