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 Member Functions | Protected Attributes | Private Attributes | List of all members
FuelCell::ApplicationCore::BlockMatrixApplication< dim > Class Template Reference

Application handling matrices and assembling linear systems of equations. More...

#include <block_matrix_application.h>

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

Public Member Functions

void serial_assemble (const FEVectors &)
 serial and PETSc assemble functions called by assemble(). More...
 
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. More...
 
virtual void solve (FEVector &dst, const FEVectors &src)
 Solve the system assembled with right hand side in FEVectors src and return the result in FEVector dst. More...
 
void serial_solve (FuelCell::ApplicationCore::FEVector system_rhs, FEVector &solution)
 
Constructors and Initialization
 BlockMatrixApplication (boost::shared_ptr< ApplicationData > data=boost::shared_ptr< ApplicationData >())
 Constructor for an object, owning its own mesh and dof handler. More...
 
 BlockMatrixApplication (DoFApplication< dim > &, bool triangulation_only)
 Constructor for an object, borrowing mesh and dof handler from another object. More...
 
void _initialize (ParameterHandler &param)
 Initialize data of this class. More...
 
virtual void declare_parameters (ParameterHandler &param)
 Declare parameters related to the linear system. More...
 
virtual void initialize (ParameterHandler &param)
 Initialize application. More...
 
ApplicationBase interface
void remesh_matrices ()
 Initialize sparsity patterns and matrices for the new mesh. More...
 
virtual void remesh ()
 Refine grid accordingly to the Refinement options under "Grid Generation" in the parameter file. More...
 
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. More...
 
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. More...

 
void residual_constraints (FEVector &dst) const
 Redefinition of residual_constraints() in DoFHandler. More...
 
Local integrators
virtual void cell_matrix (MatrixVector &cell_matrices, const typename DoFApplication< dim >::CellInfo &cell)
 Integration of local bilinear form. More...
 
virtual void bdry_matrix (MatrixVector &face_matrices, const typename DoFApplication< dim >::FaceInfo &face)
 Integration of local bilinear form. More...
 
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. More...
 
Dirichlet boundary conditions

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

virtual void dirichlet_bc (std::map< unsigned int, double > &boundary_values) const
 
- Public Member Functions inherited from 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. 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 remesh_dofs ()
 Initialize dof handler, count the dofs in each block and renumber the dofs. 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...
 
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...
 
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...
 
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 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...
 

Protected Member Functions

template<typename Vector >
void print_matrix_and_rhs (Vector &sys_rhs) const
 Internal routine to print matrix and rhs. More...
 
- Protected Member Functions inherited from FuelCell::ApplicationCore::DoFApplication< dim >
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...
 
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...
 
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

boost::shared_ptr< Quadrature
< dim > > 
quadrature_assemble_cell
 Quadrature rule for matrix assembling on cells. More...
 
boost::shared_ptr< Quadrature
< dim-1 > > 
quadrature_assemble_face
 Quadrature rule for matrix assembling on faces. More...
 
bool repair_diagonal
 Bool determining whether or not to repair diagonal before solving. More...
 
std::map< unsigned int, double > boundary_values
 Variable to store boundary values, so they only need to be computed once per mesh refinement. More...
 
System matrix and boundary conditions
BlockSparseMatrix< double > matrix
 Storage for the actual matrix. More...
 
SolverControl solver_control
 Solver control object. More...
 
Other application variables
bool assemble_numerically_flag
 Variable used to select if assembly should be done analytically or numerically. More...
 
bool mumps_additional_mem
 Variable used for configuring MUMPS to use more memory in solve() function. More...
 
bool symmetric_matrix_flag
 Variable used to specify if the matrix is symmetric to the linear solver. More...
 
bool print_debug
 Flag specifying if the matrix and rhs should be printed at each iteration. More...
 
- Protected Attributes inherited from FuelCell::ApplicationCore::DoFApplication< dim >
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...
 
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...
 
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 Attributes

BlockSparsityPattern sparsities
 Sparsity patterns. More...
 
Auxiliary data:
bool output_system_assembling_time
 If true, information about duration of the system assembling stage will be output. More...
 
Timer timer
 Object used to calculate the CPU and Run time for the system assembling stage. More...
 

Additional Inherited Members

