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

Integrate the residual using DoFApplication::integrate_1form(). More...

#include <dof_application.h>

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

Public Types

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 Member Functions

 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.
 
virtual void declare_parameters (ParameterHandler &param)
 Declare parameters related to mesh generation and finite elements.
 
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 initialize (ParameterHandler &param)
 Pure virtual initialization function.
 
virtual void init_vector (FEVector &dst) const
 Reinitialize the BlockVector dst such that it contains block_size.size() blocks.
 
virtual void remesh ()
 Refine grid accordingly to the Refinement options under "Grid Generation" in the parameter file.
 
virtual double evaluate (const FEVectors &)
 Evaluate a functional during postprocessing such as drag or lift on an aerodynamics problem.
 
virtual double estimate (const FEVectors &src)
 Estimate the error.
 
virtual double residual (FEVector &dst, const FEVectors &src, bool apply_boundaries=true)
 Loop over all cells to compute the residual.
 
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().
 
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.
 
- 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 solve (FEVector &start, const FEVectors &rhs)
 Solve the system assembled with right hand side rhs and return the result in start.
 
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.
 

Public Attributes

unsigned int verbosity
 Controls verbosity of certain functions.
 

Protected Member Functions

void constrain_boundary (FEVector &v, bool homogeneous) const
 Apply either homogeneous or inhomogeneous boundary_constraints.
 
- Protected Member Functions inherited from AppFrame::ApplicationBase
void print_caller_name (const std::string &caller_name) const
 Print caller name.
 

Protected Attributes

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.
 

Private Member Functions

virtual void sort_dofs ()
 Sort the degrees of freedom.
 
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.
 
virtual double global_from_local_errors () const
 Compute a global estimation value from cell_errors.
 

Friends

class LocalResidual< dim >
 
class LocalEstimate< dim >
 

Data output

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.
 
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)
 
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.
 

Detailed Description

template<int dim>
class AppFrame::DoFApplication< dim >

Integrate the residual using DoFApplication::integrate_1form().

Optionally, this class can use its own LocalIntegrator or use the standard residual of the DoFApplication.

Author
Guido Kanschat, 2008 Constructor initializing this object to compute the residual using application and the provided local_residual. The DoFApplication class providing the function DoFApplication::integrate_1form(). The object integrating the local residuals. Base class for all terminal applications, i.e. all applications requiring a Triangulation and a handler for degrees of freedom.

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 of a FESystem. In the latter case, the nomenclature used in the paramter file is: s et Element = FESys*tem[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]

Author
Guido Kanschat, 2006

Member Typedef Documentation

template<int dim>
typedef MeshWorker::InfoObjects::IntegrationInfo<dim, FEValuesBase<dim> > AppFrame::DoFApplication< dim >::CellInfo

This is a redefinition of CellInfo in order to call MeshWorker objects instead of the previous implementation in AppFrame.

template<int dim>
typedef MeshWorker::InfoObjects::IntegrationInfo<dim, FEFaceValuesBase<dim> > AppFrame::DoFApplication< dim >::FaceInfo

This is a redefinition of CellInfo in order to call MeshWorker objects instead of the previous implementation in AppFrame.

Constructor & Destructor Documentation

template<int dim>
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.

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

template<int dim>
AppFrame::DoFApplication< dim >::DoFApplication ( bool  multigrid)

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

template<int dim>
AppFrame::DoFApplication< dim >::DoFApplication ( DoFApplication< dim > &  ,
bool  triangulation_only,
bool  multigrid = false 
)

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.
template<int dim>
AppFrame::DoFApplication< dim >::~DoFApplication ( )

Destructor which deletes owned objects.

Member Function Documentation

template<int dim>
void AppFrame::DoFApplication< dim >::_initialize ( ParameterHandler &  param)

Initialize from parameter values.

template<int dim>
void AppFrame::DoFApplication< dim >::add_vector_for_transfer ( FEVector )

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

template<int dim>
virtual double AppFrame::DoFApplication< dim >::bdry_estimate ( const FaceInfo )
virtual

Integration of estimate on a boundary face.

The implementation in this class simply returns zero.

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

Integration of local linear form.

template<int dim>
virtual double AppFrame::DoFApplication< dim >::cell_estimate ( const CellInfo )
virtual

Integration of estimate inside a cell.

The implementation in this *class simply returns zero.

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

Integration of local linear form.

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

Apply either homogeneous or inhomogeneous boundary_constraints.

