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

Base class for all linear applications, i.e., all applications requiring Triangulation and DoFHandler classes. More...

#include <dof_application.h>

Inheritance diagram for FuelCell::ApplicationCore::DoFApplication< dim >:
Inheritance graph
[legend]
Collaboration diagram for FuelCell::ApplicationCore::DoFApplication< dim >:
Collaboration graph
[legend]

Public Types

typedef
FuelCell::ApplicationCore::IntegrationInfo
< dim, FEValuesBase< dim > > 
CellInfo
 Shortcut. More...
 
typedef
FuelCell::ApplicationCore::IntegrationInfo
< dim, FEFaceValuesBase< dim > > 
FaceInfo
 Shortcut. More...
 

Public Member Functions

Constructor, destructor and initialization:
 DoFApplication (boost::shared_ptr< ApplicationData > data=boost::shared_ptr< ApplicationData >())
 Constructor for an object owning its own mesh and dof handler. More...
 
 DoFApplication ()
 Constructor for an object owning its own mesh and dof handler and creating new ApplicationData. More...
 
 DoFApplication (DoFApplication< dim > &dof_app, bool triangulation_only)
 Constructor for an object, borrowing mesh and dof handler from another object. More...
 
 ~DoFApplication ()
 Destructor which deletes owned objects. More...
 
virtual void declare_parameters (ParameterHandler &param)
 Declare parameters related to mesh generation and finite elements. More...
 
virtual void initialize (ParameterHandler &param)
 Initialize application data based on #param object. More...
 
virtual void remesh_dofs ()
 Initialize dof handler, count the dofs in each block and renumber the dofs. More...
 
virtual void remesh ()
 Refine grid accordingly to the Refinement options under "Grid Generation" in the parameter file. More...
 
virtual void init_vector (FEVector &dst) const
 Reinitialize the BlockVector dst such that it contains block_size.size() blocks. More...
 
virtual void initialize_solution (FEVector &initial_guess, std::shared_ptr< Function< dim > > initial_function=std::shared_ptr< Function< dim > >())
 For nonlinear applications, an intial solution is needed in order to assemble the residual and the global matrix. More...
 
Other
virtual double estimate (const FEVectors &src)
 Estimate the error. More...
 
virtual double evaluate (const FEVectors &src)
 Evaluate a functional during postprocessing such as drag or lift on an aerodynamics problem. More...
 
virtual double residual (FEVector &dst, const FEVectors &src, bool apply_boundaries=true)
 Loop over all cells to compute the residual. More...
 
virtual void residual_constraints (FEVector &dst) const
 Apply boundary conditions and hanging node constraints to a residual vector after it has been computed by residual(). More...
 
void store_triangulation (Triangulation< dim > &new_tr)
 Function to copy a triangulation object for use after refinement. More...
 
void add_vector_for_transfer (FEVector *src)
 Add the vector to be transfered from one mesh to the next. More...
 
void delete_vector_for_transfer ()
 Delete the vector to be transfered from one mesh to the next. More...
 
void transfer_solution_to_coarse_mesh (Triangulation< dim > &tr_coarse, FEVector &coarse_solution, FEVector &refined_solution)
 Function to perform the transfer of a solution on a refined grid to the initial coarse grid. More...
 
unsigned int memory_consumption () const
 Compute the amount of memory needed by this object. More...
 
Output member functions:
virtual void grid_out (const std::string &basename)
 Output the grid used to solve the problem. More...
 
virtual void data_out (const std::string &basename, const FEVectors &src)
 Write data in the format specified by the ParameterHandler. More...
 
void print (const std::string &basename, const FEVector &src, const std::vector< unsigned int > &src_indices=std::vector< unsigned int >()) const
 This function prints FEVector src to a text file basename. More...
 
- Public Member Functions inherited from FuelCell::ApplicationCore::ApplicationBase
 ApplicationBase (boost::shared_ptr< ApplicationData > data=boost::shared_ptr< ApplicationData >())
 Constructor for an application. More...
 
 ApplicationBase (const ApplicationBase &other)
 Copy constructor. More...
 
virtual ~ApplicationBase ()
 Virtual destructor. More...
 
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. More...
 
virtual void start_vector (FEVector &dst, std::string) const
 Initialize vector to problem size. More...
 
virtual void solve (FEVector &dst, const FEVectors &src)=0
 Solve the system assembled with right hand side in FEVectors src and return the result in FEVector dst. More...
 
virtual void Tsolve (FEVector &, const FEVectors &)
 Solve the dual system assembled with right hand side rhs and return the result in start. More...
 
boost::shared_ptr
< ApplicationData
get_data ()
 Get access to the protected variable data. More...
 
const boost::shared_ptr
< ApplicationData
get_data () const
 Get read-only access to the protected variable data. More...
 
virtual std::string id () const
 Return a unique identification string for this application. More...
 
virtual void notify (const Event &reason)
 Add a reason for assembling. More...
 
virtual void clear ()
 All true in notifications. More...
 
virtual void clear_events ()
 All false in notifications. More...
 
virtual unsigned int get_solution_index ()
 Returns solution index. More...
 

Public Attributes

boost::shared_ptr< Boundary
< dim > > 
curved_boundary
 Curved boundary. More...
 
types::boundary_id curved_bdry_id
 Curved boundary ID. More...
 
Initial solution management
std::string filename_initial_sol
 Filename where to output the initial grid. More...
 
bool output_initial_sol
 Flag to output the initial solution used to start the solving process. More...
 
