Mesh Object

class pylagrit.MO(name, parent)

Mesh object class

addatt(attname, keyword=None, type='VDOUBLE', rank='scalar', length='nnodes', interpolate='linear', persistence='permanent', ioflag='', value=0.0)

Add a list of attributes

Parameters:
  • attnames (str) – Attribute name to add
  • keyword – Keyword used by lagrit for specific attributes
addatt_voronoi_volume(name='voronoi_volume')

Add voronoi volume attribute to mesh object

Parameters:name (str) – name of attribute in LaGriT
connect(option1='delaunay', option2=None, stride=None, big_tet_coords=[])

Connect the nodes into a Delaunay tetrahedral or triangle grid.

Parameters:
  • option1 – type of connect: delaunay, noadd, or check_interface
  • option2 – type of connect: noadd, or check_interface
  • stride (tuple(int)) – tuple of (first, last, stride) of points
connect_check_interface()

Connect the nodes into a Delaunay tetrahedral or triangle grid exhaustively checking that no edges of the mesh cross a material boundary.

connect_delaunay(option2=None, stride=None, big_tet_coords=[])

Connect the nodes into a Delaunay tetrahedral or triangle grid without adding nodes.

connect_noadd()

Connect the nodes into a Delaunay tetrahedral or triangle grid without adding nodes.

copy(name=None)

Copy mesh object

copyatt(attname_src, attname_sink=None, mo_src=None)

Add a list of attributes

Parameters:
  • attname_src (str) – Name of attribute to copy
  • attname_sink (str) – Name of sink attribute
  • mo_src (PyLaGriT Mesh Object) – Name of source mesh object
copypts(elem_type='tet', name=None)

Copy points from mesh object to new mesh object

Parameters:
  • name (str) – Name to use within lagrit for the created mesh object
  • mesh_type (str) – Mesh type for new mesh
Returns:

mesh object

create_boundary_facesets(stacked_layers=False, base_name=None)

Creates facesets for each boundary and writes associated avs faceset file :arg base_name: base name of faceset files :type base_name: str :arg stacked_layers: if mesh is created by stack_layers, user layertyp attr to determine top and bottom :type stacked_layers: bool :arg reorder_on_meds: reorder nodes on cell medians, usually needed for exodus file :type reorder_on_meds: bool :returns: Dictionary of facesets

createpts(crd, npts, mins, maxs, rz_switch=(1, 1, 1), rz_value=(1, 1, 1), connect=False)

Create and Connect Points

Parameters:
  • crd (str) – Coordinate type of either ‘xyz’ (cartesian coordinates), ‘rtz’ (cylindrical coordinates), or ‘rtp’ (spherical coordinates).
  • npts (tuple(int)) – The number of points to create in line
  • mins (tuple(int, int, int)) – The starting value for each dimension.
  • maxs (tuple(int, int, int)) – The ending value for each dimension.
  • rz_switch (tuple(int, int, int)) – Determines true or false (1 or 0) for using ratio zoning values.
createpts_brick(crd, npts, mins, maxs, ctr=(1, 1, 1), rz_switch=(1, 1, 1), rz_vls=(1, 1, 1))

Create and Connect Points

Creates a grid of points in the mesh object and connects them.

Parameters:
  • crd (str) – Coordinate type of either ‘xyz’ (cartesian coordinates), ‘rtz’ (cylindrical coordinates), or ‘rtp’ (spherical coordinates).
  • npts (tuple(int, int, int)) – The number of points to create in each dimension.
  • mins (tuple(int, int, int)) – The starting value for each dimension.
  • maxs (tuple(int, int, int)) – The ending value for each dimension.
  • ctr (tuple(int, int, int)) – Defines the center of each cell. For 0, points are placed in the middle of each cell. For 1, points are placed at the edge of each cell.
  • rz_switch (tuple(int, int, int)) – Determines true or false (1 or 0) for using ratio zmoning values.
  • rz_vls (tuple(int, int, int)) – Ratio zoning values. Each point will be multiplied by a scale of the value for that dimension.
