OpenFCST: The open-source Fuel Cell Simulation Toolbox
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Member Functions | Protected Attributes | List of all members
FuelCell::Application::AppCathode< dim > Class Template Reference

This class is used to solve a system of equations similar to the one presented in the journal article M. More...

#include <app_cathode.h>

Inheritance diagram for FuelCell::Application::AppCathode< dim >:
Inheritance graph
[legend]
Collaboration diagram for FuelCell::Application::AppCathode< dim >:
Collaboration graph
[legend]

Public Member Functions

virtual void solve (AppFrame::FEVector &start, const AppFrame::FEVectors &rhs)
 Solve the linear system of equations.
 
virtual double estimate (const AppFrame::FEVectors &sol)
 Estimate error per cell.
 
virtual double evaluate (const AppFrame::FEVectors &src)
 Post-processing.
 
virtual void data_out (const std::string &basename, const AppFrame::FEVectors &src)
 Reimplementation of the routine in the base class BaseApplication in namespace AppFrame so that the right labels are outputed and so that I can compute and output the source term.
 
void parse_init_solution (AppFrame::FEVector &src, FuelCell::OperatingConditions &OC, bool &good_solution)
 Check to see that the solution read from file can be used for the new simulation conditions.
 
bool parse_init_solution_phi_s (AppFrame::FEVector &src, double &V_cell)
 Check to see that the initial solution for the solid phase potential is usable.
 
bool parse_init_solution_x_O2 (AppFrame::FEVector &src, double &x_O2)
 Check to see that the initial solution for the oxygen mole fraction is usable.
 
Constructors, destructor, and initalization
 AppCathode (boost::shared_ptr< AppFrame::ApplicationData > data=boost::shared_ptr< AppFrame::ApplicationData >())
 Constructor.
 
 ~AppCathode ()
 Destructor.
 
virtual void declare_parameters (ParameterHandler &param)
 Declare all parameters that are needed for:
 
virtual void set_parameters (const std::vector< std::string > &name_dvar, const std::vector< double > &value_dvar, ParameterHandler &param)
 Function called by optimization loop in order to set the values in the ParameterHandler to the new design parameters.
 
void _initialize (ParameterHandler &param)
 Set up how many equations are needed and read in parameters for the parameter handler in order to initialize data.
 
virtual void initialize (ParameterHandler &param)
 Call the other initialize routines from the inherited classes.
 
void initialize_grid (ParameterHandler &param)
 Initilize and create the grid Determines whether reading from a grid file or creating new.
 
void init_solution (AppFrame::FEVector &src)
 Initialize nonlinear solution.
 
Local CG FEM based assemblers
virtual void cell_matrix (MatrixVector &cell_matrices, const typename DoFApplication< dim >::CellInfo &cell)
 Integration of local bilinear form.
 
virtual void bdry_matrix (MatrixVector &bdry_matrices, const typename DoFApplication< dim >::FaceInfo &bdry)
 Integration of local bilinear form.
 
virtual void cell_residual (AppFrame::FEVector &cell_vector, const typename DoFApplication< dim >::CellInfo &cell)
 Integration of the rhs of the equations.
 
virtual void bdry_residual (AppFrame::FEVector &bdry_vector, const typename DoFApplication< dim >::FaceInfo &bdry)
 
virtual void dirichlet_bc (std::map< unsigned int, double > &boundary_values) const
 Member function used to set dirichlet boundary conditions.
 
Optimization Member Functions
virtual void cell_dresidual_dlambda (std::vector< AppFrame::FEVector > &cell_vector, const typename DoFApplication< dim >::CellInfo &cell, std::vector< std::vector< double > > &src)
 Integration of local bilinear form.
 
virtual void check_responses ()
 This class is called by responses to make sure that all responses requested are implemented in either cell_responses, global_responses or face_responses.
 
virtual void cell_responses (std::vector< double > &resp, const typename DoFApplication< dim >::CellInfo &info, const AppFrame::FEVector &sol)
 Compute the value of all objective function and constraints.
 
virtual void global_responses (std::vector< double > &resp, const AppFrame::FEVector &sol)
 This class is used to evaluate all responses that do not require looping over cells.
 
virtual void cell_dresponses_dl (std::vector< std::vector< double > > &cell_df_dl, const typename DoFApplication< dim >::CellInfo &info, const AppFrame::FEVector &sol)
 This class is used to evaluate the derivative of all the functionals that require looping over cells with respect to the design variables.
 
virtual void global_dresponses_dl (std::vector< std::vector< double > > &df_dl, const AppFrame::FEVector &sol)
 This class is used to evaluate the sensitivities of all responses that do not require looping over cells with respect of the design variables.
 
virtual void cell_dresponses_du (std::vector< AppFrame::FEVector > &cell_df_du, const typename DoFApplication< dim >::CellInfo &info, std::vector< std::vector< double > > &src)
 This class is used to evaluate the derivative of all the functionals that require looping over cells with respect of the unknowns of the system of governing equations.
 
virtual void global_dresponses_du (std::vector< AppFrame::FEVector > &df_du, const AppFrame::FEVector &src)
 This class is used to evaluate the sensitivities of all responses that do not require looping over cells with respecto of the unknowns of the system of governing equations.
 
- Public Member Functions inherited from AppFrame::OptimizationBlockMatrixApplication< dim >
 OptimizationBlockMatrixApplication (boost::shared_ptr< AppFrame::ApplicationData > data=boost::shared_ptr< AppFrame::ApplicationData >())
 Constructor for an object, owning its own mesh and dof handler.
 
 OptimizationBlockMatrixApplication (AppFrame::DoFApplication< dim > &, bool triangulation_only)
 Constructor for an object, borrowing mesh and dof handler from another object.
 
 ~OptimizationBlockMatrixApplication ()
 Destructor.
 
virtual void responses (std::vector< double > &f, const AppFrame::FEVectors &vectors)
 Post-processing.
 
virtual void print_responses (std::vector< double > &resp)
 This function is used to print the responses to screen (deallog)
 
virtual void dresponses_dl (std::vector< std::vector< double > > &df_dl, const AppFrame::FEVectors &src)
 Post-processing.
 
virtual void dresponses_du (std::vector< AppFrame::FEVector > &dst, const AppFrame::FEVectors &src)
 Loop over all cells and assemble the vector

\[ \frac{\partial f_m}{\partial u_i} \]

that is used to solve the sensitivity equations by using the local matrix integration functions cell_dfunctional_du(), bdry_dfunctional_dlu() and edge_dfunctinal_du() provided by the derived class.

 
virtual void dresidual_dlambda (std::vector< AppFrame::FEVector > &dst, const AppFrame::FEVectors &src)
 Loop over all cells and assemble the vector

\[ \frac{\partial R_n}{\partial \lambda_k} \]

that is used to solve the sensitivity equations by using the local matrix integration functions cell_dresidual_dlambda(), bdry_dresidual_dlambda() and edge_dresidual_dlambda() provided by the derived class.

 
void solve_direct (std::vector< std::vector< double > > &df_dl, const AppFrame::FEVectors &sol)
 Solver in order to obtain the analytical sensitivities for the given objective and constraints using the direct method.
 
void solve_adjoint (std::vector< std::vector< double > > &df_dl, const AppFrame::FEVector &sol)
 Solver in order to obtain the analytical sensitivities for the given objective and constraints using the adjoint method.
 
unsigned int get_n_resp () const
 Member function that returns the number of responses.
 
unsigned int get_n_dvar () const
 Member function that returns the number of design variables.
 
std::vector< std::string > get_name_dvar () const
 Member function that returns the name of the design variables.
 
std::vector< std::string > get_name_responses () const
 Member function that returns the name of the responses.
 
void print_default_parameter_file ()
 Print the default parameter handler file.
 
bool get_bool_transfer_solution ()
 Function to see if we are transferring a solution on a refined grid to the initial coarse grid.
 
virtual void set_read_initial_solution (bool &read)
 
void set_optimization_parameters (unsigned int &n_dvar, unsigned int &n_resp, std::vector< std::string > &name_design_var, std::vector< std::string > &name_responses)
 This routine assigns the value of n_dvar, n_resp, name_design_var and name_responses to the internal variables n_dvar, n_resp, name_design_var and name_response.
 
void set_output_variables (std::vector< std::string > &dakota_name_responses)
 
- Public Member Functions inherited from AppFrame::BlockMatrixApplication< dim >
void do_assemble (const FEVectors &, typename DoFHandler< dim >::active_cell_iterator begin, typename DoFHandler< dim >::active_cell_iterator end, Threads::Mutex &mutex)
 Multithreaded assemble function called by assemble().
 
virtual void post_cell_assemble ()
 A call back function for assemble(), called after all cell matrices have been entered, but before the face matrices are computed.
 
 BlockMatrixApplication (boost::shared_ptr< ApplicationData > data=boost::shared_ptr< ApplicationData >())
 Constructor for an object, owning its own mesh and dof handler.
 
 BlockMatrixApplication (DoFApplication< dim > &, bool triangulation_only)
 Constructor for an object, borrowing mesh and dof handler from another object.
 
void _initialize (ParameterHandler &param)
 Initialize data of this class.
 
void add_vector_for_cell_matrix (std::string name, unsigned int block, unsigned int components, bool values, bool derivatives)
 Add a named vector from data to be handed to the cell_matrix() function via CellInfo::values and CellInfo::derivatives.
 
void add_vector_for_flux_matrix (std::string name, unsigned int block, unsigned int components, bool values, bool derivatives)
 Add a named vector from data to be handed to the bdry_matrix() and face_matrix() functions via FaceInfo::values and FaceInfo::derivatives.
 
void remesh_matrices ()
 Initialize sparsity patterns and matrices for the new mesh.
 
virtual void remesh ()
 Refine grid accordingly to the Refinement options under "Grid Generation" in the parameter file.
 
void assemble (const FEVectors &)
 Loop over all cells and assemble the system #matrices by using the local matrix integration functions cell_matrix(), bdry_matrix() and face_matrix() provided by the derived class.
 
void assemble_numerically (const FEVectors &src, const double delta=1e-6)
 Compute the Jacobian of the system of equations,

\[ J(i,j) = \frac{\partial R_i}{\partial u_j} \]

by using forward differences.

 
void residual_constraints (FEVector &dst) const
 Redefinition of residual_constraints() in DoFHandler.
 