template<int dim>
virtual void AppFrame::DoFApplication< dim >::data_out ( const std::string &  basename,
const FEVectors src 
)
virtual
template<int dim>
virtual void AppFrame::DoFApplication< dim >::data_out ( const std::string &  basename,
const FEVectors src,
const std::vector< std::string > &  solution_names 
)
virtual
Deprecated:
This function is deprecated.
template<int dim>
virtual void AppFrame::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 >() 
)
virtual

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 AppFrame::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.

Implements AppFrame::ApplicationBase.

Reimplemented in FuelCell::Application::AppCathode< dim >, AppFrame::MatrixApplication< dim >, FuelCell::Application::AppPemfc< dim >, FuelCell::Application::AppReadMesh< dim >, FuelCell::Application::AppLaplace< dim >, AppFrame::OptimizationBlockMatrixApplication< dim >, and AppFrame::BlockMatrixApplication< dim >.

template<int dim>
void AppFrame::DoFApplication< dim >::delete_vector_for_transfer ( )

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

template<int dim>
virtual void AppFrame::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 AppFrame::DoFApplication< dim >::estimate ( const FEVectors src)
virtual

Estimate the error.

This is the loop over all cells and faces, using the virtual local functions

While these are usually implemented for a specific problem, a default implementation measuring discontinuity of the normal derivative (and, if interior_fluxes or boundary_fluxes are true) the function itself.

Reimplemented from AppFrame::ApplicationBase.

Reimplemented in FuelCell::Application::AppCathode< dim >, FuelCell::Application::AppPemfc< dim >, and FuelCell::Application::AppLaplace< dim >.

template<int dim>
virtual double AppFrame::DoFApplication< dim >::evaluate ( const FEVectors )
virtual

Evaluate a functional during postprocessing such as drag or lift on an aerodynamics problem.

Reimplemented from AppFrame::ApplicationBase.

Reimplemented in FuelCell::Application::AppPemfc< dim >, FuelCell::Application::AppLaplace< dim >, and FuelCell::Application::AppCathode< dim >.

template<int dim>
virtual double AppFrame::DoFApplication< dim >::face_estimate ( const FaceInfo ,
const FaceInfo  
)
virtual

Integration of estimate on an interior face.

The implementation in this class simply returns the jump of the normal derivative. If interior_fluxes is set, then the jump of the function is added with appropriate weight.

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

Integration of local linear form.

template<int dim>
template<class INFO , typename TYPE >
void AppFrame::DoFApplication< dim >::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.

template<int dim>
virtual double AppFrame::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 AppFrame::DoFApplication< dim >::grid_out ( const std::string &  basename)
virtual

Output the grid used to solve the problem.

template<int dim>
virtual void AppFrame::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 AppFrame::ApplicationBase.

template<int dim>
virtual void AppFrame::DoFApplication< dim >::initialize ( ParameterHandler &  param)
virtual

Pure virtual initialization function.

Here, the parameters declared in the constructor are actually read from ParameterHandler.

This function must be overloaded in any derived class using its own values from the parameter file.

Initialization functions of derived classes must make sure to call all functions _initialize() of the base classes.

Implements AppFrame::ApplicationBase.

Reimplemented in FuelCell::Application::AppCathode< dim >, FuelCell::Application::AppPemfc< dim >, AppFrame::MatrixApplication< dim >, FuelCell::Application::AppLaplace< dim >, FuelCell::Application::AppReadMesh< dim >, AppFrame::OptimizationBlockMatrixApplication< dim >, and AppFrame::BlockMatrixApplication< dim >.

template<int dim>
template<class BOX , class WORKER >
void AppFrame::DoFApplication< dim >::integrate_1form ( BOX &  box,
WORKER &  local_forms,
FEVector dst,
const FEVectors src 
)

The function integrating 1-forms as for instance used by residual().

Since the concept of a 1-form is simply taking a set of source vectors and local integrators to produce a set of output vectors, we abstract this from the actual residual function here.

The function takes an argument of type LocalIntegrator, which is usually a derived class. LocalResidual itself will produce the residual by applying the virtual functions

template<int dim>
void AppFrame::DoFApplication< dim >::integrate_functional ( LocalEstimate< dim > &  local_forms,
FEVectors dst,
const FEVectors src 
)

The function integrating functionals using local integrators.

The result consists of two vectors containing the local contributions to the functionals of the cells and faces, respectively.

The destination argument is expected to have two FEVector objects at its front, named "cells" and "faces" for cell and face contributions, respectively.

It is possible to compute several functionals at a time in which case these two FEVector objects should be block vectors with one block for each functional. The destination vector for the LocalIntegrator will be adjusted to the correct size.

The function takes an argument of type LocalIntegrator, which is usually a derived class. LocalEstimate itself will produce the residual by applying the virtual functions