bool read_in_initial_solution
 Bool flag used to specify if the initial solution to the problem, specially important for non-linear problems, should be read from file. More...
 
bool use_predefined_solution
 Use user pre-defined initial solution. More...
 
bool output_coarse_solution
 Bool flag used to specify if the final solution should be stored in the coarse mesh in order to be used later as an initial solution to solve another problem using the flag read_in_initial_solution. More...
 
Output data:
std::vector
< DataComponentInterpretation::DataComponentInterpretation > 
solution_interpretations
 solution_interpretations identifies whether some solution_names are scalars or parts of a vector. More...
 
std::vector
< DataComponentInterpretation::DataComponentInterpretation > 
postprocessing_interpretations
 postprocessing_interpretations identifies whether some postprocessing_names are scalars or parts of a vector. More...
 
std::vector
< DataComponentInterpretation::DataComponentInterpretation > 
data_interpretation
 
bool output_materials_and_levels
 output_materials_and_levels if true then visualized, otherwise suppressed. More...
 
Vector< double > output_materials
 Vector that will be used by data_out to store the material ids. More...
 
Vector< double > output_levels
 Vector that will be used by data_out to store the refinement levels at each cell. More...
 
bool output_actual_degree
 true, if you want to output solution and postprocessing data using actual finite element fields Q_n with n >= 1. More...
 
bool print_solution
 true, if you want to print FEVector solution to a text file. More...
 
bool print_postprocessing
 true, if you want to print FEVector postprocessing to a text file. More...
 
std::vector< unsigned int > solution_printing_indices
 The indices of the FEVector solution to be printed to a text file. More...
 
std::vector< unsigned int > postprocessing_printing_indices
 The indices of the FEVector postprocessing to be printed to a text file. More...
 
bool print_blocks_instead_of_indices
 true, if the whole blocks of FEVector solution or FEVector postprocessing to be printed to a text file instead of separate indices. More...
 
bool output_matrices_and_rhs
 If true, all cell matrices and right hand sides will be output. More...
 

Protected Member Functions

virtual void data_out (const std::string &basename, const FEVector &solution, const std::vector< std::string > &solution_names, const std::vector< DataPostprocessor< dim > * > &PostProcessing)
 This routine is used to write data in the format specified by the ParameterHandler. More...
 
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. More...
 
void constrain_boundary (FEVector &v, bool homogeneous) const
 Apply either homogeneous or inhomogeneous boundary_constraints. More...
 
Initialization of application and residual
void _initialize (ParameterHandler &param)
 Initialize from parameter values. More...
 
virtual void initialize_triangulation (ParameterHandler &param)
 Function used to read in a mesh and hand it over to the boost::shared_ptr<Triangulation<dim> > tr object. More...
 
void read_init_solution (FEVector &dst, bool &good_solution) const
 Create a mesh and assign it to object tr. More...
 
Local assembly routines
virtual void cell_residual (FEVector &cell_vector, const CellInfo &cell)
 Local integration. More...
 
virtual void bdry_residual (FEVector &face_vector, const FaceInfo &face)
 Local integration. More...
 
virtual void face_residual (FEVector &face_vector1, FEVector &face_vector2, const FaceInfo &face1, const FaceInfo &face2)
 Local integration. More...
 
virtual double cell_estimate (const CellInfo &src)
 Local estimation. More...
 
virtual double bdry_estimate (const FaceInfo &src)
 Local estimation. More...
 
virtual double face_estimate (const FaceInfo &src1, const FaceInfo &src2)
 Local estimation. More...
 
- Protected Member Functions inherited from FuelCell::ApplicationCore::ApplicationBase
void print_caller_name (const std::string &caller_name) const
 Print caller name. More...
 

Protected Attributes

GridOut g_out
 The object for writing grids. More...
 
DataOut< dim, DoFHandler< dim > > d_out
 The object for writing data. More...
 
FuelCell::SystemManagement system_management
 This object knows everything about FCST equations, variables, couplings, etc. More...
 
std::map< unsigned int, double > boundary_constraints
 List of all nodes constrained by a strong boundary condition, together with a value to be assigned. More...
 
ConstraintMatrix hanging_node_constraints
 Constraint Matrix object. More...
 
boost::shared_ptr< Mapping< dim > > mapping
 The mapping used for the transformation of mesh cells. More...
 
unsigned int mapping_degree
 Degree used for polynomial mapping of arbitrary order. More...
 
boost::shared_ptr
< FiniteElement< dim > > 
element
 The finite element used in dof. More...
 
boost::shared_ptr< DoFHandler
< dim > > 
dof
 Pointer to the DoFHandler object. More...
 
BlockInfo block_info
 
Vector< float > cell_errors
 The result of error estimation by cell. More...
 
Vector< float > face_errors
 The result of error estimation by face. More...
 
std::string refinement
 Refinement parameter from parameter file. More...
 
unsigned int initial_refinement
 Initial refinement from parameter file. More...
 
double refinement_threshold
 Refinement threshold for adaptive method from parameter file. More...
 
double coarsening_threshold
 Coarsening threshold for adaptive method from parameter file. More...
 
bool sort_cuthill
 Flag for sorting with Cuthill McKee algorithm. More...
 
Point< dimsort_direction
 Direction for downstream sorting. More...
 
std::vector< FEVector * > transfer_vectors
 List of vector names to be transfered from one grid to the next. More...
 
Quadrature< dimquadrature_residual_cell
 Quadrature rule for residual computation on cells. More...
 
