|
| 1 | +"""Module for common class between Archive, and result mesh.""" |
| 2 | +from ansys.mapdl.reader import _reader, _relaxmidside |
| 3 | +from ansys.mapdl.reader.elements import ETYPE_MAP |
| 4 | +from ansys.mapdl.reader.misc import unique_rows |
| 5 | +import numpy as np |
| 6 | +import pyvista as pv |
| 7 | +from pyvista._vtk import VTK9 |
| 8 | + |
| 9 | +INVALID_ALLOWABLE_TYPES = TypeError( |
| 10 | + "`allowable_types` must be an array " "of ANSYS element types from 1 and 300" |
| 11 | +) |
| 12 | + |
| 13 | +# map MESH200 elements to a pymapdl_reader/VTK element type (see elements.py) |
| 14 | +MESH200_MAP = { |
| 15 | + 0: 2, # line |
| 16 | + 1: 2, # line |
| 17 | + 2: 2, # line |
| 18 | + 3: 2, # line |
| 19 | + 4: 3, # triangle |
| 20 | + 5: 3, # triangle |
| 21 | + 6: 3, # quadrilateral |
| 22 | + 7: 3, # quadrilateral |
| 23 | + 8: 5, # tetrahedron with 4 nodes |
| 24 | + 9: 5, # tetrahedron with 10 nodes |
| 25 | + 10: 4, # hex with 8 nodes |
| 26 | + 11: 4, |
| 27 | +} # hex with 8 nodes |
| 28 | + |
| 29 | +SHAPE_MAP = { # from ELIST definition |
| 30 | + 0: "", |
| 31 | + 1: "LINE", |
| 32 | + 2: "PARA", |
| 33 | + 3: "ARC ", |
| 34 | + 4: "CARC", |
| 35 | + 5: "", |
| 36 | + 6: "TRIA", |
| 37 | + 7: "QUAD", |
| 38 | + 8: "TRI6", |
| 39 | + 9: "QUA8", |
| 40 | + 10: "POIN", |
| 41 | + 11: "CIRC", |
| 42 | + 12: "", |
| 43 | + 13: "", |
| 44 | + 14: "CYLI", |
| 45 | + 15: "CONE", |
| 46 | + 16: "SPHE", |
| 47 | + 17: "", |
| 48 | + 18: "", |
| 49 | + 19: "PILO", |
| 50 | +} |
| 51 | +# element type to VTK conversion function call map |
| 52 | +# 0: skip |
| 53 | +# 1: Point |
| 54 | +# 2: Line (linear or quadratic) |
| 55 | +# 3: Shell |
| 56 | +# 4: 3D Solid (Hexahedral, wedge, pyramid, tetrahedral) |
| 57 | +# 5: Tetrahedral |
| 58 | +# 6: Line (always linear) |
| 59 | +TARGE170_MAP = { |
| 60 | + "TRI": 3, # 3-Node Triangle |
| 61 | + "QUAD": 3, # 4-Node Quadrilateral |
| 62 | + "CYLI": 0, # Not supported (NS) # Cylinder |
| 63 | + "CONE": 0, # NS # Cone |
| 64 | + "TRI6": 3, # 6-Node triangle |
| 65 | + "SPHE": 0, # NS # Sphere |
| 66 | + "PILO": 1, # Pilot Node |
| 67 | + "QUAD8": 3, # 8-Node Quadrilateral |
| 68 | + "LINE": 2, # Line |
| 69 | + "PARA": 2, # Parabola |
| 70 | + "POINT": 1, # Point |
| 71 | +} |
| 72 | + |
| 73 | + |
| 74 | +def _parse_vtk( |
| 75 | + mesh, |
| 76 | + allowable_types=None, |
| 77 | + force_linear=False, |
| 78 | + null_unallowed=False, |
| 79 | + fix_midside=True, |
| 80 | + additional_checking=False, |
| 81 | +): |
| 82 | + """Convert raw ANSYS nodes and elements to a VTK UnstructuredGrid |
| 83 | +
|
| 84 | + Parameters |
| 85 | + ---------- |
| 86 | + fix_midside : bool, optional |
| 87 | + Adds additional midside nodes when ``True``. When |
| 88 | + ``False``, missing ANSYS cells will simply point to the |
| 89 | + first node. |
| 90 | +
|
| 91 | + """ |
| 92 | + if not mesh._has_nodes or not mesh._has_elements: |
| 93 | + # warnings.warn('Missing nodes or elements. Unable to parse to vtk') |
| 94 | + return |
| 95 | + |
| 96 | + etype_map = ETYPE_MAP |
| 97 | + if allowable_types is not None: |
| 98 | + try: |
| 99 | + allowable_types = np.asarray(allowable_types) |
| 100 | + except: |
| 101 | + raise INVALID_ALLOWABLE_TYPES |
| 102 | + |
| 103 | + if not issubclass(allowable_types.dtype.type, np.integer): |
| 104 | + raise TypeError("Element types must be an integer array-like") |
| 105 | + |
| 106 | + if allowable_types.min() < 1 or allowable_types.max() > 300: |
| 107 | + raise INVALID_ALLOWABLE_TYPES |
| 108 | + |
| 109 | + etype_map = np.zeros_like(ETYPE_MAP) |
| 110 | + etype_map[allowable_types] = ETYPE_MAP[allowable_types] |
| 111 | + |
| 112 | + # ANSYS element type to VTK map |
| 113 | + type_ref = np.empty(2 << 16, np.int32) # 131072 |
| 114 | + type_ref[mesh._ekey[:, 0]] = etype_map[mesh._ekey[:, 1]] |
| 115 | + |
| 116 | + if allowable_types is None or 200 in allowable_types: |
| 117 | + for etype_ind, etype in mesh._ekey: |
| 118 | + |
| 119 | + # MESH200 |
| 120 | + if etype == 200 and etype_ind in mesh.key_option: |
| 121 | + # keyoption 1 contains various cell types |
| 122 | + # map them to the corresponding type (see elements.py) |
| 123 | + mapped = MESH200_MAP[mesh.key_option[etype_ind][0][1]] |
| 124 | + type_ref[etype_ind] = mapped |
| 125 | + |
| 126 | + # TARGE170 specifics |
| 127 | + if etype == 170: |
| 128 | + # edge case where missing element within the tshape_key |
| 129 | + if etype_ind not in mesh.tshape_key: # pragma: no cover |
| 130 | + continue |
| 131 | + tshape_num = mesh.tshape_key[etype_ind] |
| 132 | + if tshape_num >= 19: # weird bug when 'PILO' can be 99 instead of 19. |
| 133 | + tshape_num = 19 |
| 134 | + tshape_label = SHAPE_MAP[tshape_num] |
| 135 | + type_ref[etype_ind] = TARGE170_MAP.get(tshape_label, 0) |
| 136 | + |
| 137 | + offset, celltypes, cells = _reader.ans_vtk_convert( |
| 138 | + mesh._elem, mesh._elem_off, type_ref, mesh.nnum, True |
| 139 | + ) # for reset_midside |
| 140 | + |
| 141 | + nodes, angles, nnum = mesh.nodes, mesh.node_angles, mesh.nnum |
| 142 | + |
| 143 | + # fix missing midside |
| 144 | + if np.any(cells == -1): |
| 145 | + if fix_midside: |
| 146 | + nodes, angles, nnum = fix_missing_midside( |
| 147 | + cells, nodes, celltypes, offset, angles, nnum |
| 148 | + ) |
| 149 | + else: |
| 150 | + cells[cells == -1] = 0 |
| 151 | + |
| 152 | + if additional_checking: |
| 153 | + cells[cells < 0] = 0 |
| 154 | + # cells[cells >= nodes.shape[0]] = 0 # fails when n_nodes < 20 |
| 155 | + |
| 156 | + if VTK9: |
| 157 | + grid = pv.UnstructuredGrid(cells, celltypes, nodes, deep=True) |
| 158 | + else: |
| 159 | + grid = pv.UnstructuredGrid(offset, cells, celltypes, nodes, deep=True) |
| 160 | + |
| 161 | + # Store original ANSYS element and node information |
| 162 | + grid.point_data["ansys_node_num"] = nnum |
| 163 | + grid.cell_data["ansys_elem_num"] = mesh.enum |
| 164 | + grid.cell_data["ansys_real_constant"] = mesh.elem_real_constant |
| 165 | + grid.cell_data["ansys_material_type"] = mesh.material_type |
| 166 | + grid.cell_data["ansys_etype"] = mesh._ans_etype |
| 167 | + grid.cell_data["ansys_elem_type_num"] = mesh.etype |
| 168 | + |
| 169 | + # add components |
| 170 | + # Add element components to unstructured grid |
| 171 | + for key, item in mesh.element_components.items(): |
| 172 | + mask = np.in1d(mesh.enum, item, assume_unique=True) |
| 173 | + grid.cell_data[key] = mask |
| 174 | + |
| 175 | + # Add node components to unstructured grid |
| 176 | + for key, item in mesh.node_components.items(): |
| 177 | + mask = np.in1d(nnum, item, assume_unique=True) |
| 178 | + grid.point_data[key] = mask |
| 179 | + |
| 180 | + # store node angles |
| 181 | + if angles is not None: |
| 182 | + if angles.shape[1] == 3: |
| 183 | + grid.point_data["angles"] = angles |
| 184 | + |
| 185 | + if not null_unallowed: |
| 186 | + grid = grid.extract_cells(grid.celltypes != 0) |
| 187 | + |
| 188 | + if force_linear: |
| 189 | + # only run if the grid has points or cells |
| 190 | + if grid.n_points: |
| 191 | + grid = grid.linear_copy() |
| 192 | + |
| 193 | + # map over element types |
| 194 | + # Add tracker for original node numbering |
| 195 | + ind = np.arange(grid.n_points) |
| 196 | + grid.point_data["origid"] = ind |
| 197 | + grid.point_data["VTKorigID"] = ind |
| 198 | + return grid |
| 199 | + |
| 200 | + |
| 201 | +def fix_missing_midside(cells, nodes, celltypes, offset, angles, nnum): |
| 202 | + """Adds missing midside nodes to cells. |
| 203 | +
|
| 204 | + ANSYS sometimes does not add midside nodes, and this is denoted in |
| 205 | + the element array with a ``0``. When translated to VTK, this is |
| 206 | + saved as a ``-1``. If this is not corrected, VTK will segfault. |
| 207 | +
|
| 208 | + This function creates missing midside nodes for the quadratic |
| 209 | + elements. |
| 210 | + """ |
| 211 | + # Check for missing midside nodes |
| 212 | + mask = cells == -1 |
| 213 | + nnodes = nodes.shape[0] |
| 214 | + |
| 215 | + nextra = mask.sum() |
| 216 | + cells[mask] = np.arange(nnodes, nnodes + nextra) |
| 217 | + |
| 218 | + nodes_new = np.empty((nnodes + nextra, 3)) |
| 219 | + nodes_new[:nnodes] = nodes |
| 220 | + nodes_new[nnodes:] = 0 # otherwise, segfault disaster |
| 221 | + |
| 222 | + # Set new midside nodes directly between their edge nodes |
| 223 | + temp_nodes = nodes_new.copy() |
| 224 | + _relaxmidside.reset_midside(cells, celltypes, offset, temp_nodes) |
| 225 | + |
| 226 | + # merge midside nodes |
| 227 | + unique_nodes, idx_a, idx_b = unique_rows(temp_nodes[nnodes:]) |
| 228 | + |
| 229 | + # rewrite node numbers |
| 230 | + cells[mask] = idx_b + nnodes |
| 231 | + nextra = idx_a.shape[0] # extra unique nodes |
| 232 | + nodes_new = nodes_new[: nnodes + nextra] |
| 233 | + nodes_new[nnodes:] = unique_nodes |
| 234 | + |
| 235 | + if angles is not None: |
| 236 | + new_angles = np.empty((nnodes + nextra, 3)) |
| 237 | + new_angles[:nnodes] = angles |
| 238 | + new_angles[nnodes:] = 0 |
| 239 | + else: |
| 240 | + new_angles = None |
| 241 | + |
| 242 | + # Add extra node numbers |
| 243 | + nnum_new = np.empty(nnodes + nextra) |
| 244 | + nnum_new[:nnodes] = nnum |
| 245 | + nnum_new[nnodes:] = -1 |
| 246 | + return nodes_new, new_angles, nnum_new |
0 commit comments