Note
The function LocalIntegrator::face() receives the same result vector twice. Therefore, only one of them should be filled.
template<int dim>
unsigned int AppFrame::DoFApplication< dim >::memory_consumption ( ) const

Compute the amount of memory needed by this object.

template<int dim>
void AppFrame::DoFApplication< dim >::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<int dim>
virtual void AppFrame::DoFApplication< dim >::remesh ( )
virtual

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

Reimplemented from AppFrame::ApplicationBase.

Reimplemented in AppFrame::MatrixApplication< dim >, and AppFrame::BlockMatrixApplication< dim >.

template<int dim>
virtual void AppFrame::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 AppFrame::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

This function uses residual_adjust_boundary() right before computing the norm of the residual as return value.

Reimplemented from AppFrame::ApplicationBase.

template<int dim>
virtual void AppFrame::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 AppFrame::BlockMatrixApplication< dim >.

template<int dim>
virtual void AppFrame::DoFApplication< dim >::sort_dofs ( )
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 AppFrame::DoFApplication< dim >::store_triangulation ( Triangulation< dim > &  new_tr)
inline

Function to copy a triangulation object for use after refinement.

References AppFrame::DoFApplication< dim >::tr.

template<int dim>
void AppFrame::DoFApplication< dim >::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.

Friends And Related Function Documentation

template<int dim>
friend class LocalEstimate< dim >
friend
template<int dim>
friend class LocalResidual< dim >
friend

Member Data Documentation

template<int dim>
MeshWorker::BlockInfo AppFrame::DoFApplication< dim >::block_info
protected
template<int dim>
std::map<unsigned int, double> AppFrame::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 AppFrame::DoFApplication< dim >::boundary_fluxes
protected

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

template<int dim>
Vector<float> AppFrame::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 AppFrame::DoFApplication< dim >::coarsening_threshold
protected

Coarsening threshold for adaptive method from parameter file.

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

The object for writing data.

template<int dim>
std::vector<DataComponentInterpretation::DataComponentInterpretation> AppFrame::DoFApplication< dim >::data_interpretation
protected

Vector for adjusting output to vector data.

template<int dim>
boost::shared_ptr<DoFHandler<dim> > AppFrame::DoFApplication< dim >::dof
protected

Pointer to the DoFHandler object.

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

The finite element used in dof.

This can be a single element or a FESystem

template<int dim>
Vector<float> AppFrame::DoFApplication< dim >::face_errors
protected

The result of error estimation by face.

template<int dim>
GridOut AppFrame::DoFApplication< dim >::g_out
protected

The object for writing grids.

template<int dim>
ConstraintMatrix AppFrame::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 AppFrame::DoFApplication< dim >::initial_refinement
protected

Initial refinement from parameter file.

template<int dim>
bool AppFrame::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> > AppFrame::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 AppFrame::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<MGDoFHandler<dim> > AppFrame::DoFApplication< dim >::mg_dof
protected

Pointer to the MGDoFHandler object.

template<int dim>
Vector<double> AppFrame::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> AppFrame::DoFApplication< dim >::output_materials

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

template<int dim>
bool AppFrame::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>
std::vector< DataComponentInterpretation::DataComponentInterpretation > AppFrame::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>
Quadrature<dim-1> AppFrame::DoFApplication< dim >::quadrature_residual_bdry
protected

Quadrature rule for residual computation on boundary faces.

template<int dim>
Quadrature<dim> AppFrame::DoFApplication< dim >::quadrature_residual_cell
protected

Quadrature rule for residual computation on cells.

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

Quadrature rule for residual computation on faces.

template<int dim>
std::string AppFrame::DoFApplication< dim >::refinement
protected

Refinement parameter from parameter file.

template<int dim>
double AppFrame::DoFApplication< dim >::refinement_threshold
protected

Refinement threshold for adaptive method from parameter file.

template<int dim>
std::vector< DataComponentInterpretation::DataComponentInterpretation > AppFrame::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>
bool AppFrame::DoFApplication< dim >::sort_cuthill
protected

Flag for sorting with Cuthill McKee algorithm.

Note
Read from parameter file
template<int dim>
Point<dim> AppFrame::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>
boost::shared_ptr<Triangulation<dim> > AppFrame::DoFApplication< dim >::tr
protected
template<int dim>
std::vector<FEVector*> AppFrame::DoFApplication< dim >::transfer_vectors
protected

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

template<int dim>
unsigned int AppFrame::DoFApplication< dim >::verbosity

Controls verbosity of certain functions.

Zero is quiet and default is one.


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