Quadrature< dim-1 > quadrature_residual_bdry
 Quadrature rule for residual computation on boundary faces. More...
 
Quadrature< dim-1 > quadrature_residual_face
 Quadrature rule for residual computation on faces. More...
 
Table< 2, DoFTools::Coupling > cell_couplings
 Couplings through cell bilinear forms. More...
 
Table< 2, DoFTools::Coupling > flux_couplings
 Couplings through flux bilinear forms. More...
 
bool boundary_fluxes
 Extend the integration loops in assemble() and residual() also to boundary faces. More...
 
bool interior_fluxes
 Extend the integration loops in assemble() and residual() also to interior faces. More...
 
unsigned int verbosity
 Controls verbosity of certain functions. More...
 
Initial and boundary data information:
std::vector
< component_materialID_value_map
component_materialID_value_maps
 Each entry of this std::vector reflects the following structure (see FuelCell::InitialAndBoundaryData namespace docs): More...
 
std::vector
< component_boundaryID_value_map
component_boundaryID_value_maps
 Each entry of this std::vector reflects the following structure (see FuelCell::InitialAndBoundaryData namespace docs): More...
 
Pre-processor objects
boost::shared_ptr
< FuelCellShop::Geometry::GridBase
< dim > > 
mesh_generator
 Grid. More...
 
boost::shared_ptr
< Triangulation< dim > > 
tr
 Pointer to the Triangulation object. More...
 
- Protected Attributes inherited from FuelCell::ApplicationCore::ApplicationBase
boost::shared_ptr
< ApplicationData
data
 Object for auxiliary data. More...
 
Event notifications
 Accumulate reasons for assembling here. More...
 

Private Member Functions

virtual void sort_dofs (DoFHandler< dim > *dof_handler) const
 Sort the degrees of freedom. More...
 
virtual void distribute_face_to_cell_errors ()
 After computing the error contributions on faces and cells, this function adds half of the contribution computed by face_estimate() and the contribution computed by bdry_estimate() to the adjacent cell. More...
 
virtual double global_from_local_errors () const
 Compute a global estimation value from cell_errors. More...
 

Private Attributes

unsigned int n_ref
 Number of refinements. More...
 

Detailed Description

template<int dim>
class FuelCell::ApplicationCore::DoFApplication< dim >

Base class for all linear applications, i.e., all applications requiring Triangulation and DoFHandler classes.

The mesh as well as the dof handler may be created by this class, which is the default, or they may be provided by another object, in which case they must be specified in the constructor.

Note that in this class the received or created dof handler is associated to the finite element given by the argument "Element" on the parameter file. Therefore, this class in not responsible to generate the system of equations to be solved, only to initialize the dof handler.

Element can either be a single element or a FESystem. In the latter case, the nomenclature used in the paramter file is:

set Element = FESystem[element1_type(element1_degree)^number_of_elements1-...-elementN_type(elementN_degree)^number_of_elementsN] Example: set Element = FESystem[FE_DGQ(0)-FE_Q(1)^2]

TODO: Improve Parallelization of FEVectors. Currently all FEVectors are stored in each MPI process thereby making the calculations expensive. All FEVectors should be PETScWrappers::BlockVector so that only a part of the solution is stored. Once we do this data_out and transfer_solution_to_coarse_mesh will have to change as solution will not store everything. (Priority: low)

Author
Guido Kanschat
Valentin N. Zingan
Marc Secanell
Peter Dobson

Member Typedef Documentation

Shortcut.

Shortcut.

Constructor & Destructor Documentation

template<int dim>
FuelCell::ApplicationCore::DoFApplication< dim >::DoFApplication ( boost::shared_ptr< ApplicationData data = boost::shared_ptr< ApplicationData >())

Constructor for an object owning its own mesh and dof handler.

If data is null, a new dof handler for application data is generated.

Constructor for an object owning its own mesh and dof handler and creating new ApplicationData.

template<int dim>
FuelCell::ApplicationCore::DoFApplication< dim >::DoFApplication ( DoFApplication< dim > &  dof_app,
bool  triangulation_only 
)

Constructor for an object, borrowing mesh and dof handler from another object.

If triangulation_only is true, only the triangulation is borrowed.

Note
A pointer to DH2 must be convertible to a pointer to DH.

Destructor which deletes owned objects.

Member Function Documentation

template<int dim>
void FuelCell::ApplicationCore::DoFApplication< dim >::_initialize ( ParameterHandler &  param)
protected

Initialize from parameter values.

Called by initialize.

Deprecated:
template<int dim>
void FuelCell::ApplicationCore::DoFApplication< dim >::add_vector_for_transfer ( FEVector src)

Add the vector to be transfered from one mesh to the next.

template<int dim>
virtual double FuelCell::ApplicationCore::DoFApplication< dim >::bdry_estimate ( const FaceInfo src)
protectedvirtual

Local estimation.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::bdry_residual ( FEVector face_vector,
const FaceInfo face 
)
protectedvirtual

Local integration.

template<int dim>
virtual double FuelCell::ApplicationCore::DoFApplication< dim >::cell_estimate ( const CellInfo src)
protectedvirtual

Local estimation.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::cell_residual ( FEVector cell_vector,
const CellInfo cell 
)
protectedvirtual

Local integration.

template<int dim>
void FuelCell::ApplicationCore::DoFApplication< dim >::constrain_boundary ( FEVector v,
bool  homogeneous 
) const
protected

