Changeset 25163


Ignore:
Timestamp:
06/26/20 11:03:14 (5 years ago)
Author:
jdquinn
Message:

CHG: MATLAB -> Python translations

Location:
issm/trunk-jpl
Files:
1 added
12 edited

Legend:

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

    r25125 r25163  
    1010
    1111class slr(object):
    12     """
     12    '''
    1313    SLR class definition
    1414
    1515        Usage:
    1616          slr = slr()
    17     """
     17    '''
    1818    def __init__(self):  # {{{
    19         self.deltathickness = float('NaN')
    20         self.sealevel = float('NaN')
    21         self.spcthickness = float('NaN')
     19        self.deltathickness = np.nan
     20        self.sealevel = np.nan
     21        self.spcthickness = np.nan
    2222        self.maxiter = 0
    2323        self.reltol = 0
     
    4242        self.degacc = 0
    4343        self.horiz = 0
    44         self.Ngia = float('NaN')
    45         self.Ugia = float('NaN')
    4644        self.planetradius = planetradius('earth')
    4745        self.requested_outputs = []
     
    8886    def setdefaultparameters(self):  # {{{
    8987        #Convergence criterion: absolute, relative and residual
    90         self.reltol = 0.01  #default
    91         self.abstol = float('NaN') #1 mm of sea level rise
     88        self.reltol = 0.01 #default
     89        self.abstol = np.nan #1 mm of sea level rise
    9290        #maximum of non - linear iterations.
    9391        self.maxiter = 5
     
    110108        #numerical discretization accuracy
    111109        self.degacc = 0.01
    112         #steric:
    113         self.steric_rate = 0
     110        #hydro
    114111        self.hydro_rate = 0
    115112        #how many time steps we skip before we run SLR solver during transient
     
    148145        md = checkfield(md, 'fieldname', 'slr.maxiter', 'size', [1, 1], '>=', 1)
    149146        md = checkfield(md, 'fieldname', 'slr.geodetic_run_frequency', 'size', [1, 1], '>=', 1)
    150         md = checkfield(md, 'fieldname', 'slr.steric_rate', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
    151147        md = checkfield(md, 'fieldname', 'slr.hydro_rate', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
    152148        md = checkfield(md, 'fieldname', 'slr.degacc', 'size', [1, 1], '>=', 1e-10)
    153149        md = checkfield(md, 'fieldname', 'slr.requested_outputs', 'stringrow', 1)
    154150        md = checkfield(md, 'fieldname', 'slr.horiz', 'NaN', 1, 'Inf', 1, 'values', [0, 1])
    155         md = checkfield(md, 'fieldname', 'slr.Ngia', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
    156         md = checkfield(md, 'fieldname', 'slr.Ugia', 'NaN', 1, 'Inf', 1, 'size', [md.mesh.numberofvertices])
    157151
    158152        #check that love numbers are provided at the same level of accuracy:
     
    199193        WriteData(fid, prefix, 'object', self, 'fieldname', 'ocean_area_scaling', 'format', 'Boolean')
    200194        WriteData(fid, prefix, 'object', self, 'fieldname', 'geodetic_run_frequency', 'format', 'Integer')
    201         WriteData(fid, prefix, 'object', self, 'fieldname', 'steric_rate', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1e-3 / md.constants.yts)
    202195        WriteData(fid, prefix, 'object', self, 'fieldname', 'hydro_rate', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1e-3 / md.constants.yts)
    203         WriteData(fid, prefix, 'object', self, 'fieldname', 'Ngia', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1e-3 / md.constants.yts)
    204         WriteData(fid, prefix, 'object', self, 'fieldname', 'Ugia', 'format', 'DoubleMat', 'mattype', 1, 'scale', 1e-3 / md.constants.yts)
    205196        WriteData(fid, prefix, 'object', self, 'fieldname', 'degacc', 'format', 'Double')
    206197        WriteData(fid, prefix, 'object', self, 'fieldname', 'transitions', 'format', 'MatArray')
  • issm/trunk-jpl/src/m/coordsystems/epsg2proj.py

    r25100 r25163  
    2121
    2222    #First, get GDAL version
    23     subproc_args = "gdalsrsinfo --version | awk '{print $2}' | cut -d '.' -f1"
     23    subproc_args = shlex.split("gdalsrsinfo --version | awk '{print $2}' | cut -d '.' -f1")
    2424    subproc = subprocess.Popen(subproc_args, shell=True, stdout=subprocess.PIPE, stderr=subprocess.PIPE)
    2525    outs, errs = subproc.communicate()
     
    2929    version_major=int(outs)
    3030
    31     subproc_args = "gdalsrsinfo epsg:{} | grep PROJ.4 | tr -d '\n' | sed 's/PROJ.4 : //'".format(epsg)
     31    subproc_args = shlex.split("gdalsrsinfo epsg:{} | grep PROJ.4 | tr -d '\n' | sed 's/PROJ.4 : //'".format(epsg))
    3232    subproc = subprocess.Popen(subproc_args, shell=True, stdout=subprocess.PIPE, stderr=subprocess.PIPE)
    3333    outs, errs = subproc.communicate()
  • issm/trunk-jpl/src/m/coordsystems/gdaltransform.py

    r25125 r25163  
    3434
    3535    TODO:
    36     - Remove shlex and follow subprocess pattern implemented in src/m/coordsystems/epsg2proj.py
     36    - Follow subprocess pattern implemented in src/m/coordsystems/epsg2proj.py
    3737    '''
    3838
     
    4646    file_in.close()
    4747
    48     args = shlex.split('gdaltransform -s_srs {} -t_srs {} < {} > {}'.format(proj_in, proj_out, filename_in, filename_out))
    49     subprocess.check_call(args) # Will raise CalledProcessError if return code is not 0
     48    subproc_args = shlex.split('gdaltransform -s_srs {} -t_srs {} < {} > {}'.format(proj_in, proj_out, filename_in, filename_out))
     49    subprocess.check_call(subproc_args) # Will raise CalledProcessError if return code is not 0
    5050
    5151    A = loadvars(filename_out)
  • issm/trunk-jpl/src/m/partition/AreaAverageOntoPartition.py

    r25110 r25163  
     1import copy
     2
    13import numpy as np
    24
    35from adjacency import adjacency
    4 import copy
    56from project2d import project2d
    67from qmupart2npart import qmupart2npart
  • issm/trunk-jpl/src/m/plot/processmesh.m

    r25058 r25163  
    1515end
    1616
    17 %special case for mesg 2dvertical
     17%special case for mesh 2dvertical
    1818if strcmp(domaintype(md.mesh),'2Dvertical'),
    1919        [x y z elements is2d isplanet] = processmesh(md.mesh,options);
     
    6464
    6565%is it a 2d plot?
    66 if md.mesh.dimension()==2,
     66if md.mesh.dimension()==2 | getfieldvalue(options,'layer',0)>=1,
    6767        is2d=1;
    6868else
    69         if getfieldvalue(options,'layer',0)>=1,
    70                 is2d=1;
    71         else
    72                 is2d=0;
    73         end
     69        is2d=0;
    7470end
    7571
     
    7975                error('processmesh error message: cannot work with 3D meshes for now');
    8076        end
    81         %we modify the mesh temporarily to a 2d mesh from which the 3d mesh was extruded.
     77        %we modify the mesh temporarily to a 2d mesh from which the 3d mesh was extruded
    8278        x=x2d;
    8379        y=y2d;
  • issm/trunk-jpl/src/m/plot/processmesh.py

    r24213 r25163  
     1#import matplotlib.delaunay as delaunay
    12import numpy as np
    23
    34
    45def processmesh(md, data, options):
    5     """
    6     PROCESSMESH - process the mesh for plotting
     6    '''
     7    PROCESSMESH - process mesh to be plotted
    78
    89    Usage:
    9     x, y, z, elements, is2d = processmech(md, data, options)
     10        x, y, z, elements, is2d = processmech(md, data, options)
    1011
    1112    See also: PLOTMODEL, PROCESSDATA
    12     """
    1313
    14     # {{{ check mesh size parameters
     14    TODO:
     15    - Test application of matlplotlib.delaunay
     16    - Test that output of delaunay matches output of delaunay in MATLAB (see
     17    src/m/plot/processmesh.m)
     18    '''
     19
     20    #some checks
    1521    if md.mesh.numberofvertices == 0:
    1622        raise ValueError('processmesh error: mesh is empty')
    1723    if md.mesh.numberofvertices == md.mesh.numberofelements:
    1824        raise ValueError('processmesh error: the number of elements is the same as the number of nodes')
    19     # }}}
    20     # {{{ treating coordinates
    2125
    22     try:
     26    #special case for mesh 2dvertical
     27    if md.mesh.domaintype() == '2Dvertical':
     28        x, y, z, elements, is2d, isplanet = processmesh(md.mesh, options)
     29        return x, y, z, elements, is2d, isplanet
     30
     31    # # special case for mesh 3dsurface
     32    # if md.mesh.domaintype() == '3Dsurface':
     33    #     x, y, z, elements, is2d, isplanet = processmesh(md.mesh, options)
     34    #     if options.getfieldvalue('coord', 'xy') == 'latlon':
     35    #         x = md.mesh.long
     36    #         y = md.mesh.lat
     37    #         elements = delaunay(x, y)
     38    #         z = md.mesh.lat
     39    #         z[:] = 0
     40    #     return x, y, z, elements, is2d, isplanet
     41
     42    if hasattr(md.mesh, 'elements2d'):
     43        elements2d = md.mesh.elements2d
     44
     45    if options.exist('amr'):
     46        step = options.getfieldvalue('amr')
     47        x = md.results.TransientSolution[step].MeshX
     48        y = md.results.TransientSolution[step].MeshY
     49        elements = md.results.TransientSolution[step].MeshElements
     50    else:
     51        elements = md.mesh.elements
     52        if options.getfieldvalue('coord', 'xy') != 'latlon':
     53            x = md.mesh.x
     54            if hasattr(md.mesh, 'x2d'):
     55                x2d = md.mesh.x2d
     56            y = md.mesh.y
     57            if hasattr(md.mesh, 'y2d'):
     58                y2d = md.mesh.y2d
     59
     60    if hasattr(md.mesh, 'z'):
    2361        z = md.mesh.z
    24     except AttributeError:
     62    else:
    2563        z = np.zeros(np.shape(md.mesh.x))
    26     elements = md.mesh.elements - 1
    27 
    28     if options.getfieldvalue('layer', 0) >= 1:
    29         x = md.mesh.x2d
    30         y = md.mesh.y2d
    31         z = np.zeros(np.shape(x))
    32         elements = md.mesh.elements2d - 1
    33     elif 'latlon' in options.getfieldvalue('coord', 'xy'):
    34         x = md.mesh.long
    35         y = md.mesh.lat
    36     else:
    37         x = md.mesh.x
    38         y = md.mesh.y
     64    z = options.getfieldvalue('z', z)
     65    if isinstance(z, basestring):
     66        z = getattr(md, z)
    3967
    4068    #is it a 2D plot?
     
    4371    else:
    4472        is2d = 0
     73
     74    elements = md.mesh.elements - 1
     75
     76    #layer projection?
     77    if options.getfieldvalue('layer', 0) >= 1:
     78        if options.getfieldvalue('coord', 'xy') == 'latlon':
     79            raise Exception('processmesh error message: cannot work with 3D meshes for now')
     80        #we modify the mesh temporarily to a 2d mesh from which the 3d mesh was extruded
     81        x = x2d
     82        y = y2d
     83        z = np.zeros(np.shape(x2d))
     84        elements = md.mesh.elements2d - 1
    4585
    4686    #units
     
    5191        z = z * unit
    5292
    53     #is model a member of planet class? (workaround until planet class defined)
     93    #is model a member of planet class?
     94    #
     95    # TODO: Change this when planet class defined (see src/m/plot/processmesh.m)
     96    #
    5497    if md.__class__.__name__ != 'model':
    5598        isplanet = 1
  • issm/trunk-jpl/src/m/qmu/dakota_in_data.py

    r24213 r25163  
    99def dakota_in_data(dmeth, variables, responses, dparams, filei, *args):
    1010    '''
    11   define the data to write the dakota .in and .m files.
     11    DAKOTA_IN_DATA - define the data to write the dakota .in and .m files.
     12   
     13    Usage:
     14        [] = dakota_in_data(dmeth, variables, responses, dparams, filei, *args)
    1215
    13   [] = dakota_in_data(dmeth, variables, responses, dparams, filei, *args)
     16    where the required input is:
     17        dmeth         (dakota_method, method class object)
     18        variables     (structure array, variable class objects)
     19        responses     (structure array, response class objects)
     20        dparams       (structure array, method - independent parameters)
     21        filei         (character, name of .in and .m files)
    1422
    15   where the required input is:
    16     dmeth         (dakota_method, method class object)
    17     variables     (structure array, variable class objects)
    18     responses     (structure array, response class objects)
    19     dparams       (structure array, method - independent parameters)
    20     filei         (character, name of .in and .m files)
     23    Params may be empty, in which case defaults will be used.
    2124
    22   params may be empty, in which case defaults will be used.
     25    The optional args are passed directly through to the QmuUpdateFunctions
     26    brancher to be used by the analysis package. For example, this could be
     27    model information.
    2328
    24   the optional args are passed directly through to the
    25   QmuUpdateFunctions brancher to be used by the analysis
    26   package.  for example, this could be model information.
     29    This function defines the data to write the dakota .in and .m files. It is
     30    necessary for multiple reasons. First, it collects the parameters and
     31    applies some defaults that are unique to the environment. Second, some
     32    analysis package ariables and/or responses may be treated differently by
     33    Dakota. For example, an analysis package variable may be defined as an
     34    array, so the QmuSetupDesign brancher will create Dakota variables for each
     35    element of the array. Finally, it calls the functions to write the .in and
     36    .m files. this function is independent of the particular analysis package.
    2737
    28   this function defines the data to write the dakota .in and
    29   .m files.  it is necessary for multiple reasons.  first,
    30   it collects the parameters and applies some defaults that
    31   are unique to the environment.  second, some analysis package
    32   variables and / or responses may be treated differently by
    33   dakota.  for example, an analysis package variable may be
    34   defined as an array, so the QmuSetupDesign brancher will
    35   create dakota variables for each element of the array.
    36   finally it calls the functions to write the .in and .m files.
    37   this function is independent of the particular analysis
    38   package.
    39 
    40   this data would typically be generated by a matlab script
    41   for a specific model, using the method, variable, and
    42   response class objects.
    43 '''
     38    This data would typically be generated by a MATLAB script for a specific
     39    model, using the method, variable, and response class objects.
     40    '''
    4441
    4542    #  parameters
  • issm/trunk-jpl/src/m/qmu/postqmu.m

    r25057 r25163  
    11function md=postqmu(md)
    2 %INPUT function md=postqmu(md)
    3 %Deal with dakota output results in files.
     2%POSTQMU - Deal with Dakota output results in files
     3%
     4%   Usage:
     5%      md = postqmu(md)
    46
    57%  check to see if dakota returned errors in the err file
     
    1416           fline=fgetl(fide);
    1517       end
    16        warning(['Dakota returned error in ''' qmuerrfile ' file.  qmu directory retained.'])
     18       warning(['Dakota returned error in ''' qmuerrfile ' file. qmu directory retained.'])
    1719   end
    1820   status=fclose(fide);
     
    4648end
    4749
    48 
    4950%put dakotaresults in their right location.
    5051md.results.dakota=dakotaresults;
  • issm/trunk-jpl/src/m/qmu/postqmu.py

    r24541 r25163  
    11from os import getpid, stat
    22from os.path import isfile
    3 from subprocess import Popen
     3import shlex
     4from subprocess import call
     5
    46from dakota_out_parse import *
    57from helpers import *
     8from loadresultsfromdisk import *
    69
    710
    811def postqmu(md):
    912    '''
    10     Deal with dakota output results in files.
     13    POSTQMU - Deal with Dakota output results in files
     14   
     15        Usage:
     16            md = postqmu(md)
    1117
    12     INPUT function
    13     md = postqmu(md)
    14 
     18        TODO:
     19        - Run Dakota test to check that updates from 6/26 are working
     20        - Add checks to Popen
    1521    '''
    1622
     
    4147        dakotaresults.dresp_dat = dresp_dat
    4248
     49    if md.qmu.output:
     50        if md.qmu.method.method == 'nond_sampling':
     51            dakotaresults.modelresults = []
     52            md2 = copy.deepcopy(md)
     53            md2.qmu.isdakota = 0
     54            for i in range(md2.qmu.method.params.samples):
     55                print('reading qmu file {}.outbin.{}'.format(md2.miscellaneous.name, i))
     56                md2 = loadresultsfromdisk(md2, '{}.outbin.{}'.format(md2.miscellaneous.name, i))
     57                dakotaresults.modelresults.append(md2.results)
     58
    4359    # put dakotaresults in their right location.
    4460    md.results.dakota = dakotaresults
     
    4662    # move all the individual function evalutations into zip files
    4763    if not md.qmu.isdakota:
    48         Popen('zip - mq params.in.zip params.in.[1 - 9] * ', shell=True)
    49         Popen('zip - mq results.out.zip results.out.[1 - 9] * ', shell=True)
    50         Popen('zip - mq matlab.out.zip matlab * .out.[1 - 9] * ', shell=True)
     64        subproc_args = shlex.split('zip -mq params.in.zip params.in.[1-9]*')
     65        call(subproc_args, shell=True)
     66        subproc_args = shlex.split('zip -mq results.out.zip results.out.[1-9]*')
     67        call(subproc_args, shell=True)
     68        subproc_args = shlex.split('zip -mq matlab.out.zip matlab*.out.[1-9]*')
     69        call(subproc_args, shell=True)
    5170
    5271    return md
  • issm/trunk-jpl/src/m/shp/shpwrite.py

    r25065 r25163  
    33except ImportError:
    44    print("could not import shapefile, PyShp has not been installed, no shapefile reading capabilities enabled")
     5
    56
    67def shpwrite(shp, filename): #{{{
  • issm/trunk-jpl/test/NightlyRun/test2010.m

    r25147 r25163  
    2929mask=gmtmask(md.mesh.lat,md.mesh.long);
    3030icemask=ones(md.mesh.numberofvertices,1);
    31 pos=find(mask==0);  icemask(pos)=-1;
    32 pos=find(sum(mask(md.mesh.elements),2)<3);   icemask(md.mesh.elements(pos,:))=-1;
     31pos=find(mask==0);
     32icemask(pos)=-1;
     33pos=find(sum(mask(md.mesh.elements),2)<3);
     34icemask(md.mesh.elements(pos,:))=-1;
    3335md.mask.ice_levelset=icemask;
    3436md.mask.ocean_levelset=-icemask;
     
    6567% }}}
    6668
    67 
    6869%eustatic + rigid + elastic run:
    69 md.solidearth.settings.rigid=1; md.solidearth.settings.elastic=1; md.solidearth.settings.rotation=1;
     70md.solidearth.settings.rigid=1;
     71md.solidearth.settings.elastic=1; md.solidearth.settings.rotation=1;
    7072md.cluster=generic('name',oshostname(),'np',3);
    7173
  • issm/trunk-jpl/test/NightlyRun/test2010.py

    r25158 r25163  
    2222longe = sum(md.mesh.long[md.mesh.elements - 1], 1) / 3
    2323
    24 md.solidearth.deltathickness = np.zeros((md.mesh.numberofelements, ))
     24md.solidearth.surfaceload.icethicknesschange = np.zeros(md.mesh.numberofelements)
    2525pos = np.intersect1d(np.array(np.where(late < -75)), np.array(np.where(longe > 0)))
    2626#python does not include last element in array slices, (6:7) -> [5:7]
    27 md.solidearth.deltathickness[pos[5:7]] = -1
     27md.solidearth.surfaceload.icethicknesschange[pos[5:7]] = -1
    2828
    29 md.solidearth.sealevel = np.zeros((md.mesh.numberofvertices, ))
    30 md.dsl.global_average_thermosteric_sea_level_change=np.zeros((2, ))
    31 md.dsl.sea_surface_height_change_above_geoid=np.zeros((md.mesh.numberofvertices+1, ))
    32 md.dsl.sea_water_pressure_change_at_sea_floor=np.zeros((md.mesh.numberofvertices+1, ))
     29md.solidearth.sealevel = np.zeros(md.mesh.numberofvertices)
     30md.dsl.global_average_thermosteric_sea_level_change = np.zeros((2, 1))
     31md.dsl.sea_surface_height_change_above_geoid = np.zeros(md.mesh.numberofvertices + 1)
     32md.dsl.sea_water_pressure_change_at_sea_floor = np.zeros(md.mesh.numberofvertices + 1)
    3333
    3434md.solidearth.ocean_area_scaling = 1
     
    4040#mask:  {{{
    4141mask = gmtmask(md.mesh.lat, md.mesh.long)
    42 icemask = np.ones((md.mesh.numberofvertices, ))
     42icemask = np.ones(md.mesh.numberofvertices)
    4343pos = np.where(mask == 0)
    4444icemask[pos[0]] = -1
    45 pos = np.where(sum(mask[md.mesh.elements - 1], 1) < 3)
     45pos = np.where(np.sum(mask[md.mesh.elements - 1], 1) < 3)
    4646icemask[md.mesh.elements[pos, :] - 1] = -1
    4747md.mask.ice_levelset = icemask
    4848md.mask.ocean_levelset = -icemask
    4949
    50 #make sure that the ice level set is all inclusive:
    51 md.mask.land_levelset = np.zeros((md.mesh.numberofvertices, ))
    52 md.mask.ocean_levelset = -np.ones((md.mesh.numberofvertices, ))
    53 
    54 #make sure that the elements that have loads are fully grounded:
    55 pos = np.nonzero(md.solidearth.deltathickness)[0]
    56 md.mask.ocean_levelset[md.mesh.elements[pos, :] - 1] = 1
     50#make sure that the elements that have loads are fully grounded
     51pos = np.nonzero(md.solidearth.surfaceload.icethicknesschange)[0]
     52md.mask.ocean_levelset = -np.ones(md.mesh.numberofvertices)
    5753
    5854#make sure wherever there is an ice load, that the mask is set to ice:
     55#pos = np.nonzero(md.solidearth.surfaceload.icethicknesschange)[0] # Do we need to do this twice?
     56md.mask.ice_levelset[md.mesh.elements[pos, :] - 1] = 1
     57
     58#make sure wherever there is an ice load, that the mask is set to ice
     59md.mask.ice_levelset = np.ones(md.mesh.numberofvertices)
    5960icemask[md.mesh.elements[pos, :] - 1] = -1
    6061md.mask.ice_levelset = icemask
    6162# }}}
     63
    6264#geometry {{{
    6365di = md.materials.rho_ice / md.materials.rho_water
    64 md.geometry.thickness = np.ones((md.mesh.numberofvertices, ))
    65 md.geometry.surface = (1 - di) * np.zeros((md.mesh.numberofvertices, ))
     66md.geometry.thickness = np.ones(md.mesh.numberofvertices)
     67md.geometry.surface = (1 - di) * np.zeros(md.mesh.numberofvertices)
    6668md.geometry.base = md.geometry.surface - md.geometry.thickness
    6769md.geometry.bed = md.geometry.base
    6870# }}}
    6971#materials {{{
    70 md.initialization.temperature = 273.25 * np.ones((md.mesh.numberofvertices, ))
     72md.initialization.temperature = 273.25 * np.ones(md.mesh.numberofvertices)
    7173md.materials.rheology_B = paterson(md.initialization.temperature)
    72 md.materials.rheology_n = 3 * np.ones((md.mesh.numberofelements, ))
     74md.materials.rheology_n = 3 * np.ones(md.mesh.numberofelements)
    7375# }}}
    7476#Miscellaneous {{{
     
    7678# }}}
    7779#Solution parameters {{{
    78 md.solidearth.reltol = float('NaN')
    79 md.solidearth.abstol = 1e-3
    80 md.solidearth.geodetic = 1
     80md.solidearth.settings.reltol = np.nan
     81md.solidearth.settings.abstol = 1e-3
     82md.solidearth.settings.geodetic = 1
    8183# }}}
    8284
    83 #New stuff
    84 md.solidearth.spcthickness = np.nan * np.ones((md.mesh.numberofvertices, ))
    85 md.solidearth.Ngia = np.zeros((md.mesh.numberofvertices, ))
    86 md.solidearth.Ugia = np.zeros((md.mesh.numberofvertices, ))
    87 md.solidearth.hydro_rate = np.zeros((md.mesh.numberofvertices, ))
    88 
    8985#eustatic + rigid + elastic run:
    90 md.solidearth.rigid = 1
    91 md.solidearth.elastic = 1
    92 md.solidearth.rotation = 1
     86md.solidearth.settings.rigid = 1
     87md.solidearth.settings.elastic = 1
     88md.solidearth.settings.rotation = 1
    9389md.cluster = generic('name', gethostname(), 'np', 3)
    9490
     
    9692md = solve(md, 'Sealevelrise')
    9793eus = md.results.SealevelriseSolution.SealevelRSLEustatic
    98 solidearth = md.results.SealevelriseSolution.Sealevel
     94slr = md.results.SealevelriseSolution.Sealevel
    9995moixz = md.results.SealevelriseSolution.SealevelInertiaTensorXZ
    10096moiyz = md.results.SealevelriseSolution.SealevelInertiaTensorYZ
     
    111107
    112108#Fields and tolerances to track changes
    113 field_names = ['eus', 'solidearth', 'moixz', 'moiyz', 'moizz']
     109field_names = ['eus', 'slr', 'moixz', 'moiyz', 'moizz']
    114110field_tolerances = [1e-13, 1e-13, 1e-13, 1e-13, 1e-13]
    115111field_values = [eus, solidearth, moixz, moiyz, moizz]
Note: See TracChangeset for help on using the changeset viewer.