Getting error when trying to solve two coupled differential equations with NDSolve

Question

Mathematica newbie here. I have been trying to use NDSolve to solve two coupled differential equations. Each equation on its own (minus the coupling terms) solves fine, but when I try them both together I get an error. I apologize in advance if my question is stupid or this is not the place to ask it. Here is my code:

FIRST EQUATION:

xb = 1;
kcat = 1;
v = 1;
kp = 1;
kd = 1;
Dp = 1;
Dc = 1;
kform = 1;

NDSolve[
  {D[P[x, t], t] == Dp*D[P[x, t], x, x] - kcat*(*C[x,t]*)P[x, t] + 
     (kp/2)*(1/Sqrt[π*v])*(Exp[-((x - xb)^2)/v] + Exp[-((x + xb)^2)/v]), 
   (D[P[x, t], x] /. x -> -xb) == 0, 
   (D[P[x, t], x] /. x -> xb) == 0, 
   P[x, 0] == 0}, 
   P[x, t], {x, -xb, xb}, {t, 0, 50}]

Quit[]

SECOND EQUATION:

xb = 1;
kcat = 1;
v = 1;
kp = 1;
kd = 1;
Dp = 1;
Dc = 1;
kform = 1;

NDSolve[
  {D[C[x, t], t] == Dc*D[C[x, t], x, x] - kd*C[x, t] + 
     kform*(1/Sqrt[π*v])*Exp[-(x^2)/v], 
   (D[C[x, t], x] /. x -> -xb) == 0, 
   (D[C[x, t], x] /. x -> xb) == 0, 
   C[x, 0] == 5}, 
  C[x, t], {x, -xb, xb}, {t, 0, 50}]

Quit[]

BOTH:

xb = 1;
kcat = 1;
v = 1;
kp = 1;
kd = 1;
Dp = 1;
Dc = 1;
kform = 1;

NDSolve[
  {D[P[x, t], t] == Dp*D[P[x, t], x, x] - kcat*C[x, t]*P[x, t], 
   D[C[x, t], t] == Dc*D[C[x, t], x, x], 
   (D[P[x, t], x] /. x -> -xb) == 0, 
   (D[P[x, t], x] /. x -> xb) == 0, 
   P[x, 0] == 0, 
   (D[C[x, t], x] /. x -> -xb) == 0, 
   (D[C[x, t], x] /. x -> xb) == 0, 
   C[x, 0] == 5}, 
  {P[x, t], C[x, t]}, {x, -xb, xb}, {t, 0, 50}]

Quit[]

ERROR MESSAGE:

NDSolve::ndode: The equations are not differential equations or initial conditions in the dependent variables {P}.

Answer 1

I rewrote your equations with Q in place of C which eliminated the error message and made them reality solvable, but I don't think the solutions are much good. It is my belief that your boundary conditions are ill- formulated. Since I have no idea what system you are modeling, there is no way for me help you find correct boundary conditions. Indeed, for lack of domain knowledge, there is a good chance I could not help you even you were to describe the physical system behind your equations.

I will say that I look upon the boundary conditions given by

(D[P[x, t], x] /. x -> -xb) == 0
(D[P[x, t], x] /. x -> xb) == 0

with a great deal of suspicion.

So you can see what I found and perhaps draw some conclusions from it, I present my work.

xb = 1;
kcat = 1;
v = 1;
kp = 1;
kd = 1;
Dp = 1;
Dc = 1;
kform = 1;

p = 
  NDSolve[
   {D[P[x, t], t] == Dp*D[P[x, t], x, x] - kcat*P[x, t] + 
      (kp/2)*(1/Sqrt[π*v])*(Exp[-((x - xb)^2)/v] + Exp[-((x + xb)^2)/v]), 
    (D[P[x, t], x] /. x -> -xb) == 0, 
    (D[P[x, t], x] /. x -> xb) == 0, 
    P[x, 0] == 0}, 
   P, {x, -xb, xb}, {t, 0, 50}][[1, 1, 2]];

Plot3D[p[x, t], {x, -1., 1.}, {t, 0., 50.}, PlotRange -> All]

q = 
  NDSolve[
    {D[Q[x, t], t] == Dc*D[Q[x, t], x, x] - kd*Q[x, t] + 
       kform*(1/Sqrt[π*v])*Exp[-(x^2)/v], 
     (D[Q[x, t], x] /. x -> -xb) == 0, 
     (D[Q[x, t], x] /. x -> xb) == 0, 
     Q[x, 0] == 5},
    Q, {x, -xb, xb}, {t, 0, 50}][[1, 1, 2]];

Plot3D[q[x, t], {x, -1., 1.}, {t, 0., 50.}, PlotRange -> All]

{p2, q2} = 
   NDSolve[
     {D[P[x, t], t] == Dp*D[P[x, t], x, x] - kcat*Q[x, t]*P[x, t], 
      D[Q[x, t], t] == Dc*D[Q[x, t], x, x], 
      (D[P[x, t], x] /. x -> -xb) == 0, 
      (D[P[x, t], x] /. x -> xb) == 0, 
      P[x, 0] == 0, 
      (D[Q[x, t], x] /. x -> -xb) == 0, 
      (D[Q[x, t], x] /. x -> xb) == 0, Q[x, 0] == 5}, 
     {P, Q}, {x, -xb, xb}, {t, 0, 50}][[1, All, 2]];

Plot3D[p2[x, t], {x, -1., 1.}, {t, 0., 50.}, PlotRange -> All]

Plot3D[q2[x, t], {x, -1., 1.}, {t, 0., 50.}, PlotRange -> All]