Apply either homogeneous or inhomogeneous boundary_constraints.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::data_out ( const std::string &  basename,
const FEVectors src 
)
virtual
template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::data_out ( const std::string &  basename,
const FEVector solution,
const std::vector< std::string > &  solution_names,
const std::vector< DataPostprocessor< dim > * > &  PostProcessing 
)
protectedvirtual

This routine is used to write data in the format specified by the ParameterHandler.

This routine is usually called inside child applications wheh data_out is specified.

Usually, the solution if first processed and this routine is called last (see Usage below).

Parameters:

  • basename is a base of an output file name,
  • solution is the global solution computed in DoFs,
  • solution_names is the component names of the solution,
  • PostProcessing is a vector of DataPostprocessor objects used to evaluate other values using the solution vector. Many DataPostprocessor objects have already been developed.

Usage

The code below, from AppCathode, highlights how this routine is usually used inside an application. First, the filename, solution vector and solution_name vectors are obtained. Then, a vector of DataPostprocessor objects is used to compute any necessary additional post-processing results such as current density, relative humidity or others. Finally, this data_out routine is called:

* template<int dim>
* void NAME::AppCathode<dim>::data_out(const std::string& filename, const FuelCell::ApplicationCore::FEVectors& src){
* // --- Find solution ---
* FuelCell::ApplicationCore::FEVector solution = src.vector( src.find_vector("Solution") );
* // --- Assign solution names ---
* std::vector<std::string> solution_names;
* solution_names.push_back("oxygen_molar_fraction");
* solution_names.push_back("protonic_electrical_potential");
* solution_names.push_back("electronic_electrical_potential");
* // --- Assign solution interpretations ---
* this->solution_interpretations.clear();
* this->solution_interpretations.resize(this->element->n_blocks(), DataComponentInterpretation::component_is_scalar);
* // --- Create vector of PostProcessing objects ---
* std::vector< DataPostprocessor<dim>* > PostProcessing;
* // --- current ---
* PostProcessing.push_back(&current);
* // --- output ---
* solution,
* solution_names,
* PostProcessing);
* }
*
template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::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 >() 
)
protectedvirtual

This function outputs the results of a computation.

Parameters:

  • basename is a base of an output file name,
  • solution is the global solution computed in DoFs,
  • solution_names is the component names of the solution,
  • postprocessing is what you compute in DoFs (or in cells in the case of piece-wise constant functionals) after the solution has been computed,
  • postprocessing_names is the component names of the postprocessing.
template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::declare_parameters ( ParameterHandler &  param)
virtual

Declare parameters related to mesh generation and finite elements.

Most importantly, here we declare "Element". "Element" is declared under the subsection "Discretization". "Element" is the element that will be associated to the dof handler (by default FE_Q(1)). This can either be a single element of a FESystem. In the latter case, the nomenclature used in the paramter file is: set Element = FESystem[element1_type(element1_degree)^number_of_elements1- ...-elementN_type(elementN_degree)^number_of_elementsN] Example: set Element = FESystem[FE_DGQ(0)-FE_Q(1)^2]

Note
The function declare_parameters of derived applications should always call the one of its base class.

Reimplemented from FuelCell::ApplicationCore::ApplicationBase.

Reimplemented in FuelCell::Application::AppIncompressibleFlows< dim >, FuelCell::Application::AppCathode< dim >, FuelCell::ApplicationCore::BlockMatrixApplication< dim >, FuelCell::Application::AppPemfcNIThermal< dim >, FuelCell::Application::AppPemfcTPSaturation< dim >, FuelCell::Application::AppPemfc< dim >, FuelCell::Application::AppPemfcCapillaryTwoPhaseNIThermal< dim >, FuelCell::Application::AppReadMesh< dim >, FuelCell::Application::AppDiffusion< dim >, FuelCell::Application::AppThermalTest< dim >, FuelCell::Application::AppOhmic< dim >, FuelCell::ApplicationCore::OptimizationBlockMatrixApplication< dim >, FuelCell::Application::Thermaltesting< dim >, FuelCell::Application::AppCathodeKG< dim >, FuelCell::Application::AppCompressibleFlows< dim >, and FuelCell::Application::CapillaryTesting< dim >.

template<int dim>
void FuelCell::ApplicationCore::DoFApplication< dim >::delete_vector_for_transfer ( )

Delete the vector to be transfered from one mesh to the next.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::distribute_face_to_cell_errors ( )
privatevirtual

After computing the error contributions on faces and cells, this function adds half of the contribution computed by face_estimate() and the contribution computed by bdry_estimate() to the adjacent cell.

It deletes face_errors in the end.

This function can be overridden by derived classes to handle both estimates separately.

template<int dim>
virtual double FuelCell::ApplicationCore::DoFApplication< dim >::estimate ( const FEVectors src)
virtual

Estimate the error.

By default, the KellyErrorEstimator is used in order to estimate the error in every cell. The error for all components of the solution is added.

In general, here a loop over all cells and faces, using the virtual local functions

Reimplemented from FuelCell::ApplicationCore::ApplicationBase.

Reimplemented in FuelCell::Application::AppThermalTest< dim >.

template<int dim>
virtual double FuelCell::ApplicationCore::DoFApplication< dim >::evaluate ( const FEVectors src)
virtual
template<int dim>
virtual double FuelCell::ApplicationCore::DoFApplication< dim >::face_estimate ( const FaceInfo src1,
const FaceInfo src2 
)
protectedvirtual