virtual void face_matrix (MatrixVector &matrices11, MatrixVector &matrices12, MatrixVector &matrices21, MatrixVector &matrices22, const typename DoFApplication< dim >::FaceInfo &face1, const typename DoFApplication< dim >::FaceInfo &face2)
 Integration of local bilinear form.
 
- Public Member Functions inherited from AppFrame::DoFApplication< dim >
 DoFApplication (boost::shared_ptr< ApplicationData > data=boost::shared_ptr< ApplicationData >(), bool multigrid=false)
 Constructor for an object, owning its own mesh and dof handler.
 
 DoFApplication (bool multigrid)
 Constructor for an object owning its own mesh and dof handler and creating new ApplicationData.
 
 DoFApplication (DoFApplication< dim > &, bool triangulation_only, bool multigrid=false)
 Constructor for an object, borrowing mesh and dof handler from another object.
 
 ~DoFApplication ()
 Destructor which deletes owned objects.
 
void _initialize (ParameterHandler &param)
 Initialize from parameter values.
 
void add_vector_for_transfer (FEVector *)
 Add the vector to be transfered from one mesh to the next.
 
void delete_vector_for_transfer ()
 Delete the vector to be transfered from one mesh to the next.
 
virtual void remesh_dofs ()
 Initialize dof handler, count the dofs in each block and renumber the dofs.
 
virtual void init_vector (FEVector &dst) const
 Reinitialize the BlockVector dst such that it contains block_size.size() blocks.
 
virtual double residual (FEVector &dst, const FEVectors &src, bool apply_boundaries=true)
 Loop over all cells to compute the residual.
 
unsigned int memory_consumption () const
 Compute the amount of memory needed by this object.
 
virtual double cell_estimate (const CellInfo &)
 Integration of estimate inside a cell.
 
virtual double bdry_estimate (const FaceInfo &)
 Integration of estimate on a boundary face.
 
virtual double face_estimate (const FaceInfo &, const FaceInfo &)
 Integration of estimate on an interior face.
 
virtual void cell_residual (FEVector &cell_vector, const CellInfo &cell)
 Integration of local linear form.
 
virtual void bdry_residual (FEVector &face_vector, const FaceInfo &face)
 Integration of local linear form.
 
virtual void face_residual (FEVector &face_vector1, FEVector &face_vector2, const FaceInfo &face1, const FaceInfo &face2)
 Integration of local linear form.
 
template<class BOX , class WORKER >
void integrate_1form (BOX &box, WORKER &local_forms, FEVector &dst, const FEVectors &src)
 The function integrating 1-forms as for instance used by residual().
 
void integrate_functional (LocalEstimate< dim > &local_forms, FEVectors &dst, const FEVectors &src)
 The function integrating functionals using local integrators.
 
void store_triangulation (Triangulation< dim > &new_tr)
 Function to copy a triangulation object for use after refinement.
 
void transfer_solution_to_coarse_mesh (Triangulation< dim > &tr_coarse, AppFrame::FEVector &coarse_solution, AppFrame::FEVector &refined_solution)
 Function to perform the transfer of a solution on a refined grid to the initial coarse grid.
 
void read_init_solution (AppFrame::FEVector &dst, bool &good_solution) const
 Function to transfer a solution on a refined grid to the initial coarse grid.
 
template<class INFO , typename TYPE >
void fill_local_data (const INFO &info, const std::vector< VectorSelector > &data_vector, std::vector< std::vector< std::vector< TYPE > > > &data) const
 Fill local data vector with values of the stored finite element function vectors.
 
virtual void grid_out (const std::string &basename)
 Output the grid used to solve the problem.
 
virtual void data_out (const std::string &basename, const FEVectors &src, const std::vector< std::string > &solution_names)
 
virtual void data_out (const std::string &basename, const FEVector &solution, const std::vector< std::string > &solution_names, const FEVector &postprocessing=FEVector(), const std::vector< std::string > &postprocessing_names=std::vector< std::string >())
 This function outputs the results of a computation.
 
- Public Member Functions inherited from AppFrame::ApplicationBase
 ApplicationBase (boost::shared_ptr< ApplicationData > ex_data=boost::shared_ptr< ApplicationData >())
 Constructor for an application.
 
 ApplicationBase (const ApplicationBase &other)
 Copy constructor.
 
virtual ~ApplicationBase ()
 Virtual destructor.
 
void print_parameters_to_file (ParameterHandler &param, const std::string &file_name, const ParameterHandler::OutputStyle &style)
 Print default parameters for the application to a file.
 
virtual void start_vector (FEVector &dst, std::string caller) const
 Initialize vector to problem size.
 
virtual void Tsolve (FEVector &start, const FEVectors &rhs)
 Solve the dual system assembled with right hand side rhs and return the result in start.
 
virtual void grid_out (const std::string &filename) const
 Write the mesh in the format specified by the ParameterHandler.
 
virtual void data_out (const std::string &filename, const FEVectors &src, const std::vector< std::string >)
 
boost::shared_ptr
< ApplicationData
get_data ()
 Get access to the protected variable data.
 
const boost::shared_ptr
< ApplicationData
get_data () const
 Get read-only access to the protected variable data.
 