- Public Types inherited from FuelCell::ApplicationCore::DoFApplication< dim >
typedef
FuelCell::ApplicationCore::IntegrationInfo
< dim, FEValuesBase< dim > > 
CellInfo
 Shortcut. More...
 
typedef
FuelCell::ApplicationCore::IntegrationInfo
< dim, FEFaceValuesBase< dim > > 
FaceInfo
 Shortcut. More...
 
- Public Attributes inherited from FuelCell::ApplicationCore::DoFApplication< dim >
boost::shared_ptr< Boundary
< dim > > 
curved_boundary
 Curved boundary. More...
 
types::boundary_id curved_bdry_id
 Curved boundary ID. More...
 
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...
 
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...
 

Detailed Description

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

Application handling matrices and assembling linear systems of equations.

This class inherits the information on Triangulation and DoFHandler from its base class DoFApplication. It adds information on the linear system, in particular the sparsities.

Additionally, it holds a BlockSparseMatrix as the system matrix of the finite element problem. Generic loops over mesh cells and faces are implemented in order to assemble() the system matrices and compute the residual(). These functions use the virtual local integration functions declared here and implemented by the derived actual application class.

Note
You probably miss functions for building the right hand side of a linear system. They are missing deliberately, since the linear case is very special and can be handled in the framework of nonlinear problems. So, if you have a linear problem, code the local matrices for residual() so they ignore the BlockVector src.

Deriving your own application

Author
Guido Kanschat
Marc Secanell

Constructor & Destructor Documentation

