## Set file for 4p-r.ode on Fri Aug  9 09:48:48 2019
13   Number of equations and auxiliaries
74   Number of parameters
# Numerical stuff
1    nout
40    nullcline mesh
8   Qual RK
10000  total
0.1  DeltaT
0  T0
0  Transient
5000000000000  Bound
1e-12  DtMin
1  DtMax
0.001  Tolerance
0.001  Abs. Tolerance
0  Max Delay
100   Eigenvector iterates
0.001  Eigenvector tolerance
0.001  Newton tolerance
0  Poincare plane
1e-05  Boundary value tolerance
1e-05  Boundary value epsilon
20   Boundary value iterates
0   Poincare None
1   Poincare variable
1   Poincare sign
0   Stop on Section
0   Delay flag
10000  Current time
10000  Last Time
1   MyStart
1   INFLAG
# Delays
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
0.0
# Bndry conds
0
0
0
0
0
0
0
0
0
0
0
0
0
# Old ICs
-14.5  V1
0.195  H1
0.925925  S1
0.925925  S1_SLOW
-61.6  V2
0.199  H2
0  S2
-75  V3
0.195  H3
0  S3
-15  V4
0.195  H4
0.925925  S4
# Ending  ICs
-70.49676416791507  V1
0.9803119994014953  H1
9.089344112083678e-05  S1
0.06939266143278369  S1_SLOW
-29.6874032784189  V2
0.3496175204458359  H2
0.9144282644961006  S2
-65.7178673353848  V3
0.9907160268755857  H3
0.0004189356795631216  S3
-39.45250006278307  V4
0.2200283073605734  H4
8.375200095736193e-08  S4
# Parameters
21  c
10  gnap
2.8  gl
60  gsyni
10  gsyne
5  gnap_weak
3  gk
50  Ena
-65  El
-80  Esyn
0  EsynE
-85  Ek
-40  thetam
-6  sigmam
-47  thetah
6  sigmah
-40.25  thetah_p
-30  thetamk
-4  sigmamk
600  epsilon_i
700  epsilon_e
1500  epsilon_p
300  epsilon_p2
-35  thetatau
-35  thetasyn
-3  sigmasyn
-0.5  sigsyn_p2
-30  thetsyn_p2
-3  sigsyn_p
-35  thetsyn_p
1  pdscale
1  rtndscale
0  pd_i
0.15  pd_e
0  pd_p
0.8  pd_p2
0.3  rtnd_i
0.5  rtnd_e
0  rtnd_p
0  rtnd_p2
0.03  rtnd_ie
0  od_i
0  od_e
0  od_p
0  od_p2
0  od_ie
1  vago
1  gks_i
1  gks_e
0.3  gks_p
0.5  gks_p2
10  thr
10  conc
1  kalth_i
1  kalth_e
0  kalth_p
0  kalth_p2
1  ks
0.035  b_ei
0.001  b_ep
0.03  b_ie
0.1  b_ip
0.05  b_pp2
0.05  b_p2p
0.25  a_pe
0.15  a_ip2
0.75  dscale
1  delta_i
1  delta_e
1  delta_p
1  delta_p2
1  alpha
0.08  beta
0.005  beta_slow
# Graphics
0.7071097251469876  rm
-0.4999999999826661  rm
0.4999958366025516  rm
0.7071038372138505  rm
0.5000041633974482  rm
-0.4999999999826661  rm
0  rm
0.7071038372138505  rm
0.7071097251469876  rm
0    
1    
1    
1    
0    
0    
5    
1    
1    
9    
0    
8    
1    
1    
7    
0    
11    
1    
1    
5    
0    
1    
1    
1    
0    
0    
1    
1    
1    
0    
0    
1    
1    
1    
0    
0    
1    
1    
1    
0    
0    
1    
1    
1    
0    
0    
1    
1    
1    
0    
-1000   
1000   
0    
0   3DFlag
1   Timeflag
0   Colorflag
0   Type
1  color scale
0   minscale
10000   xmax
0   xmin
20   ymax
-80   ymin
12   zmax
-12   zmin
5000   
0.0002   
-30   
0.02   
0   
0.08333333333333333   
45   Theta
45   Phi
0    xshft
0    yshft
0    zshft
0   xlo
-80   ylo
0   
-80   
10000   xhi
20   yhi
10000   
20   
# Transpose variables etc
V1
2   n columns
1   n rows
1   row skip
1   col skip
1   row 0
# Coupling stuff for H funs
0
0
0
0
0
0
0
0
0
0
0
0
0
# Array plot stuff