virtual std::string id () const
 Return a unique identification string for this application.
 
virtual void notify (const Event &reason)
 Add a reason for assembling.
 
virtual void clear ()
 Reset the application class such that a call to initialize() is possible again and will produce the same result as if the object was fresh.
 
void clear_events ()
 Clear all notifications.
 

Protected Member Functions

For debugging purposes only
void print_cell_matrix (FullMatrix< double > matrix, std::string filename="cell_matrix.dat", bool exit=true) const
 Member function used to print a cell matrix into a file so that it can be studied.
 
void print_cell_matrix (ParameterHandler param, std::string filename="parameters_sample.prm", bool exit=true) const
 Testing file:
 
- Protected Member Functions inherited from AppFrame::OptimizationBlockMatrixApplication< dim >
void print_dresponses_dl (std::vector< std::vector< double > > pdf_pdl)
 Auxiliary routine to print the values df_dl This routine should be called once df_df is assembled.
 
void print_dresponses_du (std::vector< AppFrame::FEVector > df_du)
 Auxiliary routine to print the values of df_du This routine should be called once df_du is assembled.
 
void dfunction (std::vector< AppFrame::FEVector > &dst, const AppFrame::FEVectors &src, bool dfunctional_du, bool dresidual_dlambda)
 This is an auxiliary function called by dresidual_dlambda and dfunctional_du.
 

Protected Attributes

FuelCell::SystemManagement system_management
 Class used to store all information related to the system of equations that is being solved such as variable_names, equation_names and couplings between equations.
 
Pre-processor and operating condition classes
boost::shared_ptr
< FuelCellShop::Geometry::GridBase
< dim > > 
grid
 Create an object to generate the mesh.
 
FuelCell::OperatingConditions OC
 Operating conditions.
 
Material classes
FuelCellShop::Material::Water water
 The cathode contains water vapour, so we need to create an object water in order to compute viscosity, density, etc.
 
FuelCellShop::Material::Oxygen oxygen
 The cathode contains water vapour, so we need to create an object water in order to compute viscosity, density, etc.
 
FuelCellShop::Material::Nitrogen nitrogen
 The cathode contains water vapour, so we need to create an object water in order to compute viscosity, density, etc.
 
Layer classes
boost::shared_ptr
< FuelCellShop::Layer::GasDiffusionLayer
< dim > > 
CGDL
 The object GDL layer will contain all the information relevant to the the GDL.
 
boost::shared_ptr
< FuelCellShop::Layer::MicroPorousLayer
< dim > > 
CMPL
 The object CMPL layer will contain all the information relevant to the the micro-porous layer.
 
boost::shared_ptr
< FuelCellShop::Layer::CatalystLayer
< dim > > 
CCL
 The object CCL layer will contain all the information relevant to the the catalyst layer.
 
Physics base classes
FuelCellShop::Equation::ElectronTransportEquation
< dim
electron_transport_equation
 Equation class used to assemble the cell matrix for electron transport.
 
FuelCellShop::Equation::NewFicksTransportEquation
< dim
ficks_transport_equation
 Equation class used to assemble the cell matrix for Fick's mass transport.
 
FuelCellShop::Equation::ProtonTransportEquation
< dim
proton_transport_equation
 Equation class used to assemble the cell matrix for proton transport.
 
FuelCellShop::Equation::ReactionSourceTerms
< dim
reaction_source
 Equation class used to asemble the cell matrix and cell_residual for the reaction source terms.
 
Other internal data
std::vector< std::string > equation_names
 Structure where we store the problem we want to solve.
 
std::vector< std::string > component_names
 Structure where we store the name of each component in our problem.
 
std::vector< std::string > design_var
 Stores the design variable names so that the name can be appended to the .vtk file name.
 
std::vector< double > design_var_value
 Stores the values of the design variables so that the number can be appended to the .vtk file name.
 
- Protected Attributes inherited from AppFrame::OptimizationBlockMatrixApplication< dim >
unsigned int n_dvar
 Number of design variables.
 
unsigned int n_obj
 Number of objective functions.
 
unsigned int n_resp
 Number of responses = n_obj + n_con.
 
std::vector< std::string > name_design_var
 Member that stores the name of the design variables.
 
std::vector< std::string > name_responses
 Member that stores the name of the responses, i.e.
 
std::vector< std::string > name_output_var
 Member that stores the name of the output variables, These names will be written to name_responses if optimization is not being used and ignored otherwise.
 
bool output_coarse_solution
 Decision variable for whether the solution is to be output for transfer.
 
bool read_in_initial_solution
 Decision variable for whether the initial solution should be read from file.
 
bool optimization
 Decision variable for whether the application is being used for optimization.
 
- Protected Attributes inherited from AppFrame::BlockMatrixApplication< dim >
Table< 2, DoFTools::Coupling > cell_couplings
 Couplings through cell bilinear forms.
 
Table< 2, DoFTools::Coupling > flux_couplings
 Couplings through flux bilinear forms.
 
boost::shared_ptr< Quadrature
< dim > > 
quadrature_assemble_cell
 Quadrature rule for matrix assembling on cells.
 
boost::shared_ptr< Quadrature
< dim-1 > > 
quadrature_assemble_face
 Quadrature rule for matrix assembling on faces.
 
