NDSolve::femper for system of non-linear PDEs

Question

I am trying to solve a system of PDE in 1 dimension. They are as follows:

My system of equations looks like follows:

$$\frac{\partial \bar{c^{T}_E}}{\partial \bar{t}} = \frac{\partial^2 \bar{c^{T}_E}}{\partial \bar{x}^2} + \alpha \frac{\partial^2}{\partial \bar{x}^2}\Big(\frac{\bar{c^{T}_E} \bar{c_S}}{K_M/c_{S0} + \bar{c_S}}\Big) $$

$$\frac{\partial \bar{c_S}}{\partial \bar{t}} = \beta \frac{\partial^2 \bar{c_S}}{\partial \bar{x}^2} - k_2 \frac{c_{E0}}{c_{S0}}\frac{\bar{c^{T}_E} \bar{c_S}}{K_M/c_{S0} + \bar{c_S}}\frac{d^2}{D_E} $$

The boundary conditions that I want to impose is the following:

$$\frac{\partial c_E}{\partial x} = 0 \hspace{0.5 cm}x = 0; \hspace{0.5cm} \frac{\partial c_E}{\partial x} = 0 \hspace{0.5 cm}x = 1$$

$$c_S(t,x=0) = 1 ; \hspace{0.5cm} c_S(t,x=1) = 0$$

The initial conditions are as follows:

$$c_E(0,x) = 1; \hspace{0.5cm} c_S(0,x=0) = 1 ; c_S(0,x>0) = 0$$

When I try to implement FEM on this system, I am getting the following issue:

NDSolveValue::femper: PDE parsing error of Grad[(Ca Cs)/(2/5+Cs)]. Inconsistent equation dimensions.

Here is the code:

Needs["DifferentialEquations`InterpolatingFunctionAnatomy`"]
Needs["NDSolve`FEM`"]
(*Enzyme-substrate dyamics - hexokinase parameters*)
Clear["Global'*"];
Clear[Derivative];
Ca0 = 10^(-6);(*Initial Protein Concentration*)
d = 180*10^(-6);(*Half-width of channel*)
Da = 7.2*10^(-11);
CB0  = 10*10^(-3);
Km = 4*10^(-3);(*Michealis Menten Constant*)
beta  = 10;(*Subscript[D, s]/Subscript[D, E]*)
gamma = (Ca0/CB0)*(d^2/Da);
k2 = 8000;
alpha = -0.5;
op1 = D[Ca[t, x], t] - 
   Inactive[Div][Inactive[Grad][Ca[t, x], {x}], {x}] - 
   alpha*Inactive[Div][
     Inactive[
       Grad][(Ca[t, x]*Cs[t, x])/(Km/CB0 + Cs[t, x]), {x}], {x}];
op2 = D[Cs[t, x], t] - 
   beta*Inactive[Div][Inactive[Grad][Cs[t, x], {x}], {x}] + 
   k2*gamma*(Ca[t, x]*Cs[t, x])/(Km/CB0 + Cs[t, x]);
IC = {Ca[0, x] == 1, Cs[0, x] == Exp[-1000 x]};
Subscript[\[CapitalGamma], N1] = NeumannValue[0, x == 0];
Subscript[\[CapitalGamma], N2] = NeumannValue[0, x == 1];
Subscript[\[CapitalGamma], NS1] = 
  DirichletCondition[Cs[t, x] == 1, x == 0];
Subscript[\[CapitalGamma], NS2] = 
  DirichletCondition[Cs[t, x] == 0, x == 1];
eqn1 = {op1 == 
    Subscript[\[CapitalGamma], N1] + Subscript[\[CapitalGamma], N2]};
sol = NDSolveValue[{op1 == 
    Subscript[\[CapitalGamma], N1] + Subscript[\[CapitalGamma], N2], 
   op2 == 0, IC, Subscript[\[CapitalGamma], NS1], 
   Subscript[\[CapitalGamma], NS2]}, {Ca, Cs}, {t, 0, 1}, {x, 0, 1}]

What is the reason for this error? Can Mathematica solve this kind of PDEs? My PDE does have both Neumann and Dirichlet BC at same boundary but for different variables.

Edit1:

I simplified the 2nd term of my first equation. Here is what I wrote:

op1 = D[Ca[t, x], t] - 
   Inactive[Div][Inactive[Grad][Ca[t, x], {x}], {x}] - 
   alpha* (Inactive[Grad][Cs[t, x]/(Km/CB0 + Cs[t, x]), {x}] . 
       Inactive[Grad][
        Ca[t, x], {x}] + (Cs[t, x]/(Km/CB0 + Cs[t, x])) Inactive[Div][
        Inactive[Grad][Ca[t, x], {x}], {x}] + 
      Inactive[Grad][Ca[t, x], {x}] . 
       Inactive[Grad][Cs[t, x]/(Km/CB0 + Cs[t, x]), {x}] + 
      Ca[t, x] Inactive[Div][
        Inactive[Grad][Cs[t, x]/(Km/CB0 + Cs[t, x]), {x}], {x}]);

I now get the following error:

NDSolveValue::femper: PDE parsing error of 0.5 ((Cs Div$22707)/(2/5+Cs)+Ca Div$22709+Grad$22706.Grad$22708+Grad$22708.Grad$22706). Inconsistent equation dimensions.

What is the issue here?

Answer 1

You could try this:

Needs["DifferentialEquations`InterpolatingFunctionAnatomy`"]
Needs["NDSolve`FEM`"]
(*Enzyme-substrate dyamics-hexokinase parameters*)
Clear["Global'*"];
Clear[Derivative];
Ca0 = 10^(-6);(*Initial Protein Concentration*)d = 
 180*10^(-6);(*Half-width of channel*)Da = 7.2*10^(-11);
CB0 = 10*10^(-3);
Km = 4*10^(-3);(*Michealis Menten Constant*)beta = \
10;(*Subscript[D,s]/Subscript[D,E]*)gamma = (Ca0/CB0)*(d^2/Da);
k2 = 8000;
alpha = -0.5;
op1 = D[Ca[t, x], t] - 
   Inactive[Div][Inactive[Grad][Ca[t, x], {x}], {x}] - 
   Inactive[Div][{{alpha*(Cs[t, x])/(Km/CB0 + Cs[t, x])}} . 
     Inactive[Grad][Ca[t, x], {x}], {x}];
op2 = D[Cs[t, x], t] - 
   beta*Inactive[Div][Inactive[Grad][Cs[t, x], {x}], {x}] + 
   k2*gamma*(Ca[t, x]*Cs[t, x])/(Km/CB0 + Cs[t, x]);
IC = {Ca[0, x] == 1, Cs[0, x] == Exp[-1000 x]};
Subscript[\[CapitalGamma], N1] = NeumannValue[0, x == 0];
Subscript[\[CapitalGamma], N2] = NeumannValue[0, x == 1];
Subscript[\[CapitalGamma], NS1] = 
  DirichletCondition[Cs[t, x] == 1, x == 0];
Subscript[\[CapitalGamma], NS2] = 
  DirichletCondition[Cs[t, x] == 0, x == 1];
eqn1 = {op1 == 
    Subscript[\[CapitalGamma], N1] + Subscript[\[CapitalGamma], N2], 
   op2 == 0, IC, Subscript[\[CapitalGamma], NS1], 
   Subscript[\[CapitalGamma], NS2]};
sol = NDSolveValue[eqn1, {Ca, Cs}, {t, 0, 1}, {x, 0, 1}]

Your first equation, the third term was a bit off. Double check that I did not make a mistake transcribing your equations.