/* MemFH? Module */
/* */
/* patch file for the FitzHugh?-Nagumo membrane */
/* model, Vm goes 0 to 100 mV with AP */
/* */
/* Vrest is in mV */
/* Vt and Vp are in mV */
/* A, B, and C are in ? */
/* original paper: R FitzHugh?, "Computation of impulse initiation and */
/* saltatory conduction in a myelinated nerve fiber," Biophys J, */
/* v.2, p.445, 1961 */
/* author of this computer code: John B. Pormann */
#include "CardioWave.h"
void SetPatches?_FH( domain_t tp, vector Vm, vector Q );
int GetF_FH( vector Vm, vector Qv, vector Fv, vector Fq, vector Av );
/* FitzHugh?-Nagumo patch */
typedef struct {
real w;
} FH_patch;
typedef struct {
real Iion;
} FH_auxvar;
/* pseudo-globals */
int FH_PatchSize = sizeof(FH_patch)/sizeof(real);
int FH_AuxiliarySize = sizeof(FH_auxvar)/sizeof(real);
int FH_MemParamSize = 0;
byte FH_NodeType = ByteError?;
static char* RCSID = "$Id: MemFH.c,v 1.6 2005/02/17 20:57:04 john Exp $";
/* constants */
static real FH_Vt = 5.0;
static real FH_Vp = 100.0;
static real FH_A = 0.0;
static real FH_B = 0.2425;
static real FH_C = 1.5;
static real FH_RestVoltage = 0.0;
static FH_patch FH_RestPatch = { 0.0 };
static rword resources[] = {
{ "FH_A", 1001 },
{ "FH_B", 1002 },
{ "FH_C", 1003 },
{ "FH_Vp", 1004 },
{ "FH_Vt", 1005 },
{ "FH_IV", 1006 },
{ "FH_Node", 1100 },
{ "FH_Nodetype",1100 },
{ "FH_Patch", 1006 },
{ "FH_Type", 1100 },
{ "FH_Vp", 1004 },
{ "FH_Vr", 1007 },
{ "FH_Vrest", 1007 },
{ "FH_Vt", 1005 },
{ "showofs", 1200 },
{ "showoffset", 1200 },
{ "showoffsets",1200 },
{ "show_ofs", 1200 },
{ "show_offset",1200 },
{ "show_offsets",1200 },
{ NULL, 0 }
};
int InitMembrane?_FH( char** res ) {
int i,j,c;
int cmd;
real* iv;
real* pv;
logical show_ofs = false;
i = 0;
while( res[i] != NULL ) {
cmd = FindCommand( resources, res[i] );
switch( cmd ) {
case 1001:
FH_A = GetRealValue( res[i] );
break;
case 1002:
FH_B = GetRealValue( res[i] );
break;
case 1003:
FH_C = GetRealValue( res[i] );
break;
case 1004:
FH_Vp = GetRealValue( res[i] );
break;
case 1005:
FH_Vt = GetRealValue( res[i] );
break;
case 1006:
iv = GetRealArray( res[i] );
pv = (real*)(&FH_RestPatch);
c = GetNumValues( res[i] );
if( c > FH_PatchSize ) {
c = FH_PatchSize;
}
for(j=0;j<c;j++) {
pv[j] = iv[j];
}
break;
case 1007:
FH_RestVoltage = GetRealValue( res[i] );
break;
case 1100:
FH_NodeType = GetByteValue( res[i] );
break;
case 1200:
show_ofs = GetTFValue( res[i] );
break;
}
i++;
}
/* required!! */
if( FH_NodeType == ByteError? ) {
return( -1 );
}
PatchSize[FH_NodeType] = FH_PatchSize;
if( UseAuxvars ) {
AuxiliarySize[FH_NodeType] = FH_AuxiliarySize;
} else {
AuxiliarySize[FH_NodeType] = 0;
}
MemParamSize[FH_NodeType] = FH_MemParamSize;
FunctionTable[FH_NodeType] = GetF_FH;
SetpatchTable[FH_NodeType] = SetPatches?_FH;
if( (DebugLevel>0) and (SelfPE==0) ) {
printf("MemFH?: patchsize = 0.2lf %0.2lf %0.2lf %0.2lf %0.2lf\n",
FH_PatchSize,FH_A,FH_B,FH_C,FH_Vt,FH_Vp);
if( DebugLevel > 1 ) {
printf("MemFH?: RCSID: %s\n",RCSID);
}
}
if( show_ofs and (SelfPE==0) ) {
printf("MemFH? patch offsets:\n");
printf(" W: %i\n", OffsetOf?(FH_patch,w) );
printf("MemFH? auxvar offsets:\n");
printf(" Iion: %i\n", OffsetOf?(FH_auxvar,Iion) );
}
return( 0 );
}
void ExitMembrane?_FH( void ) {
return;
}
void SetPatches?_FH( domain_t tp, vector Vm, vector Q ) {
int i,j;
real* vm = Vm.data;
real* qp = Q.data;
real* rp = (real*)(&FH_RestPatch);
byte* nt = NodeType.data;
if( tp == Tissue ) {
for(i=0;i<TissueLocalSize;i++) {
if( nt[i] == FH_NodeType ) {
vm[i] = FH_RestVoltage;
for(j=0;j<FH_PatchSize;j++) {
qp[j] = rp[j];
}
}
qp += PatchSize[nt[i]];
}
}
return;
}
/* GetF() should fill in vectors Fv and Fq with the new currents/updates */
/* NOTE: this should be from the form: dV/dt + Fv = 0, dq/dt + Fq = 0 */
int GetF_FH( vector Vm, vector Qv, vector Fv, vector Fq, vector Av ) {
int i;
real* vmptr = Vm.data;
real* fv = Fv.data;
real* qm = Qv.data;
real* fq = Fq.data;
real* av = Av.data;
FH_patch* qp;
FH_patch* fp;
FH_auxvar* ap;
byte* nt = NodeType.data;
real vm,w,iion;
for(i=0;i<TissueLocalSize;i++) {
if( nt[i] == FH_NodeType ) {
/* map pointers to patch structures */
qp = (FH_patch*)qm;
fp = (FH_patch*)fq;
ap = (FH_auxvar*)av;
/* get local variables */
vm = vmptr[i];
w = qp->w;
iion = vm*(vm/FH_Vt-1.0)*(1.0-vm/FH_Vp) - w;
fv[i] += iion;
fp->w = FH_C*( vm + FH_A - FH_B*w );
if( UseAuxvars ) {
ap->Iion = iion;
}
} /* end-if */
/* go to next patch */
qm += PatchSize[ nt[i] ];
fq += PatchSize[ nt[i] ];
av += AuxiliarySize[ nt[i] ];
}
return( 0 );
}
