Reference documentation for deal.II version GIT relicensing-437-g81ec864850 2024-04-19 07:30:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType > Class Template Reference

#include <deal.II/base/quadrature_point_data.h>

Public Types

using CellIteratorType = typename parallel::distributed::Triangulation< dim >::cell_iterator
 

Public Member Functions

 ContinuousQuadratureDataTransfer (const FiniteElement< dim > &projection_fe, const Quadrature< dim > &mass_quadrature, const Quadrature< dim > &data_quadrature)
 
void prepare_for_coarsening_and_refinement (parallel::distributed::Triangulation< dim > &tria, CellDataStorage< CellIteratorType, DataType > &data_storage)
 
void interpolate ()
 

Private Member Functions

std::vector< char > pack_function (const typename parallel::distributed::Triangulation< dim >::cell_iterator &cell, const CellStatus status)
 
void unpack_function (const typename parallel::distributed::Triangulation< dim >::cell_iterator &cell, const CellStatus status, const boost::iterator_range< std::vector< char >::const_iterator > &data_range)
 

Private Attributes

const std::unique_ptr< const FiniteElement< dim > > projection_fe
 
std::size_t data_size_in_bytes
 
const unsigned int n_q_points
 
FullMatrix< double > project_to_fe_matrix
 
FullMatrix< double > project_to_qp_matrix
 
FullMatrix< double > matrix_dofs
 
FullMatrix< double > matrix_dofs_child
 
FullMatrix< double > matrix_quadrature
 
unsigned int handle
 
CellDataStorage< CellIteratorType, DataType > * data_storage
 
parallel::distributed::Triangulation< dim > * triangulation
 

Detailed Description

template<int dim, typename DataType>
class parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >

A class for the transfer of continuous data stored at quadrature points when performing h-adaptive refinement of parallel::distributed::Triangulation .

Implementation details

This class implements the transfer of the quadrature point data between cells in case of adaptive refinement using L2 projection. That also includes automatic shipping of information between different processors.

To that end, the constructor of the class is provided with three main objects: scalar FiniteElement projection_fe, mass_quadrature and data_quadrature Quadrature rules. First, the data located at data_quadrature of each cell is L2-projected to the continuous space defined by a single FiniteElement projection_fe . This is achieved using FETools::compute_projection_from_quadrature_points_matrix(). In doing so the mass matrix of this element is required, which will be calculated with the mass_quadrature rule . Should the cell now belong to another processor, the data is then sent to this processor. The class makes use of a feature of p4est (and parallel::distributed::Triangulation) that allows one to attach information to cells during mesh refinement and rebalancing. On receiving information on the target cell, the data is projected back to the quadrature points using the matrix calculated by FETools::compute_interpolation_to_quadrature_points_matrix() . In the case that local refinement is performed, this class first project local DoF values of the parent element to each child.

This class is templated by DataType type, however the user's DataType class has to be derived from the TransferableQuadraturePointData class. In practice that amounts to implementing the following three functions shown below for a quadrature point data with 2 scalars:

class MyQData : public TransferableQuadraturePointData
{
public:
double elasticity_parameter_lambda;
double elasticity_parameter_mu;
unsigned int number_of_values() const
{
return 2;
}
// a function to pack scalars into a vector
void pack_values(std::vector<double> &values) const
{
Assert (values.size()==2, ExcInternalError());
values[0] = elasticity_parameter_lambda;
values[1] = elasticity_parameter_mu;
}
void unpack_values(const std::vector<double> &values)
{
Assert (values.size() ==2, ExcInternalError());
elasticity_parameter_lambda = values[0];
elasticity_parameter_mu = values[1];
}
};
#define Assert(cond, exc)

Note that the order of packing and unpacking has to be the same.

This class can then be use with CellDataStorage in the following way:

data_transfer(FE_Q<dim>(2),QGauss<dim>(3),QGauss<dim>(4));
//...populate data for all active cells in data_storage
//...mark cells for refinement...
data_transfer.prepare_for_coarsening_and_refinement(triangulation,data_storage);
triangulation.execute_coarsening_and_refinement();
//...initialize quadrature point data on new cells by calling
// CellDataStorage::initialize()
data_transfer.interpolate();
Definition fe_q.h:550
parallel::distributed::Triangulation< dim > * triangulation
CellDataStorage< CellIteratorType, DataType > * data_storage

This approach can be extended to quadrature point data with Tensor objects of arbitrary order, although with a little bit more work in packing and unpacking of data inside MyQData class.

Note
Currently coarsening is not supported.
The functionality provided by this class can alternatively be achieved using parallel::distributed::SolutionTransfer. However, that would require the following steps: (i) create an auxiliary DoFHandler with a (discontinuous Galerkin) FiniteElement which has enough components to represent all data stored at the quadrature points; (ii) project the data to the FiniteElement space and thereby store results in global vectors; (iii) use parallel::distributed::SolutionTransfer to project FE vectors to the new mesh; and (iv) finally project the data back to the quadrature points on the new mesh via FEValues class. The ContinuousQuadratureDataTransfer class aims at simplifying the whole process by only requiring that the quadrature point data class is derived from the TransferableQuadraturePointData. Everything else will be done automatically.
This class is not well suited to situations where the values stored at quadrature points represent samples from a discontinuous field. An example for such a situation would be where the data stored at the quadrature points represents the elastic or plastic state of a material, i.e., a property that varies discontinuously within a solid. In such cases, trying to transfer data from the quadrature points to a finite element field that is continuous (at least within one cell) will likely yield over and undershoots that, once evaluated at a different set of quadrature points (on child or parent cells) results in values that will not make much sense.