createpts_brick_rtp(npts, mins, maxs, ctr=(1, 1, 1), rz_switch=(1, 1, 1), rz_vls=(1, 1, 1))

Create and connect spherical coordinates.

createpts_brick_rtz(npts, mins, maxs, ctr=(1, 1, 1), rz_switch=(1, 1, 1), rz_vls=(1, 1, 1))

Create and connect cylindrical coordinate points.

createpts_brick_xyz(npts, mins, maxs, ctr=(1, 1, 1), rz_switch=(1, 1, 1), rz_vls=(1, 1, 1))

Create and connect Cartesian coordinate points.

createpts_dxyz(dxyz, mins, maxs, clip='under', hard_bound='min', rz_switch=(1, 1, 1), rz_value=(1, 1, 1), connect=True)

Create and Connect Points to create an orthogonal hexahedral mesh. The vertex spacing is based on dxyz and the mins and maxs specified. mins (default, see hard_bound option) or maxs will be adhered to, while maxs (default) or mins will be modified based on the clip option to be truncated at the nearest value ‘under’ (default) or ‘over’ the range maxs-mins. clip and hard_bound options can be mixed by specifying tuples (see description below).

Parameters:
  • dxyz (tuple(float,float,float)) – The spacing between points in x, y, and z directions
  • mins (tuple(float,float,float)) – The starting value for each dimension.
  • maxs (tuple(float,float,float)) – The ending value for each dimension.
  • clip (string or tuple(string,string,string)) – How to handle bounds if range does not divide by dxyz, either clip ‘under’ or ‘over’ range
  • hard_bound (string or tuple(string,string,string)) – Whether to use the “min” or “max” as the hard constraint on dimension
  • rz_switch (tuple(int, int, int)) – Determines true or false (1 or 0) for using ratio zoning values.
  • connect (boolean) – Whether or not to connect points
createpts_line(npts, mins, maxs, rz_switch=(1, 1, 1))

Create and Connect Points in a line

Parameters:
  • npts (int) – The number of points to create in line
  • mins (tuple(int, int, int)) – The starting value for each dimension.
  • maxs (tuple(int, int, int)) – The ending value for each dimension.
  • rz_switch (tuple(int, int, int)) – Determines true or false (1 or 0) for using ratio zoning values.
delatt(attnames, force=True)

Delete a list of attributes

Parameters:
  • attnames (str or lst(str)) – Attribute names to delete
  • force (bool) – If true, delete even if the attribute permanent persistance
dump_ats_xml(filename, meshfilename, matnames={}, facenames={})

Write ats style xml file with regions :param filename: Name of xml to write :type filename: string :param meshfilename: Name of exodus file to use in xml :type meshfilename: string :param matnames: Dictionary of region names keyed by exodus material number :type matnames: dict :param facenames: Dictionary of faceset names keyed by exodus faceset number :type facenames: dict

dump_avs2(filename, points=True, elements=True, node_attr=True, element_attr=True)

Dump avs file

Parameters:
  • filename (str) – Name of avs file
  • points (bool) – Output point coordinates
  • elements (bool) – Output connectivity
  • node_attr (bool) – Output node attributes
  • element_attr (bool) – Output element attributes
dump_exo(filename, psets=False, eltsets=False, facesets=[])

Dump exo file

Parameters:
  • filename (str) – Name of exo file
  • psets (bool) – Boolean indicating that exodus will only include psets
  • eltsets (bool) – Boolean indicating that exodus will only include element sets
  • facesets (lst(FaceSet)) – Array of FaceSet objects
Example:
>>> from pylagrit import PyLaGriT
>>> l = PyLaGriT()
>>> m = l.create()
>>> m.createpts_xyz((3,3,3),(0.,0.,0.),(1.,1.,1.),rz_switch=[1,1,1],connect=True)
>>> m.status ()
>>> m.status (brief=True)
>>> fs = m.create_boundary_facesets(base_name='faceset_bounds')
>>> m.dump_exo('cube.exo',facesets=fs.values())
dump_pset(filerootname, zonetype='zone', pset=[])

