# NDSolve: there are more variables than equations, underdet

I am working on using NDSolve to solve coupled ODE. The following is the coupled equations. It works fine for sever orders. I used table to generate like 962 couple variables and equations, where the NDSolve can still work out. However my goal is to further quantize the system and generate an even larger basis, like around 5000 coupled variables. But after I did that, The NDSolve did not work out anymore with error message like "There are more variable, underdet". I am kind of sure (not completely sure) the amount of variable equals the amount of coupled equations. I am wondering is there any upper limit for the size of the coupled ODE the NDSolve can calculate? I am kind of feeling it is up to 1000 variables.

the following is the code:

\[Omega] =
2 \[Pi]*73.5*^3(*2\[Pi]*73.5*^3 Hz, the mass of lithium is low.*);
\[Delta]eff = -2 \[Pi]*200*^3(*-2\[Pi]*200*^3*);
na = 1*10^5*2.3; (*6 Default 2*10^5 *)
\[CapitalDelta]0 = -2 \[Pi]*0.2(*light shift per photon, which is \
-Subscript[g, 0]^2/Subscript[\[CapitalDelta], a]*)(* Default \
-2\[Pi]*0.2 *);
\[Kappa] = 2 \[Pi]*126*^3;(*2\[Pi]*126*^3*)
a = 2.8;(*10*)

\[CapitalDelta]t = 0.0000005;
t0 = 0.00005;

maxkc = 2;
maxkf = maxkc*0.5*1.75;
nn = 7
deltak = maxkf/nn
nc = maxkc*0.5/deltak

\[Mu] = maxkf^2 + maxkf^2; (* chemical potential, depends on maxkf *)

kT = \[Mu]/10; (* set temperature *)
coef = (* Fermi distribution with normalization process *)
Module[{coefficient},
coefficient /.
Flatten[Solve[
Sum[coefficient/(
Exp[((deltak*i)^2 + (deltak*j)^2 - \[Mu])/kT] +
1), {i, -nn - nc - nc, nn + nc + nc}, {j, -nn - nc - nc,
nn + nc + nc}] == 1, coefficient]]]
initialvalues =
Array[ini, {31,
31}, {{-nn - nc - nc, nn + nc + nc}, {-nn - nc - nc,
nn + nc + nc}}];
Do[initialvalues[[i, j]] = (coef/(
Exp[((deltak*(i - 1 - nn - nc - nc))^2 + (deltak*(j - 1 - nn - nc -
nc))^2 - \[Mu])/kT] + 1))^(1/2), {i, 1,
2*(nn + nc + nc) + 1}, {j, 1, 2*(nn + nc + nc) + 1}]
initialvalues // MatrixForm;
initialvalues[[16, 16]]
ListPlot[initialvalues[[16, 16 ;; 31]], Joined -> True, Mesh -> Full]

SetSystemOptions[
"NDSolveOptions" -> {"DefaultScanDiscontinuityTimeConstraint" ->
100., "DefaultSolveTimeConstraint" -> 1000.}]

odes = Flatten[{
Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* 12,
12 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*
Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t], {i,
nn + nc + 1, nn + nc + nc}, {j, nn + nc + 1, nn + nc + nc}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* 12,
11 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i,
nn + nc + 1, nn + nc + nc}, {j, nn + 1, nn + nc}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* 12,
10 to 12,-10 LOOP *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*
Abs[\[Beta][
t]]^2*(Subscript[\[Phi], i*deltak, j*deltak - 2][t] +
Subscript[\[Phi], i*deltak, j*deltak + 2][t]) -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i,
nn + nc + 1, nn + nc + nc}, {j, -nn, nn}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (*
12,-11 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i,
nn + nc + 1, nn + nc + nc}, {j, -nn - nc, -nn - 1}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (*
12,-12 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t], {i,
nn + nc + 1, nn + nc + nc}, {j, -nn - nc - nc, -nn - nc - 1}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* 11,
11 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i,
nn + 1, nn + nc}, {j, nn + 1, nn + nc}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* 11,
10 to 11,-10 LOOP *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*
Abs[\[Beta][
t]]^2*(Subscript[\[Phi], i*deltak, j*deltak - 2][t] +
Subscript[\[Phi], i*deltak, j*deltak + 2][t]) -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i,
nn + 1, nn + nc}, {j, -nn, nn}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (*
11,-11 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*(
Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i,
nn + 1, nn + nc}, {j, -nn - nc, -nn - 1}]],

