waves
Basic FE playground
flow::PotentialResidual Class Reference

Formulation of nonlinear potential equation residuals. More...

#include <wPotentialResidual.h>

Static Public Member Functions

static Eigen::MatrixXd buildFixed (tbox::Element const &e, std::vector< double > const &phi, Medium const &fluid)
 Build the residual matrix for a fixed-point iteration, on one element. More...
 
static Eigen::VectorXd build (tbox::Element const &e, tbox::Element const &eU, std::vector< double > const &phi, Medium const &fluid, double muC, double mCO)
 Build the residual vector, on one element. More...
 
static std::tuple< Eigen::MatrixXd, Eigen::MatrixXd > buildGradientFlow (tbox::Element const &e, tbox::Element const &eU, std::vector< double > const &phi, Medium const &fluid, double muC, double mCO)
 Build the gradient of the residual vector with respect to the flow variable (jacobian), on one element. More...
 
static std::tuple< Eigen::MatrixXd, Eigen::MatrixXd > buildGradientMesh (tbox::Element const &e, tbox::Element const &eU, std::vector< double > const &phi, Medium const &fluid, int nDim, double muC, double mCO)
 Build the gradient of the residual vector with respect to the nodes, on one element. More...
 

Detailed Description

Formulation of nonlinear potential equation residuals.

Todo:

mu as F0ElMu?

split stabilization terms from subsonic terms?

Member Function Documentation

◆ build()

Eigen::VectorXd PotentialResidual::build ( tbox::Element const &  e,
tbox::Element const &  eU,
std::vector< double > const &  phi,
Medium const &  fluid,
double  muC,
double  mCO 
)
static

Build the residual vector, on one element.

b = \int ( (1-mu)*rho + mu*rhoU ) * grad_phi * grad_psi dV

◆ buildFixed()

Eigen::MatrixXd PotentialResidual::buildFixed ( tbox::Element const &  e,
std::vector< double > const &  phi,
Medium const &  fluid 
)
static

Build the residual matrix for a fixed-point iteration, on one element.

A = \int rho * grad_phi * grad_psi dV

◆ buildGradientFlow()

std::tuple< Eigen::MatrixXd, Eigen::MatrixXd > PotentialResidual::buildGradientFlow ( tbox::Element const &  e,
tbox::Element const &  eU,
std::vector< double > const &  phi,
Medium const &  fluid,
double  muC,
double  mCO 
)
static

Build the gradient of the residual vector with respect to the flow variable (jacobian), on one element.

A = \int d( (1-mu)*rho + mu*rhoU ) * grad_phi * grad_psi dV = \int (1-mu)*drho * grad_phi * grad_psi dV

  • \int mu*drhoU * grad_phi * grad_psi dV
  • \int ( (1-mu)*rho + mu*rhoU ) * dgrad_phi * grad_psi dV
  • \int (rhoU - rho)*dmu * grad_phi * grad_psi dV

◆ buildGradientMesh()

std::tuple< Eigen::MatrixXd, Eigen::MatrixXd > PotentialResidual::buildGradientMesh ( tbox::Element const &  e,
tbox::Element const &  eU,
std::vector< double > const &  phi,
Medium const &  fluid,
int  nDim,
double  muC,
double  mCO 
)
static

Build the gradient of the residual vector with respect to the nodes, on one element.

A = d( \int ( (1-mu)*rho + mu*rhoU ) * grad_phi * grad_psi dV ) = \int (1-mu)*drho * grad_phi * grad_psi dV

  • \int mu*drhoU * grad_phi * grad_psi dV
  • \int ( (1-mu)*rho + mu*rhoU ) * dgrad_phi * grad_psi dV
  • \int ( (1-mu)*rho + mu*rhoU ) * grad_phi * dgrad_psi dV
  • \int (rhoU - rho)*dmu * grad_phi * grad_psi dV
  • \int ( (1-mu)*rho + mu*rhoU ) * grad_phi * grad_psi ddV

The documentation for this class was generated from the following files: