Problem 1

The variational formulation of the first equation that you need to specify in the file Reaction1.m is given by

    u*v*dx + k*a*du'*dv*dx

for the left-hand side and

    - k*c*w(2)*w(3)*v*dx + w(1)*v*dx

for the right-hand side.

The variational formulation of the second equation that you need to specify in the file Reaction2.m is given by

    u*v*dx + k*a*du'*dv*dx

for the left-hand side and

    k*(-c*w(2)*w(3) + f(x,d,t))*v*dx + w(1)*v*dx

for the right-hand side.

The fixed-point iteration can be implemented as follows:

while 1

  % Assemble vectors
  b1 = AssembleVector(p, e, t, 'Reaction1', [U10 U11 U21], time);
  b2 = AssembleVector(p, e, t, 'Reaction2', [U20 U11 U21], time);
    
  % Solve the linear systems
  newU1 = A1 \ b1;
  newU2 = A2 \ b2;
    
  % Check if the solution has converged
  if (norm(newU1 - U11) + norm(newU2 - U21))  < 0.01
    break;
  end
    
  % Update to the new values
  U11 = newU1;
  U21 = newU2;

end