(*10,10 to -10,-10 BIGGGGGGGGGGGGGGGGGGGGG LOOP*)
Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] == \[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*
Abs[\[Beta][
t]]^2*(Subscript[\[Phi], i*deltak, j*deltak - 2][t] +
Subscript[\[Phi], i*deltak, j*deltak + 2][t]) -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a*(Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
Subscript[\[Phi], i*deltak + 2, j*deltak][t]) +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn,
nn}, {j, -nn, nn}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] == (* -11,-11 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, i*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*(
Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn -
nc, -nn - 1}, {j, -nn - nc, -nn - 1}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* -11,
10 to -11,-10 LOOP *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*
Abs[\[Beta][
t]]^2*(Subscript[\[Phi], i*deltak, j*deltak - 2][t] +
Subscript[\[Phi], i*deltak, j*deltak + 2][t])
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*(
Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]),
{i, -nn - nc, -nn - 1}, {j, -nn, nn}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* -11,
11 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*(
Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn -
nc, -nn - 1}, {j, nn + 1, nn + nc}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] == (* -12,-12 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t], {i, -nn -
nc - nc, -nn - nc - 1}, {j, -nn - nc - nc, -nn - nc - 1}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] == (* -12,-11 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]), {i, -nn -
nc - nc, -nn - nc - 1}, {j, -nn - nc, -nn - 1}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] == (* -12,-10 to -12,
10 LOOP *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*
Abs[\[Beta][
t]]^2*(Subscript[\[Phi], i*deltak, j*deltak + 2][t] +
Subscript[\[Phi], i*deltak, j*deltak - 2][t])
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]), {i, -nn -
nc - nc, -nn - nc - 1}, {j, -nn, nn}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* -12,
11 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]), {i, -nn -
nc - nc, -nn - nc - 1}, {j, nn + 1, nn + nc}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* -12,
12 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t], {i, -nn -
nc - nc, -nn - nc - 1}, {j, nn + nc + 1, nn + nc + nc}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* 11,
12 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]), {i,
nn + 1, nn + nc}, {j, nn + nc + 1, nn + nc + nc}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* 10,
12 to -10,
12 LOOP *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a*(Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
Subscript[\[Phi], i*deltak + 2, j*deltak][t])
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak][t]), {i, -nn,
nn}, {j, nn + nc + 1, nn + nc + nc}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* -11,
12 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]), {i, -nn -
nc, -nn - 1}, {j, nn + nc + 1, nn + nc + nc}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (* 10,
11 to -10,
11 LOOP *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a*(Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
Subscript[\[Phi], i*deltak + 2, j*deltak][t])
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn,
nn}, {j, nn + 1, nn + nc}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (*
10,-11 to -10,-11 LOOP *)\[Omega]*((i*deltak)^2 + (j*
deltak)^2 - 1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a*(Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
Subscript[\[Phi], i*deltak + 2, j*deltak][t])
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn,
nn}, {j, -nn - nc, -nn - 1}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] == (* -11,-12 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i, -nn -
nc, -nn - 1}, {j, -nn - nc - nc, -nn - nc - 1}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (*
10,-12 to -10,-12 LOOP *)\[Omega]*((i*deltak)^2 + (j*
deltak)^2 - 1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a*(Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
Subscript[\[Phi], i*deltak + 2, j*deltak][t])
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i, -nn,
nn}, {j, -nn - nc - nc, -nn - nc - 1}]],

Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] == (*
11,-12 *)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t]
- 1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t]
- 1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t]
+ 1/
2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i,
nn + 1, nn + nc}, {j, -nn - nc - nc, -nn - nc - 1}]],

I*\[Beta]'[t] == (-\[Delta]eff + 1/2*na*\[CapitalDelta]0*

