:py:mod:`FELiCS.Equation.EquationCollection` ============================================ .. py:module:: FELiCS.Equation.EquationCollection Module Contents --------------- Classes ~~~~~~~ .. autoapisummary:: FELiCS.Equation.EquationCollection.EquationCollectionClass Attributes ~~~~~~~~~~ .. autoapisummary:: FELiCS.Equation.EquationCollection.logger .. py:data:: logger .. py:class:: EquationCollectionClass(param, FEMSpaces, mean, mesh) Constructs variational formulations for equations and manages boundary conditions. This class generates variational formulations for various matrices, including: - A (real and imaginary parts of the linear operator) - B (real and imaginary parts of the weight/time-derivative matrix) - B_forcing (resolvent forcing norm) - B_response (resolvent response norm) It also manages test and trial functions, boundary conditions, and mesh-based coordinate transformations. The design supports multiple physical models (momentum, mass, energy, species) and enables configuration through an external parameter object. **Initialize the EquationCollectionClass object** :param param: Parameters for the simulation, including boundary conditions, equations, and I/O settings. :type param: :py:class:`object` :param FEMSpaces: Collection of Finite Element Method spaces for the mixed formulation. :type FEMSpaces: :py:class:`object` :param mean: Mean flow fields used in variational formulations. :type mean: :py:class:`object` :param mesh: Mesh object defining the spatial domain and coordinate system. :type mesh: :py:class:`object` :ivar x: Spatial coordinates of the computational domain. :vartype x: :py:class:`ufl.SpatialCoordinate` :ivar ds: Measure for integration over boundary subdomains. :vartype ds: :py:class:`ufl.Measure` :ivar all_ds: Combined measure over all boundary subdomains. :vartype all_ds: :py:class:`ufl.Measure` :ivar n_BC: Outward-pointing unit normal on the domain boundary. :vartype n_BC: :py:class:`ufl.FacetNormal` :ivar n: Tensor form of the boundary normal vector in the chosen coordinate system. :vartype n: :py:class:`Tensor` :ivar BCs: List of Dirichlet boundary condition objects. :vartype BCs: :py:class:`list` :ivar trialFunctionsFEM: Mixed trial functions used in the variational problem. :vartype trialFunctionsFEM: :py:class:`ufl.TrialFunctions` :ivar testFunctionsFEM: Mixed test functions used in the variational problem. :vartype testFunctionsFEM: :py:class:`ufl.TestFunctions` :ivar X: List of structured test functions with tensor-awareness. :vartype X: :py:class:`list` of :py:class:`Tensor` :ivar R: Radial coordinate value for cylindrical systems or unity for Cartesian. :vartype R: :py:class:`ufl.Coefficient` :ivar equationList: List of equation objects that define the weak formulations. :vartype equationList: :py:class:`list` :ivar resolventResponseIndices: Indices selecting solution components used in the resolvent response norm. :vartype resolventResponseIndices: :py:class:`numpy.ndarray` :ivar resolventForcingIndices: Indices selecting solution components used in the resolvent forcing norm. :vartype resolventForcingIndices: :py:class:`numpy.ndarray` .. admonition:: Notes The class dynamically loads and constructs equation objects depending on the configuration provided in `param.Case.Equations`. It supports tensor-based formulations and multiple coordinate systems. .. admonition:: Examples >>> eq_collection = EquationCollectionClass(param, FEMSpaces, mean, mesh) >>> A = eq_collection.getLinearOperator(mean) >>> B = eq_collection.getWeightMatrix(mean) .. py:method:: getLinearOperator(meanFlow) Construct the linear operator matrix for the equation system. This method assembles the linear part of the weak formulation by iterating through all equations in the equation list and adding their linear expressions. :param meanFlow: Mean flow properties used in constructing the linear operator :type meanFlow: :py:class:`object` :returns: Assembled linear operator matrix with boundary conditions applied :rtype: :py:class:`petsc4py.PETSc.Mat` .. py:method:: getWeightMatrix(meanFlow) Construct the weight matrix representing the time derivative term. Assembles the matrix by summing contributions from all equations. Boundary conditions are included for resolvent analysis mode only. :param meanFlow: Mean flow properties used in constructing the weight matrix. :type meanFlow: :py:class:`object` :returns: Assembled weight matrix. :rtype: :py:class:`petsc4py.PETSc.Mat` .. admonition:: Notes - For Resolvent analysis, boundary conditions are applied - For other analysis modes, no boundary conditions are applied to avoid computational issues during eigenvalue problem solving .. py:method:: getFEMWeightMatrix() Construct the full Finite Element Method (FEM) weight matrix. This method creates a weight matrix by performing inner products of test and trial functions across all function spaces in the mixed function space. :returns: Assembled FEM weight matrix without boundary conditions :rtype: :py:class:`petsc4py.PETSc.Mat` .. admonition:: Notes - Handles both scalar and vector function spaces - Uses conjugate of test functions for matrix construction - Includes a workaround for mesh object compatibility with different versions of DOLFINx .. py:method:: getFEMDiffusionMatrix(diffusionFactor, sponge=None) Construct the Finite Element Method (FEM) diffusion matrix. This method creates a diffusion matrix by computing the weak form of the diffusion operator across all function spaces in the mixed function space. :param diffusionFactor: Coefficient multiplying the diffusion terms :type diffusionFactor: :py:class:`float` :param sponge: Additional damping term to be applied to the matrix (default is None) :type sponge: :py:class:`float` or :py:obj:`None`, *optional* :returns: Assembled FEM diffusion matrix with boundary conditions applied :rtype: :py:class:`petsc4py.PETSc.Mat` .. admonition:: Notes - Handles both scalar and vector function spaces - Computes second-order derivatives in x and y directions - Currently a quick implementation; a more comprehensive tensor framework is needed for future improvements .. py:method:: getFullRHS(func) Construct the full FEM right-hand side (RHS) vector. Forms the RHS by integrating each test function against the corresponding part of the provided function. :param func: Mixed function representing the source term or solution. :type func: :py:class:`dolfinx.fem.Function` :returns: Assembled RHS vector. :rtype: :py:class:`petsc4py.PETSc.Vec` .. py:method:: getNonlinearExpression(meanFlow, setBC=True) Construct the nonlinear vector expression for the system. Assembles the nonlinear contributions from all equations. :param meanFlow: Mean flow field used in the nonlinear formulation. :type meanFlow: :py:class:`object` :param setBC: Whether to apply boundary conditions (default is True). :type setBC: :py:class:`bool`, *optional* :returns: Assembled nonlinear vector. :rtype: :py:class:`petsc4py.PETSc.Vec` .. py:method:: getBilinearOperator(meanFlow) Construct the bilinear operator matrix. Assembles all bilinear terms across the system equations. :param meanFlow: Mean flow field used in the bilinear formulation. :type meanFlow: :py:class:`object` :returns: Assembled bilinear operator matrix with boundary conditions applied. :rtype: :py:class:`petsc4py.PETSc.Mat` .. py:method:: getForcingForInputOutput(meanFlow) Construct the forcing vector for input-output analysis. Creates ufl object with the linear equation system. :param meanFlow: Mean flow properties used to construct the forcing vector. :type meanFlow: :py:class:`object` :returns: Assembled complex forcing vector. :rtype: :py:class:`petsc4py.PETSc.Vec` .. py:method:: getResolventNorm_response(meanFlow) Construct the resolvent response norm matrix. :param meanFlow: Mean flow field. Currently not used directly here, but kept for interface consistency with other assembly routines. :type meanFlow: :py:class:`object` :returns: PETSc matrix representing the response norm, assembled from the UFL form stored in ``self.response_vf``. :rtype: :py:class:`petsc4py.PETSc.Mat` .. py:method:: getResolventNorm_forcing(meanFlow) Construct the resolvent forcing norm matrix. :param meanFlow: Mean flow field. Currently not used directly here, but kept for interface consistency with other assembly routines. :type meanFlow: :py:class:`object` :returns: PETSc matrix representing the forcing norm, assembled from the UFL form stored in ``self.forcing_vf``. :rtype: :py:class:`petsc4py.PETSc.Mat` .. py:method:: getResolventWeighting_FEM(meanFlow) Assemble the FEM weighting matrix for resolvent analysis. The underlying UFL form for the weighting is constructed beforehand by :meth:`computeResolventFEMWeights` and stored in :attr:`self.fem_weighting`. This method only assembles the corresponding PETSc matrix. :param meanFlow: Mean flow field. Currently not used directly here, but included for API compatibility with other assembly routines. :type meanFlow: :py:class:`object` :returns: PETSc matrix with FEM weights based on the previously defined variable projections. :rtype: :py:class:`petsc4py.PETSc.Mat` .. py:method:: computeResolventNorms(X, param, mean, fluc) Compute energy norms for resolvent forcing and response. :param X: Test functions with tensor awareness. :type X: :py:class:`list` of :py:class:`Tensor` :param param: Simulation configuration parameters. :type param: :py:class:`object` :param mean: Mean flow fields. :type mean: :py:class:`object` :param fluc: Fluctuation fields object. :type fluc: :py:class:`object` .. admonition:: Notes Only 'TKE' and 'Chu' norm types are currently implemented. .. py:method:: computeResolventFEMWeights(X, param, mean, fluc) Compute FEM weighting matrix for resolvent input/output scaling. :param X: Test functions. :type X: :py:class:`list` of :py:class:`Tensor` :param param: Parameter configuration. :type param: :py:class:`object` :param mean: Mean flow fields. :type mean: :py:class:`object` :param fluc: Fluctuation fields. :type fluc: :py:class:`object` .. py:method:: getRestrictorMatResponse() Construct the response restrictor matrix for the resolvent analysis. This matrix is diagonal and applies spatial restriction (if defined) to the response vector based on the `responseDomain` field in the mean flow. :returns: PETSc matrix with diagonal entries corresponding to the spatial restriction mask. :rtype: :py:class:`petsc4py.PETSc.Mat` .. admonition:: Notes - Provides a quadratic matrix, with the size of the solution space (VMixed). - Has the response restrictor values, given with the mean field, on the diagonal. .. py:method:: getRestrictorMatForcing() Construct the forcing restrictor matrix for the resolvent analysis. This matrix is diagonal and applies spatial restriction (if defined) to the forcing vector based on the ``forcingDomain`` field in the mean flow. If the forcing domain is zero everywhere, the matrix reduces to the identity. :returns: PETSc matrix with diagonal entries corresponding to the spatial restriction mask. :rtype: :py:class:`petsc4py.PETSc.Mat` .. admonition:: Notes - The matrix is quadratic with the size of the solution space (``VMixed``). - The diagonal entries are given by the forcing restrictor values provided with the mean field (or ones everywhere if no spatial restriction is prescribed). .. py:method:: getShrinkerMatResponse() Construct a shrinker matrix for the resolvent response vector. :returns: PETSc matrix projecting full response to a reduced subspace. :rtype: :py:class:`petsc4py.PETSc.Mat` .. py:method:: getShrinkerMatForcing() Creates a (possibly rectangular) matrix which serves the purpose to "shrink" the forcing vector to the requested size (e.g. only consider the velocity components when using the forcing norm "TKE"). :returns: PETSc matrix projecting full forcing vector to reduced subspace. :rtype: :py:class:`petsc4py.PETSc.Mat`