Adding solveMatrix to runSolver

pythonScripts/runSolver.py

@@ -1,4 +1,5 @@
 import ivp
+import solveMatrix
 import os
 import cModel
 import sys
@@ -7,6 +8,8 @@ import numpy
 import time
 import scipy.interpolate
 import shutil
+import importlib
+importlib.reload(solveMatrix)
 
 defaultValues={\
       'method':'LSODA',\
@@ -22,7 +25,8 @@ defaultValues={\
       'project':'Analysis/Run',\
       'schema':'lists',\
       'query':'runs',\
-      'view':'WithStrings'}
+      'view':'WithStrings',
+      'exposureScaleRef':1}
 
 def get(setup,par):
    try:
@@ -128,9 +132,39 @@ def main(parFiles,jobDir,startDir='NONE'):
          #this could be modified for different ordering of lut as well
          sIn[:,model.lutSE[x]]=s0[:,startPoint['lutSE'][x]]
    #copy t0 to setup to help get the strides right
+   exposureScale=setup['exposureScaleRef']
+   y0*=exposureScale
+   print('sIn {}'.format(sIn.shape))
+   try:
+      for x in model.lutSE:
+         sIn[:,model.lutSE[x]]*=exposureScale
+   except IndexError:
+      #if there are no initial conditions, skip this part
+      pass
    setup['t0']=t0 
    strides=getStrides(setup)
    for step in strides:
+
+      #check if solution exists
+      readable=['t','sol','se','qt','sOut']
+      l=step['label']
+      files=[os.path.join(jobDir,'{}{}.txt'.format(x,l)) for x in readable]
+      filesPresent=all( [os.path.isfile(x) for x in files])
+      if filesPresent:
+         #set t0,y0,sIn
+         latestFiles={x:'{}{}.txt'.format(x,l) for x in readable}
+         data=loadSolutionFromFiles(setup,jobDir,[latestFiles])
+         startPoint=startPointObject(data)
+         t0=startPoint['t0']
+         y0=startPoint['y0']
+         s0=startPoint['s0']
+         for x in model.lutSE:
+            sIn[:,model.lutSE[x]]=s0[:,startPoint['lutSE'][x]]
+         print('Completed step {}: {}'.format(l,filesPresent))
+         continue
+      
+      print('Step required {}'.format(l))
+
       t1=t0+step['length']
       start_time=time.time()
       if setup['mode']=='SequentialOdeint':
@@ -152,6 +186,11 @@ def main(parFiles,jobDir,startDir='NONE'):
             rtol=setup['rtol'],method=setup['method'],\
             t0=t0,y0=y0,sIn=sIn)
 
+      if setup['mode']=='solveMatrix':
+         t,sol,se,s1=solveMatrix.solveMatrix(model,t1,\
+            nt=setup['nt'],\
+            t0=t0,y0=y0,sIn=sIn,method=setup['method'])
+
       end_time=time.time()
       #interpolate on s1
       #s1 has shape ft x n x m
@@ -164,7 +203,6 @@ def main(parFiles,jobDir,startDir='NONE'):
 
       #store
       writeable={'t':t,'sol':sol,'se':se,'qt':qt,'sOut':sOut}
-      l=step['label']
       writeable={os.path.join(jobDir,'{}{}.txt'.format(x,l)):writeable[x]\
          for x in writeable}
       for x in writeable:
@@ -246,6 +284,7 @@ def parseShape(line):
       
 def read3D(fName):
    new_data=numpy.loadtxt(fName)
+#loadtxt will ignore lines starting with #, so that is where we store contextual information like shape
    #read and parse first three lines to get shape,lut and lutSE
    nLines=3
    with open(fName,'r') as f:
@@ -348,6 +387,9 @@ def loadSolutionFromDir(jobDir,loadAll=False):
    setupFile=os.path.join(jobDir,'setup.json')
    setup=parseSetup(setupFile)
    inFiles=getFiles(setup,loadAll)
+   return loadSolutionFromFiles(setup,jobDir,inFiles)
+
+def loadSolutionFromFiles(setup,jobDir,inFiles):
    isFirst=True 
    operators={x:numpy.loadtxt for x in inFiles[0]}
    operators['sOut']=read3D