template<int dim>
FuelCell::ApplicationCore::BlockMatrixApplication< dim >::BlockMatrixApplication ( 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 handler for application data is generated.

see constructor of DoFApplication.

template<int dim>
FuelCell::ApplicationCore::BlockMatrixApplication< dim >::BlockMatrixApplication ( DoFApplication< dim > &  ,
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.

see constructor of DoFApplication.

Member Function Documentation

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

Initialize data of this class.

Note
remesh_dofs() has to be called before this member function since remesh_dofs() initializes data that is used in this routine. In particular, dof and block_indices.
template<int dim>
void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::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.

template<int dim>
void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::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.

The solution is perturbed at each degree of freedom and the residual function is used to compute the residual for each perturbed solution. Forward differences are used to compute the Jacobian.

The Jacobian is stored in this->matrix

Parameters
[in]srcA vector of BlockVectors containting the solution at the previous Newton step and the residual.
[in]deltaThe step size used for numerical differenciation.
Note
This routine is used in combination with a Netwon nonlinear solver since it requires src to contain a field named "Newton iterate" with the solution at the previous Newton step.
Author
M. Secanell, 2013
template<int dim>
virtual void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::bdry_matrix ( MatrixVector face_matrices,
const typename DoFApplication< dim >::FaceInfo face 
)
virtual
template<int dim>
virtual void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::cell_matrix ( MatrixVector cell_matrices,
const typename DoFApplication< dim >::CellInfo cell 
)
virtual
template<int dim>
virtual void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::declare_parameters ( ParameterHandler &  param)
virtual

Declare parameters related to the linear system.

In particular, the following sections are declared:

*
*
* subsection Discretization
* subsection Matrix
* set Quadrature cell = -1
* set Quadrature face = -1
* end
* end
* #subsection used to select linear solver and linear solver parameters:
* subsection Linear Solver
* set Print matrix and rhs = false # print matrices and rhs for debugging purposes
* set Assemble numerically = false # Assemble Jacobian using analytical derivatives or numerically.
* set Symmetric matrix = false # If true, faster solvers will be used (only for symmetric matrices!).
*
* set Type of linear solver = MUMPS
* Options: (For parallel code (--with-petsc) MUMPS|CG|Bicgstab )
* For serial code ILU-GMRES|UMFPACK|Bicgstab )
* set Max steps = 100
* set Tolerance = 1.e-10
* set Log history = false
* set Log frequency = 1
* set Log result = true
*
* set Allocate additional memory for MUMPS = false # Configure MUMPS solver to use additional memory during solving if needed
* set Output system assembling time = false # Flag for specifying if you want to see how much time the linear system assembling takes.
* end
*
Note
The function declare_parameters of derived applications should always call the one of its base class.

Reimplemented from FuelCell::ApplicationCore::DoFApplication< dim >.

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::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>
virtual void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::dirichlet_bc ( std::map< unsigned int, double > &  boundary_values) const
virtual
template<int dim>
virtual void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::face_matrix ( MatrixVector matrices11,
MatrixVector matrices12,
MatrixVector matrices21,
MatrixVector matrices22,
const typename DoFApplication< dim >::FaceInfo face1,
const typename DoFApplication< dim >::FaceInfo face2 
)
virtual

Integration of local bilinear form.

template<int dim>
virtual void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::initialize ( ParameterHandler &  param)
virtual
template<int dim>
virtual void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::post_cell_assemble ( )
virtual

A call back function for assemble(), called after all cell matrices have been entered, but before the face matrices are computed.

template<int dim>
template<typename Vector >
void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::print_matrix_and_rhs ( Vector &  sys_rhs) const
protected

Internal routine to print matrix and rhs.

template<int dim>
virtual void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::remesh ( )
virtual

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

Reimplemented from FuelCell::ApplicationCore::DoFApplication< dim >.

template<int dim>
void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::remesh_matrices ( )

Initialize sparsity patterns and matrices for the new mesh.

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

Redefinition of residual_constraints() in DoFHandler.

In this member function, I adjust the residual in order to account for hanging nodes and I apply Dirichlet boundary conditions given by dirichlet_bc().

Reimplemented from FuelCell::ApplicationCore::DoFApplication< dim >.

template<int dim>
void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::serial_assemble ( const FEVectors )

serial and PETSc assemble functions called by assemble().

template<int dim>
void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::serial_solve ( FuelCell::ApplicationCore::FEVector  system_rhs,
FEVector solution 
)
template<int dim>
virtual void FuelCell::ApplicationCore::BlockMatrixApplication< dim >::solve ( FEVector dst,
const FEVectors src 
)
virtual

Solve the system assembled with right hand side in FEVectors src and return the result in FEVector dst.

Pure virtual.

Implements FuelCell::ApplicationCore::ApplicationBase.

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

Member Data Documentation

template<int dim>
bool FuelCell::ApplicationCore::BlockMatrixApplication< dim >::assemble_numerically_flag
protected

Variable used to select if assembly should be done analytically or numerically.

template<int dim>
std::map<unsigned int, double> FuelCell::ApplicationCore::BlockMatrixApplication< dim >::boundary_values
protected

Variable to store boundary values, so they only need to be computed once per mesh refinement.

template<int dim>
BlockSparseMatrix<double> FuelCell::ApplicationCore::BlockMatrixApplication< dim >::matrix
protected

Storage for the actual matrix.

template<int dim>
bool FuelCell::ApplicationCore::BlockMatrixApplication< dim >::mumps_additional_mem
protected

Variable used for configuring MUMPS to use more memory in solve() function.

template<int dim>
bool FuelCell::ApplicationCore::BlockMatrixApplication< dim >::output_system_assembling_time
private

If true, information about duration of the system assembling stage will be output.

template<int dim>
bool FuelCell::ApplicationCore::BlockMatrixApplication< dim >::print_debug
protected

Flag specifying if the matrix and rhs should be printed at each iteration.

template<int dim>
boost::shared_ptr<Quadrature<dim> > FuelCell::ApplicationCore::BlockMatrixApplication< dim >::quadrature_assemble_cell
protected

Quadrature rule for matrix assembling on cells.

template<int dim>
boost::shared_ptr<Quadrature<dim-1> > FuelCell::ApplicationCore::BlockMatrixApplication< dim >::quadrature_assemble_face
protected

Quadrature rule for matrix assembling on faces.

template<int dim>
bool FuelCell::ApplicationCore::BlockMatrixApplication< dim >::repair_diagonal
protected

Bool determining whether or not to repair diagonal before solving.

template<int dim>
SolverControl FuelCell::ApplicationCore::BlockMatrixApplication< dim >::solver_control
protected

Solver control object.

template<int dim>
BlockSparsityPattern FuelCell::ApplicationCore::BlockMatrixApplication< dim >::sparsities
private

Sparsity patterns.

template<int dim>
bool FuelCell::ApplicationCore::BlockMatrixApplication< dim >::symmetric_matrix_flag
protected

Variable used to specify if the matrix is symmetric to the linear solver.

template<int dim>
Timer FuelCell::ApplicationCore::BlockMatrixApplication< dim >::timer
private

Object used to calculate the CPU and Run time for the system assembling stage.


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