ogs5py.fileclasses.MSH
- class ogs5py.fileclasses.MSH(mesh_list=None, **OGS_Config)[source]
Bases:
MSHsgl
Class for a multi layer mesh file that contains multiple ‘#FEM_MSH’ Blocks.
- Parameters:
mesh_list (list of dict or None, optional) – each dictionary contains one ‘#FEM_MSH’ block of the mesh file with with the following information (sorted by keys):
- mesh_datadict
dictionary containing information about
AXISYMMETRY (bool)
CROSS_SECTION (bool)
PCS_TYPE (str)
GEO_TYPE (str)
GEO_NAME (str)
LAYER (int)
- nodesndarray
Array with all node postions
- elementsdict
contains nodelists for elements sorted by element types
- material_iddict
contains material ids for each element sorted by element types
- element_iddict
contains element ids for each element sorted by element types
task_root (str, optional) – Path to the destiny model folder. Default: cwd+”ogs5model”
task_id (str, optional) – Name for the ogs task. Default: “model”
- Attributes:
AXISYMMETRY
bool: AXISYMMETRY attribute.
CROSS_SECTION
bool: CROSS_SECTION attribute.
ELEMENTS
Get and set the ELEMENTS of the mesh.
ELEMENT_ID
Get and set the ELEMENT_IDs of the mesh.
ELEMENT_NO
int: number of ELEMENTS.
ELEMENT_TYPES
set
: ELEMENT types of the mesh.GEO_NAME
str: GEO_NAME.
GEO_TYPE
str: GEO_TYPE.
LAYER
int: LAYER.
MATERIAL_ID
Get and set the MATERIAL_IDs of the mesh.
MATERIAL_ID_flat
Get flat version of the MATERIAL_IDs of the mesh.
NODES
ndarray: (n,3) NODES of the mesh by its xyz-coordinates.
NODE_NO
int: number of NODES.
PCS_TYPE
str: PCS_TYPE.
block
int
: The actual block to access in the file.block_no
int
: The number of blocks in the file.center
Get the mesh center.
centroids
Get the centroids of the mesh.
centroids_flat
Get flat version of the centroids of the mesh.
file_name
str
: base name of the file with extension.file_path
str
: save path of the file.force_writing
bool
: state if the file is written even if empty.is_empty
State if file is empty.
name
str
: name of the file without extension.node_centroids
Get the node centroids of the mesh.
node_centroids_flat
Get flat version of the node centroids of the mesh.
volumes
Get the volumes of the mesh-elements.
volumes_flat
Get flat version of the volumes of the mesh-elements.
Methods
__call__
()Get a copy of the underlying dictionary.
add_copy_link
(path[, symlink])Add a link to copy a file instead of writing.
check
([verbose])Check if the mesh is valid.
combine_mesh
(ext_mesh, **kwargs)Combine this mesh with an external mesh.
Remove a former given link to an external file.
export_mesh
(filepath[, verbose])Export the mesh to an unstructured mesh in diffrent file-formats.
generate
([generator])Use a mesh-generator from the generator module.
Get the OGS file class name.
import_mesh
(filepath, **kwargs)Import an external unstructured mesh from diffrent file-formats.
load
(filepath, **kwargs)Load an OGS5 mesh from file.
read_file
(path[, encoding, verbose])Load an OGS5 mesh from file.
remove_dim
(remove)Remove elements by given dimensions from a mesh.
reset
()Delete every content.
rotate
(angle[, rotation_axis, rotation_point])Rotate a given mesh around a given rotation axis with a given angle.
save
(path, **kwargs)Save the mesh to an OGS5 mesh file.
set_dict
(mesh_dict)Set an mesh as returned by tools methods or generators.
set_material_id
([material_id, element_id, ...])Set material IDs by the corresponding element IDs.
shift
(vector)Shift a given mesh with a given vector.
show
([show_cell_data, show_material_id, ...])Display the mesh colored by its material ID.
swap_axis
([axis1, axis2])Swap axis of the coordinate system.
transform
(xyz_func, **kwargs)Transform a given mesh with a given function "xyz_func".
Write the actual OGS input file to the given folder.
- __call__()
Get a copy of the underlying dictionary.
- add_copy_link(path, symlink=False)
Add a link to copy a file instead of writing.
Instead of writing a file, you can give a path to an existing file, that will be copied/linked to the target folder.
- Parameters:
path (str) – path to the existing file that should be copied
symlink (bool, optional) – on UNIX systems it is possible to use a symbolic link to save time if the file is big. Default: False
- check(verbose=True)
Check if the mesh is valid.
Checked in the sence, that the contained data is consistent. Checks for correct element definitions or Node duplicates are not carried out.
- Parameters:
verbose (bool, optional) – Print information for the executed checks. Default: True
- Returns:
result – Validity of the given mesh.
- Return type:
- combine_mesh(ext_mesh, **kwargs)
Combine this mesh with an external mesh.
The node list will be updated to eliminate duplicates. Element intersections are not checked. kwargs will be forwarded to “tools.combine”
- Parameters:
ext_mesh (mesh, dict or file) – This is the mesh that should be added to the existing one.
decimals (int, optional) – Number of decimal places to round the nodes to (default: 3). This will not round the output, it is just for comparison of the node vectors.
fast (bool, optional) – If fast is True, the vector comparison is executed by a decimal comparison. If fast is False, all pairwise distances are calculated. Default: False
- del_copy_link()
Remove a former given link to an external file.
- export_mesh(filepath, verbose=False, **kwargs)
Export the mesh to an unstructured mesh in diffrent file-formats.
kwargs will be forwarded to “tools.export_mesh”
- Parameters:
filepath (string) – path to the file to export
file_format (str, optional) – Here you can specify the fileformat. If ‘None’ it will be determined by file extension. Default: None
export_material_id (bool, optional) – Here you can specify if the material_id should be exported. Default: True
export_element_id (bool, optional) – Here you can specify if the element_id should be exported. Default: True
cell_data_by_id (ndarray or dict, optional) – Here you can specify additional element data sorted by their IDs. It can be a dictionary with data-name as key and the ndarray as value. Default: None
point_data (ndarray or dict, optional) – Here you can specify additional point data sorted by their IDs. It can be a dictionary with data-name as key and the ndarray as value. Default: None
field_data (ndarray or dict, optional) – Here you can specify additional field data of the mesh. It can be a dictionary with data-name as key and the ndarray as value. Default: None
Notes
This routine calls the ‘write’ function from the meshio package and converts the input (see here: https://github.com/nschloe/meshio)
- generate(generator='rectangular', **kwargs)
Use a mesh-generator from the generator module.
See:
ogs5py.fileclasses.msh.generator
- Parameters:
generator (str) – set the generator from the generator module
**kwargs – kwargs will be forwarded to the generator in use
Notes
The following generators are available:
rectangular
([dim, mesh_origin, element_no, ...])Generate a rectangular grid in 2D or 3D.
radial
([dim, mesh_origin, angles, rad, z_arr])Generate a radial grid in 2D or 3D.
grid_adapter2D
([out_dim, in_dim, out_res, ...])Generate a grid adapter.
grid_adapter3D
([out_dim, in_dim, z_dim, ...])Generate a grid adapter.
block_adapter3D
([xy_dim, z_dim, in_res])Generate a block adapter.
gmsh
(geo_object[, import_dim])Generate mesh from gmsh code or gmsh .geo file.
- get_file_type()
Get the OGS file class name.
- import_mesh(filepath, **kwargs)
Import an external unstructured mesh from diffrent file-formats.
kwargs will be forwarded to “tools.import_mesh”
- Parameters:
filepath (string or meshio.Mesh instance) – mesh to import
file_format (str, optional) – Here you can specify the fileformat. If ‘None’ it will be determined by file extension. Default: None
ignore_unknown (bool, optional) – Unknown data in the file will be ignored. Default: False
import_dim (iterable of int, optional) – State which elements should be imported by dimensionality. Can be used to sort out unneeded elements for example from gmsh. Default: (1, 2, 3)
Notes
This routine calls the ‘read’ function from the meshio package and converts the output (see here: https://github.com/nschloe/meshio) If there is any “vertex” (0D element) in the element data, it will be removed.
- load(filepath, **kwargs)
Load an OGS5 mesh from file.
kwargs will be forwarded to “tools.load_ogs5msh”
- Parameters:
filepath (string) – path to the ‘*.msh’ OGS5 mesh file to load
verbose (bool, optional) – Print information of the reading process. Default: True
ignore_unknown (bool, optional) – Unknown data in the file will be ignored. Default: False
max_node_no (int, optional) – If you know the maximal node number per elements in the mesh file, you can optimise the reading a bit. By default the algorithm will assume hexahedrons as ‘largest’ elements in the mesh. Default: 8
encoding (str or None, optional) – encoding of the given file. If
None
is given, the system standard is used. Default:None
Notes
The $AREA keyword within the Nodes definition is NOT supported and will violate the read data if present.
- read_file(path, encoding=None, verbose=False)
Load an OGS5 mesh from file.
- Parameters:
path (str) – path to the ‘*.msh’ OGS5 mesh file to load
encoding (str or None, optional) – encoding of the given file. If
None
is given, the system standard is used. Default:None
verbose (bool, optional) – Print information of the reading process. Default: True
- remove_dim(remove)
Remove elements by given dimensions from a mesh.
- Parameters:
remove (iterable of int or single int) – State which elements should be removed by dimensionality (1, 2, 3).
- reset()
Delete every content.
- rotate(angle, rotation_axis=(0.0, 0.0, 1.0), rotation_point=(0.0, 0.0, 0.0))
Rotate a given mesh around a given rotation axis with a given angle.
- Parameters:
angle (float) – rotation angle given in radial length
rotation_axis (array_like, optional) – Array containing the vector for ratation axis. Default: (0,0,1)
rotation_point (array_like, optional) – Vector of the ratation base point. Default:(0,0,0)
- save(path, **kwargs)
Save the mesh to an OGS5 mesh file.
kwargs will be forwarded to “tools.save_ogs5msh”
- Parameters:
path (string) – path to the ‘*.msh’ OGS5 mesh file to save
verbose (bool, optional) – Print information of the writing process. Default: True
- set_dict(mesh_dict)
Set an mesh as returned by tools methods or generators.
Mesh will be checked for validity.
- Parameters:
mesh_dict (dict or None, optional) – Contains one ‘#FEM_MSH’ block of an OGS5 mesh file with the following information (sorted by keys):
- mesh_datadict
dictionary containing information about
AXISYMMETRY (bool)
CROSS_SECTION (bool)
PCS_TYPE (str)
GEO_TYPE (str)
GEO_NAME (str)
LAYER (int)
- nodesndarray
Array with all node postions
- elementsdict
contains nodelists for elements sorted by element types
- material_iddict
contains material ids for each element sorted by types
- element_iddict
contains element ids for each element sorted by types
- set_material_id(material_id=0, element_id=None, element_mask=None)
Set material IDs by the corresponding element IDs.
- Parameters:
material_id (
int
or ndarray, optional) – The new material IDs. Either one value or an array. Default: 0element_id (ndarray or None, optional) – The corresponding element IDs, where to set the material IDs. If None, all elements are assumed and the material IDs are added by their index. Default: None
element_mask (ndarray or None, optional) – Instead of the element IDs, one can specify a mask to select the element IDs. Default: None
- shift(vector)
Shift a given mesh with a given vector.
- Parameters:
vector (ndarray) – array containing the shifting vector
- show(show_cell_data=None, show_material_id=False, show_element_id=False, log_scale=False)
Display the mesh colored by its material ID.
- Parameters:
show_cell_data (ndarray or dict, optional) – Here you can specify additional element/cell data sorted by their IDs. It can be a dictionary with data-name as key and the ndarray as value. Default: None
show_material_id (bool, optional) – Here you can specify if the material_id should be shown. Default: False
show_element_id (bool, optional) – Here you can specify if the element_id should be shown. Default: False
log_scale (bool, optional) – State if the cell_data should be shown in log scale. Default: False
Notes
This routine needs “mayavi” to display the mesh. (see here: https://github.com/enthought/mayavi)
- swap_axis(axis1='y', axis2='z')
Swap axis of the coordinate system.
- transform(xyz_func, **kwargs)
Transform a given mesh with a given function “xyz_func”.
kwargs will be forwarded to “xyz_func”.
- Parameters:
xyz_func (function) – the function transforming the points:
x_new, y_new, z_new = f(x_old, y_old, z_old, **kwargs)
- write_file()
Write the actual OGS input file to the given folder.
Its path is given by “task_root+task_id+file_ext”.
- property ELEMENTS
Get and set the ELEMENTS of the mesh.
Notes
- Typedict of ndarrays
The elements are a dictionary sorted by their element-type
- “line”ndarray of shape (n_line,2)
1D element with 2 nodes
- “tri”ndarray of shape (n_tri,3)
2D element with 3 nodes
- “quad”ndarray of shape (n_quad,4)
2D element with 4 nodes
- “tet”ndarray of shape (n_tet,4)
3D element with 4 nodes
- “pyra”ndarray of shape (n_pyra,5)
3D element with 5 nodes
- “pris”ndarray of shape (n_pris,6)
3D element with 6 nodes
- “hex”ndarray of shape (n_hex,8)
3D element with 8 nodes
- property ELEMENT_ID
Get and set the ELEMENT_IDs of the mesh.
Standard element id order is given by:
“line” “tri” “quad” “tet” “pyra” “pris” “hex”
Notes
- Typedict of ndarrays
The element IDs are a dictionary containing ints sorted by their element-type
- “line”ndarray of shape (n_line,)
1D element with 2 nodes
- “tri”ndarray of shape (n_tri,)
2D element with 3 nodes
- “quad”ndarray of shape (n_quad,)
2D element with 4 nodes
- “tet”ndarray of shape (n_tet,)
3D element with 4 nodes
- “pyra”ndarray of shape (n_pyra,)
3D element with 5 nodes
- “pris”ndarray of shape (n_pris,)
3D element with 6 nodes
- “hex”ndarray of shape (n_hex,)
3D element with 8 nodes
- property MATERIAL_ID
Get and set the MATERIAL_IDs of the mesh.
Notes
- Typedict of ndarrays
The material IDs are a dictionary containing ints sorted by their element-type
- “line”ndarray of shape (n_line,)
1D element with 2 nodes
- “tri”ndarray of shape (n_tri,)
2D element with 3 nodes
- “quad”ndarray of shape (n_quad,)
2D element with 4 nodes
- “tet”ndarray of shape (n_tet,)
3D element with 4 nodes
- “pyra”ndarray of shape (n_pyra,)
3D element with 5 nodes
- “pris”ndarray of shape (n_pris,)
3D element with 6 nodes
- “hex”ndarray of shape (n_hex,)
3D element with 8 nodes
- property MATERIAL_ID_flat
Get flat version of the MATERIAL_IDs of the mesh.
See “mesh.MATERIAL_ID” This flattend MATERIAL_IDs are a stacked version of MATERIAL_ID, to get one continous array. They are stacked in order of the ELEMENT_IDs. Standard stack order is given by:
“line” “tri” “quad” “tet” “pyra” “pris” “hex”
Notes
- Typendarray
The centroids are a list containing xyz-coordiantes
- property NODES
(n,3) NODES of the mesh by its xyz-coordinates.
- Type:
ndarray
- property center
Get the mesh center.
- property centroids
Get the centroids of the mesh.
Notes
- Typedict of ndarrays
The centroids are a dictionary containing xyz-coordiantes sorted by their element-type
- “line”ndarray of shape (n_line,3)
1D element with 2 nodes
- “tri”ndarray of shape (n_tri,3)
2D element with 3 nodes
- “quad”ndarray of shape (n_quad,3)
2D element with 4 nodes
- “tet”ndarray of shape (n_tet,3)
3D element with 4 nodes
- “pyra”ndarray of shape (n_pyra,3)
3D element with 5 nodes
- “pris”ndarray of shape (n_pris,3)
3D element with 6 nodes
- “hex”ndarray of shape (n_hex,3)
3D element with 8 nodes
- property centroids_flat
Get flat version of the centroids of the mesh.
See the “mesh.get_centroids” method. This flattend centroids are a stacked version of centroids, to get one continous array. They are stacked in order of the element ids. Standard stack order is given by:
“line” “tri” “quad” “tet” “pyra” “pris” “hex”
Notes
- Typendarray
The centroids are a list containing xyz-coordiantes
- property is_empty
State if file is empty.
- property node_centroids
Get the node centroids of the mesh.
Notes
- Typedict of ndarrays
The centroids are a dictionary containing xyz-coordiantes sorted by their element-type
- “line”ndarray of shape (n_line,3)
1D element with 2 nodes
- “tri”ndarray of shape (n_tri,3)
2D element with 3 nodes
- “quad”ndarray of shape (n_quad,3)
2D element with 4 nodes
- “tet”ndarray of shape (n_tet,3)
3D element with 4 nodes
- “pyra”ndarray of shape (n_pyra,3)
3D element with 5 nodes
- “pris”ndarray of shape (n_pris,3)
3D element with 6 nodes
- “hex”ndarray of shape (n_hex,3)
3D element with 8 nodes
- property node_centroids_flat
Get flat version of the node centroids of the mesh.
See the “mesh.get_centroids” method. This flattend centroids are a stacked version of centroids, to get one continous array. They are stacked in order of the element ids. Standard stack order is given by:
“line” “tri” “quad” “tet” “pyra” “pris” “hex”
Notes
- Typendarray
The centroids are a list containing xyz-coordiantes
- property volumes
Get the volumes of the mesh-elements.
Notes
- Typedict of ndarrays
The volumes are a dictionary containing the n-dimension volumes sorted by their element-type
- “line”ndarray of shape (n_line,3)
1D element with 2 nodes
- “tri”ndarray of shape (n_tri,3)
2D element with 3 nodes
- “quad”ndarray of shape (n_quad,3)
2D element with 4 nodes
- “tet”ndarray of shape (n_tet,3)
3D element with 4 nodes
- “pyra”ndarray of shape (n_pyra,3)
3D element with 5 nodes
- “pris”ndarray of shape (n_pris,3)
3D element with 6 nodes
- “hex”ndarray of shape (n_hex,3)
3D element with 8 nodes
- property volumes_flat
Get flat version of the volumes of the mesh-elements.
This flattend volumes are a stacked version of centroids, to get one continous array. They are stacked in order of the element ids. Standard stack order is given by:
“line” “tri” “quad” “tet” “pyra” “pris” “hex”
Notes
- Typendarray
The volumes are a list containing the n-dimensional element volume