Local estimation.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::face_residual ( FEVector face_vector1,
FEVector face_vector2,
const FaceInfo face1,
const FaceInfo face2 
)
protectedvirtual

Local integration.

template<int dim>
virtual double FuelCell::ApplicationCore::DoFApplication< dim >::global_from_local_errors ( ) const
privatevirtual

Compute a global estimation value from cell_errors.

In order for this to be useful, face_errors should be transfered to cells first, for instance by distribute_face_to_cell_errors().

This function computes the sum of the cell_errors and takes the square root in the end. Reimplementation in derived classes allows for different norms.

The function is called by estimate() to compute its return value.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::grid_out ( const std::string &  basename)
virtual

Output the grid used to solve the problem.

Reimplemented from FuelCell::ApplicationCore::ApplicationBase.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::init_vector ( FEVector dst) const
virtual

Reinitialize the BlockVector dst such that it contains block_size.size() blocks.

Each block is reinitialized to dimensions block_sizes[i]. Note that block_sizes contains the number of degrees of freedom per block. So, init_vector could be used to reinitialize the solution vector.

Note
: init_vector can only be called after remesh_dofs has been called otherwise block_sizes will not have been initialized.

Reimplemented from FuelCell::ApplicationCore::ApplicationBase.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::initialize ( ParameterHandler &  param)
virtual
template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::initialize_solution ( FEVector initial_guess,
std::shared_ptr< Function< dim > >  initial_function = std::shared_ptr< Function< dim > >() 
)
virtual

For nonlinear applications, an intial solution is needed in order to assemble the residual and the global matrix.

This member function is used to either generate a solution based on initial data or to modify an existing solution in order to meet the boundary conditions imposed in the initial data file.

There are three possibilities,

  • The solution is read from a file and boundary conditions are then applied to it.
  • If solution is not read from a file and initial_function is specified, the function that is passed to initialize_solution is used to generate the initial solution and then the boundaries are applied.
  • If initial_function is not specified, then the information in the input file is used to generate a piecewise initial solution, with a constant value per material_ID.

In order to correct for boundary conditions and compute the initial solution, DoFApplication will need component_boundaryID_value_maps and component_boundaryID_value_maps respectively to be initialized in the initialization function of the application using the equation classes used in the application, for example

* // Now, initialize object that are used to setup initial solution and boundary conditions:
* component_materialID_value_maps.push_back( ficks_transport_equation.get_component_materialID_value() );
* component_materialID_value_maps.push_back( electron_transport_equation.get_component_materialID_value() );
* component_materialID_value_maps.push_back( proton_transport_equation.get_component_materialID_value() );
*
* component_boundaryID_value_maps.push_back( ficks_transport_equation.get_component_boundaryID_value() );
* component_boundaryID_value_maps.push_back( electron_transport_equation.get_component_boundaryID_value() );
* component_boundaryID_value_maps.push_back( proton_transport_equation.get_component_boundaryID_value() );
*

The value of the initial solution in each material id and each boundary id is given in the parameter file in sections like the following for, in this case, equation #FicksTransportEquation:

* subsection Equations
* subsection Ficks Transport Equation - oxygen
* subsection Initial data
* set oxygen_molar_fraction = 2: 1.0, 4: 1.0
* end
* subsection Boundary data
* set oxygen_molar_fraction = 3: 1.0
* end
* end
* end
*
*

The initial solution could also be read from a previous simulation. In order to do so, the previous simulation would have been called with the following flag in the input file: *

* subsection Adaptive refinement
* set Output solution for transfer = true
* end
*

Then, a hidden file .transfer_solution.FEVector is created and loaded afterward by setting the flag

* subsection Adaptive refinement
* set Read in initial solution from file = true
* end
*
Note
: initialize_solution can only be called after remesh_dofs has been called otherwise block_sizes will not have been initialized.

Usage

This class can be used directly or, if you would like to specify an initial solution function, you would need to reimplement the function in the child class and call it as follows:

* template <int dim>
* void NAME::AppPemfc<dim>::init_solution(FuelCell::ApplicationCore::FEVector& initial_guess) {
* std::shared_ptr< Function<dim> > initial_solution (new FuelCell::InitialSolution::AppPemfcIC<dim> (&OC, grid));
* DoFApplication<dim>::initialize_solution(initial_guess, initial_solution);
* }
*

For an example see .

Reimplemented in FuelCell::Application::AppIncompressibleFlows< dim >, FuelCell::Application::AppCathode< dim >, FuelCell::Application::AppPemfcNIThermal< dim >, FuelCell::Application::AppPemfcTPSaturation< dim >, FuelCell::Application::AppPemfc< dim >, FuelCell::Application::AppPemfcCapillaryTwoPhaseNIThermal< dim >, FuelCell::Application::AppReadMesh< dim >, FuelCell::Application::AppThermalTest< dim >, FuelCell::Application::AppDiffusion< dim >, FuelCell::Application::AppOhmic< dim >, FuelCell::Application::Thermaltesting< dim >, FuelCell::Application::AppCathodeKG< dim >, FuelCell::Application::AppCompressibleFlows< dim >, and FuelCell::Application::CapillaryTesting< dim >.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::initialize_triangulation ( ParameterHandler &  param)
protectedvirtual

Function used to read in a mesh and hand it over to the boost::shared_ptr<Triangulation<dim> > tr object.

This object stores the mesh in the application and it is used to initialize the residual and global matrix sizes.

Note
This routine is usually called by initialize or _initialize
under development

Reimplemented in FuelCell::Application::AppThermalTest< dim >.

