22
22
// State variables
23
23
//------------------------------------------------------------------------------
24
24
25
- double Y [_NB_OF_STATE_VARIABLES_ ];
25
+ extern double Y [_NB_OF_STATE_VARIABLES_ ];
26
26
double dY [_NB_OF_STATE_VARIABLES_ ];
27
+ double Ynew [_NB_OF_STATE_VARIABLES_ ];
28
+
27
29
{%- for state_var in state_vars %}
28
30
// {{loop.index0}}: {{state_var.var}} (units: {{state_var.units}}, initial value: {{state_var.initial_value}}, component: {{components[loop.index0]}})
29
31
{%- endfor %}
30
32
33
+ double Vmem ;
34
+ double time ;
35
+
31
36
char YNames [_NB_OF_STATE_VARIABLES_ ][{{stat_var_name_max_length + 1 }}];
32
37
char YUnits [_NB_OF_STATE_VARIABLES_ ][{{unit_name_max_length + 1 }}];
33
38
char YComponents [_NB_OF_STATE_VARIABLES_ ][{{component_name_max_length + 1 }}];
@@ -61,12 +66,17 @@ void init()
61
66
{%for state_var in state_vars %}
62
67
Y [{{loop .index0 }}] = {{state_var .initial_value }}; // {{state_var.var}} ({{state_var.units}}) (in {{components[loop.index0]}})
63
68
{%- endfor %}
69
+ Y [{{state_vars | length }}]; // ({{free_variable.var_name}}} (milliseconds)
70
+
64
71
{%for state_var in state_vars %}
65
72
strcpy (YNames [{{loop .index0 }}], "{{state_var.var}}" );
66
73
{%- endfor %}
74
+ strcpy (YNames [{{state_vars | length }}], "{{free_variable.var_name}}" );
75
+
67
76
{%for state_var in state_vars %}
68
77
strcpy (YUnits [{{loop .index0 }}], "{{state_var.units}}" );
69
78
{%- endfor %}
79
+ strcpy (YUnits [{{state_vars | length }}], "milliseconds" );
70
80
71
81
//------------------------------------------------------------------------------
72
82
// Constants
@@ -79,7 +89,10 @@ void init()
79
89
{% endif %}{%- endfor %}
80
90
}
81
91
82
- void compute (double {{free_variable .var_name }})
92
+ //---------------------------------------------------------------------------
93
+ // Computation
94
+ //---------------------------------------------------------------------------
95
+ void compute ()
83
96
{
84
97
// time: {{free_variable.var_name}} (millisecond)
85
98
{% for deriv in y_derivative_equations %}
@@ -96,6 +109,46 @@ void compute(double {{free_variable.var_name}})
96
109
{%- endfor %}
97
110
}
98
111
112
+ //------------------------------------------------------------------------------
113
+ // Integration & Output
114
+ //------------------------------------------------------------------------------
115
+ // Rush-Larsen method
116
+
117
+
118
+ // get tau/inf or alpha/beta
119
+ void computeTauInf ()
120
+ {% for deriv in derivative_alpha_beta %}{% if deriv .type != 'non_linear' %}double alphaOrTau_ {{loop .index0 }} = {{deriv .r_alpha_or_tau }};
121
+ double betaOrInf_ {{loop .index0 }} = {{deriv .r_beta_or_inf }};
122
+ {% endif %}{%- endfor %}
123
+ // gating variables: Exponential integration
124
+ {% for deriv in derivative_alpha_beta %}{%- if deriv .type == 'inftau' %}
125
+ Ynew [{{loop .index0 }}] = betaOrInf_ {{loop .index0 }} + (Y [{{loop .index0 }}] - betaOrInf_ {{loop .index0 }})* exp (- dt /alphaOrTau_ {{loop .index0 }});
126
+ {%- elif deriv .type != 'non_linear' %}
127
+ double tau_inv_ {{loop .index0 }} = alphaOrTau_ {{loop .index0 }} + betaOrInf_ {{loop .index0 }};
128
+ double y_inf_ {{loop .index0 }} = alphaOrTau_ {{loop .index0 }} / tau_inv_ {{loop .index0 }};
129
+ Ynew [{{loop .index0 }}] = y_inf_ {{loop .index0 }} + (Y [{{loop .index0 }}] - y_inf_ {{loop .index0 }})* exp (- dt * tau_inv_ {{loop .index0 }});
130
+ {%- endif %}
131
+ {%- endfor %}
132
+ }
133
+
134
+ // Remainder: Forward Euler
135
+ void computeRemainderForaredEuler (){
136
+ {% for deriv in derivative_alpha_beta %}
137
+ {%- if deriv .type == 'non_linear' %}
138
+ Ynew [{{loop .index0 }}] = Y [{{loop .index0 }}] + dt * {{deriv .deriv }};{%- endif %}{%- endfor %}
139
+ Ynew [{{state_vars | length }}] = Y [{{state_vars | length }}] + dt ;
140
+ }
141
+
142
+ compute ();
143
+ computeTauInf ();
144
+ computeRemainderForaredEuler ();
145
+
146
+
147
+ Vmem = Ynew [0 ];
148
+ time = Ynew [{{state_vars | length }}];
149
+
150
+
151
+
99
152
//==============================================================================
100
153
// End of file
101
154
//==============================================================================
0 commit comments