1   NCols
0   Row 1
50   NRows
8   RowSkip
0  Zmin
1  Zmax
# Torus information 
0    Torus flag 1=ON
6.283185307179586  Torus period
# Range information
c
-1   eq-range stab col
0   shoot flag 1=on
10   eq-range steps
0  eq_range low
1  eq_range high
V1
V1
20   Range steps
0   Cycle color 1=on
1   Reset data 1=on
1   Use old I.C.s 1=yes
0  Par1 low
0  Par2 low
1  Par1 high
1  Par2 high
c
0   BVP side
0   color cycle flag 1=on
10   BVP range steps
0  BVP range low
1  BVP range high
RHS etc ...
dV1/dT=(-INAP(V1,H1,GNAP)-IK(V1)-IL(V1,GL)-IAPP(V1,C_I)-ISYN(V1,S2,B_EI,KALTH_I)-IKS(V1,GKS_I)-W_I*DELTA_I*DSCALE)/C
dH1/dT=(HINF(V1,THETAH,SIGMAH)-H1)/TAUH(V1,EPSILON_I,THETATAU)
dS1/dT=ALPHA*(1-S1)*SINF(V1,THETASYN,SIGMASYN)-BETA*S1
dS1_SLOW/dT=ALPHA*(1-S1_SLOW)*SINF(V1,THETASYN,SIGMASYN)-BETA_SLOW*S1_SLOW
dV2/dT=(-INAP(V2,H2,GNAP_WEAK)-IK(V2)-IL(V2,GL)-IAPP(V2,C_E)-ISYN(V2,S1,B_IE,KALTH_E)-ISYNE(V2,S1_SLOW,A_IE*VAGO)-ISYNE(V2,S3,A_PE)-IKS(V2,GKS_E)-W_E*DELTA_E*DSCALE)/C
dH2/dT=(HINF(V2,THETAH,SIGMAH)-H2)/TAUH(V2,EPSILON_E,THETATAU)
dS2/dT=ALPHA*(1-S2)*SINF(V2,THETASYN,SIGMASYN)-BETA*S2
dV3/dT=(-INAP(V3,H3,GNAP_WEAK)-IK(V3)-IL(V3,GL)-IAPP(V3,C_P)-ISYN(V3,S2,B_EP,KALTH_P)-ISYN(V3,S1_SLOW,B_IP*VAGO,KALTH_P)-ISYN(V3,S4,B_P2P,KALTH_P)-IKS(V3,GKS_P)-W_P*DELTA_P*DSCALE)/C
dH3/dT=(HINF(V3,THETAH_P,SIGMAH)-H3)/TAUH(V3,EPSILON_P,THETATAU)
dS3/dT=ALPHA*(1-S3)*SINF(V3,THETSYN_P,SIGSYN_P)-BETA*S3
dV4/dT=(-INAP(V4,H4,GNAP_WEAK)-IK(V4)-IL(V4,GL)-IAPP(V4,C_P2)-ISYN(V4,S3,B_PP2,KALTH_P2)-ISYNE(V1,S1,A_IP2)-IKS(V4,GKS_P2)-W_P2*DELTA_P2*DSCALE)/C
dH4/dT=(HINF(V4,THETAH,SIGMAH)-H4)/TAUH(V4,EPSILON_P2,THETATAU)
dS4/dT=ALPHA*(1-S4)*SINF(V4,THETSYN_P2,SIGSYN_P2)-BETA*S4

where ...
C_I = PDSCALE*PD_I+RTNDSCALE*RTND_I+OD_I
   
C_E = PDSCALE*PD_E+RTNDSCALE*RTND_E+OD_E
   
C_P = PDSCALE*PD_P+RTNDSCALE*RTND_P+OD_P
   
C_P2 = PDSCALE*PD_P2+RTNDSCALE*RTND_P2+OD_P2
   
A_IE = (RTND_IE*RTNDSCALE)
   
MKS = 1/(1+EXP(-(CONC-THR)/2))
   

User-defined functions:
HINF(V,THET,SIG) = 1/(1+EXP((V-THET)/SIG))
  
SINF(V,THETA,SIGMA) = 1/(1+EXP((V-THETA)/SIGMA))
  
TAUH(V,EP,THETA) = EP/(COSH((V-THETA)/SIGMAH/2))
  
MINF(V) = 1/(1+EXP((V-THETAM)/SIGMAM))
  
MK(V) = 1/(1+EXP((V-THETAMK)/SIGMAMK))
  
ISYN(V,S,B,KALTH) = GSYNI*B*S*(1+KS*KALTH)*(V-ESYN)
  
ISYNE(V,S,A) = GSYNE*A*S*(V-ESYNE)
  
INAP(V,H,G) = G*MINF(V)*H*(V-ENA)
  
IL(V,G) = G*(V-EL)
  
IAPP(V,C) = C*V
  
IKS(V,G) = G*MKS*(V-EK)
  
IK(V) = GK*MK(V)^4*(V-EK)