template<int dim>
unsigned int FuelCell::ApplicationCore::DoFApplication< dim >::memory_consumption ( ) const

Compute the amount of memory needed by this object.

template<int dim>
void FuelCell::ApplicationCore::DoFApplication< dim >::print ( const std::string &  basename,
const FEVector src,
const std::vector< unsigned int > &  src_indices = std::vector< unsigned int >() 
) const

This function prints FEVector src to a text file basename.

Note
If src_indices is empty, the whole FEVector src will be printed to a text file basename.
template<int dim>
void FuelCell::ApplicationCore::DoFApplication< dim >::read_init_solution ( FEVector dst,
bool &  good_solution 
) const
protected

Create a mesh and assign it to object tr.

This member function is usually called by initialize Function to transfer a solution on a refined grid to the initial coarse grid

Set read to true if you would like the application to read the initial solution from a previous simulation. The initial solution from a previous iteration is stored in the hiddern file .transfer_solution.FEVector. In order to generate this file, you need to setup the following parameter in the input file to true:

* subsection Adaptive refinement
* set Output solution for transfer = true
* end
*
Note
This function is called by initialize_solution
This function should be private
template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::remesh ( )
virtual

Refine grid accordingly to the Refinement options under "Grid Generation" in the parameter file.

Reimplemented from FuelCell::ApplicationCore::ApplicationBase.

Reimplemented in FuelCell::ApplicationCore::BlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::remesh_dofs ( )
virtual

Initialize dof handler, count the dofs in each block and renumber the dofs.

Different numbering schemes can be easily implemented by reimplementing this function in a derived class and do the sorting after having called this function. Remember though to always sort by component in the end.

template<int dim>
virtual double FuelCell::ApplicationCore::DoFApplication< dim >::residual ( FEVector dst,
const FEVectors src,
bool  apply_boundaries = true 
)
virtual

Loop over all cells to compute the residual.

Uses the virtual local functions

Reimplemented from FuelCell::ApplicationCore::ApplicationBase.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::residual_constraints ( FEVector dst) const
virtual

Apply boundary conditions and hanging node constraints to a residual vector after it has been computed by residual().

This function is used in residual() and its default implementation is void.

Note
Reimplement this function if you want to constrain the residual. A typical implementation would be
hanging_node_constraints.condense(dst);
constrain_boundary(v, true);

Reimplemented in FuelCell::ApplicationCore::BlockMatrixApplication< dim >.

template<int dim>
virtual void FuelCell::ApplicationCore::DoFApplication< dim >::sort_dofs ( DoFHandler< dim > *  dof_handler) const
privatevirtual

Sort the degrees of freedom.

In a derived class, sorting order and schemes can be changed by overloading this function.

template<int dim>
void FuelCell::ApplicationCore::DoFApplication< dim >::store_triangulation ( Triangulation< dim > &  new_tr)

Function to copy a triangulation object for use after refinement.

template<int dim>
void FuelCell::ApplicationCore::DoFApplication< dim >::transfer_solution_to_coarse_mesh ( Triangulation< dim > &  tr_coarse,
FEVector coarse_solution,
FEVector refined_solution 
)

Function to perform the transfer of a solution on a refined grid to the initial coarse grid.

Member Data Documentation

template<int dim>
BlockInfo FuelCell::ApplicationCore::DoFApplication< dim >::block_info
protected
template<int dim>
std::map<unsigned int, double> FuelCell::ApplicationCore::DoFApplication< dim >::boundary_constraints
protected

List of all nodes constrained by a strong boundary condition, together with a value to be assigned.

This map must be filled by the application and is used in constrain_boundary().

template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::boundary_fluxes
protected

Extend the integration loops in assemble() and residual() also to boundary faces.

template<int dim>
Table<2, DoFTools::Coupling> FuelCell::ApplicationCore::DoFApplication< dim >::cell_couplings
protected

Couplings through cell bilinear forms.

cell_coupling allows to specify which variables couple with which equations. This is used by DoFTools to generate a sparsity pattern that does not contain elements unless specified in the cell_coupings, therefore, reducing the amount of memory needed.

Note
cell_coupings is a two dimensional table (i.e. a matrix) of DoFTools::Coupling objects. DoFTools::Coulings takes the values: always (components couple), zero (no coupling) and nonzero

Example: For the 2D Stokes equations

\begin{eqnarray*} \frac{d^2u}{dx^2} +\frac{d^2u}{dy^2} + \frac{dp}{dx} = 0 \\ \frac{d^2v}{dx^2} +\frac{d^2v}{dy^2} + \frac{dp}{dy} = 0 \\ \frac{du}{dx} +\frac{dv}{dy} = 0 \end{eqnarray*}

we need to initialize cell_coupling to:

* cell_coupling.reinit(3,3);
* cell_coupling(0,0) = DoFTools::always;
* cell_coupling(0,1) = DoFTools::zero;
* cell_couplings(0,2) = DoFTools::always;
* cell_couplings(1,0) = DoFTools::zero;
* ...
* 
template<int dim>
Vector<float> FuelCell::ApplicationCore::DoFApplication< dim >::cell_errors
protected

The result of error estimation by cell.

This variable will be used for adaptive mesh refinement.

template<int dim>
double FuelCell::ApplicationCore::DoFApplication< dim >::coarsening_threshold
protected

Coarsening threshold for adaptive method from parameter file.

template<int dim>
std::vector< component_boundaryID_value_map > FuelCell::ApplicationCore::DoFApplication< dim >::component_boundaryID_value_maps
protected