Dump zone file of psets :arg filerootname: rootname of files to create, pset name will be added to name :type filerootname: string :arg zonetype: Type of zone file to dump, ‘zone’ or ‘zonn’ :type zonetype: string :arg pset: list of psets to dump, all psets dumped if empty list :type pset: list[strings]

eltset_bool(eset_list, boolstr='union', name=None)

Create element set from boolean operation of set of element sets

Parameters:
  • eset_list (lst(PyLaGriT element set)) – List of elements to perform boolean operation on
  • boolstr (str) – type of boolean operation to perform on element sets, one of [union,inter,not]
  • name (str) – The name to be assigned to the EltSet within LaGriT
Returns:

PyLaGriT element set object

eltset_object(mo, name=None)

Create element set from the intersecting elements with another mesh object

extrude(offset, offset_type='const', return_type='volume', direction=[], name=None)

Extrude mesh object to new mesh object This command takes the current mesh object (topologically 1d or 2d mesh (a line, a set of line segments, or a planar or non-planar surface)) and extrudes it into three dimensions along either the normal to the curve or surface (default), along a user defined vector, or to a set of points that the user has specified. If the extrusion was along the normal of the surface or along a user defined vector, the command can optionally find the external surface of the volume created and return that to the user. Refer to http://lagrit.lanl.gov/docs/commands/extrude.html for more details on arguments.

Parameters:
  • name (str) – Name to use within lagrit for the created mesh object
  • offset (float) – Distance to extrude
  • offset_type (str) – either const or min (interp will be handled in the PSET class in the future)
  • return_type (str) – either volume for entire mesh or bubble for just the external surface
  • direction (lst[float,float,float]) – Direction to extrude in, defaults to normal of the object
Returns:

mesh object

grid2grid(ioption, name=None)

Convert a mesh with one element type to a mesh with another

Parameters:
  • ioption – type of conversion: quadtotri2 quad to 2 triangles, no new points. prismtotet3 prism to 3 tets, no new points. quadtotri4 quad to 4 triangles, with one new point. pyrtotet4 pyramid to 4 tets, with one new point. hextotet5 hex to 5 tets, no new points. hextotet6 hex to 6 tets, no new points. prismtotet14 prism to 14 tets, four new points (1 + 3 faces). prismtotet18 prism to 18 tets, six new points (1 + 5 faces). hextotet24 hex to 24 tets, seven new points (1 + 6 faces). tree_to_fe quadtree or octree grid to grid with no parent-type elements.
  • name (str) – Internal Lagrit name of new mesh object, automatically created if None

Returns MO object

grid2grid_hextotet24(name=None)

Hex to 24 tets, seven new points (1 + 6 faces) :arg name: Internal Lagrit name of new mesh object, automatically created if None :type name: str

Returns MO object

grid2grid_hextotet5(name=None)

Hex to 5 tets, no new points :arg name: Internal Lagrit name of new mesh object, automatically created if None :type name: str

Returns MO object

grid2grid_hextotet6(name=None)

Hex to 6 tets, no new points :arg name: Internal Lagrit name of new mesh object, automatically created if None :type name: str

Returns MO object

grid2grid_prismtotet14(name=None)

Prism to 14 tets, four new points (1 + 3 faces) :arg name: Internal Lagrit name of new mesh object, automatically created if None :type name: str

Returns MO object

grid2grid_prismtotet18(name=None)

Prism to 18 tets, four new points (1 + 3 faces) :arg name: Internal Lagrit name of new mesh object, automatically created if None :type name: str

Returns MO object

grid2grid_prismtotet3(name=None)

Quad to 2 triangles, no new points. Prism to 3 tets, no new points. :arg name: Internal Lagrit name of new mesh object, automatically created if None :type name: str

Returns MO object

grid2grid_pyrtotet4(name=None)

Pyramid to 4 tets, with one new point :arg name: Internal Lagrit name of new mesh object, automatically created if None :type name: str

Returns MO object

grid2grid_quadtotri2(name=None)

