Source code for prismo.solvers.fem_solver

"""
Finite Element Method (FEM) solver for electromagnetic simulations.

This module provides a FEM solver wrapper using FEniCS/dolfinx for frequency-domain
and eigenvalue electromagnetic problems.
"""

from typing import Any, Optional, Union

import numpy as np

from prismo.backends import Backend, get_backend
from prismo.core.fields import ElectromagneticFields
from prismo.core.grid import YeeGrid
from prismo.solvers.base import FrequencyDomainSolver

# Try to import FEniCS
try:
    import dolfinx
    from dolfinx import fem, mesh
    from mpi4py import MPI

    FENICS_AVAILABLE = True
except ImportError:
    try:
        import fenics

        FENICS_AVAILABLE = True
    except ImportError:
        FENICS_AVAILABLE = False
        dolfinx = None
        fenics = None


[docs] class FEMSolver(FrequencyDomainSolver): """ Finite Element Method solver for electromagnetic simulations. This solver uses FEniCS/dolfinx for frequency-domain and eigenvalue problems. It supports PMC and magnetic PML boundary conditions. Parameters ---------- grid : YeeGrid The simulation grid. geometry : Any, optional Geometry definition (to be integrated with Prismo geometry system). materials : dict, optional Material properties dictionary. boundary_conditions : dict, optional Boundary condition specification, e.g., {'x_min': 'pmc', 'x_max': 'pml'}. backend : Backend, optional Computational backend to use. """
[docs] def __init__( self, grid: YeeGrid, geometry: Optional[Any] = None, materials: Optional[dict] = None, boundary_conditions: Optional[dict] = None, backend: Optional[Union[str, Backend]] = None, ): if not FENICS_AVAILABLE: raise ImportError( "FEniCS is not available. FEniCS is not on PyPI.\n" "Install with: conda install -c conda-forge fenics-dolfinx\n" "or from source: https://fenicsproject.org/download/" ) super().__init__(grid) # Initialize backend if isinstance(backend, str): self.backend = get_backend(backend) elif isinstance(backend, Backend): self.backend = backend elif backend is None: self.backend = get_backend() else: raise TypeError("backend must be a Backend instance or string name") self.geometry = geometry self.materials = materials or {} self.boundary_conditions = boundary_conditions or {} # Physical constants self.eps0 = 8.854187817e-12 # Vacuum permittivity (F/m) self.mu0 = 4 * np.pi * 1e-7 # Vacuum permeability (H/m) self.c = 299792458.0 # Speed of light (m/s) # Initialize mesh (will be created from grid/geometry) self.mesh = None self._setup_mesh() # Initialize fields self.fields = ElectromagneticFields(grid, backend=self.backend) # FEM function space and solution self.function_space = None self.solution = None
def _setup_mesh(self) -> None: """Set up FEM mesh from grid or geometry.""" # This is a placeholder - full implementation would: # 1. Create mesh from grid dimensions # 2. Or use provided geometry to generate mesh # 3. Tag boundaries for boundary conditions if dolfinx is not None: # Use dolfinx for mesh creation # For now, create a simple box mesh nx, ny, nz = self.grid.dimensions dx, dy, dz = self.grid.spacing # Create mesh (simplified - would need proper 3D mesh generation) # This is a placeholder implementation pass else: # Fallback to fenics pass
[docs] def solve(self, frequency: float) -> ElectromagneticFields: """ Solve for fields at a specific frequency. Parameters ---------- frequency : float Frequency in Hz. Returns ------- ElectromagneticFields Solution fields at the specified frequency. """ if not FENICS_AVAILABLE: raise ImportError("FEniCS is required for FEM solver") # Angular frequency omega = 2 * np.pi * frequency k0 = omega / self.c # Free-space wavenumber # This is a placeholder - full implementation would: # 1. Set up weak form of Maxwell's equations # 2. Apply boundary conditions (PMC, magnetic PML, etc.) # 3. Assemble and solve linear system # 4. Extract solution to Prismo field format # For now, return zero fields as placeholder self.fields.zero_fields() return self.fields
[docs] def solve_eigenvalue(self, num_modes: int = 1) -> tuple[np.ndarray, list]: """ Solve eigenvalue problem (e.g., for waveguide modes). Parameters ---------- num_modes : int Number of modes to compute. Returns ------- tuple (eigenvalues, eigenmodes) where eigenmodes is a list of fields. """ if not FENICS_AVAILABLE: raise ImportError("FEniCS is required for FEM solver") # This is a placeholder - full implementation would: # 1. Set up eigenvalue problem for waveguide modes # 2. Solve using FEniCS eigenvalue solver # 3. Extract modes and effective indices # 4. Convert to Prismo field format # Placeholder return eigenvalues = np.zeros(num_modes, dtype=complex) eigenmodes = [ElectromagneticFields(self.grid) for _ in range(num_modes)] return eigenvalues, eigenmodes
[docs] def step(self, fields: Optional[ElectromagneticFields] = None) -> None: """ FEM is frequency-domain, so step() is not applicable. Use solve() or solve_eigenvalue() instead. """ raise NotImplementedError( "FEM solver is frequency-domain. Use solve() or solve_eigenvalue() instead." )
[docs] def get_fields(self) -> ElectromagneticFields: """Get the current electromagnetic fields.""" return self.fields
[docs] def reset(self) -> None: """Reset simulation to initial state.""" self.fields.zero_fields() self.time = 0.0 self.step_count = 0 self.solution = None
[docs] def get_time_step(self) -> float: """ FEM is frequency-domain, so time step is not applicable. Returns ------- float Returns 0.0 as FEM doesn't use time stepping. """ return 0.0
[docs] def __repr__(self) -> str: """String representation.""" return f"FEMSolver(grid={self.grid.dimensions}, FEniCS={FENICS_AVAILABLE})"