Each entry of this std::vector reflects the following structure (see FuelCell::InitialAndBoundaryData namespace docs):

  • first argument : name of the solution component,
  • second argument : boundary id,
  • third argument : value of the solution component.
template<int dim>
std::vector< component_materialID_value_map > FuelCell::ApplicationCore::DoFApplication< dim >::component_materialID_value_maps
protected

Each entry of this std::vector reflects the following structure (see FuelCell::InitialAndBoundaryData namespace docs):

  • first argument : name of the solution component,
  • second argument : material id,
  • third argument : value of the solution component.
template<int dim>
types::boundary_id FuelCell::ApplicationCore::DoFApplication< dim >::curved_bdry_id

Curved boundary ID.

template<int dim>
boost::shared_ptr< Boundary<dim> > FuelCell::ApplicationCore::DoFApplication< dim >::curved_boundary

Curved boundary.

template<int dim>
DataOut<dim, DoFHandler<dim> > FuelCell::ApplicationCore::DoFApplication< dim >::d_out
protected

The object for writing data.

template<int dim>
std::vector<DataComponentInterpretation::DataComponentInterpretation> FuelCell::ApplicationCore::DoFApplication< dim >::data_interpretation
Deprecated:
Vector for adjusting output to vector data.
template<int dim>
boost::shared_ptr<DoFHandler<dim> > FuelCell::ApplicationCore::DoFApplication< dim >::dof
protected

Pointer to the DoFHandler object.

template<int dim>
boost::shared_ptr<FiniteElement<dim> > FuelCell::ApplicationCore::DoFApplication< dim >::element
protected

The finite element used in dof.

This can be a single element or a FESystem

template<int dim>
Vector<float> FuelCell::ApplicationCore::DoFApplication< dim >::face_errors
protected

The result of error estimation by face.

template<int dim>
std::string FuelCell::ApplicationCore::DoFApplication< dim >::filename_initial_sol

Filename where to output the initial grid.

This value is initialized via the input file using:

* subsection Adaptive refinement
* set Initial solution filename = initial_sol # Set flag to true if you want to output the initial solution to file
* end
*
template<int dim>
Table<2, DoFTools::Coupling> FuelCell::ApplicationCore::DoFApplication< dim >::flux_couplings
protected

Couplings through flux bilinear forms.

template<int dim>
GridOut FuelCell::ApplicationCore::DoFApplication< dim >::g_out
protected

The object for writing grids.

template<int dim>
ConstraintMatrix FuelCell::ApplicationCore::DoFApplication< dim >::hanging_node_constraints
protected

Constraint Matrix object.

This object contains a list of the constraints from the hanging nodes. It needs to be used to remove hanging nodes from the system matrix and rhs before the system of equations is solved using

* hanging_node_constraints.condense(system_matrix);
* hanging_node_constraints.condense(system_rhs);
* 

and it needs to be used to add the hanging nodes to the solution once the system is solved using

* hanging_node_constraints.distribute(solution);
* 
template<int dim>
unsigned int FuelCell::ApplicationCore::DoFApplication< dim >::initial_refinement
protected

Initial refinement from parameter file.

template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::interior_fluxes
protected

Extend the integration loops in assemble() and residual() also to interior faces.

template<int dim>
boost::shared_ptr<Mapping<dim> > FuelCell::ApplicationCore::DoFApplication< dim >::mapping
protected

The mapping used for the transformation of mesh cells.

Note
Should be set by initialize() in a derived class before _initialize() is called.
template<int dim>
unsigned int FuelCell::ApplicationCore::DoFApplication< dim >::mapping_degree
protected

Degree used for polynomial mapping of arbitrary order.

By mapping we mean the transformation between the unit cell (i.e. the unit line, square, or cube) to the cells in real space. Before we have implicitly used linear or d-linear mappings and have not noticed this at all, since this is what happens if we do not do anything special. However, if the domain has curved boundaries, there are cases where the piecewise linear approximation of the boundary (i.e. by straight line segments) is not sufficient, and we may want that the computational domain is an approximation to the real domain using curved boundaries as well.

If the boundary approximation uses piecewise quadratic parabolas to approximate the true boundary, then we say that this is a quadratic or $ Q_2 $ approximation. If we use piecewise graphs of cubic polynomials, then this is a $ Q_3 $ approximation, and so on.

For some differential equations, it is known that piecewise linear approximations of the boundary, i.e. $ Q_1 $ mappings, are not sufficient if the boundary of the domain is curved. Examples are the biharmonic equation using $ C_1 $ elements, or the Euler equation on domains with curved reflective boundaries. In these cases, it is necessary to compute the integrals using a higher order mapping. The reason, of course, is that if we do not use a higher order mapping, the order of approximation of the boundary dominates the order of convergence of the entire numerical scheme, irrespective of the order of convergence of the discretization in the interior of the domain.

template<int dim>
boost::shared_ptr< FuelCellShop::Geometry::GridBase<dim> > FuelCell::ApplicationCore::DoFApplication< dim >::mesh_generator
protected

Grid.

template<int dim>
unsigned int FuelCell::ApplicationCore::DoFApplication< dim >::n_ref
private

Number of refinements.

Note
Normally, n_ref is stored in the data object. But if adaptive refinement class is not used, we need to initialize this variable in this class so that the mesh refinement saturation does not happen (see the source file of this class, search for n_ref).
template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::output_actual_degree

true, if you want to output solution and postprocessing data using actual finite element fields Q_n with n >= 1.