Quad to 2 triangles, no new points. :arg name: Internal Lagrit name of new mesh object, automatically created if None :type name: str

Returns MO object

grid2grid_quadtotri4(name=None)

Quad to 4 triangles, with one new point :arg name: Internal Lagrit name of new mesh object, automatically created if None :type name: str

Returns MO object

grid2grid_tree_to_fe(name=None)

Quadtree or octree grid to grid with no parent-type elements. :arg name: Internal Lagrit name of new mesh object, automatically created if None :type name: str

Returns MO object

interpolate(method, attsink, cmosrc, attsrc, stride=[1, 0, 0], tie_option=None, flag_option=None, keep_option=None, interp_function=None)

Interpolate values from attribute attsrc from mesh object cmosrc to current mesh object

intersect_elements(mo, attr_name='attr00')

This command takes two meshes and creates an element-based attribute in mesh1 that contains the number of elements in mesh2 that intersected the respective element in mesh1. We define intersection as two elements sharing any common point.

Parameters:
  • mo (PyLaGriT mesh object) – Mesh object to intersect with current mesh object to determine where to refine
  • attr_name (str) – Name to give created attribute
Returns:

attr_name

pset_attribute(attribute, value, comparison='eq', stride=(1, 0, 0), name=None)

Define PSet by attribute

Parameters:
  • attribute (str) – Nodes defined by attribute ID.
  • value (integer) – attribute ID value.
  • comparison (can use default without specifiy anything, or list[lt|le|gt|ge|eq|ne]) – attribute comparison, default is eq.
  • stride (list[int, int, int]) – Nodes defined by ifirst, ilast, and istride.
  • name (str) – The name to be assigned to the PSet created.

Returns: PSet object

pset_bool(pset_list, boolean='union', name=None)

Return PSet from boolean operation on list of psets

Defines and returns a PSet from points that are not inside the PSet, ps.

pset_geom(mins, maxs, ctr=(0, 0, 0), geom='xyz', stride=(1, 0, 0), name=None)

Define PSet by Geometry

Selects points from geomoetry specified by string geom and returns a PSet.

Parameters:
  • mins (tuple(int, int, int)) – Coordinate of one of the shape’s defining points. xyz (Cartesian): (x1, y1, z1); rtz (Cylindrical): (radius1, theta1, z1); rtp (Spherical): (radius1, theta1, phi1);
  • maxs (tuple(int, int, int)) – Coordinate of one of the shape’s defining points. xyz (Cartesian): (x2, y2, z2); rtz (Cylindrical): (radius2, theta2, z2); rtp (Spherical): (radius2, theta2, phi2);
  • ctr (tuple(int, int, int)) – Coordinate of the relative center.
  • geom (str) – Type of geometric shape: ‘xyz’ (spherical), ‘rtz’ (cylindrical), ‘rtp’ (spherical)
  • stride (list[int, int, int]) – Nodes defined by ifirst, ilast, and istride.
  • name (str) – The name to be assigned to the PSet created.

Returns: PSet object

pset_geom_rtp(mins, maxs, ctr=(0, 0, 0), stride=(1, 0, 0), name=None)

Forms a pset of nodes within the sphere, sperical shell or sperical section given by radius1 to radius2, and angles theta1 to theta2 (0 - 180) and angles phi1 to phi2 (0 - 360). Refer to http://lagrit.lanl.gov/docs/conventions.html for an explanation of angles

Parameters:
  • mins (tuple(int, int, int)) – Defines radius1, theta1, and phi1.
  • maxs (tuple(int, int, int)) – Defines radius2, theta2, and phi2.
  • stride (list[int, int, int]) – Nodes defined by ifirst, ilast, and istride.
  • name (str) – The name to be assigned to the PSet created.
  • ctr (tuple(int, int, int)) – Coordinate of the relative center.
  • stride – Nodes defined by ifirst, ilast, and istride.
  • name – The name to be assigned to the PSet created.

Returns: PSet object

pset_geom_rtz(mins, maxs, ctr=(0, 0, 0), stride=(1, 0, 0), name=None)

