mesh2vec.mesh2vec_cae
.Mesh2VecCae
- class mesh2vec.mesh2vec_cae.Mesh2VecCae(distance, mesh, mesh_info, calc_strategy='dfs')
Class to use finite element meshes as hypergraphs. Provide methods to add features from common CAE tools (for now: ANSA meshes, LSDYNA binout and d3plot information from lasso-cae) Inherits from Mesh2VecBase
All vertices in this class are assumed to be elements in CAE meshes (restricted to LSDYNA so far).
Methods
__init__
(distance, mesh, mesh_info[, ...])Create a Mesh2VecCae instance
add_feature_from_d3plot
(feature, d3plot[, ...])Map feature from a d3plot to the CAE elements which are the vertices of the hg.
add_features_from_ansa
(features[, ansafile, ...])Add values derived or calculated from ANSA shell elements (currently restricted to LSDYNA models) for each element.
add_features_from_csv
(csv_file[, ...])Map the content of a CSV file to the vertices of the hypergraph.
add_features_from_dataframe
(df)Map the content of a Pandas dataframe to the vertices of the hypergraph.
aggregate
(feature, dist, aggr[, aggr_name, ...])Aggregate features from neighborhoods for each distance in
dist
aggregate_angle_diff
(dist[, aggr, ...])Aggregate angle differences
aggregate_categorical
(feature, dist[, ...])For categorical features, aggregate the numbers of occurrences of each categorical value.
available_aggregated_features
()returns a list the names of all aggregated features
available_features
()returns a list the names of all features
features
()Returns a Pandas dataframe with all feature columns.
from_ansa_shell
(distance, ansafile[, ...])Read the given ANSA file and use the shell elements corresponding with
partid
to generate a hypergraph, using CAE nodes as hyperedges, and adjacent elements as hypervertices.from_d3plot_shell
(distance, d3plot[, ...])Read the given d3plot file and use the shell elements corresponding with
partid
to generate a hypergraph, using CAE nodes as hyperedges, and adjacent elements as hypervertices.from_file
(hg_file, distance[, calc_strategy])Read a hypergraph (hg) from a text file.
from_keyfile_shell
(distance, keyfile[, ...])Read the given keyfile and use the shell elements to generate a hypergraph, using mesh nodes as hyperedges, and adjacent elements as hypervertices.
from_stl_shell
(distance, stlfile)Read the given stl file and use the shell elements to generate a hypergraph, using mesh nodes as hyperedges, and adjacent elements as hypervertices.
Return a Pandas dataframe containing a row for each element in the hypergraph vertices with
get_feature_from_d3plot
(feature, d3plot_data)Map a single feature from a d3plot to the CAE elements which are the vertices of the hg.
get_max_distance
()returns the distance value used to generate the hypergraph neighborhood
get_nbh
(vtx, dist)Get a list of neighbors with the exact distance
dist
of a given vertexvtx
get_visualization_plotly
(feature)visualize an aggregated feature on mesh Useage in Notebook import plotly.io as pio pio.renderers.default = 'sphinx_gallery' fig = m2v.get_visualization_plotly(name) fig
get_visualization_trimesh
(feature)Get a trimesh with face colored by feature values Use trimesh_mesh.show(smooth=False, flags={"cull": False}) to visualize.
load
(path)Load the Mesh2Vec object from a file with joblib
mesh
()return the FE mesh
save
(path)Save the Mesh2Vec object to a file with joblib
to_array
([vertices])Returns a numpy array with all the beforehand aggregated feature columns.
to_dataframe
([vertices])Returns a Pandas dataframe with all the beforehand aggregated feature columns.
vtx_ids
()returns a list the ids of all hyper vertices
- Parameters:
distance (int) –
mesh (CaeShellMesh) –
mesh_info (DataFrame) –
calc_strategy (str) –
- __init__(distance, mesh, mesh_info, calc_strategy='dfs')
Create a Mesh2VecCae instance
- Parameters:
distance (int) – the maximum distance for neighborhood generation and feature aggregation
mesh (CaeShellMesh) – points, point_ids/uids, connectivity and element_ids/uids
mesh_info (DataFrame) –
additional info about the elements in mesh (same order is required) columns “part_name”, “part_id”, “file_path”, “element_id” are required calc_strategy: choose the algorithm to calculate adjacencies
”dfs”: depth first search (defaultl fast)
”bfs”: breadth first search (low memory consumption)
”matmul”: matrix multiplication (deprecated, for compatibility only)
calc_strategy (str) –
- Return type:
None
Example
>>> import numpy as np >>> import pandas as pd >>> from mesh2vec.mesh2vec_cae import Mesh2VecCae >>> from mesh2vec.mesh_features import CaeShellMesh >>> point_coordinates = np.array([[v, v, v] for v in range(6)]) >>> pnt_ids = np.array(["0", "1", "2", "3", "4", "5"]) >>> elem_ids = np.array(["0", "1", "2", "3"]) >>> elem_node_idxs = np.array([[0, 1, 2, 2], [1, 2, 3, 3], [2, 3, 4, 4], [3, 4, 5, 5]]) >>> mesh = CaeShellMesh(point_coordinates, pnt_ids, elem_ids, elem_node_idxs) >>> mesh_info = pd.DataFrame({"element_id": elem_ids}) >>> mesh_info["part_name"] = "part_name" >>> mesh_info["part_id"] = "part_id" >>> mesh_info["file_path"] = "file_path" >>> m2v = Mesh2VecCae(2, mesh, mesh_info) >>> m2v._hyper_edges OrderedDict([('0', ['0']), ('1', ['0', '1']), ('2', ['0', '1', '2']), ('3', ['1', '2', '3']), ('4', ['2', '3']), ('5', ['3'])])
- static from_stl_shell(distance, stlfile)
Read the given stl file and use the shell elements to generate a hypergraph, using mesh nodes as hyperedges, and adjacent elements as hypervertices.
- Parameters:
distance (int) –
stlfile (Path) –
- Return type:
- static from_ansa_shell(distance, ansafile, partid='', json_mesh_file=None, ansa_executable=None, ansa_script=None, verbose=False, calc_strategy='dfs')
Read the given ANSA file and use the shell elements corresponding with
partid
to generate a hypergraph, using CAE nodes as hyperedges, and adjacent elements as hypervertices. This is a simple wrapper around an ANSA Python script that generates a json_mesh_file, ifjson_mesh_file
is a valid file path. Otherwise, a temporary file is generated and deleted after importing.json_mesh_file
is not overwritten if it already exists . For now, this function is only guaranteed to work with shell meshes for LSDYNA. Each element gets a unique internal ID consisting of its element ID (which may be not unique), and a unique hash value.- Parameters:
distance (int) – the maximum distance for neighborhood generation and feature aggregation
ansafile (Path) – path to ansa file
partid (str) – part id to use for hypergraph generation
json_mesh_file (Optional[Path]) – path to json mesh file, if exists, it will be loaded instead of the ansafile else it will be generated.
ansa_executable (Optional[Path]) – path to ansa executable Path to ANSA executable can also be provided in environment var: ANSA_EXECUTABLE
ansa_script (Optional[Path]) – path to ansa script You can use a customized script to include more features ansa_script (see Customize Ansa script)
verbose (bool) – print additional information
calc_strategy (str) –
choose the algorithm to calculate adjacencies
”dfs”: depth first search (defaultl fast)
”bfs”: breadth first search (low memory consumption)
”matmul”: matrix multiplication (deprecated, for compatibility only)
- Return type:
Example
>>> from pathlib import Path >>> from mesh2vec.mesh2vec_cae import Mesh2VecCae >>> m2v = Mesh2VecCae.from_ansa_shell(4, ... Path("data/hat/Hatprofile.k"), ... json_mesh_file=Path("data/hat/cached_hat_key.json")) >>> len(m2v._hyper_edges) 6666
- static from_d3plot_shell(distance, d3plot, partid=None, calc_strategy='dfs')
Read the given d3plot file and use the shell elements corresponding with
partid
to generate a hypergraph, using CAE nodes as hyperedges, and adjacent elements as hypervertices.- Parameters:
distance (int) – the maximum distance for neighborhood generation and feature aggregation
d3plot (Path) – path to d3plot file
partid (Optional[str]) – part id to use
calc_strategy –
choose the algorithm to calculate adjacencies
”dfs”: depth first search (defaultl fast)
”bfs”: breadth first search (low memory consumption)
”matmul”: matrix multiplication (deprecated, for compatibility only)
- Return type:
Example
>>> from pathlib import Path >>> from mesh2vec.mesh2vec_cae import Mesh2VecCae >>> m2v = Mesh2VecCae.from_d3plot_shell(3, Path("data/hat/HAT.d3plot")) >>> len(m2v._hyper_edges) 6666
- static from_keyfile_shell(distance, keyfile, partid='', calc_strategy='bfs')
Read the given keyfile and use the shell elements to generate a hypergraph, using mesh nodes as hyperedges, and adjacent elements as hypervertices.
- Parameters:
distance (int) – the maximum distance for neighborhood generation and feature aggregation
keyfile (Path) – path to keyfile
partid – part id to use for hypergraph generation (default empty string, use all shell parts)
calc_strategy –
choose the algorithm to calculate adjacencies
”dfs”: depth first search (defaultl fast)
”bfs”: breadth first search (low memory consumption)
”matmul”: matrix multiplication (deprecated, for compatibility only)
- Return type:
Example
>>> from pathlib import Path >>> from mesh2vec.mesh2vec_cae import Mesh2VecCae >>> m2v = Mesh2VecCae.from_keyfile_shell(4, Path("data/hat/Hatprofile.k")) >>> len(m2v._hyper_edges) 6666
- get_elements_info()
Return a Pandas dataframe containing a row for each element in the hypergraph vertices with
its internal ID (vtx_id)
element ID (element_id)
part ID (part_id)
Part name, or None if not available (part_name)
File path, or None if not available (file_path)
- Return type:
DataFrame
- add_features_from_ansa(features, ansafile=None, json_mesh_file=None, ansa_executable=None, ansa_script=None, verbose=False, allow_additional_ansa_features=False)
Add values derived or calculated from ANSA shell elements (currently restricted to LSDYNA models) for each element.
Path to ANSA executable can also be provided in environment var: ANSA_EXECUTABLE
You can use a customized script to include more features ansa_script (see Customize Ansa script)
features
is a subset of:aspect: The aspect ratio of each element (ansafile is required)
warpage: (ansafile is required)
num_border: number of border edges
is_tria
midpoint x,y,z
normal vector x,y,z (ansafile is required)
area (ansafile is required)
- custom features from your customized ansa_script (scalar values only,
allow_additional_ansa_features must be True)
Example
>>> from pathlib import Path >>> from mesh2vec.mesh2vec_cae import Mesh2VecCae >>> m2v = Mesh2VecCae.from_ansa_shell( ... 4, ... Path("data/hat/Hatprofile.k"), ... json_mesh_file=Path("data/hat/cached_hat_key.json")) >>> m2v.add_features_from_ansa( ... ["aspect", "warpage"], ... Path("data/hat/Hatprofile.k"), ... json_mesh_file=Path("data/hat/cached_hat_key.json")) ['aspect', 'warpage'] >>> print(f'{m2v._features["warpage"][14]:.4f}') 0.0188
- Parameters:
features (List[str]) –
ansafile (Optional[Path]) –
json_mesh_file (Optional[Path]) –
ansa_executable (Optional[Path]) –
ansa_script (Optional[Path]) –
verbose (bool) –
allow_additional_ansa_features (bool) –
- Return type:
List[str]
- get_feature_from_d3plot(feature, d3plot_data, timestep=None, shell_layer=None, history_var_index=None)
Map a single feature from a d3plot to the CAE elements which are the vertices of the hg. Restricted to arrays available from lasso-cae.
Example
>>> from pathlib import Path >>> from lasso.dyna import ArrayType, D3plot >>> from mesh2vec.mesh2vec_cae import Mesh2VecCae >>> m2v = Mesh2VecCae.from_d3plot_shell(3, Path("data/hat/HAT.d3plot")) >>> names, values, ids = m2v.get_feature_from_d3plot( ... ArrayType.element_shell_strain, ... D3plot(Path("data/hat/HAT.d3plot").as_posix()), ... timestep=1, shell_layer=0) >>> print([f'{v:.4f}' for v in values[42]]) ['0.0010', '-0.0003', '-0.0000', '-0.0012', '-0.0000', '-0.0003']
- Parameters:
feature (str) –
d3plot_data (D3plot) –
timestep (Optional[int]) –
shell_layer (Optional[Union[int, Callable]]) –
history_var_index (Optional[int]) –
- Return type:
Tuple[str, List[Any]]
- add_feature_from_d3plot(feature, d3plot, timestep=None, shell_layer=None, history_var_index=None)
Map feature from a d3plot to the CAE elements which are the vertices of the hg. Restricted to arrays available from lasso-cae.
- Parameters:
feature (str) – name of the feature to add (a shell_array name of lasso.dyna.ArrayType)
d3plot (Union[Path, D3plot]) – path to d3plot file or loaded d3plot data
timestep (Optional[int]) – timestep to extract (required for time dependend arrays, ignored otherwise)
shell_layer (Optional[Union[int, Callable]]) – integration point index or callabe to accumulate over all integration points (required for layerd arrays, ignored otherwise)
history_var_index (Optional[int]) – index of the history variable to extract (required for element_shell_history_vars, ignored otherwise)
- Return type:
str
Example
>>> from pathlib import Path >>> from lasso.dyna import ArrayType >>> from mesh2vec.mesh2vec_cae import Mesh2VecCae >>> m2v = Mesh2VecCae.from_d3plot_shell(3, Path("data/hat/HAT.d3plot")) >>> name = m2v.add_feature_from_d3plot( ... ArrayType.element_shell_strain, ... Path("data/hat/HAT.d3plot"), ... timestep=1, shell_layer=0) >>> # print eps_xx, eps_yy, eps_zz, eps_xy, eps_yz, eps_xz of layer 0 >>> print([f'{v:.4f}' for v in m2v.features()[name][42]]) ['0.0010', '-0.0003', '-0.0000', '-0.0012', '-0.0000', '-0.0003']
Example
>>> from pathlib import Path >>> from lasso.dyna import ArrayType >>> from mesh2vec.mesh2vec_cae import Mesh2VecCae >>> m2v = Mesh2VecCae.from_d3plot_shell(3, Path("data/hat/HAT.d3plot")) >>> def mean_over_components_all_layers(v): ... return np.mean(v, axis=-1) >>> name = m2v.add_feature_from_d3plot( ... ArrayType.element_shell_strain, ... Path("data/hat/HAT.d3plot"), ... timestep=-1, shell_layer=mean_over_components_all_layers) >>> # print mean of all components for each layer (0,1) >>> print([f'{v:.4f}' for v in m2v.features()[name][42]]) ['0.0002', '0.0001']
- aggregate_angle_diff(dist, aggr=None, agg_add_ref=True, default_value=0.0, skip_arcos=False)
Aggregate angle differences
Aggregate a new feature calculated from the angle difference of each element’s normal vector to the reference vector given by the center element (in radian).
- Parameters:
dist (Union[List[int], int]) –
either
distance of the maximum neighborhood of vertex for aggregation, 0 <=
dist
<=self.distance
, ora list of distances, e.g.
[0, 2, 5]
.
aggr (Optional[Callable]) – aggregation function, default is np.mean
agg_add_ref (bool) – the aggregation callable needs the feature value of the center element as reference as 2nd argument. (default is True)
default_value (float) – value to use in aggregation when a feature is missing for a neighbor or no neighor with the given dist exist.
skip_arcos (bool) – if True, the angle difference is not calculated, but the dot product of the normal vectors is used.
- Return type:
Union[str, List[str]]
Example
>>> from pathlib import Path >>> from lasso.dyna import ArrayType >>> from mesh2vec.mesh2vec_cae import Mesh2VecCae >>> m2v = Mesh2VecCae.from_ansa_shell( ... 4, ... Path("data/hat/Hatprofile.k"), ... json_mesh_file=Path("data/hat/cached_hat_key.json")) >>> m2v.add_features_from_ansa( ... ["normal"], ... Path("data/hat/Hatprofile.k"), ... json_mesh_file=Path("data/hat/cached_hat_key.json")) ['normal'] >>> name = m2v.aggregate_angle_diff(2) >>> print(f'{ m2v._aggregated_features[name][14]:.4f}') 0.6275
- get_visualization_trimesh(feature)
Get a trimesh with face colored by feature values Use trimesh_mesh.show(smooth=False, flags={“cull”: False}) to visualize.
- Parameters:
feature (str) –
- Return type:
Trimesh
- get_visualization_plotly(feature)
visualize an aggregated feature on mesh Useage in Notebook import plotly.io as pio pio.renderers.default = ‘sphinx_gallery’ fig = m2v.get_visualization_plotly(name) fig
- Parameters:
feature (str) –
- Return type:
Figure
- mesh()
return the FE mesh
- Return type:
CaeShellMesh