2
$\begingroup$

I would like to solve the following reaction-diffusion problem in Mathematica using NDSolve:

DiffConst = 150;
kon = 0.14/60;
koff = 0.18/60;
leftBoundary = 0.;
rightBoundary = 10.;
eqs = {D[L[x, t], t] == DiffConst*D[L[x, t], {x, 2}] + koff*LR[x, t] - 
    kon*L[x, t]*R[x, t] + NeumannValue[0, x == leftBoundary] + 
    NeumannValue[0, x == rightBoundary],
  D[LR[x, t], t] == -koff*LR[x, t] + kon*L[x, t]*R[x, t] + NeumannValue[0, x == leftBoundary] + 
    NeumannValue[0, x == rightBoundary],
  D[R[x, t], t] == koff*LR[x, t] - kon*L[x, t]*R[x, t] + NeumannValue[0, x == leftBoundary] + 
    NeumannValue[0, x == rightBoundary],
  L[x, 0] == 10,
  R[0, 0] == 10,
  LR[x, 0] == 0}
NDSolve[eqs, {L, R, LR}, {t, 0, 100}, {x, 0, 10}]

The idea is to apply no-flux boundary conditions (with the NeumannValue functions) and apply the initial conditions in the last three lines (LR is zero all over, L is 10 everywhere, and R is 10 at x=0). I get the "NDSolve::bcedge: Boundary condition R[0,0]==10 is not specified on a single edge of the boundary of the computational domain" error message. R[0,0] is supposed to be an initial condition, but Mathematica believes it is a boundary condition. What's the problem? Thanks for any help.

$\endgroup$
2
  • $\begingroup$ Perhaps it should be R[0, t] == 10 ? $\endgroup$ Sep 16, 2023 at 16:36
  • $\begingroup$ Or rather R[x,0] == 10, since it's an initial condition. Setting R[0,0]=10 clashes with your Neumann conditions for the left boundary at x=0, since they are independent of t. You might also consider using the predicate in the Neumann condition to excise the point {0,0} to avoid any further conflict using an inequality instead of the default True everywhere. $\endgroup$
    – user87932
    Sep 16, 2023 at 17:52

1 Answer 1

1
$\begingroup$

Here is a way to solve it. Note that there is no need for the NeumannValue, as the 0 NeumannValue is the default. I have changed the initial condition on R and gave a method option to use the finite element method:

DiffConst = 150;
kon = 0.14/60;
koff = 0.18/60;
leftBoundary = 0.;
rightBoundary = 10.;
eqs = {D[L[x, t], t] == 
    DiffConst*D[L[x, t], {x, 2}] + koff*LR[x, t] - 
     kon*L[x, t]*R[x, t](*+NeumannValue[0,x==leftBoundary]+
   NeumannValue[0,x==rightBoundary]*), 
   D[LR[x, t], t] == -koff*LR[x, t] + 
     kon*L[x, t]*R[x, t](*+NeumannValue[0,x==leftBoundary]+
   NeumannValue[0,x==rightBoundary]*), 
   D[R[x, t], t] == 
    koff*LR[x, t] - kon*L[x, t]*R[x, t](*+NeumannValue[0,x==
   leftBoundary]+NeumannValue[0,x==rightBoundary]*),
   L[x, 0] == 10, R[x, 0] == 10, LR[x, 0] == 0};
res = NDSolveValue[eqs, {L, R, LR}, {t, 0, 100}, {x, 0, 10}, 
  Method -> {"PDEDiscretization" -> {"MethodOfLines", 
      "SpatialDiscretization" -> "FiniteElement"}}]
$\endgroup$
1
  • $\begingroup$ Thanks. This is syntactically correct, but not what I wanted to achieve. I wanted R to be 0 everywhere at t=0, except at a certain location. This is the code I replaced R[x,0]==10 with: R[x, 0] == If[x == 0, 10, 0] or R[x, 0] == If[ 0 <= x <= 0.1, 10, 0] $\endgroup$
    – Peter
    Sep 18, 2023 at 15:03

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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