Forms a pset of nodes within the cylinder or cylindrical shell section given by radius1 to radius2, and angles theta1 to theta2 and height z1 to z2. Refer to http://lagrit.lanl.gov/docs/conventions.html for an explanation of angles

Parameters:
  • mins (tuple(int, int, int)) – Defines radius1, theta1, and z1.
  • maxs (tuple(int, int, int)) – Defines radius2, theta2, and z2.
  • stride (list[int, int, int]) – Nodes defined by ifirst, ilast, and istride.
  • name (str) – The name to be assigned to the PSet created.
  • ctr (tuple(int, int, int)) – Coordinate of the relative center.
  • stride – Nodes defined by ifirst, ilast, and istride.
  • name – The name to be assigned to the PSet created.

Returns: PSet object

pset_geom_xyz(mins, maxs, ctr=(0, 0, 0), stride=(1, 0, 0), name=None)

Define PSet by Tetrahedral Geometry

Selects points from a Tetrahedral region.

Parameters:
  • mins (tuple(int, int, int)) – Coordinate point of 1 of the tetrahedral’s corners.
  • maxs (tuple(int, int, int)) – Coordinate point of 1 of the tetrahedral’s corners.
  • ctr (tuple(int, int, int)) – Coordinate of the relative center.
  • stride (list[int, int, int]) – Nodes defined by ifirst, ilast, and istride.
  • name (str) – The name to be assigned to the PSet created.

Returns: PSet object

pset_region(region, stride=(1, 0, 0), name=None)

Define PSet by region

Parameters:
  • region – region to create pset
  • stride (list[int, int, int]) – Nodes defined by ifirst, ilast, and istride.
  • name (str) – The name to be assigned to the PSet created.

Returns: PSet object

pset_surface(surface, stride=(1, 0, 0), name=None)

Define PSet by surface

Parameters:
  • surface – surface to create pset
  • stride (list[int, int, int]) – Nodes defined by ifirst, ilast, and istride.
  • name (str) – The name to be assigned to the PSet created.

Returns: PSet object

refine_to_object(mo, level=None, imt=None, prd_choice=None)

Refine mesh at locations that intersect another mesh object

Parameters:
  • mo (PyLaGriT mesh object) – Mesh object to intersect with current mesh object to determine where to refine
  • level (int) – max level of refinement
  • imt (int) – Value to assign to imt (LaGriT material type attribute)
  • prd_choice (int) – directions of refinement
resetpts_itp()

set node type from connectivity of mesh

rmpoint_compress(filter_bool=False, resetpts_itp=True)

remove all marked nodes and correct the itet array

Parameters:resetpts_itp (bool) – set node type from connectivity of mesh
rotateln(coord1, coord2, theta, center=[0, 0, 0], copy=False, stride=(1, 0, 0))

Rotates a point distribution (specified by ifirst,ilast,istride) about a line. The copy option allows the user to make a copy of the original points as well as the rotated points, while copy=False just keeps the rotated points themselves. The line of rotation defined by coord1 and coord2 needs to be defined such that the endpoints extend beyond the point distribution being rotated. theta (in degrees) is the angle of rotation whose positive direction is determined by the right-hand-rule, that is, if the thumb of your right hand points in the direction of the line (1 to 2), then your fingers will curl in the direction of rotation. center is the point where the line can be shifted to before rotation takes place. If the copy option is chosen, the new points will have only coordinate values (xic, yic, zic); no values will be set for any other mesh object attribute for these points. Note: The end points of the line segment must extend well beyond the point set being rotated.

