1
$\begingroup$

I need to solve 4 coupled PDEs. Thanks to @xzczd's help, a few days ago I can obtain the solution of the first 2 equations with some proper boundary conditions, see this post. Now, these two equations have been coupled with the other two, in which I encountered some new issues.

L = 30; tmax = 20; c = 1/10; \[Epsilon] = 1/10;

eqs = {D[Ar[x, t], t] + 1/4 (Ar[x, t]^3 + Ai[x, t]^2 Ar[x, t]) - D[Ar[x, t], {x, 2}] - 2 Ar[x, t] == 0, 
   D[Ai[x, t], t] + 1/4 (Ai[x, t]^3 + Ar[x, t]^2 Ai[x, t]) - D[Ai[x, t], {x, 2}] - 2 Ai[x, t] == 0,
   D[Br[x, t], t] - D[Br[x, t], {x, 2}] + (Ar[x, t]^2 - Ai[x, t]^2)*(1/4*Br[x, t] + 10*x*D[Ar[x, t], x]
+ (50*x^2 - 1/4)* D[Ai[x, t], x]) + (Ar[x, t]^2 + Ai[x, t]^2)*(1/2*Br[x, t]
+ 10*x*D[Ar[x, t], x] - (50*x^2 + 3/4)*D[Ai[x, t], x])
+ (1/2*Bi[x, t] + 20*x*D[Ai[x, t], x] - (100*x^2 - 1/2)*D[Ar[x, t], x])*Ar[x, t]*Ai[x, t]
- 2*Br[x, t] - 2*D[Ai[x, t], {x, 3}] - 6*D[Ai[x, t], x] - 2*Ar[x, t] == 0,
   D[Bi[x, t], t] - D[Bi[x, t], {x, 2}] + (Ar[x, t]^2 - Ai[x, t]^2)*(-(1/4)*Bi[x, t] - 10*x*D[Ai[x, t], x]
+ (50*x^2 - 1/4)*D[Ar[x, t], x]) + (Ar[x, t]^2 + Ai[x, t]^2)*(1/2*Bi[x, t]
+ 10*x*D[Ai[x, t], x] + (50*x^2 + 3/4)*D[Ar[x, t], x])
+ (1/2*Br[x, t] + 20*x*D[Ar[x, t], x] + (100*x^2 - 1/2)*D[Ai[x, t], x])*Ar[x, t]*Ai[x, t]
- 2*Bi[x, t] + 2*D[Ar[x, t], {x, 3}] + 6*D[Ar[x, t], x] - 2*Ai[x, t] == 0};

subject to the boundary conditions:

bc = {Ar[-L, t] == Ar[L, t], Ai[-L, t] == -Ai[L, t], 
   Br[-L, t] == Br[L, t], Bi[-L, t] == -Bi[L, t],
   \[Epsilon]*Derivative[1, 0][Ar][-L, t] - 
     Ai[-L, t] == \[Epsilon]*Derivative[1, 0][Ar][L, t] - Ai[L, t],
   \[Epsilon]*Derivative[1, 0][Ai][-L, t] + 
     Ar[-L, t] == -\[Epsilon]*Derivative[1, 0][Ai][L, t] - Ar[L, t],
   \[Epsilon]*Derivative[1, 0][Br][-L, t] - 
     Bi[-L, t] == \[Epsilon]*Derivative[1, 0][Br][L, t] - Bi[L, t],
   \[Epsilon]*Derivative[1, 0][Bi][-L, t] + 
     Br[-L, t] == -\[Epsilon]*Derivative[1, 0][Bi][L, t] - Br[L, t]};

The initial conditions are set randomly:

ic = {Ar[x, 0] == iniAr[x], Ai[x, 0] == iniAi[x], Br[x, 0] == iniBr[x], Bi[x, 0] == iniBi[x]};

iniAr = ListInterpolation[RandomReal[{-c, c}, 20], {{-L, L}}];
iniAi = ListInterpolation[RandomReal[{-c, c}, 20], {{-L, L}}];
iniBr = ListInterpolation[RandomReal[{-c, c}, 20], {{-L, L}}];
iniBi = ListInterpolation[RandomReal[{-c, c}, 20], {{-L, L}}];

As mentioned in the previous post, when solving it with NDSolve

sol = NDSolveValue[{eqs, ic, bc}, {Ar[x, t], Ai[x, t], Br[x, t], Bi[x, t]}, {x, -L, L}, {t, 0, tmax}]

it returns an error

NDSolveValue::bcedge: Boundary condition Ai[-30,t]==-Ai[30,t] is not specified on a single edge of the boundary of the computational domain.

So I tried @xzczd's pdetoode

points = 200; domain = {-L, L}; difforder = 4;
grid = Array[# &, points, domain];
ptoofunc = pdetoode[{Ar, Ai, Br, Bi}[x, t], t, grid, difforder];
odebc = Map[ptoofunc, bc, {2}]

del = #[[2 ;; -2]] &;
odeic = del /@ ptoofunc@ic;
ode = del /@ ptoofunc@eqs;
sollst = NDSolveValue[{ode, odeic, odebc}, Table[v[x], {v, {Ar, Ai, Br, Bi}}, {x, grid}], {t, 0, tmax}];

However, I got another error:

NDSolveValue::mconly: For the method NDSolve`IDA, only machine real code is available. Unable to continue with complex values or beyond floating-point exceptions.

What is machine real code? Could anyone help with solving the problem and explaining the error? Thank you in advance!

$\endgroup$
0
$\begingroup$

Strongly related:

Why does NDSolve fail to solve the PDEs and spit out mconly warning?

odebc = With[{sf = 1}, Map[ptoofunc[sf # + D[#, t]] &, bc, {2}]];

odeic = ptoofunc@ic;

var = Table[v[x], {v, {Ar, Ai, Br, Bi}}, {x, grid}];

lhs = D[Flatten[var][t] // Through, t];

{barray, marray} = 
   CoefficientArrays[Flatten@{ode, odebc}, lhs]; // AbsoluteTiming
(* {0.118406, Null} *)

rhs = LinearSolve[marray, -barray]; // AbsoluteTiming
(* {2.87564, Null} *)

sollst = NDSolveValue[{lhs == rhs, odeic}, var, {t, 0, tmax}]; // AbsoluteTiming
(* {3.15959, Null} *)

Necessary explanation has been included in the post above so I'd like not to repeat it here. (Don't miss those links in that post. )

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.