No subject

Reading back the dumped output and save it to chkpt format won't be
difficult. However, I won't have time to work on it right now. I welcome
anybody who can contribute the code.


Cheers,
Eh




-- 
Eh Tan
Staff Scientist
Computational Infrastructure for Geodynamics
California Institute of Technology, 158-79
Pasadena, CA 91125
(626) 395-1693
http://www.geodynamics.org


--------------080309050805090500020908
Content-Type: text/x-python;
 name="chkptreader.py"
Content-Transfer-Encoding: 7bit
Content-Disposition: inline;
 filename="chkptreader.py"

#!/usr/bin/python

import sys
import numpy as np

class ChkptReader(object):

        def __init__(self, filename, tracer=False, composition=False):
                f = open(filename, 'rb')
                self.read_general_checkpoint(f)
                self.read_energy_checkpoint(f)
                self.read_momentum_checkpoint(f)

                if tracer:
                        self.read_tracer_checkpoint(f)
                        if composition:
                                self.read_composition_checkpoint(f)


        def read_sentinel(self, f):
                tmp = np.fromfile(f, count=4, dtype=np.int)
                for i in tmp:
                        if i != 0: raise ValueError('Error in reading checkpoint file: wrong sentinel')
                return


        def read_general_checkpoint(self, f):
                tmp = np.fromfile(f, count=8, dtype=np.int)
                self.nodex = tmp[0]
                self.nodey = tmp[1]
                self.nodez = tmp[2]
                self.step = tmp[7]
                self.nno = self.nodex * self.nodey * self.nodez
                self.nel = (self.nodex-1) * (self.nodey-1) * (self.nodez-1)

                tmp = np.fromfile(f, count=3, dtype=np.single)
                self.time = tmp[0]

                #print self.nodex, self.nodey, self.nodez, self.step, self.time
                return


        def read_energy_checkpoint(self, f):
                self.read_sentinel(f)

                self.t = np.fromfile(f, count=self.nno+1, dtype=np.double)[1:]
                self.t.shape = (self.nodey, self.nodex, self.nodez)
                self.temperature = self.temp = self.t

                self.tdot = np.fromfile(f, count=self.nno+1, dtype=np.double)[1:]
                self.tdot.shape = (self.nodey, self.nodex, self.nodez)
                #print self.t[1,1,:]
                return

        def read_momentum_checkpoint(self, f):
                self.read_sentinel(f)
                tmp = np.fromfile(f, count=2, dtype=np.single)
                self.P = np.fromfile(f, count=self.nel+1, dtype=np.double)[1:]
                self.P.shape = (self.nodey-1, self.nodex-1, self.nodez-1)

                self.v = np.fromfile(f, count=3*self.nno, dtype=np.double)
                self.v.shape = (self.nodey, self.nodex, self.nodez, 3)
                self.velocity = self.velo = self.v
                #print self.v[1,1,:,2]
                return

        def read_tracer_checkpoint(self, f):
                self.read_sentinel(f)
                tmp = np.fromfile(f, count=5, dtype=np.int)
                nbasicq = tmp[0]
                nextraq = tmp[1]
                ntracers = tmp[4]

                tmp = np.fromfile(f, count=(6+nextraq)*(ntracers+1), dtype=np.double)
                tmp.shape = (6+nextraq, ntracers+1)
                self.tracer_tfr = tmp[:3,1:]
                self.tracer_xyz = tmp[3:6,1:]
                self.extraq = tmp[6:,1:]

                self.tracer_elem = np.fromfile(f, count=ntracers+1, dtype=np.int)[1:]
                #print self.tracer_tfr
                return

        def read_composition_checkpoint(self, f):
                self.read_sentinel(f)
                self.ncomp = np.fromfile(f, count=1, dtype=np.int)
                tmp = np.fromfile(f, count=2*self.ncomp, dtype=np.double)

                tmp = np.fromfile(f, count=self.ncomp*(self.nel+1), dtype=np.double)
                tmp.shape = (self.ncomp, self.nel+1)
                self.comp_el = tmp[:,1:]
                self.comp_el.shape = (self.ncomp, self.nodey-1, self.nodex-1, self.nodez-1)
                #print self.comp_el[:,0,0,:]
                return


if __name__ == '__main__':

        if len(sys.argv) < 2:
                print '''usuage: read_citcoms_chkpt.py chkpt_file [tracer_flag [composition_flag] ]'''
                sys.exit(1)

        chkpt = sys.argv[1]
        tracer_flag = composition_flag = False
        if len(sys.argv) > 2:
                tracer_flag = True
                if len(sys.argv) > 3:
                        composition_flag = True
        ChkptReader(chkpt, tracer_flag, composition_flag)



--------------080309050805090500020908--