Example 1:
>>> from pylagrit import PyLaGriT
>>> import numpy
>>> x = numpy.arange(0,10.1,1)
>>> y = x
>>> z = [0,1]
>>> lg = PyLaGriT()
>>> mqua = lg.gridder(x,y,z,elem_type='hex',connect=True)
>>> mqua.rotateln([mqua.xmin-0.1,0,0],[mqua.xmax+0.1,0,0],25)
>>> mqua.dump_exo('rotated.exo')
>>> mqua.dump_ats_xml('rotated.xml','rotated.exo')
>>> mqua.paraview()
Example 2:
>>> from pylagrit import PyLaGriT
>>> import numpy
>>> x = numpy.arange(0,10.1,1)
>>> y = [0,1]
>>> #z = [0,1]
>>> lg = PyLaGriT()
>>> layer = lg.gridder(x=x,y=y,elem_type='quad',connect=True)
>>> layer.rotateln([0,layer.ymin-0.10,0],[0,layer.ymax+0.1,0],25)
>>> layer.dump('tmp_lay_top.inp')
>>> # Layer depths?
>>> #           1   2   3    4    5    6    7   8    9   10
>>> layers = [ .1, 1.]
>>> addnum = [  4, 2]
>>> #matnum = [2]*len(layers)
>>> matnum = [2, 1]
>>> layer_interfaces = numpy.cumsum(layers)
>>> mtop = layer.copy()
>>> stack_files = ['tmp_lay_top.inp 1,9']
>>> #stack_files.append('tmp_lay_peat_bot.inp 1,33')
>>> i = 1
>>> for li,m,a in zip(layer_interfaces,matnum,addnum):
>>>     layer.math('sub',li,'zic',cmosrc=mtop)
>>>     stack_files.append('tmp_lay'+str(i)+'.inp '+str(int(m))+', '+str(a))
>>>     layer.dump('tmp_lay'+str(i)+'.inp')
>>>     i += 1
>>> layer.math('sub',2,'zic',cmosrc=mtop)
>>> #layer.setatt('zic',-2.)
>>> layer.dump('tmp_lay_bot.inp')
>>> stack_files.append('tmp_lay_bot.inp 2')
>>> stack_files.reverse()
>>> # Create stacked layer mesh and fill
>>> stack = lg.create()
>>> stack.stack_layers('avs',stack_files,flip_opt=True)
>>> stack_hex = stack.stack_fill()
>>> stack_hex.dump_exo('rotated.exo')
>>> stack_hex.dump_ats_xml('rotated.xml','rotated.exo')
>>> stack_hex.paraview()
subset(mins, maxs, geom='xyz')

Return Mesh Object Subset

Creates a new mesh object that contains only a geometric subset defined by mins and maxs.

Parameters:
  • mins – Coordinate of one of the shape’s defining points. xyz (Cartesian): (x1, y1, z1); rtz (Cylindrical): (radius1, theta1, z1); rtp (Spherical): (radius1, theta1, phi1);
  • maxs (tuple(int, int, int)) – Coordinate of one of the shape’s defining points. xyz (Cartesian): (x2, y2, z2); rtz (Cylindrical): (radius2, theta2, z2); rtp (Spherical): (radius2, theta2, phi2);
  • geom (str) – Type of geometric shape: ‘xyz’ (spherical), ‘rtz’ (cylindrical), ‘rtp’ (spherical)
Typep mins:

tuple(int, int, int)

Returns: MO object

Example:
>>> #To use pylagrit, import the module.
>>> import pylagrit
>>> #Start the lagrit session.
>>> lg = pylagrit.PyLaGriT()
>>> #Create a mesh object.
>>> mo = lg.create()
>>> mo.createpts_brick_xyz((5,5,5), (0,0,0), (5,5,5))
>>> #Take the subset from (3,3,3)
>>> mo.subset((3,3,3),(5,5,5))
subset_rtp(mins, maxs)

Return Spherical MO Subset

Creates a new mesh object that contains only a spherical subset defined by mins and maxs.

Parameters:
  • mins (tuple(int, int, int)) – Defines radius1, theta1, and phi1.
  • maxs (tuple(int, int, int)) – Defines radius2, theta2, and phi2.

Returns: MO object

subset_rtz(mins, maxs)

Return Cylindrical MO Subset

Creates a new mesh object that contains only a cylindrical subset defined by mins and maxs.

Parameters:
  • mins (tuple(int, int, int)) – Defines radius1, theta1, and z1.
  • maxs (tuple(int, int, int)) – Defines radius2, theta2, and z2.

Returns: MO object