BlockSparseMatrix< double > matrix
 Storage for the actual matrix.
 
FilteredMatrix< BlockVector
< double > > 
filtered_matrix
 Wrapper for the actual matrix that contains all the boundary conditions.
 
- Protected Attributes inherited from AppFrame::DoFApplication< dim >
std::map< unsigned int, double > boundary_constraints
 List of all nodes constrained by a strong boundary condition, together with a value to be assigned.
 
ConstraintMatrix hanging_node_constraints
 Constraint Matrix object.
 
boost::shared_ptr< Mapping< dim > > mapping
 The mapping used for the transformation of mesh cells.
 
unsigned int mapping_degree
 Degree used for polynomial mapping of arbitrary order.
 
boost::shared_ptr
< Triangulation< dim > > 
tr
 Pointer to the Triangulation object.
 
boost::shared_ptr
< FiniteElement< dim > > 
element
 The finite element used in dof.
 
boost::shared_ptr< DoFHandler
< dim > > 
dof
 Pointer to the DoFHandler object.
 
boost::shared_ptr
< MGDoFHandler< dim > > 
mg_dof
 Pointer to the MGDoFHandler object.
 
MeshWorker::BlockInfo block_info
 
Vector< float > cell_errors
 The result of error estimation by cell.
 
Vector< float > face_errors
 The result of error estimation by face.
 
GridOut g_out
 The object for writing grids.
 
DataOut< dim, DoFHandler< dim > > d_out
 The object for writing data.
 
std::vector
< DataComponentInterpretation::DataComponentInterpretation > 
data_interpretation
 Vector for adjusting output to vector data.
 
std::string refinement
 Refinement parameter from parameter file.
 
unsigned int initial_refinement
 Initial refinement from parameter file.
 
double refinement_threshold
 Refinement threshold for adaptive method from parameter file.
 
double coarsening_threshold
 Coarsening threshold for adaptive method from parameter file.
 
bool sort_cuthill
 Flag for sorting with Cuthill McKee algorithm.
 
Point< dimsort_direction
 Direction for downstream sorting.
 
std::vector< FEVector * > transfer_vectors
 List of vector names to be transfered from one grid to the next.
 
Quadrature< dimquadrature_residual_cell
 Quadrature rule for residual computation on cells.
 
Quadrature< dim-1 > quadrature_residual_bdry
 Quadrature rule for residual computation on boundary faces.
 
Quadrature< dim-1 > quadrature_residual_face
 Quadrature rule for residual computation on faces.
 
bool boundary_fluxes
 Extend the integration loops in assemble() and residual() also to boundary faces.
 
bool interior_fluxes
 Extend the integration loops in assemble() and residual() also to interior faces.
 
- Protected Attributes inherited from AppFrame::ApplicationBase
boost::shared_ptr
< ApplicationData
data
 Object for auxiliary data.
 
Event notifications
 Accumulate reasons for reassembling here.
 

Additional Inherited Members

- Public Types inherited from AppFrame::DoFApplication< dim >
typedef
MeshWorker::InfoObjects::IntegrationInfo
< dim, FEValuesBase< dim > > 
CellInfo
 This is a redefinition of CellInfo in order to call MeshWorker objects instead of the previous implementation in AppFrame.
 
typedef
MeshWorker::InfoObjects::IntegrationInfo
< dim, FEFaceValuesBase< dim > > 
FaceInfo
 This is a redefinition of CellInfo in order to call MeshWorker objects instead of the previous implementation in AppFrame.
 
- Public Attributes inherited from AppFrame::DoFApplication< dim >
unsigned int verbosity
 Controls verbosity of certain functions.
 
std::vector
< DataComponentInterpretation::DataComponentInterpretation > 
solution_interpretations
 solution_interpretations identifies whether some solution_names are scalars or parts of a vector.
 
std::vector
< DataComponentInterpretation::DataComponentInterpretation > 
postprocessing_interpretations
 postprocessing_interpretations identifies whether some postprocessing_names are scalars or parts of a vector.
 
bool output_materials_and_levels
 output_materials_and_levels if true then visualized, otherwise suppressed.
 
Vector< double > output_materials
 Vector that will be used by data_out to store the material ids.
 
Vector< double > output_levels
 Vector that will be used by data_out to store the refinement levels at each cell.
 
- Static Public Attributes inherited from AppFrame::OptimizationBlockMatrixApplication< dim >
static const AppFrame::Event sensitivity_analysis
 The Event set by OptimizationMatrixApplication if if we want to compute the sensitivities and a new matrix should be assembled.
 
- Protected Types inherited from AppFrame::OptimizationBlockMatrixApplication< dim >
typedef void(* CELL_Dvalues )(std::vector< AppFrame::FEVector > &, const typename DoFApplication< dim >::CellInfo &, std::vector< std::vector< double > > &)
 Definition of a pointer to a function that is called in a loop over cells to compute the derivatives of either residual or responses with respect to design variables.
 

Detailed Description

template<int dim>
class FuelCell::Application::AppCathode< dim >

This class is used to solve a system of equations similar to the one presented in the journal article M.