Sum[Re[Subscript[\[Phi], i*deltak, j*deltak][t]*
Conjugate[
Subscript[\[Phi], i*deltak, j*deltak + 2][t]]], {i, -nn -
nc - nc, nn + nc + nc}, {j, -nn - nc - nc, nn}] -
I*\[Kappa])*\[Beta][t] -
I*1/8*na*\[CapitalDelta]0*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5
*(
Sum[Subscript[\[Phi], i*deltak, j*deltak][t](*-12 to 11,12*)
*Conjugate[Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]]
+ Conjugate[Subscript[\[Phi], i*deltak, j*deltak][t]]
*
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t], {i, -nn -
nc - nc, nn + nc}, {j, nn + nc + 1, nn + nc + nc}]

+
Sum[Subscript[\[Phi], i*deltak, j*deltak][
t]*(Conjugate[
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]] +
Conjugate[
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][
t]])(*-12 to 11, -11 to 11*)
+
Conjugate[
Subscript[\[Phi], i*deltak, j*deltak][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][
t]), {i, -nn - nc - nc, nn + nc}, {j, -nn - nc, nn + nc}]

+
Sum[Subscript[\[Phi], i*deltak, j*deltak][
t](*-12 to 11, -12*)
*Conjugate[Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]]
+ Conjugate[Subscript[\[Phi], i*deltak, j*deltak][t]]
*
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t], {i, -nn -
nc - nc, nn + nc}, {j, -nn - nc - nc, -nn - nc - 1}])}]
Text[Style[
"There are " <> ToString[Dimensions[odes]] <>
" coupled differential equations", Red, 18, Italic]]

ics = Flatten[{Table[
Subscript[\[Phi], i*deltak, j*deltak][0] ==
initialvalues[[i + nn + nc + nc + 1,
j + nn + nc + nc + 1]], {i, -nn - nc - nc,
nn + nc + nc}, {j, -nn - nc - nc, nn + nc + nc}],

(*Table[Subscript[\[Phi],i*deltak,j*deltak][0]\[Equal]0.3,{i,-nn,
nn},{j,-nn,nn}],Table[Subscript[\[Phi],i*deltak,j*deltak][
0]\[Equal]0.02,{i,nn+1,nn+nc+nc},{j,-nn-nc-nc,nn+nc+nc}],
Table[Subscript[\[Phi],i*deltak,j*deltak][0]\[Equal]0.02,{i,-nn-
nc-nc,-nn-1},{j,-nn-nc-nc,nn+nc+nc}],
Table[Subscript[\[Phi],i*deltak,j*deltak][0]\[Equal]0.02,{i,-nn,
nn},{j,nn+1,nn+nc+nc}],
Table[Subscript[\[Phi],i*deltak,j*deltak][0]\[Equal]0.02,{i,-nn,
nn},{j,-nn-nc-nc,-nn-1}],*)

\[Beta][0] == 0}];
Text[Style[
"There are " <> ToString[Dimensions[ics]] <>
" initial condition equations", Red, 18, Italic]]

variable =
Flatten[ {\[Beta],
Flatten[Table[
Subscript[\[Phi], i*deltak, j*deltak], {i, -nn - nc - nc,
nn + nc + nc}, {j, -nn - nc - nc, nn + nc + nc}]]}];
Text[Style["There are " <> ToString[Dimensions[ics]] <> " variables",
Red, 18, Italic]]

solve = NDSolve[{odes, ics}, variable, {t, 0, 0.0001}];



I generated coupled equations, initial conditions and variables separately before using NDSolve. nn and nc are used to control the order of the couple equations.

deltak needs to be either 0.1,0.2,0.25,0.5 or 1. (it has to be something periodic within 1) maxkc needs to be 2. maxkf needs to be 1.5 (although here I used 1.75).

Thus there are some combinations between nn and nc. Current the limit is nn=6, thus nc =4, and the order is (2*(6+4+4)+1)^2 = (29)^2 =841;

But if I use nn=15, thus nc =10, the order is (2*(15+10+10)+1)^2 = (71)^2 = 5041, the NDSolve gave error: more varibles and underdet.

Can someone have a look at this problem?

Best

There are several typos in the code. First, initial data array probably should be in a form

initialvalues =  Array[ini, {2 (nn + nc + nc) + 1,  2 (nn + nc + nc) + 1}, {{-nn - nc - nc,      nn + nc + nc}, {-nn - nc - nc, nn + nc + nc}}]