subset_xyz(mins, maxs)

Return Tetrehedral MO Subset

Creates a new mesh object that contains only a tetrehedral subset defined by mins and maxs.

Parameters:
  • mins (tuple(int, int, int)) – Coordinate point of 1 of the tetrahedral’s corners.
  • maxs (tuple(int, int, int)) – Coordinate point of 1 of the tetrahedral’s corners.

Returns: MO object

trans(xold, xnew, stride=(1, 0, 0))

Translate mesh according to old coordinates “xold” to new coordinates “xnew”

Parameters:
  • xold (tuple(float,float,float)) – old position
  • xnew (tuple(float,float,float)) – new position
  • stride (tuple(int,int,int)) – tuple of (first, last, stride) of points
tri_mesh_output_prep()

Prepare tri mesh for output, remove dudded points, ensure delaunay volumes, etc. Combination of lagrit commands: filter/1 0 0 rmpoint/compress recon/1 resetpts/itp

triangulate(order='clockwise')

triangulate will take an ordered set of nodes in the current 2d mesh object that define a perimeter of a polygon and create a trangulation of the polygon. The nodes are assumed to lie in the xy plane; the z coordinate is ignored. No checks are performed to verify that the nodes define a legal perimeter (i.e. that segments of the perimeter do not cross). The code will connect the last node to the first node to complete the perimeter.

This code support triangulation of self-intersecting polygons (polygon with holes), assuming that the order of the nodes are correct. Moreover the connectivity of the polyline must also be defined correctly. No checks are made.

One disadvantage of the algorithm for triangulating self-intersecting polygons is that it does not always work. For example, if the holes have complicated shapes, with many concave vertices, the code might fail. In this case, the user may try to rotate the order of the nodes:
NODE_ID:
1 -> 2 2 -> 3 ... N -> 1
param order:direction of point ordering
type order:string
Example:
>>> from pylagrit import PyLaGriT
>>> # Create pylagrit object
>>> lg = PyLaGriT()
>>> # Define polygon points in clockwise direction
>>> # and create tri mesh object
>>> coords = [[0.0, 0.0, 0.0],
>>>           [0.0, 1000.0, 0.0],
>>>           [2200.0, 200.0, 0.0],
>>>           [2200.0, 0.0, 0.0]]
>>> motri = lg.tri_mo_from_polyline(coords)
>>> # Triangulate polygon
>>> motri.triangulate()
>>> motri.setatt('imt',1)
>>> motri.setatt('itetclr',1)
>>> # refine mesh with successively smaller edge length constraints
>>> edge_length = [1000,500,250,125,75,40,20,15]
>>> for i,l in enumerate(edge_length):
>>>     motri.resetpts_itp()
>>>     motri.refine(refine_option='rivara',refine_type='edge',values=[l],inclusive_flag='inclusive')
>>>     motri.smooth()
>>>     motri.recon(0)
>>> # provide additional smoothing after the last refine
>>> for i in range(5):
>>>     motri.smooth()
>>>     motri.recon(0)
>>> # create delaunay mesh and clean up
>>> motri.tri_mesh_output_prep()
>>> # dump fehm files
>>> motri.dump_fehm('nk_mesh00')
>>> # view results
>>> motri.paraview() 
upscale(method, attsink, cmosrc, attsrc=None, stride=(1, 0, 0), boundary_choice=None, keepatt=False, set_id=False)

The upscale command is used to interpolate attribute values from nodes of a fine source mesh to node attributes of a coarse sink mesh. The subroutine finds nodes of the fine source mesh within the Voronoi cell of every node in the coarser sink mesh. Nodes on cell boundaries are assigned to two or more sink nodes. Then the attributes of all the source nodes within a source node’s cell are upscaled into a single value based on the chosen method. Mesh elements and connectivity are ignored and only node values are used to upscale values on to the sink mesh nodes.

