Wrapping MAPDL as a Tesseract¶

This example wraps Ansys MAPDL as a differentiable Tesseract that can e.g. be used for SIMP (Solid Isotropic Material with Penalization) topology optimization.

The Tesseract computes structural compliance using finite element analysis and provides analytical sensitivities for gradient-based optimization.

Prerequisites¶

This example requires a running Ansys MAPDL server accessible via gRPC, e.g. via:

# Linux
/usr/ansys_inc/v241/ansys/bin/ansys241 -grpc -port 50052

# Windows
"C:/Program Files/ANSYS Inc/v241/ansys/bin/winx64/ANSYS241.exe" -grpc -port 50052

Tesseract definition¶

Core functionality — schemas and `apply` function¶

The Tesseract accepts a hexahedral mesh, density field, and boundary conditions as inputs. The density field rho is marked as Differentiable to enable gradient computation:

class HexMesh(BaseModel):
    """Hexahedral mesh representation with 8-node elements."""

    points: Array[(None, 3), Float32] = Field(description="Array of vertex positions.")
    elems: Array[(None, 8), Int32] = Field(
        description="Array of hexahedral elems defined by indices into the points array."
    )
    n_points: Int32 = Field(
        default=0, description="Number of valid points in the points array."
    )
    n_elems: Int32 = Field(
        default=0, description="Number of valid elems in the elems array."
    )

class InputSchema(BaseModel):
    """Input specification for SIMP compliance computation via ANSYS MAPDL.

    Defines all inputs required to run a static structural analysis with SIMP
    material interpolation. Requires a running ANSYS MAPDL server accessible
    via gRPC protocol.

    Example MAPDL server startup:
        Linux: /usr/ansys_inc/v241/ansys/bin/ansys241 -grpc -port 50052
        Windows: "C:/Program Files/ANSYS Inc/v241/ansys/bin/winx64/ANSYS241.exe" -grpc -port 50052

    The density field (rho) is marked as Differentiable, enabling gradient
    computation via vector-Jacobian products for topology optimization.
    """

    host: str = Field(description="The IP of your MAPDL grpc server.")
    port: str = Field(description="The port of your MAPDL grpc server.")
    rho: Differentiable[
        Array[
            (
                None,
                None,
            ),
            Float32,
        ]
    ] = Field(description="2D density field for topology optimization")
    von_neumann_mask: Array[(None,), Int32] = Field(
        description="Mask for von Neumann boundary conditions",
    )
    von_neumann_values: Array[(None, None), Float32] = Field(
        description="Values for von Neumann boundary conditions",
    )
    dirichlet_mask: Array[(None,), Int32] = Field(
        description="Mask for Dirichlet boundary conditions",
    )
    dirichlet_values: Array[(None,), Float32] = Field(
        description="Values for Dirichlet boundary conditions",
    )
    hex_mesh: HexMesh = Field(
        description="Hexahedral mesh representation of the geometry",
    )

    E0: float = Field(
        default=1.0,
        description="Base Young's modulus in Pa for SIMP material interpolation.",
    )

    rho_min: float = Field(
        default=1e-6,
        description="Minimum density value to avoid singular stiffness matrix in SIMP.",
    )

    p: float = Field(
        default=3.0,
        description="SIMP penalty parameter for material interpolation (default: 3.0).",
    )

    log_level: str = Field(
        default="WARNING",
        description="Logging level for output messages (DEBUG, INFO, WARNING, ERROR).",
    )

    vtk_output: str | None = Field(
        default=None,
        description="The path to write the results in VTK format.",
    )

The outputs include compliance (the optimization objective), element-wise strain energy, and sensitivity values:

class OutputSchema(BaseModel):
    """Output specification from SIMP compliance analysis.

    Contains the results of a static structural analysis including compliance
    (objective function for minimum compliance optimization), element-wise
    strain energy, and sensitivity information for gradient-based optimization.

    The compliance is marked as Differentiable to support automatic differentiation
    through the ANSYS solver using cached sensitivities from the forward pass.
    """

    compliance: Differentiable[
        Array[
            (),
            Float32,
        ]
    ] = Field(
        description="Compliance of the structure (scalar), a measure of structural flexibility (inverse of stiffness)"
    )
    strain_energy: Array[(None,), Float32] = Field(
        description="Element-wise strain energy array of shape (n_elements,)"
    )
    sensitivity: Array[(None,), Float32] = Field(
        description="Derivative of compliance with respect to each density variable: shape (n_elements,)"
    )

The apply function connects to MAPDL, creates the mesh, assigns SIMP-interpolated materials based on density, applies boundary conditions, solves the static analysis, and extracts results:

def apply(inputs: InputSchema) -> OutputSchema:
    """Run the pymapdl Tesseract for SIMP-Elasticity analysis."""
    # Initialize MAPDL
    mapdl = Mapdl(inputs.host, port=inputs.port)
    mapdl.clear()

    # Create solver instance and run analysis
    solver = SIMPElasticity(inputs, mapdl)
    return solver.solve()

SIMP material interpolation¶

The solver creates a unique material for each element with properties scaled by the density field using the SIMP formula:

    @log_timing
    def _define_simp_materials(self) -> None:
        """Define materials using SIMP approach with batch commands for efficiency."""
        # Flatten rho array and ensure it matches the number of elements
        rho_flat = np.array(self.rho).flatten()

        if len(rho_flat) != self.n_elements:
            raise ValueError(
                f"Density field size {len(rho_flat)} does not match "
                f"number of elements {self.n_elements}"
            )

        # Apply minimum density constraint to avoid singular stiffness matrix
        # SIMP formula: E = E0 * rho^p
        E_values = self.E0 * (self.rho_min + (1 - self.rho_min) * (rho_flat**self.p))
        dens_values = self.rho_min + (1 - self.rho_min) * rho_flat

        # Build all commands as a list for batch execution
        commands = []
        for i in range(len(rho_flat)):
            mat_id = i + 1
            commands.append(f"MP,EX,{mat_id},{E_values[i]}")
            commands.append(f"MP,DENS,{mat_id},{dens_values[i]}")
            commands.append(f"MP,NUXY,{mat_id},0.3")
        self.mapdl.input_strings(commands)

This interpolation scheme penalizes intermediate densities, encouraging binary (solid/void) designs in topology optimization.

Vector-Jacobian product endpoint¶

Abstract evaluation¶

The abstract_eval function computes output shapes based on input shapes without running the Ansys solver, enabling static analysis and memory pre-allocation:

def abstract_eval(abstract_inputs: Any) -> dict:
    """Calculate output shapes and dtypes without executing the forward pass."""
    n_elems = abstract_inputs.hex_mesh.elems.shape[0]
    return {
        "compliance": ShapeDType(shape=(), dtype="float32"),
        "strain_energy": ShapeDType(shape=(n_elems,), dtype="float32"),
        "sensitivity": ShapeDType(shape=(n_elems,), dtype="float32"),
    }

Wrapping MAPDL as a Tesseract¶

Prerequisites¶

Tesseract definition¶

Core functionality — schemas and apply function¶

SIMP material interpolation¶

Vector-Jacobian product endpoint¶

Abstract evaluation¶

Demo script¶

Core functionality — schemas and `apply` function¶