It could be better to use ArrayPlot[initialvalues] for initial data visualization, then we can see all initial values. Second, don't use Subscript[\[Phi], i*deltak, j*deltak] with real induces. It could be better to use it with rational induces. For example, let put deltak = 1/10 then we have

Clear["Global*"]

\[Omega] =
2 \[Pi]*73.5*^3(*2\[Pi]*73.5*^3 Hz,the mass of lithium is low.*);
\[Delta]eff = -2 \[Pi]*200*^3(*-2\[Pi]*200*^3*);
na = 1*10^5*2.3;(*6 Default 2*10^5*)\[CapitalDelta]0 = -2 \
\[Pi]*0.2(*light shift per photon,which \
is-Subscript[g,0]^2/Subscript[\[CapitalDelta],a]*)(*Default-2\[Pi]*0.\
2*);
\[Kappa] =
2 \[Pi]*126*^3;(*2\[Pi]*126*^3*)a = 2.8;(*10*)\[CapitalDelta]t = \
0.0000005;
t0 = 0.00005;

maxkc = 2;
maxkf = maxkc*0.5*1.75;
nn = 15;
deltak = 1/10;
nc = maxkc*1/2/deltak;

\[Mu] = maxkf^2 +
maxkf^2;(*chemical potential,depends on maxkf*)kT = \[Mu]/
10;(*set temperature*)coef =(*Fermi distribution with normalization \
process*)Module[{coefficient},
coefficient /.
Flatten[Solve[
Sum[coefficient/(Exp[((deltak*i)^2 + (deltak*j)^2 - \[Mu])/kT] +
1), {i, -nn - nc - nc, nn + nc + nc}, {j, -nn - nc - nc,
nn + nc + nc}] == 1, coefficient]]];
initialvalues =
Array[ini, {2 (nn + nc + nc) + 1,
2 (nn + nc + nc) + 1}, {{-nn - nc - nc,
nn + nc + nc}, {-nn - nc - nc, nn + nc + nc}}];
Do[initialvalues[[i,
j]] = (coef/(Exp[((deltak*(i - 1 - nn - nc -
nc))^2 + (deltak*(j - 1 - nn - nc - nc))^2 - \[Mu])/
kT] + 1))^(1/2), {i, 1, 2*(nn + nc + nc) + 1}, {j, 1,
2*(nn + nc + nc) + 1}];
initialvalues // MatrixForm;
ArrayPlot[initialvalues, Frame -> False, ColorFunction -> "Rainbow"]


SetSystemOptions[
"NDSolveOptions" -> {"DefaultScanDiscontinuityTimeConstraint" ->
100., "DefaultSolveTimeConstraint" -> 1000.}];

