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.

get_elements_info()

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 vertex vtx

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:

Mesh2VecCae

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, if json_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:

Mesh2VecCae

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:

Mesh2VecCae

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:

Mesh2VecCae

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, or

    • a 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