Secanell et al. "Numerical Optimization of Proton Exchange Membrane Fuel Cell Cathode Electrodes", Electrochimica Acta, 52(7):2668-2682, February 2007 and also presented at the ECCM 2006 conference under the abstract titled M. Secanell et al, "A PEM fuel cell electrode model for gradient-based optimization", June 2006.

The cathode catalyst layer and the gas diffusion layer are modelled using a set of equations found in the journal article:

The simulation is isothermal, isotropic, models the electrochemical reactions using tafel kinetics, oxygen diffusion using Ficks Law and electronic and protonic transport using Ohms Law.

The equations solved are written as follows:

\[ R_1(\vec{u}) = \nabla \cdot (c_{total}D^{eff}_{O_2} \nabla x_{O_2} ) = \frac{1}{nF}\nabla\cdot i \]

\[ R_2(\vec{u}) = \nabla \cdot (\sigma^{eff}_{m}\nabla\phi_m) = \nabla \cdot i \]

\[ R_3(\vec{u}) = \nabla \cdot (\sigma^{eff}_{s}\nabla\phi_s) = -\nabla \cdot i \]

where, for the case of a macro-homogeneous catalyst layer,

\[ \nabla \cdot i = A_v i^{ref}_0 \left( \frac{c_{O_2}} {c^{ref}_{O_2}} \right) \mbox{exp} \left( \frac{\alpha F}{RT}(\phi_m - \phi_s) \right) \]

The solution variables, \( \vec{u} \) are the protonic potential, \(\phi_m\), the electronic potential, \(\phi_s\) and, instead of the oxygen concentration, we use the oxygen molar fraction, \(x_{O_2}\), that also accounts for the oxygen dissolving in the ionomer by using Henrys Law.

The govering equations above are nonlinear and therefore they cannot be solved directly.

In OpenFCST, we have decided to solve the system of equations using a nonlinear Newton solver. Therefore, instead of implemeting the equations above, cell_matrix and cell_residual implement a linearization of the governing equations above, i.e.

\[ \frac{\partial R_i(\vec{u^{n-1}})}{\partial u_j} u_j = R_i(\vec{u^{n-1}}) \]

For more information on the governing equations, discretization and solution methodology please read the following references:

Also, the folder /data/cathode/Secanell_EA07_Numerical_Optimization_PEMFC_Cathode_Electrodes contains a sample datafile to run this application with a macro-homogeneous model.

The application can be used to solve a cathode with and without an MPL. You can select the appropriate geometry by changing the following in the input file:

subsection Grid generation
set Type of mesh = Cathode OR CathodeMPL
subsection Internal mesh generator parameters
subsection Dimensions
set Cathode current collector width [cm] = 0.1 # Thickness of the rib of the bipolar plates (BPP)
set Cathode channel width [cm] = 0.1 # Thickness of the channel on the BPP
set Cathode CL thickness [cm] = 1E-3 # Thickness of the cathode catalyst layer [cm]
set Cathode GDL thickness [cm] = 2E-2 # Thickness of the cathode gas diffusion layer [cm]
set Cathode MPL thickness [cm] = 2E-3 # Thickness of the cathode microporous layer [cm]
end
end
end

The model can also be used to solve an agglomerate catalyst layer model. The governing equations are similar to the ones outlined above, however, the volumetric current density source, i.e. \( \nabla \cdot i \) is obtained as specified in the following article:

A sample dataset used to solve an agglomerate cathode catalyst layer model is given in data/cathode/Secanell_EA07_MultiVariable_Optimization_PEMFC_Cathodes_Agglomerate_Model

If using this function please cite the articles above as references.

Author
Marc Secanell, 2006-13

Constructor & Destructor Documentation

template<int dim>
FuelCell::Application::AppCathode< dim >::AppCathode ( boost::shared_ptr< AppFrame::ApplicationData data = boost::shared_ptr< AppFrame::ApplicationData >())

Constructor.

Note
the pointer data is initialized to boost::shared_ptr<> (), this means that the pointer is empty and when we do data.get() it will return 0. This is good because at ApplicationBase constructor an ApplicationData will be constructed.
template<int dim>
FuelCell::Application::AppCathode< dim >::~AppCathode ( )

Destructor.

Member Function Documentation

template<int dim>
void FuelCell::Application::AppCathode< dim >::_initialize ( ParameterHandler &  param)

Set up how many equations are needed and read in parameters for the parameter handler in order to initialize data.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::bdry_matrix ( MatrixVector face_matrices,
const typename DoFApplication< dim >::FaceInfo face 
)
virtual

Integration of local bilinear form.

Reimplemented from AppFrame::BlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::bdry_residual ( AppFrame::FEVector bdry_vector,
const typename DoFApplication< dim >::FaceInfo bdry 
)
virtual
template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::cell_dresidual_dlambda ( std::vector< AppFrame::FEVector > &  cell_vector,
const typename DoFApplication< dim >::CellInfo cell,
std::vector< std::vector< double > > &  src 
)
virtual

Integration of local bilinear form.

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::cell_dresponses_dl ( std::vector< std::vector< double > > &  cell_df_dl,
const typename DoFApplication< dim >::CellInfo info,
const AppFrame::FEVector sol 
)
virtual

This class is used to evaluate the derivative of all the functionals that require looping over cells with respect to the design variables.

