Simulate this system of $n$ nonlinear DEs?

Questions
, , ,
1

The following code solves a system of n nonlinear DEs in one dimension:

p[i_, j_] = P[j][t] - P[i][t];
A@p[j, i] = (a*p[j, i])/Exp[b*p[j, i]];
R@p[j, i] = c/p[j, i];
summand1 = R@p[j, i] - A@p[j, i];
A@p[i, j] = (a*p[i, j])/Exp[b*p[i, j]];
R@p[i, j] = c/p[i, j];
summand2 = A@p[i, j] - R@p[i, j];
n = 12;

sys = Table[P[i]'[t] == Sum[If[j < i, summand1, 0], {j, 1, n}] + 
Sum[If[i < j, summand2, 0], {j, 1, n}], {i, 1, n}];

sol = ParametricNDSolveValue[{sys, P[1][0] == 1, P[2][0] == 1.01, P[3][0] == 1.02, 
P[4][0] == 5, P[5][0] == 6.8, P[6][0] == 6.9, P[7][0] == 7, P[8][0] == 7.2, P[9][0] == 8.9, 
P[10][0] == 10, P[11][0] == 10.25, P[12][0] == 11}, P /@ Range[12], {t, 0, 100}, {a, b, c}];

Manipulate[Plot[sol[a, b, c][t] // Through // Evaluate, {t, 0, T}, PlotRange -> Automatic, 
ImageSize -> Large, AspectRatio -> 1], {{a, 1}, .1, 10}, {{b, 1}, .1, 10}, 
{{c, .1}, .1, 100}, {{T, 1}, .1, 100}]

However, trying the two-dimensional case, which I tried as follows:

ClearAll["Global`*"]
P[i_][t] = {X[i][t], Y[i][t]}
p[j_, i_] = Norm[P[i][t] - P[j][t]]
A@p[j, i] = (a*(P[i][t] - P[j][t]))/Exp[b*p[j, i]];
R@p[j, i] = (c*(P[i][t] - P[j][t]))/p[j, i]^2;
summand1 = R@p[j, i] - A@p[j, i];
A@p[i, j] = (a*(P[j][t] - P[i][t]))/Exp[b*p[i, j]];    
R@p[i, j] = (c*(P[j][t] - P[i][t]))/p[i, j]^2;    
summand2 = A@p[i, j] - R@p[i, j];
n = 5;

sys = Table[P[i]'[t] == Sum[If[j < i, summand1, 0], {j, 1, n}] + 
 Sum[If[i < j, summand2, 0], {j, 1, n}], {i, 1, n}];

sol = ParametricNDSolveValue[{sys, P[1][0] == {1, 1}, P[2][0] == {1, 6}, P[3][0] == {1, -1}, 
P[4][0] == {5, 3}, P[5][0] == {6, 2}}, P /@ Range[5], {t, 0, 5}, {a, b, c}];

Manipulate[ParametricPlot[sol[a, b, c][t] // Through // Evaluate, {t, 0, T}, 
PlotRange -> Automatic, ImageSize -> Large, AspectRatio -> 1], {{a, 1}, .1, 5}, 
{{b, 1}, .1, 5}, {{c, .1}, .1, 5}, {{T, 1}, 0, 5}]

gives the error of the system being underdetermined.

I think that it has to do with the way I defined the array

P[i_][t]

I don’t see how to improve the code. Can someone help me out of this?

2

It looks so there is no equation for X, Y. Try working code that can be developed further

ClearAll["Global`*"]
P[i_][t] = {X[i][t], Y[i][t]};
p[j_, i_] = Norm[P[i][t] - P[j][t]];
A@p[j, i] = (a*(P[i][t] - P[j][t]))/Exp[b*p[j, i]];
R@p[j, i] = (c*(P[i][t] - P[j][t]))/p[j, i]^2;
summand1 = R@p[j, i] - A@p[j, i];
A@p[i, j] = (a*(P[j][t] - P[i][t]))/Exp[b*p[i, j]];
R@p[i, j] = (c*(P[j][t] - P[i][t]))/p[i, j]^2;
summand2 = A@p[i, j] - R@p[i, j];
n = 5;
sys = Table[{X[i]'[t], Y[i]'[t]} == 
    Sum[If[j < i, summand1, 0], {j, 1, n}] + 
     Sum[If[i < j, summand2, 0], {j, 1, n}], {i, 1, n}];
Xsol = Table[
  ParametricNDSolveValue[{sys, {X[1][0], Y[1][0]} == {1, 1}, {X[2][0],
       Y[2][0]} == {1, 6}, {X[3][0], Y[3][0]} == {1, -1}, {X[4][0], 
      Y[4][0]} == {5, 3}, {X[5][0], Y[5][0]} == {6, 2}}, 
   X[i], {t, 0, 5}, {a, b, c}], {i, 1, n}]; Ysol = 
 Table[ParametricNDSolveValue[{sys, {X[1][0], Y[1][0]} == {1, 
      1}, {X[2][0], Y[2][0]} == {1, 6}, {X[3][0], 
      Y[3][0]} == {1, -1}, {X[4][0], Y[4][0]} == {5, 3}, {X[5][0], 
      Y[5][0]} == {6, 2}}, Y[i], {t, 0, 5}, {a, b, c}], {i, 1, n}];

Manipulate[
 ParametricPlot[{Xsol[[1]][a, b, c][t], Ysol[[1]][a, b, c][t]}, {t, 0,
    5}, PlotRange -> Automatic, ImageSize -> Large, 
  AspectRatio -> 1], {{a, 1}, .1, 5}, {{b, 1}, .1, 5}, {{c, .1}, .1, 
  5}]