odes = Flatten[{Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*12,
12*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*
Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t], {i,
nn + nc + 1, nn + nc + nc}, {j, nn + nc + 1, nn + nc + nc}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*12,
11*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i,
nn + nc + 1, nn + nc + nc}, {j, nn + 1, nn + nc}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*12,
10 to 12,-10 LOOP*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*
Abs[\[Beta][t]]^2*(Subscript[\[Phi], i*deltak, j*deltak - 2][
t] + Subscript[\[Phi], i*deltak, j*deltak + 2][t]) -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i,
nn + nc + 1, nn + nc + nc}, {j, -nn, nn}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] ==(*12,-11*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i,
nn + nc + 1, nn + nc + nc}, {j, -nn - nc, -nn - 1}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] ==(*12,-12*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t], {i,
nn + nc + 1, nn + nc + nc}, {j, -nn - nc - nc, -nn - nc - 1}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*11,
11*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i,
nn + 1, nn + nc}, {j, nn + 1, nn + nc}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*11,
10 to 11,-10 LOOP*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*
Abs[\[Beta][t]]^2*(Subscript[\[Phi], i*deltak, j*deltak - 2][
t] + Subscript[\[Phi], i*deltak, j*deltak + 2][t]) -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i,
nn + 1, nn + nc}, {j, -nn, nn}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] ==(*11,-11*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i,
nn + 1, nn + nc}, {j, -nn - nc, -nn - 1}]],(*10,
10 to-10,-10 BIGGGGGGGGGGGGGGGGGGGGG LOOP*)
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] == \[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*
Abs[\[Beta][t]]^2*(Subscript[\[Phi], i*deltak, j*deltak - 2][
t] + Subscript[\[Phi], i*deltak, j*deltak + 2][t]) -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a*(Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
Subscript[\[Phi], i*deltak + 2, j*deltak][t]) +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn,
nn}, {j, -nn, nn}]],
Flatten[
Table[I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] ==(*-11,-11*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, i*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn -
nc, -nn - 1}, {j, -nn - nc, -nn - 1}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*-11,
10 to-11,-10 LOOP*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*
Abs[\[Beta][t]]^2*(Subscript[\[Phi], i*deltak, j*deltak - 2][
t] + Subscript[\[Phi], i*deltak, j*deltak + 2][t]) -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn -
nc, -nn - 1}, {j, -nn, nn}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*-11,
11*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn -
nc, -nn - 1}, {j, nn + 1, nn + nc}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] ==(*-12,-12*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t], {i, -nn -
nc - nc, -nn - nc - 1}, {j, -nn - nc - nc, -nn - nc - 1}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] ==(*-12,-11*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]), {i, -nn -
nc - nc, -nn - nc - 1}, {j, -nn - nc, -nn - 1}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*-12,-10 to-12,
10 LOOP*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*
Abs[\[Beta][t]]^2*(Subscript[\[Phi], i*deltak, j*deltak + 2][
t] + Subscript[\[Phi], i*deltak, j*deltak - 2][t]) -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]), {i, -nn -
nc - nc, -nn - nc - 1}, {j, -nn, nn}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*-12,
11*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]), {i, -nn -
nc - nc, -nn - nc - 1}, {j, nn + 1, nn + nc}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*-12,
12*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][t]]*
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t], {i, -nn -
nc - nc, -nn - nc - 1}, {j, nn + nc + 1, nn + nc + nc}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*11,
12*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]), {i,
nn + 1, nn + nc}, {j, nn + nc + 1, nn + nc + nc}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*10,12 to-10,
12 LOOP*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a*(Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
Subscript[\[Phi], i*deltak + 2, j*deltak][t]) +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak][t]), {i, -nn,
nn}, {j, nn + nc + 1, nn + nc + nc}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*-11,
12*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]), {i, -nn -
nc, -nn - 1}, {j, nn + nc + 1, nn + nc + nc}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[t] ==(*10,11 to-10,
11 LOOP*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak - 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a*(Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
Subscript[\[Phi], i*deltak + 2, j*deltak][t]) +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn,
nn}, {j, nn + 1, nn + nc}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] ==(*10,-11 to-10,-11 LOOP*)\[Omega]*((i*deltak)^2 + (j*
deltak)^2 - 1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a*(Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
Subscript[\[Phi], i*deltak + 2, j*deltak][t]) +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak - 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]), {i, -nn,
nn}, {j, -nn - nc, -nn - 1}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] ==(*-11,-12*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak + 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i, -nn -
nc, -nn - 1}, {j, -nn - nc - nc, -nn - nc - 1}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] ==(*10,-12 to-10,-12 LOOP*)\[Omega]*((i*deltak)^2 + (j*
deltak)^2 - 1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a*(Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
Subscript[\[Phi], i*deltak + 2, j*deltak][t]) +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i, -nn,
nn}, {j, -nn - nc - nc, -nn - nc - 1}]],
Flatten[Table[
I*Subscript[\[Phi], i*deltak, j*deltak]'[
t] ==(*11,-12*)\[Omega]*((i*deltak)^2 + (j*deltak)^2 -
1/2*Abs[\[Beta][t]]^2 -
1/2*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a))*
Subscript[\[Phi], i*deltak, j*deltak][t] -
1/4*\[Omega]*Abs[\[Beta][t]]^2*
Subscript[\[Phi], i*deltak, j*deltak + 2][t] -
1/4*\[Omega]*(Tanh[(t - t0)/\[CapitalDelta]t] + 1)*a*
Subscript[\[Phi], i*deltak - 2, j*deltak][t] +
1/2*\[Omega]*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5 Im[\[Beta][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak - 1, j*deltak + 1][t]), {i,
nn + 1, nn + nc}, {j, -nn - nc - nc, -nn - nc - 1}]],
I*\[Beta]'[
t] == (-\[Delta]eff +
1/2*na*\[CapitalDelta]0*
Sum[Re[Subscript[\[Phi], i*deltak, j*deltak][t]*

Conjugate[
Subscript[\[Phi], i*deltak, j*deltak + 2][t]]], {i, -nn -
nc - nc, nn + nc + nc}, {j, -nn - nc - nc, nn}] -
I*\[Kappa])*\[Beta][t] -
I*1/8*na*\[CapitalDelta]0*((Tanh[(t - t0)/\[CapitalDelta]t] + 1)*
a)^0.5*(Sum[
Subscript[\[Phi], i*deltak, j*deltak][t](*-12 to 11,12*)*
Conjugate[
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t]] +
Conjugate[Subscript[\[Phi], i*deltak, j*deltak][t]]*
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][t], {i, -nn -
nc - nc, nn + nc}, {j, nn + nc + 1, nn + nc + nc}] +
Sum[Subscript[\[Phi], i*deltak, j*deltak][
t]*(Conjugate[
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]] +
Conjugate[
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][
t]])(*-12 to 11,-11 to 11*)+
Conjugate[
Subscript[\[Phi], i*deltak, j*deltak][
t]]*(Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t] +
Subscript[\[Phi], i*deltak + 1, j*deltak - 1][
t]), {i, -nn - nc - nc, nn + nc}, {j, -nn - nc,
nn + nc}] +
Sum[Subscript[\[Phi], i*deltak, j*deltak][t](*-12 to 11,-12*)*
Conjugate[
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t]] +
Conjugate[Subscript[\[Phi], i*deltak, j*deltak][t]]*
Subscript[\[Phi], i*deltak + 1, j*deltak + 1][t], {i, -nn -
nc - nc, nn + nc}, {j, -nn - nc - nc, -nn - nc - 1}])}]
Text[Style[
"There are " <> ToString[Dimensions[odes]] <>
" coupled differential equations", Red, 18, Italic]]