This class is called by responses to evaluate the response at each cell.

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::cell_dresponses_du ( std::vector< AppFrame::FEVector > &  cell_df_du,
const typename DoFApplication< dim >::CellInfo info,
std::vector< std::vector< double > > &  src 
)
virtual

This class is used to evaluate the derivative of all the functionals that require looping over cells with respect of the unknowns of the system of governing equations.

This class is called by responses to evaluate the response at each cell.

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::cell_matrix ( MatrixVector cell_matrices,
const typename DoFApplication< dim >::CellInfo cell 
)
virtual

Integration of local bilinear form.

Here we loop over the quadrature points and over degrees of freedom in order to compute the matrix for the cell This routine depends on the problem at hand and is called by assemble() in DoF_Handler class The matrix to be assembled is:

\[ \begin{array}{l} M(i,j).block(0) = \int_{\Omega} a \nabla \phi_i \nabla \phi_j d\Omega + \int_{\Omega} \phi_i \frac{\partial f}{\partial u_0}\Big|_n \phi_j d\Omega \\ M(i,j).block(1) = \int_{\Omega} \phi_i \frac{\partial f}{\partial u_1}\Big|_n \phi_j d\Omega \\ M(i,j).block(2) = \int_{\Omega} \phi_i \frac{\partial f}{\partial u_2}\Big|_n \phi_j d\Omega \end{array} \]

Reimplemented from AppFrame::BlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::cell_residual ( AppFrame::FEVector cell_vector,
const typename DoFApplication< dim >::CellInfo cell 
)
virtual

Integration of the rhs of the equations.

Here we loop over the quadrature points and over degrees of freedom in order to compute the right hand side for each cell This routine depends on the problem at hand and is called by residual() in DoF_Handler class

Note
This function is called residual because in the case of nonlinear systems the rhs is equivalent to the residual
template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::cell_responses ( std::vector< double > &  resp,
const typename DoFApplication< dim >::CellInfo info,
const AppFrame::FEVector sol 
)
virtual

Compute the value of all objective function and constraints.

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::check_responses ( )
virtual

This class is called by responses to make sure that all responses requested are implemented in either cell_responses, global_responses or face_responses.

Note
Every time we add a new response that we can calculate we need to update this file.

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::data_out ( const std::string &  basename,
const AppFrame::FEVectors src 
)
virtual

Reimplementation of the routine in the base class BaseApplication in namespace AppFrame so that the right labels are outputed and so that I can compute and output the source term.

Reimplemented from AppFrame::DoFApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::declare_parameters ( ParameterHandler &  param)
virtual

Declare all parameters that are needed for:

  • the computation of the equation coefficients
  • the control of the linear system solution
  • ...

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::dirichlet_bc ( std::map< unsigned int, double > &  boundary_values) const
virtual

Member function used to set dirichlet boundary conditions.

This function is application specific and it only computes the boundary_value values that are used to constraint the linear system of equations that is being solved

Reimplemented from AppFrame::BlockMatrixApplication< dim >.

template<int dim>
virtual double FuelCell::Application::AppCathode< dim >::estimate ( const AppFrame::FEVectors sol)
virtual

Estimate error per cell.

Reimplemented from AppFrame::DoFApplication< dim >.

template<int dim>
virtual double FuelCell::Application::AppCathode< dim >::evaluate ( const AppFrame::FEVectors src)
virtual

Post-processing.

Evaluate a functional such as the objective function of an optimization problem

Reimplemented from AppFrame::DoFApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::global_dresponses_dl ( std::vector< std::vector< double > > &  df_dl,
const AppFrame::FEVector sol 
)
virtual

This class is used to evaluate the sensitivities of all responses that do not require looping over cells with respect of the design variables.

An example of one of this types of constraints is the solid volume fraction.

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::global_dresponses_du ( std::vector< AppFrame::FEVector > &  df_du,
const AppFrame::FEVector src 
)
virtual

This class is used to evaluate the sensitivities of all responses that do not require looping over cells with respecto of the unknowns of the system of governing equations.

An example of one of this types of constraints is the solid volume fraction.

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::global_responses ( std::vector< double > &  resp,
const AppFrame::FEVector sol 
)
virtual

This class is used to evaluate all responses that do not require looping over cells.

An example of one of this types of constraints is the solid volume fraction.

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
void FuelCell::Application::AppCathode< dim >::init_solution ( AppFrame::FEVector src)
virtual

Initialize nonlinear solution.

Note
This routine is not working well yet. In fact, I don't know how to do it...

Implements AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::initialize ( ParameterHandler &  param)
virtual

Call the other initialize routines from the inherited classes.

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
void FuelCell::Application::AppCathode< dim >::initialize_grid ( ParameterHandler &  param)

Initilize and create the grid Determines whether reading from a grid file or creating new.

template<int dim>
void FuelCell::Application::AppCathode< dim >::parse_init_solution ( AppFrame::FEVector src,
FuelCell::OperatingConditions OC,
bool &  good_solution 
)

Check to see that the solution read from file can be used for the new simulation conditions.

Uses member functions to compare the initial conditions in the parameter file to those in the loaded solution.

template<int dim>
bool FuelCell::Application::AppCathode< dim >::parse_init_solution_phi_s ( AppFrame::FEVector src,
double &  V_cell 
)