Definition at line 439 of file quadrature_point_data.h.

Member Typedef Documentation

◆ CellIteratorType

template<int dim, typename DataType >
using parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::CellIteratorType = typename parallel::distributed::Triangulation<dim>::cell_iterator

An alias for a cell.

Definition at line 449 of file quadrature_point_data.h.

Constructor & Destructor Documentation

◆ ContinuousQuadratureDataTransfer()

template<int dim, typename DataType >
parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::ContinuousQuadratureDataTransfer ( const FiniteElement< dim > &  projection_fe,
const Quadrature< dim > &  mass_quadrature,
const Quadrature< dim > &  data_quadrature 
)

Constructor which takes the FiniteElement projection_fe , the quadrature rule mass_quadrature used to integrate its local mass matrix and finally the quadrature rule data_quadrature which is used to store DataType.

Precondition
projection_fe has to be scalar-valued.
Note
Since this class does projection on cell-by-cell basis, projection_fe is only required to be continuous within the cell.

Member Function Documentation

◆ prepare_for_coarsening_and_refinement()

template<int dim, typename DataType >
void parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::prepare_for_coarsening_and_refinement ( parallel::distributed::Triangulation< dim > &  tria,
CellDataStorage< CellIteratorType, DataType > &  data_storage 
)

Prepare for coarsening and refinement of a triangulation tria . data_storage represents the cell data which should be transferred and it should be initialized for each locally owned active cell.

Note
Although CellDataStorage class allows storing on different cells different objects derived from the base class, here we assume that data_storage contains objects of the same type, more specifically they pack/unpack the same data.

◆ interpolate()

template<int dim, typename DataType >
void parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::interpolate ( )

Interpolate the data previously stored in this object before the mesh was refined or coarsened onto the quadrature points of the currently active set of cells.

Note
Before calling this function the user is expected to populate the data stored in the data_storage object provided to prepare_for_coarsening_and_refinement() at new cells using CellDataStorage::initialize(). If that is not the case, an exception will be thrown in debug mode.

◆ pack_function()

template<int dim, typename DataType >
std::vector< char > parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::pack_function ( const typename parallel::distributed::Triangulation< dim >::cell_iterator &  cell,
const CellStatus  status 
)
private

A callback function used to pack the data on the current mesh into objects that can later be retrieved after refinement, coarsening and repartitioning.

◆ unpack_function()

template<int dim, typename DataType >
void parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::unpack_function ( const typename parallel::distributed::Triangulation< dim >::cell_iterator &  cell,
const CellStatus  status,
const boost::iterator_range< std::vector< char >::const_iterator > &  data_range 
)
private

A callback function used to unpack the data on the current mesh that has been packed up previously on the mesh before refinement, coarsening and repartitioning.

Member Data Documentation

◆ projection_fe

template<int dim, typename DataType >
const std::unique_ptr<const FiniteElement<dim> > parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::projection_fe
private

FiniteElement used to project data from and to quadrature points.

Definition at line 522 of file quadrature_point_data.h.

◆ data_size_in_bytes

template<int dim, typename DataType >
std::size_t parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::data_size_in_bytes
private

The size of the data that will be sent, which depends on the DataType class.

Definition at line 528 of file quadrature_point_data.h.

◆ n_q_points

template<int dim, typename DataType >
const unsigned int parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::n_q_points
private

Number of quadrature points at which DataType is stored.

Definition at line 533 of file quadrature_point_data.h.

◆ project_to_fe_matrix

template<int dim, typename DataType >
FullMatrix<double> parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::project_to_fe_matrix
private

Projection matrix from the quadrature points to local DoFs for a single scalar.

Definition at line 539 of file quadrature_point_data.h.

◆ project_to_qp_matrix

template<int dim, typename DataType >
FullMatrix<double> parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::project_to_qp_matrix
private

Projection matrix from the local DoFs to quadrature points for a single scalar.

Definition at line 545 of file quadrature_point_data.h.

◆ matrix_dofs

template<int dim, typename DataType >
FullMatrix<double> parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::matrix_dofs
private

Auxiliary matrix which represents projection of each internal value stored at the quadrature point (second index) to the local DoFs of the projection_fe (first index).

Definition at line 553 of file quadrature_point_data.h.

◆ matrix_dofs_child

template<int dim, typename DataType >
FullMatrix<double> parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::matrix_dofs_child
private

Projection of matrix_dofs to each child cell in case of adaptive refinement.

Definition at line 558 of file quadrature_point_data.h.

◆ matrix_quadrature

template<int dim, typename DataType >
FullMatrix<double> parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::matrix_quadrature
private

Auxiliary matrix which represents data (second index) stored at each quadrature point (first index).

Definition at line 564 of file quadrature_point_data.h.

◆ handle

template<int dim, typename DataType >
unsigned int parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::handle
private

The handle that the parallel::distributed::Triangulation has assigned to this object while registering the pack_callback function.

Definition at line 570 of file quadrature_point_data.h.

◆ data_storage

template<int dim, typename DataType >
CellDataStorage<CellIteratorType, DataType>* parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::data_storage
private

A pointer to the CellDataStorage class whose data will be transferred.

Definition at line 575 of file quadrature_point_data.h.

◆ triangulation

template<int dim, typename DataType >
parallel::distributed::Triangulation<dim>* parallel::distributed::ContinuousQuadratureDataTransfer< dim, DataType >::triangulation
private

A pointer to the distributed triangulation to which cell data is attached.

Definition at line 581 of file quadrature_point_data.h.


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