ics = Flatten[{Table[
Subscript[\[Phi], i*deltak, j*deltak][0] ==
initialvalues[[i + nn + nc + nc + 1,
j + nn + nc + nc + 1]], {i, -nn - nc - nc,
nn + nc + nc}, {j, -nn - nc - nc, nn + nc + nc}],(*Table[
Subscript[\[Phi],i*deltak,j*deltak][0]\[Equal]0.3,{i,-nn,
nn},{j,-nn,nn}],Table[Subscript[\[Phi],i*deltak,j*deltak][
0]\[Equal]0.02,{i,nn+1,nn+nc+nc},{j,-nn-nc-nc,nn+nc+nc}],Table[
Subscript[\[Phi],i*deltak,j*deltak][0]\[Equal]0.02,{i,-nn-nc-
nc,-nn-1},{j,-nn-nc-nc,nn+nc+nc}],Table[Subscript[\[Phi],i*deltak,
j*deltak][0]\[Equal]0.02,{i,-nn,nn},{j,nn+1,nn+nc+nc}],Table[
Subscript[\[Phi],i*deltak,j*deltak][0]\[Equal]0.02,{i,-nn,
nn},{j,-nn-nc-nc,-nn-1}],*)\[Beta][0] == 0}];
Text[Style[
"There are " <> ToString[Dimensions[ics]] <>
" initial condition equations", Red, 18, Italic]]

variable =
Flatten[{\[Beta],
Flatten[Table[
Subscript[\[Phi], i*deltak, j*deltak], {i, -nn - nc - nc,
nn + nc + nc}, {j, -nn - nc - nc, nn + nc + nc}]]}];
Text[Style["There are " <> ToString[Dimensions[ics]] <> " variables",
Red, 18, Italic]]
solve = NDSolve[{odes, ics}, variable, {t, 0, 0.0001}];

Plot[Abs[\[Beta][t]] /. solve, {t, 0, 0.0001}, PlotRange -> All]
`