FITZ: Parallel integration for FitzHugh-Nagumo ODEs

Purpose

• The package "Fitz" solves the FitzHugh-Nagumo equation using the parareal method. To construct serial solvers with coarse and fine resolutions, the first-order forward Euler scheme was used.

Specifications

• Name: Fitz.
• Author: Feng Chen.
• Finishing date: 07/06/2013.
• Languages: MATLAB.

Simple Example

• Equation: \begin{equation} \begin{aligned} & v' = v-v^3-w+I_{\text{ext}},\\ & \tau w' = v + a - bw. \end{aligned} \end{equation}
• Parameters: \begin{equation} I_{\text{ext}}=0.5, \, a=0.7,\, b=0.8, \,\tau=12.5,\, T=200, \,\Delta T =0.1, \, \Delta t = 0.05, \, \delta t=0.001. \end{equation}
• Initial conditions: \begin{equation} v(0)=1, \quad w(0)=0. \end{equation}

Quick Start

• Compiling and running:

  matlab Fitz_Driver_Parareal
  

• Output:

  iteration 0, 0.000000e+00, 7.040824e-01 
  iteration 1, 1.235809e-02, 6.917243e-01 
  iteration 2, 8.904816e-04, 6.908338e-01 
  iteration 3, 6.172277e-05, 6.907721e-01 
  iteration 4, 2.008283e-06, 6.907701e-01 
  iteration 5, 9.682426e-08, 6.907702e-01 
  iteration 6, 2.131734e-08, 6.907702e-01 
  iteration 7, 2.023203e-09, 6.907702e-01 
  iteration 8, 1.402405e-10, 6.907702e-01 
  iteration 9, 7.806200e-12, 6.907702e-01 
  iteration 10, 3.443912e-13, 6.907702e-01 
  sucessful return! 
  ...
  

• Graphics:
  
• CPU: Intel(R) Xeon(R) CPU X5550 @2.67GHz.
• OS: CentOS release 6.4 (Final).
• Release: MATLAB R2012a.

References

• Lions, J. L and Maday, Y. and Turinici, G. A "parareal" in time discretization of PDE's, Comptes Rendus de l'Academie des Sciences Series I Mathematics, Vol 332, (2001).
• The FitzHugh-Nagumo model. Link.

Code Highlights

% parareal iteration
for k = 2:imax+1
 % fine solver
 for i = k:n+1
  ytF = Fitz_Time(sF, yp(:,i-1), l*m); yF(:,i) = ytF(:,end);
 end
 
 % coarse solver
 yc(:,k) = yF(:,k); 
 for i = k+1:n+1
  ytG = Fitz_Time(sG, yc(:,i-1), m); yGc(:,i) = ytG(:,end);
  yc(:,i) = yF(:,i) + yGc(:,i) - yGp(:,i); 
 end
  
  % check errors
 ed = abs(yc(1,end)-yp(1,end));  % relative error
 et = abs(ye(1,end)-yc(1,end));  % true error
 fprintf('iteration %d, %e, %e \n', k-1, ed, et);  

 % check convergence
 if (ed < tol)
  fprintf('sucessful return! \n');
  break
 end

 % update solution
 yGp(:,k:n+1) = yGc(:,k:n+1); 
 yp(:,k:n+1) = yc(:,k:n+1); 
end