Numerically solve a system with 64 equations

Question

I'm trying to numerically solve a system of 64 equations (see picture below for the equation, where Erf is the error function and Q is a 64 by 64 matrix). But typing every single equation seems very tedious. So I was wondering (more like hoping) if there is some way to implement this formula in mathematica without explicitly typing all the equations, i.e. using the notion below somehow with e.g. NDSolve?

You have to have Q defined beforehand; then you can use Table to construct a List of such equations. — corey979, Nov 18 '16 at 18:26
note the issue that can arise with large numbers of equations: http://mathematica.stackexchange.com/q/131411/2079 — george2079, Nov 18 '16 at 19:16

score 6 · Accepted Answer · edited Apr 13 '17 at 12:55

6

NDSolve can work with this in vector form .. see Vector form using NDSolve

n = 64;
lam = .2;
mu = 3;
Q = RandomReal[1, {n, n}];
vf[x : {_?NumberQ ..}] := 
 Table[Sum[Q[[i, j]] x[[j]] + mu Q[[i, j]] x[[i]] x[[j]] ,  {j, n}] - 
   lam x[[i]], {i, 64}]
ic = Table[RandomReal[{-1, 1}], {n}];

vsoln = NDSolveValue[{x'[t] == Erf[vf[x[t]]], x[0] == ic}, x[t], {t, 0, 2}];

Plot[vsoln, {t, 0, 2}]

aside, I'm a bit annoyed I cant figure how to seperately style the output of the resulting vector valued interpolation function..except like this:

ListPlot[Transpose[Table[vsoln, {t, 0, 2, .01}]], Joined -> True, 
 DataRange -> {0, 2}]

edited Apr 13 '17 at 12:55

Community

1

answered Nov 18 '16 at 19:13

george2079

38,913
1
43
110

Elegant solutions, thanks very much!! – holistic Nov 18 '16 at 19:33
Try ListLinePlot@Head@vsoln. – Michael E2 Nov 19 '16 at 00:52

corey979 · Answer 2 · 2016-11-18T18:46:22.680

4

SeedRandom[1];
n = 2;
μ = RandomReal[];
λ = RandomReal[];
Q = RandomReal[1, {n, n}];

eq = Table[
  x[i]'[t] == 
   Erf @ Sum[Q[[i, j]] x[j][t] + μ Q[[i, j]] x[j][t] x[i][t], {j, 1, n}] - λ x[i][t], {i, 1, n}
      ]

ic = Table[x[i][0] == RandomReal[], {i, 1, n}]

sol = NDSolve[eq~Join~ic, Table[x[i][t], {i, 1, n}], {t, 0, 1}]

func = Table[x[i][t], {i, 1, n}] /. sol

Plot[func, {t, 0, 1}]

Works fine when n = 64:

In your real problem you need to, of course, set the values of parameters of interest (i.e., change all RandomReals to what you want/need).

edited Nov 18 '16 at 18:46

answered Nov 18 '16 at 18:39

corey979

23,947
7
58
101

Nice solution, thank you! – holistic Nov 18 '16 at 19:33

score 1 · Answer 3 · answered Nov 19 '16 at 01:00

Here's another vector form, without any reference to components, and using Dot instead of Sum:

n = 20;
λ = 1/2;
μ = 3;
SeedRandom[0];
Q = RandomReal[1, {n, n}];
ic = Table[RandomReal[{-1, 1}], {n}];

vsoln = NDSolveValue[
   {x'[t] == Erf[Q.x[t] + μ Q.x[t] x[t] - λ x[t]], x[0] == ic}, 
   x, {t, 0, 2}];

ListLinePlot[vsoln, PlotRange -> All]

Mathematica graphics

Numerically solve a system with 64 equations

3 Answers3