Changeset 22267


Ignore:
Timestamp:
11/18/17 15:19:15 (7 years ago)
Author:
kruegern
Message:

ADD: added python versions and dependecies for tests: 243,260,261,293,340,342,343,344,350,430,435,436,437,438,439,441,442,460,461,462,463,464,465,540,701,808,2424,2425

Location:
issm/trunk-jpl
Files:
43 added
18 edited

Legend:

Unmodified
Added
Removed
  • issm/trunk-jpl/src/m/classes/SMBgemb.m

    r21232 r22267  
    358358           
    359359                        WriteData(fid,prefix,'data',dt,'name','md.smb.dt','format','Double','scale',yts);
    360            
     360
    361361                        % Check if smb_dt goes evenly into transient core time step
    362362                        if (mod(md.timestepping.time_step,dt) >= 1e-10)
  • issm/trunk-jpl/src/m/classes/SMBgradientsela.py

    r21469 r22267  
    99
    1010           Usage:
    11               SMBgradientsela=SMBgradientsela();
     11              SMBgradientsela=SMBgradientsela()
    1212        """
    1313
     
    1919                self.b_min   = float('NaN')
    2020                self.requested_outputs      = []
     21                self.setdefaultparameters()
    2122                #}}}
    2223        def __repr__(self): # {{{
    23                 string="   surface forcings parameters:"
     24                string = "   surface forcings parameters:"
     25                string+= '\n   SMB gradients ela parameters:'
    2426
    25                 string="%s\n%s"%(string,fielddisplay(self,'issmbgradientsela','is smb gradients ela method activated (0 or 1, default is 0)'))
    2627                string="%s\n%s"%(string,fielddisplay(self,'ela',' equilibrium line altitude from which deviation is used to calculate smb using the smb gradients ela method [m a.s.l.]'))
    2728                string="%s\n%s"%(string,fielddisplay(self,'b_pos',' vertical smb gradient (dB/dz) above ela'))
     
    4748                return self
    4849        #}}}
     50        def setdefaultparameters(self): # {{{
     51                self.b_max=9999.
     52                self.b_min=-9999.
     53                return self
     54        #}}}
    4955        def checkconsistency(self,md,solution,analyses):    # {{{
    5056
     
    5662                        md = checkfield(md,'fieldname','smb.b_min','timeseries',1,'NaN',1,'Inf',1)
    5763
    58                 md = checkfield(md,'fieldname','masstransport.requested_outputs','stringrow',1)
     64                md = checkfield(md,'fieldname','smb.requested_outputs','stringrow',1)
    5965                return md
    6066        # }}}
  • issm/trunk-jpl/src/m/classes/amr.py

    r22241 r22267  
    2828        self.deviatoricerror_threshold = 0.
    2929        self.deviatoricerror_maximum    = 0.
     30        self.restart=0.
    3031        #set defaults
    3132        self.setdefaultparameters()
     
    4849        string="%s\n%s"%(string,fielddisplay(self,"deviatoricerror_threshold","maximum threshold deviatoricstress error permitted"))
    4950        string="%s\n%s"%(string,fielddisplay(self,"deviatoricerror_maximum","maximum deviatoricstress error permitted"))
     51        string="%s\n%s"%(string,fielddisplay(self,'restart',['indicates if ReMesh() will call before first time step']))
    5052        return string
    5153    #}}}
     
    6769        self.deviatoricerror_threshold = 0
    6870        self.deviatoricerror_maximum    = 0
     71        self.restart = 0.
    6972        return self
    7073    #}}}
     
    8487        md = checkfield(md,'fieldname','amr.deviatoricerror_threshold','numel',[1],'>=',0,'<=',1,'NaN',1);       
    8588        md = checkfield(md,'fieldname','amr.deviatoricerror_maximum','numel',[1],'>=',0,'NaN',1,'Inf',1);
     89        md = checkfield(md,'fieldname','amr.restart','numel',[1],'>=',0,'<=',1,'NaN',1)
    8690        return md
    8791    # }}}
     
    104108        WriteData(fid,prefix,'object',self,'fieldname','deviatoricerror_threshold','format','Double');
    105109        WriteData(fid,prefix,'object',self,'fieldname','deviatoricerror_maximum','format','Double');
     110        WriteData(fid,prefix,'object',self,'class','amr','fieldname','restart','format','Integer')
    106111    # }}}
  • issm/trunk-jpl/src/m/classes/basalforcings.py

    r21303 r22267  
    4545                        self.floatingice_melting_rate=np.zeros((md.mesh.numberofvertices))
    4646                        print "      no basalforcings.floatingice_melting_rate specified: values set as zero"
     47                #if np.all(np.isnan(self.geothermalflux)):
     48                        #self.geothermalflux=np.zeros((md.mesh.numberofvertices))
     49                        #print "      no basalforcings.geothermalflux specified: values set as zero"
    4750
    4851                return self
  • issm/trunk-jpl/src/m/classes/constants.py

    r22122 r22267  
    1212
    1313        def __init__(self): # {{{
    14                 self.g                    = 0
    15                 self.omega                = 0
    16                 self.yts                  = 0
    17                 self.referencetemperature = 0
     14                self.g                    = 0.
     15                self.omega                = 0.
     16                self.yts                  = 0.
     17                self.referencetemperature = 0.
    1818               
    1919                #set defaults
     
    4949        def checkconsistency(self,md,solution,analyses):    # {{{
    5050
    51                 md = checkfield(md,'fieldname','constants.g','>',0,'size',[1])
    52                 md = checkfield(md,'fieldname','constants.omega','>=',0,'size',[1])
    53                 md = checkfield(md,'fieldname','constants.yts','>',0,'size',[1])
    54                 md = checkfield(md,'fieldname','constants.referencetemperature','size',[1])
     51                md = checkfield(md,'fieldname','constants.g','>=',0,'size',[1,1])
     52                md = checkfield(md,'fieldname','constants.omega','>=',0,'size',[1,1])
     53                md = checkfield(md,'fieldname','constants.yts','>',0,'size',[1,1])
     54                md = checkfield(md,'fieldname','constants.referencetemperature','size',[1,1])
    5555
    5656                return md
  • issm/trunk-jpl/src/m/classes/flowequation.py

    r21303 r22267  
    9292                md = checkfield(md,'fieldname','flowequation.fe_SSA','values',['P1','P1bubble','P1bubblecondensed','P2','P2bubble'])
    9393                md = checkfield(md,'fieldname','flowequation.fe_HO' ,'values',['P1','P1bubble','P1bubblecondensed','P1xP2','P2xP1','P2','P2bubble','P1xP3','P2xP4'])
    94                 md = checkfield(md,'fieldname','flowequation.fe_FS' ,'values',['P1P1','P1P1GLS','MINIcondensed','MINI','TaylorHood','XTaylorHood','OneLayerP4z','CrouzeixRaviart'])
     94                md = checkfield(md,'fieldname','flowequation.fe_FS' ,'values',['P1P1','P1P1GLS','MINIcondensed','MINI','TaylorHood','LATaylorHood','XTaylorHood','OneLayerP4z','CrouzeixRaviart','LACrouzeixRaviart'])
    9595                md = checkfield(md,'fieldname','flowequation.borderSSA','size',[md.mesh.numberofvertices],'values',[0,1])
    9696                md = checkfield(md,'fieldname','flowequation.borderHO','size',[md.mesh.numberofvertices],'values',[0,1])
     
    104104                        md = checkfield(md,'fieldname','flowequation.vertex_equation','size',[md.mesh.numberofvertices],'values',[1,2])
    105105                        md = checkfield(md,'fieldname','flowequation.element_equation','size',[md.mesh.numberofelements],'values',[1,2])
     106                elif m.strcmp(md.mesh.domaintype(),'2Dvertical'):
     107                        md = checkfield(md,'fieldname','flowequation.vertex_equation','size',[md.mesh.numberofvertices],'values',[2,4,5])
     108                        md = checkfield(md,'fieldname','flowequation.element_equation','size',[md.mesh.numberofelements],'values',[2,4,5])
    106109                elif m.strcmp(md.mesh.domaintype(),'3D'):
    107110                        md = checkfield(md,'fieldname','flowequation.vertex_equation','size',[md.mesh.numberofvertices],'values',np.arange(0,8+1))
  • issm/trunk-jpl/src/m/classes/mismipbasalforcings.py

    r21354 r22267  
    4040    def initialize(self,md): # {{{
    4141        if np.all(np.isnan(self.groundedice_melting_rate)):
    42             self.groundedice_melting_rate=np.zeros(md.mesh.numberofvertices)
     42            self.groundedice_melting_rate=np.zeros((md.mesh.numberofvertices))
    4343            print ' no basalforcings.groundedice_melting_rate specified: values set as zero'
     44        if np.all(np.isnan(self.geothermalflux)):
     45                        self.geothermalflux=np.zeros((md.mesh.numberofvertices))
     46                        print "      no basalforcings.geothermalflux specified: values set as zero"
    4447        return self
    4548    #}}}
     
    8487
    8588        floatingice_melting_rate = np.zeros((md.mesh.numberofvertices))
    86         floatingice_melting_rate = md.basalforcings.meltrate_factor*np.tanh((md.geometry.base-md.geometry.bed)/md.basalforcings.threshold_thickness)*np.amax(md.basalforcings.upperdepth_melt-md.geometry.base,0)
     89        floatingice_melting_rate = md.basalforcings.meltrate_factor*np.tanh((md.geometry.base-md.geometry.bed)/md.basalforcings.threshold_thickness)*(md.basalforcings.upperdepth_melt-md.geometry.base)
    8790
    8891        WriteData(fid,prefix,'name','md.basalforcings.model','data',3,'format','Integer')
  • issm/trunk-jpl/src/m/classes/taoinversion.py

    r21303 r22267  
    66from IssmConfig import IssmConfig
    77from marshallcostfunctions import marshallcostfunctions
     8from supportedcontrols import *
     9from supportedcostfunctions import *
    810
    9 
    10 class taoinversion:
     11class taoinversion(object):
    1112        def __init__(self):
    12                 iscontrol                   = 0
    13                 incomplete_adjoint          = 0
    14                 control_parameters          = float('NaN')
    15                 maxsteps                    = 0
    16                 maxiter                     = 0
    17                 fatol                       = 0
    18                 frtol                       = 0
    19                 gatol                       = 0
    20                 grtol                       = 0
    21                 gttol                       = 0
    22                 algorithm                   = ''
    23                 cost_functions              = float('NaN')
    24                 cost_functions_coefficients = float('NaN')
    25                 min_parameters              = float('NaN')
    26                 max_parameters              = float('NaN')
    27                 vx_obs                      = float('NaN')
    28                 vy_obs                      = float('NaN')
    29                 vz_obs                      = float('NaN')
    30                 vel_obs                     = float('NaN')
    31                 thickness_obs               = float('NaN')
    32                 surface_obs                 = float('NaN')
     13                self.iscontrol                   = 0
     14                self.incomplete_adjoint          = 0
     15                self.control_parameters          = float('NaN')
     16                self.maxsteps                    = 0
     17                self.maxiter                     = 0
     18                self.fatol                       = 0
     19                self.frtol                       = 0
     20                self.gatol                       = 0
     21                self.grtol                       = 0
     22                self.gttol                       = 0
     23                self.algorithm                   = ''
     24                self.cost_functions              = float('NaN')
     25                self.cost_functions_coefficients = float('NaN')
     26                self.min_parameters              = float('NaN')
     27                self.max_parameters              = float('NaN')
     28                self.vx_obs                      = float('NaN')
     29                self.vy_obs                      = float('NaN')
     30                self.vz_obs                      = float('NaN')
     31                self.vel_obs                     = float('NaN')
     32                self.thickness_obs               = float('NaN')
     33                self.surface_obs                 = float('NaN')
     34                self.setdefaultparameters()
    3335
    3436        def __repr__(self):
     
    98100                #several responses can be used:
    99101                self.cost_functions=101;
    100 
    101102                return self
    102103
     
    119120
    120121        def checkconsistency(self,md,solution,analyses):
    121                 if not self.control:
     122                if not self.iscontrol:
    122123                        return md
    123124                if not IssmConfig('_HAVE_TAO_')[0]:
     
    125126
    126127
    127                 num_controls= np.numel(md.inversion.control_parameters)
    128                 num_costfunc= np.size(md.inversion.cost_functions,2)
     128                num_controls= np.size(md.inversion.control_parameters)
     129                num_costfunc= np.size(md.inversion.cost_functions)
    129130
    130131                md = checkfield(md,'fieldname','inversion.iscontrol','values',[0, 1])
    131                 md = checkfield(md,'fieldname','inversion.incomplete_adjoint','values',[0, 1])
     132                md = checkfield(md,'fieldname','inversion.incomplete_adjoint','values',[0,1])
    132133                md = checkfield(md,'fieldname','inversion.control_parameters','cell',1,'values',supportedcontrols())
    133134                md = checkfield(md,'fieldname','inversion.maxsteps','numel',1,'>=',0)
     
    143144                PETSCMINOR = IssmConfig('_PETSC_MINOR_')[0]
    144145                if(PETSCMAJOR>3 or (PETSCMAJOR==3 and PETSCMINOR>=5)):
    145                         md = checkfield(md,'fieldname','inversion.algorithm','values',{'blmvm','cg','lmvm'})
     146                        md = checkfield(md,'fieldname','inversion.algorithm','values',['blmvm','cg','lmvm'])
    146147                else:
    147                         md = checkfield(md,'fieldname','inversion.algorithm','values',{'tao_blmvm','tao_cg','tao_lmvm'})
     148                        md = checkfield(md,'fieldname','inversion.algorithm','values',['tao_blmvm','tao_cg','tao_lmvm'])
    148149
    149150
    150                 md = checkfield(md,'fieldname','inversion.cost_functions','size',[1, num_costfunc],'values',supportedcostfunctions())
     151                md = checkfield(md,'fieldname','inversion.cost_functions','size', [num_costfunc],'values',supportedcostfunctions())
    151152                md = checkfield(md,'fieldname','inversion.cost_functions_coefficients','size',[md.mesh.numberofvertices, num_costfunc],'>=',0)
    152153                md = checkfield(md,'fieldname','inversion.min_parameters','size',[md.mesh.numberofvertices, num_controls])
     
    162163                        md = checkfield(md,'fieldname','inversion.vy_obs','size',[md.mesh.numberofvertices],'NaN',1,'Inf',1)
    163164
    164                 def marshall(self, md, fid):
     165        def marshall(self,prefix,md,fid):
    165166
    166                         yts=md.constants.yts;
    167                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','iscontrol','format','Boolean')
    168                         WriteData(fid,prefix,'name','md.inversion.type','data',1,'format','Integer')
    169                         if not self.iscontrol:
    170                                 return
    171                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','incomplete_adjoint','format','Boolean')
    172                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','maxsteps','format','Integer')
    173                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','maxiter','format','Integer')
    174                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','fatol','format','Double')
    175                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','frtol','format','Double')
    176                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','gatol','format','Double')
    177                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','grtol','format','Double')
    178                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','gttol','format','Double')
    179                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','algorithm','format','String')
    180                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','cost_functions_coefficients','format','DoubleMat','mattype',1)
    181                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','min_parameters','format','DoubleMat','mattype',3)
    182                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','max_parameters','format','DoubleMat','mattype',3)
    183                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','vx_obs','format','DoubleMat','mattype',1,'scale',1./yts)
    184                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','vy_obs','format','DoubleMat','mattype',1,'scale',1./yts)
    185                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','vz_obs','format','DoubleMat','mattype',1,'scale',1./yts)
    186                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','thickness_obs','format','DoubleMat','mattype',1)
    187                         WriteData(fid,prefix,'object',self,'class','inversion','fieldname','surface_obs','format','DoubleMat','mattype',1)
     167                yts=md.constants.yts;
     168                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','iscontrol','format','Boolean')
     169                WriteData(fid,prefix,'name','md.inversion.type','data',1,'format','Integer')
     170                if not self.iscontrol:
     171                        return
     172                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','incomplete_adjoint','format','Boolean')
     173                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','maxsteps','format','Integer')
     174                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','maxiter','format','Integer')
     175                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','fatol','format','Double')
     176                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','frtol','format','Double')
     177                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','gatol','format','Double')
     178                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','grtol','format','Double')
     179                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','gttol','format','Double')
     180                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','algorithm','format','String')
     181                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','cost_functions_coefficients','format','DoubleMat','mattype',1)
     182                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','min_parameters','format','DoubleMat','mattype',3)
     183                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','max_parameters','format','DoubleMat','mattype',3)
     184                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','vx_obs','format','DoubleMat','mattype',1,'scale',1./yts)
     185                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','vy_obs','format','DoubleMat','mattype',1,'scale',1./yts)
     186                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','vz_obs','format','DoubleMat','mattype',1,'scale',1./yts)
     187                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','thickness_obs','format','DoubleMat','mattype',1)
     188                WriteData(fid,prefix,'object',self,'class','inversion','fieldname','surface_obs','format','DoubleMat','mattype',1)
    188189
    189                         #process control parameters
    190                         num_control_parameters = np.numel(self.control_parameters)
    191                         WriteData(fid,prefix,'object',self,'fieldname','control_parameters','format','StringArray')
    192                         WriteData(fid,prefix,'data',num_control_parameters,'name','md.inversion.num_control_parameters','format','Integer')
     190                #process control parameters
     191                num_control_parameters = np.size(self.control_parameters)
     192                WriteData(fid,prefix,'object',self,'fieldname','control_parameters','format','StringArray')
     193                WriteData(fid,prefix,'data',num_control_parameters,'name','md.inversion.num_control_parameters','format','Integer')
    193194
    194                         #process cost functions
    195                         num_cost_functions = np.size(self.cost_functions,2)
    196                         data= marshallcostfunctions(self.cost_functions)
    197                         WriteData(fid,prefix,'data',data,'name','md.inversion.cost_functions','format','StringArray')
    198                         WriteData(fid,prefix,'data',num_cost_functions,'name','md.inversion.num_cost_functions','format','Integer')
     195                #process cost functions
     196                num_cost_functions = np.size(self.cost_functions)
     197                data= marshallcostfunctions(self.cost_functions)
     198                WriteData(fid,prefix,'data',data,'name','md.inversion.cost_functions','format','StringArray')
     199                WriteData(fid,prefix,'data',num_cost_functions,'name','md.inversion.num_cost_functions','format','Integer')
  • issm/trunk-jpl/src/m/classes/thermal.py

    r21542 r22267  
    101101
    102102                if 'EnthalpyAnalysis' in analyses and md.thermal.isenthalpy and md.mesh.dimension()==3:
    103                         pos=np.where(~np.isnan(md.thermal.spctemperature[0:md.mesh.numberofvertices]))
     103                        TEMP = md.thermal.spctemperature[:-1].flatten(-1)
     104                        pos=np.where(~np.isnan(TEMP))
    104105                        try:
    105106                                spccol=np.size(md.thermal.spctemperature,1)
    106107                        except IndexError:
    107108                                spccol=1
    108                         replicate=np.tile(md.geometry.surface-md.mesh.z,(spccol))
    109                         control=md.materials.meltingpoint-md.materials.beta*md.materials.rho_ice*md.constants.g*replicate
    110                         md = checkfield(md,'fieldname','thermal.spctemperature','field',md.thermal.spctemperature[pos],'<=',control[pos],'message',"spctemperature should be below the adjusted melting point")
     109
     110                        replicate=np.tile(md.geometry.surface-md.mesh.z,(spccol)).flatten(-1)
     111
     112                        control=md.materials.meltingpoint-md.materials.beta*md.materials.rho_ice*md.constants.g*replicate+10**-5
     113
     114                        md = checkfield(md,'fieldname','thermal.spctemperature','field',md.thermal.spctemperature.flatten(-1)[pos],'<=',control[pos],'message',"spctemperature should be below the adjusted melting point")
    111115                        md = checkfield(md,'fieldname','thermal.isenthalpy','numel',[1],'values',[0,1])
    112116                        md = checkfield(md,'fieldname','thermal.isdynamicbasalspc','numel',[1],'values',[0,1]);
  • issm/trunk-jpl/src/m/classes/transient.py

    r21802 r22267  
    1212
    1313        def __init__(self): # {{{
    14                 self.issmb   = False
     14                self.issmb             = False
    1515                self.ismasstransport   = False
    1616                self.isstressbalance   = False
     
    2323                self.ishydrology       = False
    2424                self.isslr             = False
     25                self.iscoupler         = False
     26                self.amr_frequency     = 0
    2527                self.isoceancoupling   = False
    26                 self.iscoupler         = False
    27                 amr_frequency                     = 0
    2828                self.requested_outputs = []
    2929
     
    7676                self.iscoupler         = False
    7777                self.amr_frequency        = 0
     78
     79                #default output
     80                self.requested_outputs=[]
     81                return self
     82        #}}}
     83        def deactivateall(self):#{{{
     84                self.issmb             = False
     85                self.ismasstransport   = False
     86                self.isstressbalance   = False
     87                self.isthermal         = False
     88                self.isgroundingline   = False
     89                self.isgia             = False
     90                self.isesa             = False
     91                self.isdamageevolution = False
     92                self.ismovingfront     = False
     93                self.ishydrology       = False
     94                self.isslr             = False
     95                self.isoceancoupling   = False
     96                self.iscoupler         = False
     97                self.amr_frequency     = 0
    7898
    7999                #default output
  • issm/trunk-jpl/src/m/consistency/checkfield.py

    r22150 r22267  
    8484        if options.exist('numel'):
    8585                fieldnumel=options.getfieldvalue('numel')
    86                 if np.size(field) not in fieldnumel:
     86                if (type(fieldnumel) == int and np.size(field) != fieldnumel) or (type(fieldnumel) == list and np.size(field) not in fieldnumel):
    8787                        if   len(fieldnumel)==1:
    8888                                md = md.checkmessage(options.getfieldvalue('message',\
     
    101101                                "NaN values found in field '%s'" % fieldname))
    102102
     103
    103104        #check Inf
    104105        if options.getfieldvalue('Inf',0):
     
    106107                        md = md.checkmessage(options.getfieldvalue('message',\
    107108                                "Inf values found in field '%s'" % fieldname))
     109
    108110
    109111        #check cell
     
    129131        #check greater
    130132        if options.exist('>='):
    131                 lowerbound=options.getfieldvalue('>=')
    132                 if np.any(field<lowerbound):
     133                lowerbound = options.getfieldvalue('>=')
     134                field2 = np.reshape(field,(np.prod(np.shape(field)),1),order='F').copy()
     135                if options.getfieldvalue('timeseries',0):
     136                        field2 = np.reshape(field[:-1],np.prod(np.shape(field[:-1])),1)
     137
     138                if options.getfieldvalue('singletimeseries',0):
     139                        field2 = np.reshape(field[0],np.prod(np.shape(field[0])),1)
     140
     141                if np.any(field2<lowerbound):
    133142                        md = md.checkmessage(options.getfieldvalue('message',\
    134143                                "field '%s' should have values above %d" % (fieldname,lowerbound)))
     144
    135145        if options.exist('>'):
    136146                lowerbound=options.getfieldvalue('>')
    137                 if np.any(field<=lowerbound):
     147                field2 = np.reshape(field,(np.prod(np.shape(field)),1),order='F').copy()
     148                if options.getfieldvalue('timeseries',0):
     149                        field2 = np.reshape(field[:-1],np.prod(np.shape(field[:-1])),1)
     150
     151                if options.getfieldvalue('singletimeseries',0):
     152                        field2 = np.reshape(field[0],np.prod(np.shape(field[0])),1)
     153
     154                if np.any(field2<=lowerbound):
    138155                        md = md.checkmessage(options.getfieldvalue('message',\
    139156                                "field '%s' should have values above %d" % (fieldname,lowerbound)))
     
    142159        if options.exist('<='):
    143160                upperbound=options.getfieldvalue('<=')
    144                 if np.any(field>upperbound):
     161                field2 = np.reshape(field,(np.prod(np.shape(field)),1),order='F').copy()
     162                if options.getfieldvalue('timeseries',0):
     163                        field2 = np.reshape(field[:-1],np.prod(np.shape(field[:-1])),1)
     164
     165                if options.getfieldvalue('singletimeseries',0):
     166                        field2 = np.reshape(field[0],np.prod(np.shape(field[0])),1)
     167
     168                if np.any(field2>upperbound):
    145169                        md = md.checkmessage(options.getfieldvalue('message',\
    146170                                "field '%s' should have values below %d" % (fieldname,upperbound)))
    147171        if options.exist('<'):
    148172                upperbound=options.getfieldvalue('<')
    149                 if np.any(field>=upperbound):
     173                field2 = np.reshape(field,(np.prod(np.shape(field)),1),order='F').copy()
     174                if options.getfieldvalue('timeseries',0):
     175                        field2 = np.reshape(field[:-1],np.prod(np.shape(field[:-1])),1)
     176
     177                if options.getfieldvalue('singletimeseries',0):
     178                        field2 = np.reshape(field[0],np.prod(np.shape(field[0])),1)
     179
     180                if np.any(field2>=upperbound):
    150181                        md = md.checkmessage(options.getfieldvalue('message',\
    151182                                "field '%s' should have values below %d" % (fieldname,upperbound)))
     
    164195        #Check forcings (size and times)
    165196        if options.getfieldvalue('timeseries',0):
    166                 if   np.size(field,0)==md.mesh.numberofvertices:
     197                if np.size(field,0)==md.mesh.numberofvertices or np.size(field,0)==md.mesh.numberofelements:
    167198                        if np.ndim(field)>1 and not np.size(field,1)==1:
    168199                                md = md.checkmessage(options.getfieldvalue('message',\
    169200                                        "field '%s' should have only one column as there are md.mesh.numberofvertices lines" % fieldname))
    170                 elif np.size(field,0)==md.mesh.numberofvertices+1 or np.size(field,0)==2:
    171                         if not all(field[-1,:]==np.sort(field[-1,:])):
     201                elif np.size(field,0)==md.mesh.numberofvertices+1 or np.size(field,0)==md.mesh.numberofelements+1:
     202                        if np.ndim(field) > 1 and not all(field[-1,:]==np.sort(field[-1,:])):
    172203                                md = md.checkmessage(options.getfieldvalue('message',\
    173204                                        "field '%s' columns should be sorted chronologically" % fieldname))
    174                         if any(field[-1,0:-1]==field[-1,1:]):
     205                        if np.ndim(field) > 1 and any(field[-1,0:-1]==field[-1,1:]):
    175206                                md = md.checkmessage(options.getfieldvalue('message',\
    176207                                        "field '%s' columns must not contain duplicate timesteps" % fieldname))
     
    198229        return md
    199230
     231
  • issm/trunk-jpl/src/m/dev/devpath.py

    r20679 r22267  
    3131#Manual imports for commonly used functions
    3232from plotmodel import plotmodel
     33from runme import runme
    3334
    3435#c = get_ipython().config
  • issm/trunk-jpl/src/m/mesh/bamg.py

    r22216 r22267  
    11import os.path
    22import numpy as  np
    3 from mesh2d import mesh2d
     3from mesh2dvertical import *
    44from collections import OrderedDict
    55from pairoptions import pairoptions
     
    8080                #Check that file exists
    8181                domainfile=options.getfieldvalue('domain')
    82                 if not os.path.exists(domainfile):
    83                         raise IOError("bamg error message: file '%s' not found" % domainfile)
    84                 domain=expread(domainfile)
     82                if type(domainfile) == str:
     83                        if not os.path.exists(domainfile):
     84                                raise IOError("bamg error message: file '%s' not found" % domainfile)
     85                        domain=expread(domainfile)
     86                else:
     87                        domain=domainfile
    8588
    8689                #Build geometry
     
    8992
    9093                        #Check that the domain is closed
    91                         if (domaini['x'][0]!=domaini['x'][-1] or domaini['y'][0]!=domaini['y'][-1]):
     94                        if (domaini.x[0]!=domaini.x[-1] or domaini.y[0]!=domaini.y[-1]):
    9295                                raise RuntimeError("bamg error message: all contours provided in ''domain'' should be closed")
    9396
    9497                        #Checks that all holes are INSIDE the principle domain outline
    9598                        if i:
    96                                 flags=ContourToNodes(domaini['x'],domaini['y'],domainfile,0)[0]
     99                                flags=ContourToNodes(domaini.x,domaini.y,domainfile,0)[0]
    97100                                if np.any(np.logical_not(flags)):
    98101                                        raise RuntimeError("bamg error message: All holes should be strictly inside the principal domain")
    99102
    100103                        #Add all points to bamg_geometry
    101                         nods=domaini['nods']-1    #the domain are closed 0=end
    102                         bamg_geometry.Vertices=np.vstack((bamg_geometry.Vertices,np.vstack((domaini['x'][0:nods],domaini['y'][0:nods],np.ones((nods)))).T))
     104                        nods=domaini.nods-1    #the domain are closed 0=end
     105                        bamg_geometry.Vertices=np.vstack((bamg_geometry.Vertices,np.vstack((domaini.x[0:nods],domaini.y[0:nods],np.ones((nods)))).T))
    103106                        bamg_geometry.Edges   =np.vstack((bamg_geometry.Edges,np.vstack((np.arange(count+1,count+nods+1),np.hstack((np.arange(count+2,count+nods+1),count+1)),1.*np.ones((nods)))).T))
    104107                        if i:
     
    115118                        #update counter
    116119                        count+=nods
     120
     121                if options.getfieldvalue('vertical',0):
     122                        if np.size(options.getfieldvalue('Markers',[]))!=np.size(bamg_geometry.Edges,0):
     123                                raise RuntimeError('for 2d vertical mesh, ''Markers'' option is required, and should be of size ' + str(np.size(bamg_geometry.Edges,0)))
     124                        if np.size(options.getfieldvalue('Markers',[]))==np.size(bamg_geometry.Edges,0):
     125                                bamg_geometry.Edges[:,2]=options.getfieldvalue('Markers')
    117126
    118127                #take care of rifts
     
    323332
    324333        # plug results onto model
     334        if options.getfieldvalue('vertical',0):
     335                md.mesh=mesh2dvertical()
     336                md.mesh.x=bamgmesh_out['Vertices'][:,0].copy()
     337                md.mesh.y=bamgmesh_out['Vertices'][:,1].copy()
     338                md.mesh.elements=bamgmesh_out['Triangles'][:,0:3].astype(int)
     339                md.mesh.edges=bamgmesh_out['IssmEdges'].astype(int)
     340                md.mesh.segments=bamgmesh_out['IssmSegments'][:,0:3].astype(int)
     341                md.mesh.segmentmarkers=bamgmesh_out['IssmSegments'][:,3].astype(int)
     342
     343                #Fill in rest of fields:
     344                md.mesh.numberofelements=np.size(md.mesh.elements,axis=0)
     345                md.mesh.numberofvertices=np.size(md.mesh.x)
     346                md.mesh.numberofedges=np.size(md.mesh.edges,axis=0)
     347                md.mesh.vertexonboundary=np.zeros(md.mesh.numberofvertices,bool)
     348                md.mesh.vertexonboundary[md.mesh.segments[:,0:2]-1]=True
     349
     350                #Now, build the connectivity tables for this mesh. Doubled in matlab for some reason
     351                md.mesh.vertexonboundary=np.zeros(md.mesh.numberofvertices,)
     352                md.mesh.vertexonboundary[md.mesh.segments[:,0:2]-1]=1
     353
     354        elif options.getfieldvalue('3dsurface',0):
     355                md.mesh=mesh3dsurface()
     356                md.mesh.x=bamgmesh_out['Vertices'][:,0].copy()
     357                md.mesh.y=bamgmesh_out['Vertices'][:,1].copy()
     358                md.mesh.z=md.mesh.x
     359                md.mesh.z[:]=0
     360                md.mesh.elements=bamgmesh_out['Triangles'][:,0:3].astype(int)
     361                md.mesh.edges=bamgmesh_out['IssmEdges'].astype(int)
     362                md.mesh.segments=bamgmesh_out['IssmSegments'][:,0:3].astype(int)
     363                md.mesh.segmentmarkers=bamgmesh_out['IssmSegments'][:,3].astype(int)
     364
     365                #Fill in rest of fields:
     366                md.mesh.numberofelements=np.size(md.mesh.elements,axis=0)
     367                md.mesh.numberofvertices=np.size(md.mesh.x)
     368                md.mesh.numberofedges=np.size(md.mesh.edges,axis=0)
     369                md.mesh.vertexonboundary=np.zeros(md.mesh.numberofvertices,bool)
     370                md.mesh.vertexonboundary[md.mesh.segments[:,0:2]-1]=True
     371
     372        else:
     373                md.mesh=mesh2d()
     374                md.mesh.x=bamgmesh_out['Vertices'][:,0].copy()
     375                md.mesh.y=bamgmesh_out['Vertices'][:,1].copy()
     376                md.mesh.elements=bamgmesh_out['Triangles'][:,0:3].astype(int)
     377                md.mesh.edges=bamgmesh_out['IssmEdges'].astype(int)
     378                md.mesh.segments=bamgmesh_out['IssmSegments'][:,0:3].astype(int)
     379                md.mesh.segmentmarkers=bamgmesh_out['IssmSegments'][:,3].astype(int)
     380
     381                #Fill in rest of fields:
     382                md.mesh.numberofelements=np.size(md.mesh.elements,axis=0)
     383                md.mesh.numberofvertices=np.size(md.mesh.x)
     384                md.mesh.numberofedges=np.size(md.mesh.edges,axis=0)
     385                md.mesh.vertexonboundary=np.zeros(md.mesh.numberofvertices,bool)
     386                md.mesh.vertexonboundary[md.mesh.segments[:,0:2]-1]=True
     387
     388        #Bamg private fields
    325389        md.private.bamg=OrderedDict()
    326390        md.private.bamg['mesh']=bamgmesh(bamgmesh_out)
    327391        md.private.bamg['geometry']=bamggeom(bamggeom_out)
    328         md.mesh = mesh2d()
    329         md.mesh.x=bamgmesh_out['Vertices'][:,0].copy()
    330         md.mesh.y=bamgmesh_out['Vertices'][:,1].copy()
    331         md.mesh.elements=bamgmesh_out['Triangles'][:,0:3].astype(int)
    332         md.mesh.edges=bamgmesh_out['IssmEdges'].astype(int)
    333         md.mesh.segments=bamgmesh_out['IssmSegments'][:,0:3].astype(int)
    334         md.mesh.segmentmarkers=bamgmesh_out['IssmSegments'][:,3].astype(int)
    335 
    336         #Fill in rest of fields:
    337         md.mesh.numberofelements=np.size(md.mesh.elements,axis=0)
    338         md.mesh.numberofvertices=np.size(md.mesh.x)
    339         md.mesh.numberofedges=np.size(md.mesh.edges,axis=0)
    340         md.mesh.vertexonboundary=np.zeros(md.mesh.numberofvertices,bool)
    341         md.mesh.vertexonboundary[md.mesh.segments[:,0:2]-1]=True
    342392        md.mesh.elementconnectivity=md.private.bamg['mesh'].ElementConnectivity
    343393        md.mesh.elementconnectivity[np.nonzero(np.isnan(md.mesh.elementconnectivity))]=0
  • issm/trunk-jpl/src/m/solve/WriteData.py

    r22204 r22267  
    5151                yts = options.getfieldvalue('yts')
    5252                if np.ndim(data) > 1:
    53                         data[-1,:] = data[-1,:]*yts
    54                 else:
    55                         data[-1] = data[-1]*yts
     53                        data[-1,:] = yts*data[-1,:]
     54                else:
     55                        data[-1] = yts*data[-1]
    5656
    5757        #Step 1: write the enum to identify this record uniquely
  • issm/trunk-jpl/src/m/solve/parseresultsfromdisk.py

    r21674 r22267  
    175175                if fieldname=='BalancethicknessThickeningRate':
    176176                        field = field*yts
    177                 elif fieldname=='Time':
    178                         field = field/yts
    179177                elif fieldname=='HydrologyWaterVx':
    180178                        field = field*yts
     
    191189                elif fieldname=='BasalforcingsGroundediceMeltingRate':
    192190                        field = field*yts
     191                elif fieldname=='BasalforcingsFloatingiceMeltingRate':
     192                        field = field*yts
    193193                elif fieldname=='TotalFloatingBmb':
    194                         field = field/10.**12.*yts #(GigaTon/year)
     194                        field = field/10.**12*yts #(GigaTon/year)
    195195                elif fieldname=='TotalGroundedBmb':
    196                         field = field/10.**12.*yts #(GigaTon/year)
     196                        field = field/10.**12*yts #(GigaTon/year)
    197197                elif fieldname=='TotalSmb':
    198                         field = field/10.**12.*yts #(GigaTon/year)
     198                        field = field/10.**12*yts #(GigaTon/year)
    199199                elif fieldname=='SmbMassBalance':
    200200                        field = field*yts
     201                elif fieldname=='SmbPrecipitation':
     202                        field = field*yts
     203                elif fieldname=='SmbRunoff':
     204                        field = field*yts
     205                elif fieldname=='SmbEC':
     206                        field = field*yts
     207                elif fieldname=='SmbAccumulation':
     208                        field = field*yts
     209                elif fieldname=='SmbMelt':
     210                        field = field*yts
     211                elif fieldname=='SmbDz_add':
     212                        field = field*yts
     213                elif fieldname=='SmbM_add':
     214                        field = field*yts
    201215                elif fieldname=='CalvingCalvingrate':
    202216                        field = field*yts
     217
     218                if time !=-9999:
     219                        time = time/yts
    203220
    204221                saveres=OrderedDict()
  • issm/trunk-jpl/test/Par/ISMIPE.par

    r21168 r22267  
    1414        md.geometry.surface(i)=(1.-coeff)*data(point1,3)+coeff*data(point2,3);
    1515end
     16
    1617md.geometry.thickness=md.geometry.surface-md.geometry.base;
    1718md.geometry.thickness(find(~md.geometry.thickness))=0.01;
  • issm/trunk-jpl/test/Par/ISMIPE.py

    r21409 r22267  
    1313        point1=numpy.floor(y/100.)
    1414        point2=numpy.minimum(point1+1,50)
    15         coeff=(y-(point1-1.)*100.)/100.
     15        coeff=(y-(point1)*100.)/100.
    1616        md.geometry.base[i]=(1.-coeff)*data[point1,1]+coeff*data[point2,1]
    1717        md.geometry.surface[i]=(1.-coeff)*data[point1,2]+coeff*data[point2,2]
     18
    1819md.geometry.thickness=md.geometry.surface-md.geometry.base
    1920md.geometry.thickness[numpy.nonzero(numpy.logical_not(md.geometry.thickness))]=0.01
  • issm/trunk-jpl/test/Par/SquareThermal.py

    r21409 r22267  
    2626print "      creating temperatures"
    2727md.initialization.temperature=(273.-20.)*numpy.ones((md.mesh.numberofvertices))
     28md.initialization.pressure=numpy.zeros((md.mesh.numberofvertices,))
     29md.initialization.waterfraction=numpy.zeros((md.mesh.numberofvertices,))
     30md.initialization.watercolumn=numpy.zeros((md.mesh.numberofvertices,))
    2831
    2932print "      creating flow law parameter"
     
    3336print "      creating surface mass balance"
    3437md.smb.mass_balance=numpy.ones((md.mesh.numberofvertices))/md.constants.yts    #1m/a
    35 md.basalforcings.melting_rate=0.*numpy.ones((md.mesh.numberofvertices))/md.constants.yts    #1m/a
     38#md.basalforcings.melting_rate=0.*numpy.ones((md.mesh.numberofvertices))/md.constants.yts    #1m/a
     39md.basalforcings.groundedice_melting_rate=0.*numpy.ones((md.mesh.numberofvertices))/md.constants.yts    #1m/a
     40md.basalforcings.floatingice_melting_rate=0.*numpy.ones((md.mesh.numberofvertices))/md.constants.yts    #1m/a
    3641
    3742#Deal with boundary conditions:
     
    4146
    4247print "      boundary conditions for thermal model"
    43 md.thermal.spctemperature=md.initialization.temperature
     48md.thermal.spctemperature[:]=md.initialization.temperature
    4449md.basalforcings.geothermalflux=numpy.zeros((md.mesh.numberofvertices))
    4550md.basalforcings.geothermalflux[numpy.nonzero(md.mask.groundedice_levelset>0.)[0]]=1.*10**-3    #1 mW/m^2
Note: See TracChangeset for help on using the changeset viewer.