Parameters:
  • method (str) – Type of upscaling: sum, min, max, and averages ariave, harave, geoave
  • attsink (str) – attribute sink
  • cmosrc (PyLaGriT Mesh Object) – PyLaGriT mesh object source
  • attsrc (str) – attribute src, defaults to name of attsink
  • stride (tuple(int)) – tuple of (first, last, stride) of points
  • boundary_choice (str) – method of choice when source nodes are found on the boundary of multiple Voronoi volumes of sink nodes: single, divide, or multiple
upscale_ariave(attsink, cmosrc, attsrc=None, stride=(1, 0, 0), boundary_choice=None, keepatt=False, set_id=False)

Upscale using arithmetic average of cmosrc points within Voronoi volumes of current mesh

Parameters:
  • attsink (str) – attribute sink
  • cmosrc (PyLaGriT Mesh Object) – PyLaGriT mesh object source
  • attsrc (str) – attribute src
  • stride (tuple(int)) – tuple of (first, last, stride) of points
  • boundary_choice (str) – method of choice when source nodes are found on the boundary of multiple Voronoi volumes of sink nodes: single, divide, or multiple
upscale_geoave(attsink, cmosrc, attsrc=None, stride=(1, 0, 0), boundary_choice=None, keepatt=False, set_id=False)

Upscale using geometric average of cmosrc points within Voronoi volumes of current mesh

Parameters:
  • attsink (str) – attribute sink
  • cmosrc (PyLaGriT Mesh Object) – PyLaGriT mesh object source
  • attsrc (str) – attribute src
  • stride (tuple(int)) – tuple of (first, last, stride) of points
  • boundary_choice (str) – method of choice when source nodes are found on the boundary of multiple Voronoi volumes of sink nodes: single, divide, or multiple
upscale_harave(attsink, cmosrc, attsrc=None, stride=(1, 0, 0), boundary_choice=None, keepatt=False, set_id=False)

Upscale using harmonic average of cmosrc points within Voronoi volumes of current mesh

Parameters:
  • attsink (str) – attribute sink
  • cmosrc (PyLaGriT Mesh Object) – PyLaGriT mesh object source
  • attsrc (str) – attribute src
  • stride (tuple(int)) – tuple of (first, last, stride) of points
  • boundary_choice (str) – method of choice when source nodes are found on the boundary of multiple Voronoi volumes of sink nodes: single, divide, or multiple
upscale_max(attsink, cmosrc, attsrc=None, stride=(1, 0, 0), boundary_choice=None, keepatt=False, set_id=False)

Upscale using maximum of cmosrc points within Voronoi volumes of current mesh

Parameters:
  • attsink (str) – attribute sink
  • cmosrc (PyLaGriT Mesh Object) – PyLaGriT mesh object source
  • attsrc (str) – attribute src
  • stride (tuple(int)) – tuple of (first, last, stride) of points
  • boundary_choice (str) – method of choice when source nodes are found on the boundary of multiple Voronoi volumes of sink nodes: single, divide, or multiple
upscale_min(attsink, cmosrc, attsrc=None, stride=(1, 0, 0), boundary_choice=None, keepatt=False, set_id=False)

Upscale using minimum of cmosrc points within Voronoi volumes of current mesh

Parameters:
  • attsink (str) – attribute sink
  • cmosrc (PyLaGriT Mesh Object) – PyLaGriT mesh object source
  • attsrc (str) – attribute src
  • stride (tuple(int)) – tuple of (first, last, stride) of points
  • boundary_choice (str) – method of choice when source nodes are found on the boundary of multiple Voronoi volumes of sink nodes: single, divide, or multiple
upscale_sum(attsink, cmosrc, attsrc=None, stride=(1, 0, 0), boundary_choice=None, keepatt=False, set_id=False)

Upscale using sum of cmosrc points within Voronoi volumes of current mesh

Parameters:
  • attsink (str) – attribute sink
  • cmosrc (PyLaGriT Mesh Object) – PyLaGriT mesh object source
  • attsrc (str) – attribute src
  • stride (tuple(int)) – tuple of (first, last, stride) of points
  • boundary_choice (str) – method of choice when source nodes are found on the boundary of multiple Voronoi volumes of sink nodes: single, divide, or multiple