|
|
@@ -12,6 +12,7 @@ class model:
|
|
|
def __init__(self):
|
|
|
self.compartments={}
|
|
|
self.seJ={}
|
|
|
+ self.scaled=[]
|
|
|
|
|
|
def add_source(self,compartmentName,formula):
|
|
|
self.compartments[compartmentName]['source']=formula
|
|
|
@@ -233,7 +234,8 @@ class model:
|
|
|
pcParName=self.mod['partitionCoefficients'][pc]
|
|
|
pcPar=pars[pcParName]
|
|
|
print('\t{}:{} [{}]'.format(pc,pcParName,self.get(pcParName)))
|
|
|
-
|
|
|
+
|
|
|
+ def inspectSE(self):
|
|
|
|
|
|
print('SE parameters')
|
|
|
for p in self.seJ:
|
|
|
@@ -253,8 +255,12 @@ class model:
|
|
|
with open(parameterFile,'r') as f:
|
|
|
self.parSetup=json.load(f)
|
|
|
|
|
|
+ self.mod['compartments'].append('total')
|
|
|
for m in self.mod['compartments']:
|
|
|
self.add_compartment(m)
|
|
|
+
|
|
|
+ for m in self.mod['scaled']:
|
|
|
+ self.scaled.append(m)
|
|
|
|
|
|
self.add_default_parameters()
|
|
|
#standard parameters such as one,zero etc.
|
|
|
@@ -322,15 +328,44 @@ class model:
|
|
|
pars['zero']={'value':0}
|
|
|
|
|
|
|
|
|
- def M(self,t):
|
|
|
+ def M(self,t,y=numpy.array([])):
|
|
|
for i in self.dM:
|
|
|
for j in self.dM[i]:
|
|
|
self.fM[i,j]=self.dM[i][j](t)
|
|
|
#create an array and fill it with outputs of function at t
|
|
|
+ if (y.size==0):
|
|
|
+ return self.fM
|
|
|
+ self.set_scaledM(t,y)
|
|
|
return self.fM
|
|
|
|
|
|
+ def set_scaledM(self,t,y):
|
|
|
+ #prevent zero division
|
|
|
+ eps=1e-8
|
|
|
+ for c in self.scaled:
|
|
|
+ i=self.lut[c]
|
|
|
+ it=self.lut['total']
|
|
|
+ try:
|
|
|
+ k=numpy.copy(self.originalK[i])
|
|
|
+ except AttributeError:
|
|
|
+ k=numpy.copy(self.fM[i,:])
|
|
|
+ self.originalK={}
|
|
|
+ self.originalK[i]=k
|
|
|
+ #make another copy
|
|
|
+ k=numpy.copy(self.originalK[i])
|
|
|
+ except KeyError:
|
|
|
+ k=numpy.copy(self.fM[i,:])
|
|
|
+ self.originalK[i]=k
|
|
|
+ #make another copy
|
|
|
+ k=numpy.copy(self.originalK[i])
|
|
|
+
|
|
|
+ k[i]=k[i]-self.u(t)[it]
|
|
|
+ #scale all inputs by total input mass
|
|
|
+ for j in range(self.n):
|
|
|
+ self.fM[i,j]=k[j]/(y[it]+eps)
|
|
|
+
|
|
|
def u(self,t):
|
|
|
ub=[f(t) for f in self.fu]
|
|
|
+ ub[self.lut['total']]=sum(ub)
|
|
|
return numpy.array(ub)
|
|
|
|
|
|
def jacobiFull(self,t):
|
|
|
@@ -344,14 +379,36 @@ class model:
|
|
|
|
|
|
|
|
|
|
|
|
- def fSS(self,t):
|
|
|
+ def fSS(self,t,y=numpy.array([])):
|
|
|
for k in self.qSS:
|
|
|
for i in self.qSS[k]:
|
|
|
for j in self.qSS[k][i]:
|
|
|
#print('[{},{},{}] {}'.format(k,i,j,self.qSS[k][i][j]))
|
|
|
self.SS[k,i,j]=(self.qSS[k][i][j])(t)
|
|
|
+ if y.size==0:
|
|
|
+ return self.SS
|
|
|
+ self.set_scaledSS(t,y)
|
|
|
return self.SS
|
|
|
|
|
|
+ def set_scaledSS(self,t,y):
|
|
|
+ #prevent zero division
|
|
|
+ eps=1e-8
|
|
|
+ for c in self.scaled:
|
|
|
+ it=self.lut['total']
|
|
|
+ i=self.lut[c]
|
|
|
+ try:
|
|
|
+ dkdp=numpy.copy(self.originalSS[i])
|
|
|
+ except AttributeError:
|
|
|
+ dkdp=numpy.copy(self.SS[:,i,:])
|
|
|
+ self.originalSS={}
|
|
|
+ self.originalSS[i]=dkdp
|
|
|
+ dkdp=numpy.copy(self.originalSS[i])
|
|
|
+ except KeyError:
|
|
|
+ dkdp=numpy.copy(self.SS[:,i,:])
|
|
|
+ self.originalSS[i]=dkdp
|
|
|
+ dkdp=numpy.copy(self.originalSS[i])
|
|
|
+ self.SS[:,i,:]=dkdp/(y[it]+eps)
|
|
|
+ #should add error on u!
|
|
|
|
|
|
def fSY(self,y,t):
|
|
|
#M number of sensitivity parameters
|
|
|
@@ -360,7 +417,7 @@ class model:
|
|
|
|
|
|
#assume a tabulated solution y(t) at t spaced intervals
|
|
|
|
|
|
- qS=self.fSS(t).dot(y)
|
|
|
+ qS=self.fSS(t,y).dot(y)
|
|
|
#qS is MxN
|
|
|
#but NxM is expected, so do a transpose
|
|
|
|
