# NDSolve for non-analytical forms

I came across a problem of coupled differential equations of a non-analytical results. It follows that NDsolve requires an evaluated form of equations to be fed into it for a differential equation to process.But, my functions can't be evaluated analytically. I want to know if there is any possible way to tackle the problem. I have created a sample problem here:

First I create a Complex system just to confirm I will get non-analytical solutions:

Mat[n_, x_, y_] := SparseArray[{Band[{1, 1}, {n, n}] -> {x^3, x I + 5 y^2 + 4, Sqrt[x]},
Band[{1, 2}, {n, n}] -> {y^3, Sqrt[x + I x^3 - y^2 + 4],
Sqrt[x - y^2]}, Band[{3, 1}, {n, n}] -> {I x^3, x + y x^2 + 4, Sqrt[x + y^2]}}]

eval[n_, x_, y_] := Eigensystem[Mat[n, x, y]][[1]]

evec[n_, x_, y_] := Eigensystem[Mat[n, x, y]][[2]]


Secondly, I will form some time-dependent functions so as to cover the complexity of my problem:

val = D[Mat[8, x, y], x];
x2[x1_, t_] := x1 + t^2
y2[y1_, t_] := y1 + t
T1[x1_, y1_] := {I Conjugate[#], #} &@(Conjugate[evec[8, x1, y1][[6]]].SparseArray[
ArrayRules[val] /. {x -> x1, y -> y1}, Dimensions[val]].evec[8, x1, y1][[7]] // N)
T2[t_] := {t, t^2}
T3[x1_, y1_, t_] := T1[x2[x1, t], y2[y1, t]]
T4[x1_, y1_, t_] := Re[eval[8, x2[x1, t], y2[x1, t]][[7]]] - Re[eval[8, x2[x1, t], y2[x1, t]][[8]]]//N


Last, I will list them in the required form of equation and initial conditons for processing:

t0=-5;
eqns[x1_, y1_, t_] := {A1'[t] == (T2[t].T3[x1, y1, t]) A2[t], A2'[t] == A1[t] (T2[t].Conjugate[T3[x1, y1, t]]),
i'[t] == T4[x1, y1, t],A1[t0] == 0, A2[t0] == 1, i[t0] == 0}
sol1 = ParametricNDSolve[eqns[x1, y1, t], {A1[t], A2[t]}, {t, t0, 5}, {x1, y1}]


As you see, say, eqns[x1, y1, t] can't be evaluated unless you provide numerical values of all the parameters. How do we solve the equations in that case. I would be grateful for your help.

(Note: this is a sample just to reflect my problem, feel free to make reasonable changes)

• eqns also has a function i. Have you defined it or should it be solved for? Jun 19, 2020 at 11:05
• It’s a way of writing if there is an integral in your calculations . I want to evaluate integral of ‘T4’ alongside. So, Yes, ‘i’ should be solved for. I have provided initial conditions for it. Jun 19, 2020 at 11:17
• I am really not sure what you are trying to do here. Can you simplify your example to something much simpler? Are you looking for the NumericQ functionality, as described here? Jun 19, 2020 at 15:25
• Sure, Let's just solve for a parameter containing T4. Much more simplified, how can i solve this: ParametricNDSolve[ i'[t] == T4[x1, y1, t], i[t0] == 0, {i[t]}, {t, t0, 5}, {x1, y1}]. Jun 19, 2020 at 17:18
• @MarcoB All I am trying is to solve a differential equation using NDsolve. The equations can only be numerically determined. Above, In the comment I have presented simplified version of it. Please do mention if you want any further info Jun 19, 2020 at 17:55

Your equations are of the form x' = f(x), x = {x1,...,xn}, flow part f(x) can be defined implicitly, i.e. as a black box outside of NDSolve.

This explicit form:

sol = NDSolve[{x'[t]==-y[t]-x[t]^2,y'[t]==2x[t]-y[t]^3,x[0]==y[0]==1},{x,y},{t,20}]
ParametricPlot[Evaluate[{x[t],y[t]}/.sol],{t,0,20}]


can be converted to:

ClearAll[flow] ;
flow[x_,y_] := flow[x,y] = {-y-x^2,2 x-y^3} ;

ClearAll[fx,fy] ;
fx[arg__?NumericQ] := flow[arg][[1]] ;
fy[arg__?NumericQ] := flow[arg][[2]] ;

sol = NDSolve[{x'[t] == fx[x[t],y[t]] ,y'[t]==fy[x[t],y[t]],x[0]==y[0]==1},{x,y},{t,20}] ;
ParametricPlot[Evaluate[{x[t],y[t]}/.sol],{t,0,20}]


## Edit

ClearAll[flow] ;
flow[x1_, y1_, t_,A1_,A2_,i_] :=  flow[x1,y1,t,A1,A2,i] = {
(T2[t].T3[x1, y1, t]) A2,
A1 (T2[t].Conjugate[T3[x1, y1, t]]),
T4[x1, y1, t]
} ;

ClearAll[f1,f2,f3] ;
f1[arg__?NumericQ] := flow[arg][[1]] ;
f2[arg__?NumericQ] := flow[arg][[2]] ;
f3[arg__?NumericQ] := flow[arg][[3]] ;

t0=-5;
(* return only i *)
sol = ParametricNDSolveValue[
{
A1'[t] == f1[x1,y1,t,A1[t],A2[t],i[t]],
A2'[t] == f2[x1,y1,t,A1[t],A2[t],i[t]],
i'[t] == f3[x1,y1,t,A1[t],A2[t],i[t]],
A1[t0] == 0, A2[t0] == 1, i[t0] == 0
},
i,
{t, t0, 5},
{x1, y1}
]
sol[1,1]

• Thank you. When you defined flow[x,y] = {-y-x^2,2 x-y^3} that still becomes an analytical representation. Whereas, if you see in mine it has a set delayed. Could you just possibly obtain ParametricNDSolve[ i'[t] == T4[x1, y1, t], i[t0] == 0, {i[t]}, {t, t0, 5}, {x1, y1}] from my example. Jun 21, 2020 at 6:21
• @Rupesh, you can define flow[x_,y_] := flow[x,y] = whatever, flow[x_,y_] := flow[x,y] part is a kind of memorization, since both fx and fy evaluate flow with the same argument
– I.M.
Jun 21, 2020 at 6:33
• I understand what you are saying but since there are many intermediary steps from defining a flow to NDSolve in my case, how do you go on about it. Plus, I am trying to find a time evolution of a separate function but not x an y itself. I am afraid that I can't make a resemblance between your suggestion and my problem. I would be glad if you could take my own example and at least find an evolution for i[t]. You can change the time intervals or I expect the solution parametric in terms of x and y Jun 21, 2020 at 6:54
• @Rupesh, see the edit.
– I.M.
Jun 21, 2020 at 9:12
• @ I.M I want to thank you a lot and appreciate your time and patience. I was totally unaware of this approach and now it does make sense to me. I tried fitting your first example to my case and messed up, it was bad on my part. Sorry if I behaved impolitely. Jun 21, 2020 at 12:38