Check to see that the initial solution for the solid phase potential is usable.

If the variable varies by more than 30% at the boundary, the solution is considered unusable, and a default solution will be generated.

template<int dim>
bool FuelCell::Application::AppCathode< dim >::parse_init_solution_x_O2 ( AppFrame::FEVector src,
double &  x_O2 
)

Check to see that the initial solution for the oxygen mole fraction is usable.

If the variable varies by more than 30% at the boundary, the solution is considered unusable, and a default solution will be generated.

template<int dim>
void FuelCell::Application::AppCathode< dim >::print_cell_matrix ( FullMatrix< double >  matrix,
std::string  filename = "cell_matrix.dat",
bool  exit = true 
) const
inlineprotected

Member function used to print a cell matrix into a file so that it can be studied.

This routine can be called as follows:

print_cell_matrix(cell_matrices[0].matrix, "marc.dat");

This will store the data in the first block in the MatrixVector in a file named marc.dat.

template<int dim>
void FuelCell::Application::AppCathode< dim >::print_cell_matrix ( ParameterHandler  param,
std::string  filename = "parameters_sample.prm",
bool  exit = true 
) const
inlineprotected

Testing file:

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::set_parameters ( const std::vector< std::string > &  name_dvar,
const std::vector< double > &  value_dvar,
ParameterHandler &  param 
)
virtual

Function called by optimization loop in order to set the values in the ParameterHandler to the new design parameters.

Since ParameterHandler depends on the problem we are solving, set_parameters() is set at the most inner loop of the application.

Reimplemented from AppFrame::OptimizationBlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::Application::AppCathode< dim >::solve ( AppFrame::FEVector start,
const AppFrame::FEVectors rhs 
)
virtual

Solve the linear system of equations.

Reimplemented from AppFrame::ApplicationBase.

Member Data Documentation

template<int dim>
boost::shared_ptr<FuelCellShop::Layer::CatalystLayer<dim> > FuelCell::Application::AppCathode< dim >::CCL
protected

The object CCL layer will contain all the information relevant to the the catalyst layer.

We can request any effective property from this class

template<int dim>
boost::shared_ptr<FuelCellShop::Layer::GasDiffusionLayer<dim> > FuelCell::Application::AppCathode< dim >::CGDL
protected

The object GDL layer will contain all the information relevant to the the GDL.

We can request any effective property from this class

template<int dim>
boost::shared_ptr<FuelCellShop::Layer::MicroPorousLayer<dim> > FuelCell::Application::AppCathode< dim >::CMPL
protected

The object CMPL layer will contain all the information relevant to the the micro-porous layer.

We can request any effective property from this class

template<int dim>
std::vector<std::string> FuelCell::Application::AppCathode< dim >::component_names
protected

Structure where we store the name of each component in our problem.

The component names are stored in the same way as they are stored in the solution.

template<int dim>
std::vector<std::string> FuelCell::Application::AppCathode< dim >::design_var
protected

Stores the design variable names so that the name can be appended to the .vtk file name.

template<int dim>
std::vector<double> FuelCell::Application::AppCathode< dim >::design_var_value
protected

Stores the values of the design variables so that the number can be appended to the .vtk file name.

template<int dim>
FuelCellShop::Equation::ElectronTransportEquation<dim> FuelCell::Application::AppCathode< dim >::electron_transport_equation
protected

Equation class used to assemble the cell matrix for electron transport.

template<int dim>
std::vector<std::string> FuelCell::Application::AppCathode< dim >::equation_names
protected

Structure where we store the problem we want to solve.

Each vector component contains a string with the name of the equation we want to solve Then, the number of components is equation_names.size()

template<int dim>
FuelCellShop::Equation::NewFicksTransportEquation<dim> FuelCell::Application::AppCathode< dim >::ficks_transport_equation
protected

Equation class used to assemble the cell matrix for Fick's mass transport.

Note
The order in the .h file is important here. Objects that need to be used in the initialization of classes need to be declared before this class. In this case FicksTransportEquation takes in two material classes. Therefore, they need to be defined before FicksTransportEquation.
template<int dim>
boost::shared_ptr< FuelCellShop::Geometry::GridBase<dim> > FuelCell::Application::AppCathode< dim >::grid
protected

Create an object to generate the mesh.

The cathode contains water vapour, so we need to create an object water in order to compute viscosity, density, etc.

for water

Operating conditions.

The cathode contains water vapour, so we need to create an object water in order to compute viscosity, density, etc.

for water

template<int dim>
FuelCellShop::Equation::ProtonTransportEquation<dim> FuelCell::Application::AppCathode< dim >::proton_transport_equation
protected

Equation class used to assemble the cell matrix for proton transport.

template<int dim>
FuelCellShop::Equation::ReactionSourceTerms< dim > FuelCell::Application::AppCathode< dim >::reaction_source
protected

Equation class used to asemble the cell matrix and cell_residual for the reaction source terms.

template<int dim>
FuelCell::SystemManagement FuelCell::Application::AppCathode< dim >::system_management
protected

Class used to store all information related to the system of equations that is being solved such as variable_names, equation_names and couplings between equations.

The cathode contains water vapour, so we need to create an object water in order to compute viscosity, density, etc.

for water


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