false, if only Q_1 is used even for higher degree data outputs.

template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::output_coarse_solution

Bool flag used to specify if the final solution should be stored in the coarse mesh in order to be used later as an initial solution to solve another problem using the flag read_in_initial_solution.

* subsection Adaptive refinement
* set Output solution for transfer = false # Check whether a solution on a refined grid should be output to file on a coarse mesh
* end
*
template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::output_initial_sol

Flag to output the initial solution used to start the solving process.

This flag is set to true to make sure the initial solution satisfies initial and boundary conditions. For nonlinear problems, the final solution might also be dependent on the initial guess, so this flag can be used to asses if that is the case.

This flag is initialized via the input file using:

* subsection Adaptive refinement
* set Output initial solution = false # Set flag to true if you want to output the initial solution to file
* end
*
template<int dim>
Vector<double> FuelCell::ApplicationCore::DoFApplication< dim >::output_levels

Vector that will be used by data_out to store the refinement levels at each cell.

template<int dim>
Vector<double> FuelCell::ApplicationCore::DoFApplication< dim >::output_materials

Vector that will be used by data_out to store the material ids.

template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::output_materials_and_levels

output_materials_and_levels if true then visualized, otherwise suppressed.

Initialized as true in the constructors.

template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::output_matrices_and_rhs

If true, all cell matrices and right hand sides will be output.

Warning
Might lead to huge output depending on number of degrees of freedom.
template<int dim>
std::vector< DataComponentInterpretation::DataComponentInterpretation > FuelCell::ApplicationCore::DoFApplication< dim >::postprocessing_interpretations

postprocessing_interpretations identifies whether some postprocessing_names are scalars or parts of a vector.

Note that if postprocessing_interpretations is empty, all postprocessing_names are treated as scalars.

template<int dim>
std::vector<unsigned int> FuelCell::ApplicationCore::DoFApplication< dim >::postprocessing_printing_indices

The indices of the FEVector postprocessing to be printed to a text file.

Note
If print_postprocessing=true and postprocessing_printing_indices is not specified in the parameter file, the whole FEVector postprocessing will be printed to a text file.
template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::print_blocks_instead_of_indices

true, if the whole blocks of FEVector solution or FEVector postprocessing to be printed to a text file instead of separate indices.

false, otherwise.

template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::print_postprocessing

true, if you want to print FEVector postprocessing to a text file.

false, otherwise.

template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::print_solution

true, if you want to print FEVector solution to a text file.

false, otherwise.

template<int dim>
Quadrature<dim-1> FuelCell::ApplicationCore::DoFApplication< dim >::quadrature_residual_bdry
protected

Quadrature rule for residual computation on boundary faces.

template<int dim>
Quadrature<dim> FuelCell::ApplicationCore::DoFApplication< dim >::quadrature_residual_cell
protected

Quadrature rule for residual computation on cells.

template<int dim>
Quadrature<dim-1> FuelCell::ApplicationCore::DoFApplication< dim >::quadrature_residual_face
protected

Quadrature rule for residual computation on faces.

template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::read_in_initial_solution

Bool flag used to specify if the initial solution to the problem, specially important for non-linear problems, should be read from file.

This flag is specified in the input file using:

* subsection Adaptive refinement
* set Read in initial solution from file = false # Check whether a stored solution should be read from file and applied to the grid
* end
*
template<int dim>
std::string FuelCell::ApplicationCore::DoFApplication< dim >::refinement
protected

Refinement parameter from parameter file.

template<int dim>
double FuelCell::ApplicationCore::DoFApplication< dim >::refinement_threshold
protected

Refinement threshold for adaptive method from parameter file.

template<int dim>
std::vector< DataComponentInterpretation::DataComponentInterpretation > FuelCell::ApplicationCore::DoFApplication< dim >::solution_interpretations

solution_interpretations identifies whether some solution_names are scalars or parts of a vector.

Note that if solution_interpretations is empty, all solution_names are treated as scalars.

template<int dim>
std::vector<unsigned int> FuelCell::ApplicationCore::DoFApplication< dim >::solution_printing_indices

The indices of the FEVector solution to be printed to a text file.

Note
If print_solution=true and solution_printing_indices is not specified in the parameter file, the whole FEVector solution will be printed to a text file.
template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::sort_cuthill
protected

Flag for sorting with Cuthill McKee algorithm.

Note
Read from parameter file
template<int dim>
Point<dim> FuelCell::ApplicationCore::DoFApplication< dim >::sort_direction
protected

Direction for downstream sorting.

No downstream sorting if this vector is zero.

Note
Read from parameter file
template<int dim>
FuelCell::SystemManagement FuelCell::ApplicationCore::DoFApplication< dim >::system_management
protected

This object knows everything about FCST equations, variables, couplings, etc.

template<int dim>
boost::shared_ptr<Triangulation<dim> > FuelCell::ApplicationCore::DoFApplication< dim >::tr
protected

Pointer to the Triangulation object.

template<int dim>
std::vector<FEVector*> FuelCell::ApplicationCore::DoFApplication< dim >::transfer_vectors
protected

List of vector names to be transfered from one grid to the next.

template<int dim>
bool FuelCell::ApplicationCore::DoFApplication< dim >::use_predefined_solution

Use user pre-defined initial solution.

template<int dim>
unsigned int FuelCell::ApplicationCore::DoFApplication< dim >::verbosity
protected

Controls verbosity of certain functions.

Zero is quiet and default is one.


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