I am trying to build an Action potential simulation (a graph of voltage relative to time, where the voltage values are obtained through numerically solving a couple of differential equiations). I have used the code from an online simulation of this in javaScript and am trying to translate it into python. However, due to the different ways in which python and java handle floats, I get an Overflowerror: math range error. Does anyone know how I can circumvent this in this specific situation? I am posting everything i've written and my error output:
import numpy as np
import math
import matplotlib.pyplot as plt
import random
class HH():
def __init__(self):
self.dt = 0.025
self.Capacitance = 1
self.GKMax = 36
self.GNaMax = 120
self.Gm = 0.3
self.EK = 12
self.ENa = 115
self.VRest = 10.613
self.V = 0
self.n = 0.32
self.m = 0.05
self.h = 0.60
self.INa = 0
self.IK = 0
self.Im = 0
self.Iinj = 0
self.tStop = 200
self.tInjStart = 25
self.tInjStop = 175
self.IDC = 0
self.IRand = 35
self.Itau = 1
@classmethod
def alphaN(cls,V):
if(V == 10):
return alphaN(V+0.001)
else:
a = (10-V)/(100 * (math.exp((10-V)/10) -1))
print("alphaN: ", a)
return a
@classmethod
def betaN(cls,V):
a = 0.125*math.exp(-V/80)
print("betaN: ", a)
return a
@classmethod
def alphaM(cls, V):
if (V == 25):
return alphaM(V+0.0001)
else:
a = (25 - V)/10 * (math.exp( (25-V)/10) - 1)
print("alphaM", a)
return (a)
@classmethod
def betaM(cls, V):
return 4 * math.exp(-V/18)
@classmethod
def alphaH(cls, V):
return 0.07 * math.exp(-V/20)
@classmethod
def betaH(cls, V):
return 1/(math.exp((30-V/10)+1))
def iteration(self):
aN = self.alphaN(self.V)
bN = self.betaN(self.V)
aM = self.alphaM(self.V)
bM = self.betaM(self.V)
aH = self.alphaH(self.V)
bH = self.betaH(self.V)
tauN = 1/(aN + bN)
print("tauN: ", tauN)
tauM = 1/(aM + bM)
print("tauM", tauM)
tauH = 1/(aH + bH)
print("tauH: ", tauH)
nInf = aN * tauN
print("nInf: ", nInf)
mInf = aM * tauM
print("mInf: ", mInf)
hInf = aH * tauH
print("hInf: ", hInf)
self.n += self.dt/ tauN * (nInf - self.n)
print("n: ", self.n)
self.m += self.dt/tauM * (mInf - self.m)
print("m: ", self.m)
self.h += self.dt/tauH * (hInf - self.h)
print("h: ", self.h)
self.IK = self.GKMax * self.n * self.n * self.n * self.n * (self.VRest - self.EK)
print("IK: ", self.IK)
self.INa = self.GNaMax * self.m * self.m * self.m * self.h * (self.VRest - self.ENa)
print("INa: ", self.INa)
self.Im = self.Gm * (self.VRest * self.V)
print("Im: ", self.Im)
self.V += self.dt/self.Capacitance * (self.INa + self.IK + self.Im + self.Iinj)
print("V: ", self.V)
print("nnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnn")
return self.V
if __name__ == '__main__':
hodgkinHuxley = HH()
V_vector = np.zeros(math.floor(hodgkinHuxley.tStop/hodgkinHuxley.dt)) #v_vector
Iinj_vector = np.zeros(math.floor(hodgkinHuxley.tStop/hodgkinHuxley.dt)) #stim_vector
n_vector = np.zeros(math.floor(hodgkinHuxley.tStop/hodgkinHuxley.dt))
m_vector = np.zeros(math.floor(hodgkinHuxley.tStop/hodgkinHuxley.dt))
h_vector = np.zeros(math.floor(hodgkinHuxley.tStop/hodgkinHuxley.dt))
plotSampleRate = math.ceil(hodgkinHuxley.tStop/2000/hodgkinHuxley.dt) #t_vector
print("plotsamplerate: ", plotSampleRate)
i = 0
t_vector = []
for t in range(0, hodgkinHuxley.tStop+1):
t= t+hodgkinHuxley.dt
if(math.floor(t)>hodgkinHuxley.tInjStart & math.ceil(t)<hodgkinHuxley.tInjStop):
rawInj = hodgkinHuxley.IDC + hodgkinHuxley.IRand * 2 * (random.random()-0.5)
else:
rawInj = 0
hodgkinHuxley.Iinj += hodgkinHuxley.dt/hodgkinHuxley.Itau * (rawInj - hodgkinHuxley.Iinj)
hodgkinHuxley.iteration()
counter = 0
if(i == plotSampleRate):
V_vector[counter] = hodgkinHuxley.V
Iinj_vector[counter] = hodgkinHuxley.Iinj
n_vector[counter] = hodgkinHuxley.n
m_vector[counter] = hodgkinHuxley.m
h_vector[counter] = hodgkinHuxley.h
i=0;
counter+=1
i+=1
Error:
plotsamplerate: 4
alphaN: 0.05819767068693265
betaN: 0.125
alphaM 27.956234901758684
tauN: 5.458584687514421
tauM 0.03129279788667988
tauH: 14.285714285707257
nInf: 0.3176769140606974
mInf: 0.8748288084532805
hInf: 0.9999999999995081
n: 0.31998936040167786
m: 0.7089605789167688
h: 0.6006999999999995
IK: -0.5235053389507994
INa: -2681.337959108097
Im: 0.0
V: -67.0465366111762
nnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnn
alphaN: 0.00034742442109860416
betaN: 0.2889909708647963
alphaM 91515.37543280824
tauN: 3.456160731837547
tauM 1.0907358211913938e-05
tauH: 0.5000401950825147
nInf: 0.0012007546414823879
mInf: 0.9981909817434281
hInf: 1.0
n: 0.31768341581102577
m: 663.6339264765652
h: 0.6206633951393696
IK: -0.5085774931025618
INa: -2272320232559.0464
Im: -213.46946791632385
V: -56808005886.37157
nnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnn
Traceback (most recent call last):
File "C:\Users\Huinqk 2.0\Documents\projects\python\tariq_neural networks\plotting function.py", line 131, in <module>
hodgkinHuxley.iteration()
File "C:\Users\Huinqk 2.0\Documents\projects\python\tariq_neural networks\plotting function.py", line 71, in iteration
aN = self.alphaN(self.V)
File "C:\Users\Huinqk 2.0\Documents\projects\python\tariq_neural networks\plotting function.py", line 38, in alphaN
a = (10-V)/(100 * (math.exp((10-V)/10) -1))
OverflowError: math range error
[Finished in 0.8s]
