SERiF 0.0.1a
3+1D Stellar Structure and Evolution
Loading...
Searching...
No Matches
serif::polytrope::formBundle Struct Reference

Structure to hold the various bilinear and nonlinear forms for the polytrope problem. More...

#include <polySolver.h>

Public Attributes

std::unique_ptr< mfem::MixedBilinearForm > M
 Mixed bilinear form $ M(\psi^\theta, N^\phi) = \int_\Omega \nabla\psi^\theta \cdot N^\phi \,dV $. This form arises from the term $\nabla \theta$ in the first equation, tested with a vector test function $N^\phi$. It couples the scalar field $\theta$ (potential) with the vector field $\boldsymbol{\phi}$ (flux).
 
std::unique_ptr< mfem::MixedBilinearForm > Q
 Mixed bilinear form $ Q(\psi^\phi, N^\theta) = \int_\Omega \psi^\phi \cdot \nabla N^\theta \,dV $. This form arises from the term $\boldsymbol{\phi}$ in the first equation, tested with a scalar test function $N^\theta$, after integration by parts of the $\nabla \theta$ term. It can also be seen as related to the transpose of a discrete gradient operator.
 
std::unique_ptr< mfem::BilinearForm > D
 Bilinear form $ D(\psi^\phi, N^\phi) = \int_\Omega \psi^\phi \cdot N^\phi \,dV $. This is a mass matrix for the vector field $\boldsymbol{\phi}$. It arises from the $\boldsymbol{\phi}$ term in the first equation when tested with a vector test function $N^\phi$.
 
std::unique_ptr< mfem::BilinearForm > S
 Bilinear form $ S(\psi^\theta, N^\theta) = \int_\Omega \nabla\psi^\theta \cdot \nabla N^\theta \,dV $. This is a stiffness matrix (Laplacian) for the scalar field $\theta$. It is used for stabilization terms or alternative formulations.
 
std::unique_ptr< mfem::NonlinearForm > f
 Nonlinear form $ f(\theta)(\psi^\theta) = \int_\Omega \psi^\theta \cdot \theta^n \,dV $. This form arises from the nonlinear source term $\theta^n$ in the second equation, tested with a scalar test function $\psi^\theta$.
 

Detailed Description

Structure to hold the various bilinear and nonlinear forms for the polytrope problem.

This structure bundles the MFEM forms that represent the discretized weak formulation of the mixed variable polytropic equations. The system being solved is typically:

\[\begin{aligned}
  \boldsymbol{\phi} + \nabla \theta &= \mathbf{0} \\
  \nabla \cdot \boldsymbol{\phi} - \theta^n &= 0
\end{aligned}
\]

After integration by parts and discretization, these lead to matrix equations involving the forms stored here.

Examples
/Users/tboudreaux/Programming/SERiF/src/polytrope/solver/public/polySolver.h.

Definition at line 137 of file polySolver.h.

Member Data Documentation

◆ D

std::unique_ptr<mfem::BilinearForm> serif::polytrope::formBundle::D

Bilinear form $ D(\psi^\phi, N^\phi) = \int_\Omega \psi^\phi \cdot N^\phi \,dV $. This is a mass matrix for the vector field $\boldsymbol{\phi}$. It arises from the $\boldsymbol{\phi}$ term in the first equation when tested with a vector test function $N^\phi$.

Examples
/Users/tboudreaux/Programming/SERiF/src/polytrope/solver/public/polySolver.h.

Definition at line 158 of file polySolver.h.

◆ f

std::unique_ptr<mfem::NonlinearForm> serif::polytrope::formBundle::f

Nonlinear form $ f(\theta)(\psi^\theta) = \int_\Omega \psi^\theta \cdot \theta^n \,dV $. This form arises from the nonlinear source term $\theta^n$ in the second equation, tested with a scalar test function $\psi^\theta$.

Examples
/Users/tboudreaux/Programming/SERiF/src/polytrope/solver/public/polySolver.h.

Definition at line 172 of file polySolver.h.

◆ M

std::unique_ptr<mfem::MixedBilinearForm> serif::polytrope::formBundle::M

Mixed bilinear form $ M(\psi^\theta, N^\phi) = \int_\Omega \nabla\psi^\theta \cdot N^\phi \,dV $. This form arises from the term $\nabla \theta$ in the first equation, tested with a vector test function $N^\phi$. It couples the scalar field $\theta$ (potential) with the vector field $\boldsymbol{\phi}$ (flux).

Examples
/Users/tboudreaux/Programming/SERiF/src/polytrope/solver/public/polySolver.h.

Definition at line 143 of file polySolver.h.

◆ Q

std::unique_ptr<mfem::MixedBilinearForm> serif::polytrope::formBundle::Q

Mixed bilinear form $ Q(\psi^\phi, N^\theta) = \int_\Omega \psi^\phi \cdot \nabla N^\theta \,dV $. This form arises from the term $\boldsymbol{\phi}$ in the first equation, tested with a scalar test function $N^\theta$, after integration by parts of the $\nabla \theta$ term. It can also be seen as related to the transpose of a discrete gradient operator.

Examples
/Users/tboudreaux/Programming/SERiF/src/polytrope/solver/public/polySolver.h.

Definition at line 151 of file polySolver.h.

◆ S

std::unique_ptr<mfem::BilinearForm> serif::polytrope::formBundle::S

Bilinear form $ S(\psi^\theta, N^\theta) = \int_\Omega \nabla\psi^\theta \cdot \nabla N^\theta \,dV $. This is a stiffness matrix (Laplacian) for the scalar field $\theta$. It is used for stabilization terms or alternative formulations.

Examples
/Users/tboudreaux/Programming/SERiF/src/polytrope/solver/public/polySolver.h.

Definition at line 165 of file polySolver.h.


The documentation for this struct was generated from the following file: