1
$\begingroup$

I want to solve a system of ODEs in sweeping a parameter $\omega$, and in that final conditions need to be carried for the initial condition for the next value of parameter $\omega$. Also, I want to export the data to Excel sheet. The given answer is useful, but exporting the values seems very difficult. Here, I have written my code.

tmax = 30; X0 = {1, 0.9, 2, 0.5, -1}; X0d = {0, 0, 0, 0, 0}; ω = {0.5, 0.6, 0.7, 0.8, 0.9, 1}; mt = ConstantArray[0, {1000, 6}]; K = 2*IdentityMatrix[5]; i = 1; While[i < 5, K[[i, i + 1]] = -1; i++]; i = 1; While[i < 5, K[[i + 1, i]] = -1; i++] Table[ X[t_] := Table[Subscript[x, i][t], {i, 1, 5}]; Xb[t_] := Subscript[x, 5][t]*Sin[ω [[j]] t]; {s} = NDSolve[ {D[X[t], t, t] + K.X[t] == {Xb[t], 1, Xb[t], 1, Xb[t]} * Sin[t], X[0] == X0, X'[0] == X0d}, X[t], {t, 0, tmax} ]; (*Using Subscript[x,1][t] /. s /. t -> tmax and Updating X0 && X0d*) i = 1; While[ i < 6, X0[[i]] = Subscript[x, i][t] /. s /. t -> tmax; X0d[[i]] = D[Subscript[x, i][t] /. s, t] /. t -> tmax; i++ ]; (* To get the list of all values for perticular Subscript[x,i][t] *) Needs["DifferentialEquations`InterpolatingFunctionAnatomy`"]; xvals = InterpolatingFunctionValuesOnGrid[First[First[s]]]; (* But even for x1 also it is not working *) mt[[1 ;; Length[xvals], j]] = xvals;, {j, 1, 6} ] Export["narhari.xls", {"trials" -> mt}] 

Here, I am not able to export the list for $x_i$ values and also D[X[t],t]/.t->tmax also not possible.

I am confused among different types of usage of NDSolve.

$\endgroup$

1 Answer 1

2
$\begingroup$

The problem with this code is to track output NDSolve[] must exactly match input InterpolatingFunctionCoordinates[]. To do this, separate the variables X[t] and X.

tmax = 30; X0 = {1, 0.9, 2, 0.5, -1}; X0d = {0, 0, 0, 0, 0}; \[Omega] = {0.5, 0.6, 0.7, 0.8, 0.9, 1}; mt = ConstantArray[0, {1000, 6}]; M = SparseArray[{{i_, i_} -> 2, {i_, j_} /; Abs[i - j] == 1 -> -1}, {5, 5}]; X = Table[Subscript[x, i][t], {i, 1, 5}]; X1 = Table[Subscript[x, i]'[t], {i, 1, 5}]; Xb[t_, j_] := Subscript[x, 5][t]*Sin[\[Omega][[j]] t]; f[t_, j_] := {Xb[t, j], 1, Xb[t, j], 1, Xb[t, j]}*Sin[t] eq[j_] := {D[X, t, t] == -M.X + f[t, j], (X /. t -> 0) == X0, (X1 /. t -> 0) == X0d}; var = Table[Subscript[x, i], {i, 1, 5}]; var1 = Table[Subscript[x, i]', {i, 1, 5}]; T = {t, 0, tmax}; sol = Table[NDSolve[eq[j], var, T], {j, 6}]; sol1 = Table[NDSolve[eq[j], var1, T], {j, 6}]; fun = X /. sol;fun1 = X1 /. sol1; 

Let's see what we got

Table[Plot[Evaluate[fun[[i]]], T, PlotLegends -> Automatic, PlotLabel -> \[Omega][[i]]], {i, 6}] Table[Plot[Evaluate[fun1[[i]]], T, PlotLegends -> Automatic, PlotLabel -> \[Omega][[i]]], {i, 6}] 

Figure 1

Now we compose a list of variable $x_1$ for export

ifun = var /. sol; ifun1 = First[First[First[ifun]]]; Needs["DifferentialEquations`InterpolatingFunctionAnatomy`"]; coords = First[InterpolatingFunctionCoordinates[ifun1]]; Table[mt[[i]] = Table[First[ifun[[j]]][[1]][#], {j, 6}] &@coords[[i]], {i, Length[coords]}]; Export["C:\\…\\var.xls", mt] 

Check var.xls and compare the function $x_1(t)and its values on thecoords`

Show[Plot[First[First[fun[[1]]]], T], ListPlot[Transpose[{coords, ifun1[coords]}], PlotStyle -> Red]] 

Figure 2

$\endgroup$
2
  • $\begingroup$ Thanks, It is working and It is so compact coding. $\endgroup$ Commented Nov 13, 2019 at 6:36
  • $\begingroup$ @Hari You're welcome! Code can be made even more compact using Verlet algorithm from mathematica.stackexchange.com/questions/208590/… $\endgroup$ Commented Nov 13, 2019 at 11:16

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.