"""Load MAPDL CDB input decks into an femorph-solver :class:`~femorph_solver.Model`.
Thin wrapper around :mod:`mapdl_archive`. Only the file-format reader
is involved — MAPDL itself is never invoked.
The CDB's ``ET`` table (MAPDL element-type declarations such as
``ET, 1, 186``) is translated on load so the returned Model already
knows each ET id's kernel — callers don't have to re-declare them.
References
----------
* CDB format definition (MAPDL input deck, produced by the
``CDWRITE`` command): Ansys, Inc., *Ansys Mechanical APDL
Command Reference*, Release 2022R2, entry ``CDWRITE``; and
*Ansys Mechanical APDL Basic Analysis Guide*, chapter on the
database file.
* Open-source parser: ``mapdl-archive``
(https://github.com/akaszynski/mapdl-archive, MIT) — the
production parser we delegate to. Does not invoke MAPDL, only
reads the text grammar.
"""
from __future__ import annotations
import re
from pathlib import Path
from typing import TYPE_CHECKING, Any
import numpy as np
from femorph_solver.labels import UnitSystem, unit_system_from_mapdl
if TYPE_CHECKING:
from femorph_solver.model import Model
#: MAPDL element-type-number → registered kernel name. Covers every
#: type femorph-solver ships a kernel for; a deck that uses an
#: unrecognised type silently skips the declaration (the caller gets
#: a clean "element type N not declared" error the first time they
#: try to use it). Keys match the numeric codes stored in
#: ``mapdl_archive.Archive.ekey``.
_MAPDL_NUM_TO_NAME: dict[int, str] = {
185: "HEX8",
186: "HEX20",
187: "TET10",
188: "BEAM2",
181: "QUAD4_SHELL",
182: "QUAD4_PLANE",
180: "TRUSS2",
14: "SPRING",
21: "POINT_MASS",
}
#: VTK cell-type → kernel name for the quadratic-3D family. SOLID186
#: in MAPDL is a single declaration (ET = HEX20) but degenerates to
#: wedge / pyramid / tet at element instantiation time. ``from_cdb``
#: splits one ET-id-per-shape so the assembler picks the right kernel.
#:
#: Keys come from :mod:`vtkmodules.util.numpy_support`'s VTK constants
#: (24 = QUADRATIC_TETRA, 25 = QUADRATIC_HEXAHEDRON, 26 =
#: QUADRATIC_WEDGE, 27 = QUADRATIC_PYRAMID).
_VTK_CELLTYPE_TO_KERNEL: dict[int, str] = {
24: "TET10",
25: "HEX20",
26: "WEDGE15",
27: "PYR13",
# Linear hexahedra reuse the HEX8 kernel (SOLID185).
12: "HEX8",
}
# ``/UNITS,<token>`` or ``/UNITS,<token>,...`` — whitespace allowed
# around the comma; token is the leading comma-separated field. CDB
# parsers (including MAPDL itself) only look at the first token.
_UNITS_RE = re.compile(r"^\s*/UNITS\s*,\s*([A-Za-z0-9_]+)", re.IGNORECASE)
def _detect_unit_system(path: str | Path) -> UnitSystem:
"""Scan the CDB for a ``/UNITS,<token>`` command and return its :class:`UnitSystem`.
Returns :attr:`UnitSystem.UNSPECIFIED` when no ``/UNITS`` command
appears in the first 8 KiB of the file (the header region). MAPDL
decks place ``/UNITS`` near the top; scanning further down would
just slow this call without improving accuracy.
"""
try:
with open(path, encoding="utf-8", errors="ignore") as fh:
head = fh.read(8192)
except OSError: # pragma: no cover - defensive; from_cdb opens again
return UnitSystem.UNSPECIFIED
for line in head.splitlines():
m = _UNITS_RE.match(line)
if m:
return unit_system_from_mapdl(m.group(1))
return UnitSystem.UNSPECIFIED
[docs]
def from_cdb(path: str, **kwargs: Any) -> Model:
"""Load a MAPDL CDB file as a femorph-solver :class:`Model`.
Reads a MAPDL ASCII archive deck via :mod:`mapdl_archive`,
parses the deck's ``/UNITS`` directive (if present) onto
:attr:`Model.unit_system`, and registers the CDB's element-type
(``ET``) table so the resulting :class:`Model` is ready for a
follow-up :meth:`Model.assign` call (or, if the deck already
carries materials and BCs, ready to :meth:`Model.solve_static`).
Parameters
----------
path : str or pathlib.Path
Path to a CDB archive file.
**kwargs
Forwarded to :class:`mapdl_archive.Archive`. See its docs
for the supported keys (e.g. ``parse_vtk=True``).
Returns
-------
femorph_solver.Model
Native :class:`Model` whose grid carries the CDB's nodes,
elements, element-type table, and (when ``/UNITS`` is
declared) :attr:`UnitSystem` stamp.
Raises
------
FileNotFoundError
If ``path`` does not point at an existing file.
ImportError
If the ``mapdl`` extra is not installed.
See Also
--------
femorph_solver.interop.nastran.from_bdf : NASTRAN BDF loader.
femorph_solver.interop.abaqus.from_inp : Abaqus INP loader.
femorph_solver.interop.optistruct.from_fem : OptiStruct FEM loader.
Examples
--------
The bundled `mapdl_archive` test fixtures ship a small CDB; round-
trip it through :func:`from_cdb` to confirm the call shape.
>>> import femorph_solver as fs # doctest: +SKIP
>>> import mapdl_archive # doctest: +SKIP
>>> from femorph_solver.mapdl_api.cdb import from_cdb # doctest: +SKIP
>>> path = mapdl_archive.examples.sector_archive_file # doctest: +SKIP
>>> model = from_cdb(path) # doctest: +SKIP
>>> model.n_nodes > 0 # doctest: +SKIP
True
"""
try:
import mapdl_archive
except ImportError as exc: # pragma: no cover - import-guard path
raise ImportError(
"Reading MAPDL CDB decks requires the 'mapdl' extra. "
"Install with: pip install femorph_solver[mapdl]"
) from exc
from femorph_solver.model import Model
# ``mapdl_archive.Archive`` raises an opaque ``RuntimeError: Error
# opening file`` from the C++ reader when the path doesn't exist.
# Surface the missing file as a ``FileNotFoundError`` instead so
# callers get the conventional Python error type with the path
# echoed in the message.
p = Path(path)
if not p.is_file():
raise FileNotFoundError(f"CDB file not found: {p}")
archive = mapdl_archive.Archive(path, **kwargs)
grid = archive.grid
# ``mapdl_archive`` populates ``cell_data["ansys_elem_type_num"]``
# with the MAPDL element-type *number* (e.g. 186 for SOLID186),
# but :class:`Model` keys ``_etypes`` on the deck's ET *id* (1, 2,
# …) and the assembly path looks up ``ansys_elem_type_num`` as an
# ET id. ``ansys_etype`` already carries the per-cell ET id so we
# remap the field before handing the grid to the Model.
if "ansys_etype" in grid.cell_data:
grid.cell_data["ansys_elem_type_num"] = np.asarray(
grid.cell_data["ansys_etype"], dtype=np.int32
)
# SOLID186 (and its lower-order kin) is declared once per ET id
# but instantiates as a mix of HEX20 / TET10 / WEDGE15 / PYR13 by
# element-time degeneration. The assembler dispatches on a single
# kernel per ET id, so split each affected ET into one ET-id-per-
# VTK-celltype and remap ``ansys_elem_type_num`` cell-by-cell.
ekey = getattr(archive, "ekey", None)
et_kernel: dict[int, str] = {} # final ET id → kernel name
if ekey is not None:
existing_ids = {int(row[0]) for row in ekey}
next_synth_id = (max(existing_ids) if existing_ids else 0) + 1
celltypes = np.asarray(grid.celltypes, dtype=np.int32)
et_per_cell = np.asarray(grid.cell_data["ansys_elem_type_num"], dtype=np.int32)
ekey_kernel: dict[int, str] = {}
for row in ekey:
mapdl_num = int(row[1])
name = _MAPDL_NUM_TO_NAME.get(mapdl_num)
if name is not None:
ekey_kernel[int(row[0])] = name
# Map (orig_et_id, vtk_celltype) → (assigned_et_id, kernel).
synth: dict[tuple[int, int], tuple[int, str]] = {}
new_et_per_cell = et_per_cell.copy()
for i in range(et_per_cell.size):
orig = int(et_per_cell[i])
ct = int(celltypes[i])
kernel = _VTK_CELLTYPE_TO_KERNEL.get(ct)
base = ekey_kernel.get(orig)
if kernel is None:
# Fall back to the ET's declared kernel — the assembler
# will raise a clear error on shape mismatch if the
# caller's kernel doesn't actually match the cell.
if base is not None:
et_kernel[orig] = base
continue
if base is not None and base == kernel:
# ET's declared kernel already matches this shape.
et_kernel[orig] = kernel
continue
key = (orig, ct)
assigned = synth.get(key)
if assigned is None:
if base is None:
# ET id not in our table — claim it for this kernel.
assigned_id = orig
else:
assigned_id = next_synth_id
next_synth_id += 1
while assigned_id in existing_ids:
assigned_id = next_synth_id
next_synth_id += 1
existing_ids.add(assigned_id)
synth[key] = (assigned_id, kernel)
et_kernel[assigned_id] = kernel
assigned = (assigned_id, kernel)
new_et_per_cell[i] = assigned[0]
grid.cell_data["ansys_elem_type_num"] = new_et_per_cell
model = Model.from_grid(grid)
model.set_unit_system(_detect_unit_system(path))
# Register the (possibly expanded) ET table. ``_et_impl`` uses
# the private path to skip the ``Model.et()`` deprecation warning
# on every ``from_cdb`` call.
for et_id, name in et_kernel.items():
model._et_impl(et_id, name)
# Preserve the CDB's named-component table on the Model so
# downstream callers can ``model.fix(component="_FIXEDSU")``
# rather than re-opening the archive to pull the node-id list.
# Component dicts are ``{name: int32 array of ids}`` — match
# the layout :meth:`Model.fix` expects.
nc = getattr(archive, "node_components", None)
if nc:
for name, ids in nc.items():
model._node_components[str(name)] = np.asarray(ids, dtype=np.int64)
ec = getattr(archive, "element_components", None)
if ec:
for name, ids in ec.items():
model._element_components[str(name)] = np.asarray(ids, dtype=np.int64)
return model