OpenFCST: The open-source Fuel Cell Simulation Toolbox
|
This class deals with Fick's Transport Equation. More...
#include <ficks_transport_equation.h>
Public Member Functions | |
Constructors, destructor, and initalization | |
FicksTransportEquation (FuelCell::SystemManagement &system_management, FuelCellShop::Material::PureGas *solute, FuelCellShop::Material::PureGas *solvent, boost::shared_ptr< FuelCell::ApplicationCore::ApplicationData > data=boost::shared_ptr< FuelCell::ApplicationCore::ApplicationData >()) | |
Constructor. More... | |
FicksTransportEquation (FuelCell::SystemManagement &system_management, std::string &name_section, boost::shared_ptr< FuelCell::ApplicationCore::ApplicationData > data=boost::shared_ptr< FuelCell::ApplicationCore::ApplicationData >()) | |
Constructor. More... | |
virtual | ~FicksTransportEquation () |
Destructor. More... | |
virtual void | declare_parameters (ParameterHandler ¶m) |
Declare parameters. More... | |
virtual void | initialize (ParameterHandler ¶m) |
Initialize parameters. More... | |
void | set_solute_and_solvent (FuelCellShop::Material::PureGas *solute, FuelCellShop::Material::PureGas *solvent, ParameterHandler ¶m) |
Method to set solute and solve if other constructor (not passing solute and solvent in the constructor) is being used. More... | |
Local CG FEM based assemblers | |
virtual void | assemble_cell_residual (FuelCell::ApplicationCore::FEVector &cell_rhs, const typename FuelCell::ApplicationCore::DoFApplication< dim >::CellInfo &cell_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble local cell residual for nonlinear problems. More... | |
virtual void | assemble_bdry_matrix (FuelCell::ApplicationCore::MatrixVector &bdry_matrices, const typename FuelCell::ApplicationCore::DoFApplication< dim >::FaceInfo &bdry_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble local boundary matrix. More... | |
virtual void | assemble_bdry_residual (FuelCell::ApplicationCore::FEVector &bdry_rhs, const typename FuelCell::ApplicationCore::DoFApplication< dim >::FaceInfo &bdry_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble local boundary residual. More... | |
Public Member Functions inherited from FuelCellShop::Equation::EquationBase< dim > | |
virtual void | assemble_cell_matrix (FuelCell::ApplicationCore::MatrixVector &cell_matrices, const typename FuelCell::ApplicationCore::DoFApplication< dim >::CellInfo &cell_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble local cell matrix. More... | |
const couplings_map & | get_internal_cell_couplings () const |
This function returns internal_cell_couplings of a derived equation class. More... | |
const couplings_map & | get_internal_flux_couplings () const |
This function returns internal_flux_couplings (DG FEM only) of a derived equation class. More... | |
const component_materialID_value_map & | get_component_materialID_value () const |
This function returns component_materialID_value of a derived equation class. More... | |
const component_boundaryID_value_map & | get_component_boundaryID_value () const |
This function returns component_boundaryID_value of a derived equation class. More... | |
const std::vector< BoundaryType > & | get_boundary_types () const |
This function returns boundary_types of a derived equation class. More... | |
const std::vector< std::vector < BoundaryType > > & | get_multi_boundary_types () const |
This function returns multi_boundary_types of a derived equation class. More... | |
const std::vector< OutputType > & | get_output_types () const |
This function returns output_types of a derived equation class. More... | |
const std::vector< std::vector < OutputType > > & | get_multi_output_types () const |
This function returns multi_output_types of a derived equation class. More... | |
const std::string & | get_equation_name () const |
This function returns equation_name of a derived equation class. More... | |
const std::vector< unsigned int > & | get_matrix_block_indices () const |
This function returns matrix_block_indices of a derived equation class. More... | |
const std::vector< unsigned int > & | get_residual_indices () const |
This function returns residual_indices of a derived equation class. More... | |
Accessors and info | |
virtual void | print_equation_info () const |
The function printing out the equations info. More... | |
void | class_test () |
Member function used to test the functionality of the class. More... | |
double | compute_Knudsen_diffusivity (int cell_index) |
Member function used to compute the effective Knudsen diffusivity using the Knudsen radius entered using the field data in VTK mesh. More... | |
Tensor< 2, dim, double > | effective_diffusion_coefficient (Tensor< 2, dim, double > &D_bulk, double D_Knud) |
Member function used to compute effective diffusion coefficient using the Bosanquet approximation by combining the Knudsen and Bulk diffusion coefficients. More... | |
double | compute_current (double x_gas) |
Computes the current using oxygen dissolution in the ionomer film by solving a Newton loop. More... | |
Additional Inherited Members | |
Public Attributes inherited from FuelCellShop::Equation::EquationBase< dim > | |
bool | variable_initial_data |
true , if variable initial data is prescribed on a part of the domain. More... | |
bool | variable_boundary_data |
true , if variable Dirichlet boundary conditions are prescribed on a part of the boundary. More... | |
Protected Member Functions inherited from FuelCellShop::Equation::EquationBase< dim > | |
EquationBase (FuelCell::SystemManagement &sys_management, boost::shared_ptr< FuelCell::ApplicationCore::ApplicationData > data=boost::shared_ptr< FuelCell::ApplicationCore::ApplicationData >()) | |
Constructor. More... | |
virtual | ~EquationBase () |
Destructor. More... | |
virtual void | declare_parameters (ParameterHandler ¶m) const |
Declare parameters. More... | |
virtual void | set_parameters (const std::vector< std::string > &name_dvar, const std::vector< double > &value_dvar, ParameterHandler ¶m) |
Set parameters using the parameter file, in order to run parametric/optimization studies. More... | |
virtual void | make_assemblers_generic_constant_data () |
Function used to initialize variable information that will be needed to assemble matrix and residual and that will remain constant throughout the simulation. More... | |
virtual void | make_assemblers_cell_constant_data (const typename FuelCell::ApplicationCore::DoFApplication< dim >::CellInfo &cell_info) |
Function used to initialize cell speciific information that remains constant regardless of the cell being computed. More... | |
virtual void | make_assemblers_bdry_constant_data (const typename FuelCell::ApplicationCore::DoFApplication< dim >::FaceInfo &bdry_info) |
virtual void | make_assemblers_cell_variable_data (const typename FuelCell::ApplicationCore::DoFApplication< dim >::CellInfo &cell_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Function used to compute cell specific information such as shape functions, shape function gradients, solution at each quadrature point. More... | |
virtual void | make_assemblers_bdry_variable_data (const typename FuelCell::ApplicationCore::DoFApplication< dim >::FaceInfo &bdry_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
void | select_cell_assemblers (const typename FuelCell::ApplicationCore::DoFApplication< dim >::CellInfo &cell_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
This routine is used to select the make_assembly routines that need to be called inside assemble_cell_matrix to compute. More... | |
virtual void | make_internal_cell_couplings () |
This function fills out internal_cell_couplings of a derived equation class. More... | |
virtual void | make_internal_flux_couplings () |
This function fills out internal_flux_couplings (DG FEM only) of a derived equation class. More... | |
virtual void | make_component_materialID_value () |
This function fills out component_materialID_value of a derived equation class. More... | |
virtual void | make_component_boundaryID_value () |
This function fills out component_boundaryID_value of a derived equation class. More... | |
virtual void | make_boundary_types () |
This function fills out boundary_types of a derived equation class. More... | |
virtual void | make_multi_boundary_types () |
This function fills out multi_boundary_types of a derived equation class. More... | |
virtual void | make_output_types () |
This function fills out output_types of a derived equation class. More... | |
virtual void | make_multi_output_types () |
This function fills out multi_output_types of a derived equation class. More... | |
virtual void | make_matrix_block_indices () |
This function is only needed to provide the last argument to dealII_to_appframe. More... | |
virtual void | make_residual_indices () |
This function is only needed to provide the last argument to dealII_to_appframe. More... | |
virtual void | assemble_cell_Jacobian_matrix (FuelCell::ApplicationCore::MatrixVector &cell_matrices, const typename FuelCell::ApplicationCore::DoFApplication< dim >::CellInfo &cell_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble the local Jacobian Matrix for Non-Linear problems. More... | |
virtual void | assemble_bdry_Jacobian_matrix (FuelCell::ApplicationCore::MatrixVector &bdry_matrices, const typename FuelCell::ApplicationCore::DoFApplication< dim >::FaceInfo &bdry_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble local Jacobian boundary matrix for Non-Linear problems. More... | |
virtual void | assemble_cell_linear_matrix (FuelCell::ApplicationCore::MatrixVector &cell_matrices, const typename FuelCell::ApplicationCore::DoFApplication< dim >::CellInfo &cell_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble the local cell matrix for Linear problems. More... | |
virtual void | assemble_cell_residual_rhs (FuelCell::ApplicationCore::FEVector &cell_residual, const typename FuelCell::ApplicationCore::DoFApplication< dim >::CellInfo &cell_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble local cell RHS for nonlinear problems. More... | |
virtual void | assemble_cell_linear_rhs (FuelCell::ApplicationCore::FEVector &cell_residual, const typename FuelCell::ApplicationCore::DoFApplication< dim >::CellInfo &cell_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble local cell RHS for Linear problems. More... | |
virtual void | assemble_bdry_linear_matrix (FuelCell::ApplicationCore::MatrixVector &bdry_matrices, const typename FuelCell::ApplicationCore::DoFApplication< dim >::FaceInfo &bdry_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble local boundary matrix for linear problems. More... | |
virtual void | assemble_bdry_linear_rhs (FuelCell::ApplicationCore::FEVector &bdry_residual, const typename FuelCell::ApplicationCore::DoFApplication< dim >::FaceInfo &bdry_info, FuelCellShop::Layer::BaseLayer< dim > *const layer) |
Assemble local boundary RHS for linear problems. More... | |
void | standard_to_block_wise (FullMatrix< double > &target) const |
This function changes the order of dealii::FullMatrix<double> target from standard to block-wise. More... | |
void | standard_to_block_wise (Vector< double > &target) const |
This function changes the order of dealii::Vector<double> target from standard to block-wise. More... | |
void | dealII_to_appframe (FuelCell::ApplicationCore::MatrixVector &dst, const FullMatrix< double > &src, const std::vector< unsigned int > &matrix_block_indices) const |
This function converts the standard ordered structure dealii::FullMatrix<double> src into the block-wise ordered structure FuelCell::ApplicationCore::MatrixVector dst . More... | |
void | dealII_to_appframe (FuelCell::ApplicationCore::FEVector &dst, const Vector< double > &src, const std::vector< unsigned int > &residual_indices) const |
This function converts the standard ordered structure dealii::Vector<double> src into the block-wise ordered structure FuelCell::ApplicationCore::FEVector dst . More... | |
bool | belongs_to_boundary (const unsigned int &tria_boundary_id, const unsigned int ¶m_boundary_id) const |
This function returns true if a boundary indicator of an external face on the triangulation coincides with a boundary indicator defined in the parameters file of a derived equation class. More... | |
void | print_caller_name (const std::string &caller_name) const |
This function is used to print out the name of another function that has been declared in the scope of this class, but not yet been implemented. More... | |
Protected Attributes inherited from FuelCellShop::Equation::EquationBase< dim > | |
unsigned int | dofs_per_cell |
Number of degrees of freedom per cell. More... | |
unsigned int | n_q_points_cell |
Number of quadrature points per cell. More... | |
unsigned int | n_q_points_bdry |
Number of quadrature points per boundary. More... | |
DoFHandler< dim > ::active_cell_iterator | cell |
Currently active DoFHandler<dim> active cell iterator. More... | |
DoFHandler< dim > ::active_face_iterator | bdry |
Currently active DoFHandler<dim> active boundary iterator. More... | |
std::vector< double > | JxW_cell |
Jacobian of mapping by Weight in the quadrature points of a cell. More... | |
std::vector< double > | JxW_bdry |
Jacobian of mapping by Weight in the quadrature points of a boundary. More... | |
std::vector< Point< dim > > | normal_vectors |
Normal vectors in the quadrature points of a boundary. More... | |
std::vector< std::vector < Point< dim > > > | tangential_vectors |
Tangential vectors in the quadrature points of a boundary. More... | |
FuelCell::SystemManagement * | system_management |
Pointer to the external YourApplication<dim>::system_management object. More... | |
couplings_map | internal_cell_couplings |
This object contains the info on how the equations and solution variables of a derived equation class are coupled. More... | |
couplings_map | internal_flux_couplings |
This object contains the info on how the "X" and "Y" of a derived equation class are coupled (DG FEM only). More... | |
component_materialID_value_map | component_materialID_value |
This object reflects the following structure (see FuelCell::InitialAndBoundaryData namespace docs): More... | |
component_boundaryID_value_map | component_boundaryID_value |
This object reflects the following structure (see FuelCell::InitialAndBoundaryData namespace docs): More... | |
std::vector< BoundaryType > | boundary_types |
The list of boundary types of a derived equation class. More... | |
std::vector< std::vector < BoundaryType > > | multi_boundary_types |
The list of multiple boundary types of a derived equation class. More... | |
std::vector< OutputType > | output_types |
The list of output types of a derived equation class. More... | |
std::vector< std::vector < OutputType > > | multi_output_types |
The list of multiple output types of a derived equation class. More... | |
std::string | equation_name |
The name of a derived equation class. More... | |
std::string | name_base_variable |
Const std::string member storing name of the base solution variable corresponding to the equation represented by this class. More... | |
std::vector< unsigned int > | matrix_block_indices |
The system matrix block indices (a derived equation class) drawn from the global structure (a derived equation class + other active equation classes included into the computation). More... | |
std::vector< unsigned int > | residual_indices |
The residual indices (a derived equation class) drawn from the global structure (a derived equation class + other active equation classes included into the computation). More... | |
std::vector< bool > | counter |
This vector contains the collection of internal "counters" used by the derived equation classes. More... | |
EquationFlags | assemble_flags |
This vector contains a collection of internal flags to tell derived equation classes what needs to be re-computed. More... | |
boost::shared_ptr < FuelCell::ApplicationCore::ApplicationData > | data |
Data object for the application data to be passed to the equation classes. More... | |
std::string | solution_vector_name |
The name of the solution vector in FEVectors. More... | |
std::string | residual_vector_name |
The name of the residual vector name in FEVectors. More... | |
This class deals with Fick's Transport Equation.
This equation class solves for Fick's law of diffusion inside the porous layers. The solute gas and solvent gas are normally passed inside the constructor of this equation class.
It is solved with respect to:
(solutegas_molar_fraction
)This equation can be written as:
cm^2/s
], which can be function of
(temperature_of_REV
) and
(liquid_water_saturation
).mol/cm^3
].To be well-posed, these equations are equipped with the appropriate boundary conditions. All the boundary conditions can be described by boundary_id
(s
) and
boundary_type
. Besides, this some boundary types require additional information, which can also be provided by the parameter file. We consider following types of boundary conditions:
The boundary_ids
are specified in the parameter file under subsection "Dirichlet Boundary Indicators", as a list of comma-separated values.
e.g.
Neumann
boundary condition.No
gas
flux
or Symmetric
boundary conditions, as FEM formulation automatically implies a particular boundary is one of these cases, by default.set_gases
method.We solve the whole problem by linearizing the governing equation at each Newton iteration with subsequent CG FEM discretization in space. The class contains the necessary class members to add the necessary contributions to cell_matrix and cell_residual to the governing equations used to analyze gas transport via ficks diffusion model,
ReactionSourceTerms
class. Please read the documentation of ReactionSourceTerms class, for additional methods to be implemented in the application.adjust_internal_cell_couplings
member function of ReactionSourceTerms class, before using make_cell_couplings
of SystemManagement at the application level.TODO Old Boundary conditions including Dirichlet Boundary Indicators are supposed to be replaced with the new subsections, see TO BE REMOVED comments in .cc file.
FuelCellShop::Equation::FicksTransportEquation< dim >::FicksTransportEquation | ( | FuelCell::SystemManagement & | system_management, |
FuelCellShop::Material::PureGas * | solute, | ||
FuelCellShop::Material::PureGas * | solvent, | ||
boost::shared_ptr< FuelCell::ApplicationCore::ApplicationData > | data = boost::shared_ptr< FuelCell::ApplicationCore::ApplicationData >() |
||
) |
Constructor.
FuelCellShop::Equation::FicksTransportEquation< dim >::FicksTransportEquation | ( | FuelCell::SystemManagement & | system_management, |
std::string & | name_section, | ||
boost::shared_ptr< FuelCell::ApplicationCore::ApplicationData > | data = boost::shared_ptr< FuelCell::ApplicationCore::ApplicationData >() |
||
) |
Constructor.
|
virtual |
Destructor.
|
virtual |
Assemble local boundary matrix.
Reimplemented from FuelCellShop::Equation::EquationBase< dim >.
|
virtual |
Assemble local boundary residual.
Reimplemented from FuelCellShop::Equation::EquationBase< dim >.
|
virtual |
Assemble local cell residual for nonlinear problems.
Reimplemented from FuelCellShop::Equation::EquationBase< dim >.
void FuelCellShop::Equation::FicksTransportEquation< dim >::class_test | ( | ) |
Member function used to test the functionality of the class.
It create an object of this class and test functionality.
|
protected |
Computes the current using oxygen dissolution in the ionomer film by solving a Newton loop.
The ICCP model is used to solve the current for the microstructure. This function runs the Newton loop to solve for the oxygen concentration at the ionomer-solid interface. The function returns the computed current density. The coupled equations solved are:
where is the ionomer dissolution rate, is the diffusion rate coefficient for oxygen in nafion, is the ionomer film thickness, is the oxygen concentration at gas-ionomer interface, is the equilibrium oxygen concentration computed using Henry's Law, is the reactant concentration at the ionomer-solid interface and is the current density. The residual for the Newton loops is implemented in the function get_ICCP_residual
. The derivative of the residual wrt to is approximated using the central difference method.
double FuelCellShop::Equation::FicksTransportEquation< dim >::compute_Knudsen_diffusivity | ( | int | cell_index | ) |
Member function used to compute the effective Knudsen diffusivity using the Knudsen radius entered using the field data in VTK mesh.
Cell index is used as the key to reference the local cell radius which is used to compute the Knudsen diffusion coefficient.
|
virtual |
Declare parameters.
Tensor<2,dim,double> FuelCellShop::Equation::FicksTransportEquation< dim >::effective_diffusion_coefficient | ( | Tensor< 2, dim, double > & | D_bulk, |
double | D_Knud | ||
) |
Member function used to compute effective diffusion coefficient using the Bosanquet approximation by combining the Knudsen and Bulk diffusion coefficients.
|
virtual |
Initialize parameters.
Reimplemented from FuelCellShop::Equation::EquationBase< dim >.
|
virtual |
The function printing out the equations info.
Reimplemented from FuelCellShop::Equation::EquationBase< dim >.
|
inline |
Method to set solute and solve if other constructor (not passing solute and solvent in the constructor) is being used.
It will also setup name_equation and name_solution
References FuelCellShop::Equation::EquationBase< dim >::equation_name, FuelCellShop::Equation::EquationBase< dim >::name_base_variable, and FuelCellShop::Material::BaseMaterial::name_material().