};
-/**
- * A worker class for integration_loop(), separating assembling and
- * local integration. Objects of this class rely on the base class
- * provided by the template argument ASSEMBLER for assembling local
- * data into global. The integration of local data is delegated to the
- * other template parameter class INTEGRATOR.
- *
- * In order to use this class, it will be necessary to create an
- * INTEGRATOR class following this template:
- *
- * @code
- * template <int dim>
- * class MyLocalOperator : public Subscriptor
- * {
- * public:
- * void cell(typename IntegrationWorker<dim>::CellInfo& info) const;
- * void bdry(typename IntegrationWorker<dim>::FaceInfo& info) const;
- * void face(typename IntegrationWorker<dim>::FaceInfo& info1,
- * typename IntegrationWorker<dim>::FaceInfo& info2) const;
- * };
- * @endcode
- *
- * This class will do whatever your problem requires locally on each
- * cell and/or face. Once this class is defined, you choose a suitable
- * assembler for your problem from the Assembler namespace and set up
- * objects:
- *
- * @code
- * MyLocalOperator<dim> myop;
- *
- * AssemblingIntegrator<dim, Assembler::MyAssembler, MyLocalOperator<dim> >
- * integrator(myop);
- * @endcode
- *
- * You do the necessary initializations of this @p integrator and then
- * you have a worker object suitable for integration_loop().
- *
- * @ingroup MeshWorker
- * @author Guido Kanschat, 2009
- */
- template <int dim, class ASSEMBLER>
- class AssemblingIntegrator : public IntegrationWorker<dim>,
- public ASSEMBLER
- {
- public:
- typedef
- typename IntegrationWorker<dim>::CellInfo
- CellInfo;
-
- typedef
- typename IntegrationWorker<dim>::FaceInfo
- FaceInfo;
-
- /**
- * Constructor, initializing
- * the local data.
- *
- * The three arguments are
- * objects that that will do
- * the local computations on
- * cells, boundary and interior
- * faces. They may be specified
- * as pointers to global or
- * static functions, pointers
- * to objects with an
- * operator() that takes the
- * required number of
- * arguments, or functions with
- * bound arguments.
- */
- AssemblingIntegrator(const std_cxx1x::function<void (CellInfo &)> &cell_worker,
- const std_cxx1x::function<void (FaceInfo &)> &boundary_worker,
- const std_cxx1x::function<void (FaceInfo &, FaceInfo &)> &face_worker);
-
- /**
- * The cell operator called by
- * integration_loop().
- */
- void cell(CellInfo& info);
-
- /**
- * The local boundary operator
- * called by
- * integration_loop().
- */
- void bdry(FaceInfo& info);
-
- /**
- * The interior face operator
- * called by
- * integration_loop().
- */
- void face(FaceInfo& info1,
- FaceInfo& info2);
-
- private:
- /**
- * Pointers to the three
- * functions that will do the
- * local computations on cells,
- * boundary and interior faces.
- */
- const std_cxx1x::function<void (CellInfo &)> cell_worker;
- const std_cxx1x::function<void (FaceInfo &)> boundary_worker;
- const std_cxx1x::function<void (FaceInfo &, FaceInfo &)> face_worker;
- };
-
//----------------------------------------------------------------------//
template <int dim>
inline
:
interior_fluxes(true), boundary_fluxes(true)
{}
-
-//----------------------------------------------------------------------//
-
- template <int dim, class ASSEMBLER>
- inline
- AssemblingIntegrator<dim,ASSEMBLER>::
- AssemblingIntegrator(const std_cxx1x::function<void (CellInfo &)> &cell_worker,
- const std_cxx1x::function<void (FaceInfo &)> &boundary_worker,
- const std_cxx1x::function<void (FaceInfo &, FaceInfo &)> &face_worker)
- :
- cell_worker (cell_worker),
- boundary_worker (boundary_worker),
- face_worker (face_worker)
- {}
-
-
- template <int dim, class ASSEMBLER>
- inline void
- AssemblingIntegrator<dim,ASSEMBLER>::cell(CellInfo& info)
- {
- cell_worker(info);
- ASSEMBLER::assemble(info);
- }
-
-
- template <int dim, class ASSEMBLER>
- inline void
- AssemblingIntegrator<dim,ASSEMBLER>::bdry(FaceInfo& info)
- {
- boundary_worker(info);
- ASSEMBLER::assemble(info);
- }
-
-
- template <int dim, class ASSEMBLER>
- inline void
- AssemblingIntegrator<dim,ASSEMBLER>::face(FaceInfo& info1,
- FaceInfo& info2)
- {
- face_worker(info1, info2);
- ASSEMBLER::assemble(info1, info2);
- }
}
+
DEAL_II_NAMESPACE_CLOSE
#endif
void initialize_local(MatrixBlock<FullMatrix<number> >& M,
const unsigned int row,
const unsigned int col);
-
+
public:
void initialize_numbers(const unsigned int n);
void initialize_vectors(const unsigned int n);
* provide row and column info.
*/
void initialize_matrices(unsigned int n, bool both);
-
+
/**
* Allocate a local matrix
* for each of the global
template <class MatrixPtr>
void initialize_matrices(const std::vector<MatrixPtr>& matrices,
bool both);
-
+
/**
* Reinitialize matrices for
* new cell. Resizes the
* them to zero.
*/
void reinit(const BlockIndices& local_sizes);
-
+
/**
* The local numbers,
* computed on a cell or on a
* face.
*/
std::vector<number> J;
-
+
/**
* The local vectors. This
* field is public, so that
* write to it.
*/
std::vector<BlockVector<number> > R;
-
+
/**
* The local matrices
* coupling degrees of
* first cell on a face.
*/
std::vector<MatrixBlock<FullMatrix<number> > > M1;
-
+
/**
* The local matrices
* coupling test functions on
std::vector<MatrixBlock<FullMatrix<number> > > M2;
};
-
+
/**
* Basic info class only containing information on geometry and
* degrees of freedom of the mesh object.
template<int dim, int spacedim = dim>
class DoFInfo : public LocalResults<double>
{
- public:
+ public:
/// The current cell
typename Triangulation<dim>::cell_iterator cell;
-
+
/// The current face
typename Triangulation<dim>::face_iterator face;
-
+
/**
* The number of the current
* face on the current cell.
* if the info object was
* initialized with a cell.
*/
-
+
unsigned int face_number;
/**
* The number of the current
/// The DoF indices of the current cell
std::vector<unsigned int> indices;
-
+
/**
* Constructor setting the
* #block_info pointer.
*/
template <class DH>
DoFInfo(const DH& dof_handler);
-
+
/**
* Set the current cell and
* fill #indices.
void reinit(const DHCellIterator& c,
const DHFaceIterator& f,
const unsigned int n);
-
+
/**
* Set the current subface
* and fill #indices if the
const unsigned int s);
const BlockIndices& local_indices() const;
-
-
+
+
/// The block structure of the system
SmartPointer<const BlockInfo,DoFInfo<dim,spacedim> > block_info;
-
+
private:
/// Fill index vector
void get_indices(const typename DoFHandler<dim, spacedim>::cell_iterator c);
-
+
/// Fill index vector with level indices
void get_indices(const typename MGDoFHandler<dim, spacedim>::cell_iterator c);
-
+
/// Auxiliary vector
std::vector<unsigned int> indices_org;
-
+
/**
* An auxiliary local
* BlockIndices object created
*/
BlockIndices aux_local_indices;
};
-
+
/**
* Class for objects handed to local integration functions.
*
* information to DoFInfo.
*/
IntegrationInfo(const BlockInfo& block_info);
-
+
/**
* Constructor forwarding
* information to DoFInfo.
*/
template <class DH>
IntegrationInfo(const DH& dof_handler);
-
+
/**
* Build all internal
* structures, in particular
* selector.
*/
void initialize_data(const boost::shared_ptr<VectorDataBase<dim,spacedim> > data);
-
+
/**
* Delete the data created by initialize().
*/
void clear();
-
+
/// This is true if we are assembling for multigrid
bool multigrid;
/// Access to finite element
* vector of elements.
*/
const FEVALUESBASE& fe() const;
-
+
/// Access to finite elements
/**
* This access function must
* @see DGBlockSplitApplication
*/
const FEVALUESBASE& fe(const unsigned int i) const;
-
+
/**
* The vector containing the
* values of finite element
* point.
*/
std::vector<std::vector<std::vector<double> > > values;
-
+
/**
* The vector containing the
* derivatives of finite
* point.
*/
std::vector<std::vector<std::vector<Tensor<1,dim> > > > gradients;
-
+
/**
* The vector containing the
* second derivatives of finite
* point.
*/
std::vector<std::vector<std::vector<Tensor<2,dim> > > > hessians;
-
+
/**
* Reinitialize internal data
* structures for use on a cell.
void reinit(const DHCellIterator& c,
const DHFaceIterator& f,
const unsigned int fn);
-
+
/**
* Reinitialize internal data
* structures for use on a subface.
*/
template <typename T>
IntegrationInfoBox(const T&);
-
- template <class WORKER>
+
+ template <class WORKER, class ASSEMBLER>
void initialize(const WORKER&,
+ ASSEMBLER &assembler,
const FiniteElement<dim, spacedim>& el,
const Mapping<dim, spacedim>& mapping);
const FiniteElement<dim, spacedim>& el,
const Mapping<dim, spacedim>& mapping,
const NamedData<VECTOR*>& data);
-
+
// private:
boost::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > cell_data;
boost::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > bdry_data;
boost::shared_ptr<MeshWorker::VectorDataBase<dim, spacedim> > face_data;
-
+
CellInfo cell_info;
FaceInfo bdry_info;
FaceInfo face_info;
FaceInfo subface_info;
FaceInfo neighbor_info;
};
-
+
//----------------------------------------------------------------------//
J.resize(n);
}
-
+
template <typename number>
inline void
LocalResults<number>::initialize_vectors(const unsigned int n)
R.resize(n);
}
-
+
template <typename number>
template <class MatrixPtr>
inline void
}
}
}
-
-
+
+
template <typename number>
inline void
LocalResults<number>::initialize_matrices(const unsigned int n,
}
}
}
-
-
+
+
template <typename number>
inline void
LocalResults<number>::reinit(const BlockIndices& bi)
for (unsigned int i=0;i<M2.size();++i)
M2[i].matrix.reinit(bi.block_size(M2[i].row),
bi.block_size(M2[i].column));
- }
-
+ }
+
//----------------------------------------------------------------------//
-
+
template <int dim, int spacedim>
template <class DH>
DoFInfo<dim,spacedim>::DoFInfo(const DH& dof_handler)
aux[0] = dof_handler.get_fe().dofs_per_cell;
aux_local_indices.reinit(aux);
}
-
-
+
+
template <int dim, int spacedim>
template <class DHCellIterator>
inline void
if (block_info)
LocalResults<double>::reinit(block_info->local());
else
- LocalResults<double>::reinit(aux_local_indices);
+ LocalResults<double>::reinit(aux_local_indices);
}
-
-
+
+
template<int dim, int spacedim>
template <class DHCellIterator, class DHFaceIterator>
inline void
DoFInfo<dim, spacedim>::reinit(const DHCellIterator& c,
const DHFaceIterator& f,
unsigned int n)
- {
+ {
if ((cell.state() != IteratorState::valid)
|| cell != static_cast<typename Triangulation<dim>::cell_iterator> (c))
get_indices(c);
if (block_info)
LocalResults<double>::reinit(block_info->local());
else
- LocalResults<double>::reinit(aux_local_indices);
+ LocalResults<double>::reinit(aux_local_indices);
}
-
-
+
+
template<int dim, int spacedim>
template <class DHCellIterator, class DHFaceIterator>
inline void
if (block_info)
LocalResults<double>::reinit(block_info->local());
else
- LocalResults<double>::reinit(aux_local_indices);
+ LocalResults<double>::reinit(aux_local_indices);
}
-
-
+
+
template<int dim, int spacedim>
inline const BlockIndices&
DoFInfo<dim, spacedim>::local_indices() const
return block_info->local();
return aux_local_indices;
}
-
+
//----------------------------------------------------------------------//
-
+
template <int dim, class FVB, int spacedim>
template <class DH>
IntegrationInfo<dim,FVB,spacedim>::IntegrationInfo(const DH& dof_handler)
global_data(boost::shared_ptr<VectorDataBase<dim, spacedim> >(new VectorDataBase<dim, spacedim>))
{}
-
+
template <int dim, class FVB, int spacedim>
inline const FVB&
IntegrationInfo<dim,FVB,spacedim>::fe() const
const bool split_fevalues = this->block_info != 0;
fill_local_data(split_fevalues);
}
-
-
+
+
template <int dim, class FVB, int spacedim>
template <class DHCellIterator, class DHFaceIterator>
inline void
FEFaceValues<dim>& fe = dynamic_cast<FEFaceValues<dim>&> (febase);
fe.reinit(this->cell, fn);
}
-
+
const bool split_fevalues = this->block_info != 0;
fill_local_data(split_fevalues);
}
-
-
+
+
template <int dim, class FVB, int spacedim>
template <class DHCellIterator, class DHFaceIterator>
inline void
FESubfaceValues<dim>& fe = dynamic_cast<FESubfaceValues<dim>&> (febase);
fe.reinit(this->cell, fn, sn);
}
-
+
const bool split_fevalues = this->block_info != 0;
fill_local_data(split_fevalues);
}
subface_info(t),
neighbor_info(t)
{}
-
-
+
+
template <int dim, int sdim>
- template <class WORKER>
+ template <class WORKER, class ASSEMBLER>
void
- IntegrationInfoBox<dim,sdim>::initialize(
- const WORKER& integrator,
- const FiniteElement<dim,sdim>& el,
- const Mapping<dim,sdim>& mapping)
+ IntegrationInfoBox<dim,sdim>::
+ initialize(const WORKER& integrator,
+ ASSEMBLER &assembler,
+ const FiniteElement<dim,sdim>& el,
+ const Mapping<dim,sdim>& mapping)
{
- integrator.initialize_info(cell_info, false);
- integrator.initialize_info(bdry_info, false);
- integrator.initialize_info(face_info, true);
- integrator.initialize_info(subface_info, true);
- integrator.initialize_info(neighbor_info, true);
-
+ assembler.initialize_info(cell_info, false);
+ assembler.initialize_info(bdry_info, false);
+ assembler.initialize_info(face_info, true);
+ assembler.initialize_info(subface_info, true);
+ assembler.initialize_info(neighbor_info, true);
+
cell_info.initialize<FEValues<dim,sdim> >(el, mapping, integrator.cell_quadrature,
integrator.cell_flags);
bdry_info.initialize<FEFaceValues<dim,sdim> >(el, mapping, integrator.bdry_quadrature,
{
initialize(integrator, el, mapping);
boost::shared_ptr<VectorData<VECTOR, dim, sdim> > p;
-
+
p = boost::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (integrator.cell_selector));
p->initialize(data);
cell_data = p;
cell_info.initialize_data(p);
-
+
p = boost::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (integrator.bdry_selector));
p->initialize(data);
bdry_data = p;
bdry_info.initialize_data(p);
-
+
p = boost::shared_ptr<VectorData<VECTOR, dim, sdim> >(new VectorData<VECTOR, dim, sdim> (integrator.face_selector));
p->initialize(data);
face_data = p;
* @ingroup MeshWorker
* @author Guido Kanschat, 2009
*/
- template<class CELLINFO, class FACEINFO, class ITERATOR>
+ template<class CELLINFO, class FACEINFO, class ITERATOR, typename ASSEMBLER>
void loop(ITERATOR begin,
typename identity<ITERATOR>::type end,
CELLINFO& cell_info, FACEINFO& boundary_info,
const std_cxx1x::function<void (CELLINFO &)> &cell_worker,
const std_cxx1x::function<void (FACEINFO &)> &boundary_worker,
const std_cxx1x::function<void (FACEINFO &, FACEINFO &)> &face_worker,
+ ASSEMBLER &assembler,
bool cells_first = true)
{
const bool integrate_cell = (cell_worker != 0);
{
cell_info.reinit(cell);
cell_worker(cell_info);
+ assembler.assemble (cell_info);
}
if (integrate_interior_face || integrate_boundary)
{
boundary_info.reinit(cell, face, face_no);
boundary_worker(boundary_info);
+ assembler.assemble (boundary_info);
}
}
else if (integrate_interior_face)
// conform to
// old version
face_worker(face_info, subface_info);
+ assembler.assemble (face_info, subface_info);
}
else
{
neighbor_face_no);
face_worker(face_info, neighbor_info);
+ assembler.assemble (face_info, neighbor_info);
neighbor->face(neighbor_face_no)->set_user_flag ();
}
-
}
} // faces
// Execute this, if faces
{
cell_info.reinit(cell);
cell_worker(cell_info);
+ assembler.assemble (cell_info);
}
}
}
* @ingroup MeshWorker
* @author Guido Kanschat, 2009
*/
- template<class CELLINFO, class FACEINFO, int dim, class ITERATOR>
+ template<class CELLINFO, class FACEINFO, int dim, class ITERATOR, typename ASSEMBLER>
void integration_loop(ITERATOR begin,
typename identity<ITERATOR>::type end,
IntegrationInfoBox<dim, dim>& box,
const std_cxx1x::function<void (CELLINFO &)> &cell_worker,
const std_cxx1x::function<void (FACEINFO &)> &boundary_worker,
const std_cxx1x::function<void (FACEINFO &, FACEINFO &)> &face_worker,
+ ASSEMBLER &assembler,
bool cells_first = true)
{
loop<CELLINFO,FACEINFO>
cell_worker,
boundary_worker,
face_worker,
+ assembler,
cells_first);
}
}
// result of std::bind if the local
// integrators were, for example,
// non-static member functions.
- typedef
- MeshWorker::AssemblingIntegrator
- <dim,
- MeshWorker::Assembler::SystemSimple<SparseMatrix<double>,
- Vector<double> > >
- Integrator;
- Integrator integrator(&DGMethod<dim>::integrate_cell_term,
- &DGMethod<dim>::integrate_boundary_term,
- &DGMethod<dim>::integrate_face_term);
+ MeshWorker::IntegrationWorker<dim> integration_worker;
// First, we initialize the
// quadrature formulae and the
// independently, we have to hand
// over this value three times.
const unsigned int n_gauss_points = dof_handler.get_fe().degree+1;
- integrator.initialize_gauss_quadrature(n_gauss_points,
- n_gauss_points,
- n_gauss_points);
+ integration_worker.initialize_gauss_quadrature(n_gauss_points,
+ n_gauss_points,
+ n_gauss_points);
// These are the types of values we
// need for integrating our
UpdateFlags update_flags = update_quadrature_points |
update_values |
update_gradients;
- integrator.add_update_flags(update_flags, true, true, true, true);
+ integration_worker.add_update_flags(update_flags, true, true, true, true);
// Finally, we have to tell the
// assembler base class where to
// put the local data. These will
// be our system matrix and the
// right hand side.
- integrator.initialize(system_matrix, right_hand_side);
+ MeshWorker::Assembler::SystemSimple<SparseMatrix<double>, Vector<double> >
+ assembler;
+ assembler.initialize(system_matrix, right_hand_side);
// We are now ready to get to the
// integration loop. @p info_box is
// receives all the stuff we
// created so far.
MeshWorker::IntegrationInfoBox<dim> info_box(dof_handler);
- info_box.initialize(integrator, fe, mapping);
+ info_box.initialize(integration_worker, assembler, fe, mapping);
// Finally, the integration loop
// over all active cells
MeshWorker::integration_loop<CellInfo,FaceInfo>
(dof_handler.begin_active(), dof_handler.end(),
info_box,
- std_cxx1x::bind (&Integrator::cell, integrator, _1),
- std_cxx1x::bind (&Integrator::bdry, integrator, _1),
- std_cxx1x::bind (&Integrator::face, integrator, _1, _2));
+ &DGMethod<dim>::integrate_cell_term,
+ &DGMethod<dim>::integrate_boundary_term,
+ &DGMethod<dim>::integrate_face_term,
+ assembler);
}