std::vector<MatrixBlock<FullMatrix<number> > > M2;
};
+ template <int dim, int spacedim> class DoFInfoBox;
/**
* A class containing information on geometry and degrees of freedom
* the #aux_local_indices.
*/
template <class DH>
- DoFInfo(const DH& dof_handler);
-
+ DoFInfo (const DH& dof_handler);
+
/**
* Set the current cell and
* fill #indices.
SmartPointer<const BlockInfo,DoFInfo<dim,spacedim> > block_info;
private:
- /// Fill index vector with active indices
+ /**
+ * Standard constructor, not
+ * setting any block
+ * indices. Use of this
+ * constructor is not
+ * recommended, but it is
+ * needed for the arrays in
+ * DoFInfoBox.
+ */
+ DoFInfo ();
+
+ /// Fill index vector with active indices
void get_indices(const typename DoFHandler<dim, spacedim>::cell_iterator& c);
/// Fill index vector with level indices
* degrees of freedom per cell.
*/
BlockIndices aux_local_indices;
+
+ friend class DoFInfoBox<dim, spacedim>;
+ };
+
+
+ /**
+ * A class bundling the MeshWorker::DoFInfo objects used on a cell.
+ *
+ * @todo Currently, we are storing an object for the cells and two for
+ * each face. We could gather all face data pertaining to the cell
+ * itself in one object, saving a bit of memory and a few operations,
+ * but sacrificing some cleanliness.
+ *
+ * @author Guido Kanschat, 2010
+ */
+ template < int dim, int spacedim=dim>
+ class DoFInfoBox
+ {
+ public:
+ /**
+ * Constructor copying the seed
+ * into all other objects.
+ */
+ DoFInfoBox(const MeshWorker::DoFInfo<dim, spacedim>& seed);
+
+ /**
+ * Reset all the availability flags.
+ */
+ void reset();
+
+ /**
+ * The data for the cell.
+ */
+ MeshWorker::DoFInfo<dim> cell;
+ /**
+ * The data for the faces from inside.
+ */
+ MeshWorker::DoFInfo<dim> interior[GeometryInfo<dim>::faces_per_cell];
+ /**
+ * The data for the faces from outside.
+ */
+ MeshWorker::DoFInfo<dim> exterior[GeometryInfo<dim>::faces_per_cell];
+
+ /**
+ * A set of flags, indicating
+ * whether data on an interior
+ * face is available.
+ */
+ bool interior_face_available[GeometryInfo<dim>::faces_per_cell];
+ /**
+ * A set of flags, indicating
+ * whether data on an exterior
+ * face is available.
+ */
+ bool exterior_face_available[GeometryInfo<dim>::faces_per_cell];
};
+
+
+
/**
* Class for objects handed to local integration functions.
*
return aux_local_indices;
}
+//----------------------------------------------------------------------//
+
+ template < int dim, int spacedim>
+ inline
+ DoFInfoBox<dim, spacedim>::DoFInfoBox(const DoFInfo<dim, spacedim>& seed)
+ :
+ cell(seed)
+ {
+ for (unsigned int i=0;i<GeometryInfo<dim>::faces_per_cell;++i)
+ {
+ exterior[i] = seed;
+ interior[i] = seed;
+ interior_face_available[i] = false;
+ exterior_face_available[i] = false;
+ }
+ }
+
+ template < int dim, int spacedim>
+ inline void
+ DoFInfoBox<dim, spacedim>::reset ()
+ {
+ for (unsigned int i=0;i<GeometryInfo<dim>::faces_per_cell;++i)
+ {
+ interior_face_available[i] = false;
+ exterior_face_available[i] = false;
+ }
+ }
+
//----------------------------------------------------------------------//
template <int dim, class FVB, int spacedim>
* @ingroup MeshWorker
* @author Guido Kanschat, 2009
*/
- template<class DOFINFO, class INFOBOX, class ITERATOR, typename ASSEMBLER>
+ template<class INFOBOX, int dim, int spacedim, typename ASSEMBLER, class ITERATOR>
void loop(ITERATOR begin,
typename identity<ITERATOR>::type end,
- DOFINFO cell_dof_info,
+ DoFInfo<dim, spacedim> dinfo,
INFOBOX& info,
- const std_cxx1x::function<void (DOFINFO&, typename INFOBOX::CellInfo &)> &cell_worker,
- const std_cxx1x::function<void (DOFINFO&, typename INFOBOX::FaceInfo &)> &boundary_worker,
- const std_cxx1x::function<void (DOFINFO&, DOFINFO&,
+ const std_cxx1x::function<void (DoFInfo<dim, spacedim>&, typename INFOBOX::CellInfo &)> &cell_worker,
+ const std_cxx1x::function<void (DoFInfo<dim, spacedim>&, typename INFOBOX::FaceInfo &)> &boundary_worker,
+ const std_cxx1x::function<void (DoFInfo<dim, spacedim>&, DoFInfo<dim, spacedim>&,
typename INFOBOX::FaceInfo &,
typename INFOBOX::FaceInfo &)> &face_worker,
ASSEMBLER &assembler,
const bool integrate_boundary = (boundary_worker != 0);
const bool integrate_interior_face = (face_worker != 0);
- ;
- DOFINFO face_dof_info = cell_dof_info;
- DOFINFO neighbor_dof_info = cell_dof_info;
- assembler.initialize_info(cell_dof_info, false);
- assembler.initialize_info(face_dof_info, true);
- assembler.initialize_info(neighbor_dof_info, true);
+ DoFInfoBox<dim, spacedim> dof_info(dinfo);
+
+ assembler.initialize_info(dof_info.cell, false);
+ for (unsigned int i=0;i<GeometryInfo<dim>::faces_per_cell;++i)
+ {
+ assembler.initialize_info(dof_info.interior[i], true);
+ assembler.initialize_info(dof_info.exterior[i], true);
+ }
// Loop over all cells
for (ITERATOR cell = begin; cell != end; ++cell)
// before faces
if (integrate_cell && cells_first)
{
- cell_dof_info.reinit(cell);
- info.cell.reinit(cell_dof_info);
- cell_worker(cell_dof_info, info.cell);
- assembler.assemble(cell_dof_info);
+ dof_info.cell.reinit(cell);
+ info.cell.reinit(dof_info.cell);
+ cell_worker(dof_info.cell, info.cell);
+ assembler.assemble(dof_info.cell);
}
if (integrate_interior_face || integrate_boundary)
{
if (integrate_boundary)
{
- cell_dof_info.reinit(cell, face, face_no);
- info.boundary.reinit(cell_dof_info);
- boundary_worker(cell_dof_info, info.boundary);
- assembler.assemble(cell_dof_info);
+ dof_info.interior_face_available[face_no] = true;
+ dof_info.interior[face_no].reinit(cell, face, face_no);
+ info.boundary.reinit(dof_info.interior[face_no]);
+ boundary_worker(dof_info.interior[face_no], info.boundary);
+ assembler.assemble(dof_info.interior[face_no]);
}
}
else if (integrate_interior_face)
typename ITERATOR::AccessorType::Container::face_iterator nface
= neighbor->face(neighbor_face_no.first);
- face_dof_info.reinit(cell, face, face_no);
- info.face.reinit(face_dof_info);
- neighbor_dof_info.reinit(neighbor, nface,
- neighbor_face_no.first, neighbor_face_no.second);
- info.subface.reinit(neighbor_dof_info);
+ dof_info.interior_face_available[face_no] = true;
+ dof_info.exterior_face_available[face_no] = true;
+ dof_info.interior[face_no].reinit(cell, face, face_no);
+ info.face.reinit(dof_info.interior[face_no]);
+ dof_info.exterior[face_no].reinit(
+ neighbor, nface, neighbor_face_no.first, neighbor_face_no.second);
+ info.subface.reinit(dof_info.exterior[face_no]);
+
// Neighbor
// first to
// conform to
// old version
- face_worker(face_dof_info, neighbor_dof_info, info.face, info.subface);
- assembler.assemble (face_dof_info, neighbor_dof_info);
+ face_worker(dof_info.interior[face_no], dof_info.exterior[face_no],
+ info.face, info.subface);
+ assembler.assemble (dof_info.interior[face_no], dof_info.exterior[face_no]);
}
else
{
unsigned int neighbor_face_no = cell->neighbor_of_neighbor(face_no);
Assert (neighbor->face(neighbor_face_no) == face, ExcInternalError());
// Regular interior face
- face_dof_info.reinit(cell, face, face_no);
- info.face.reinit(face_dof_info);
- neighbor_dof_info.reinit(neighbor, neighbor->face(neighbor_face_no),
- neighbor_face_no);
- info.neighbor.reinit(neighbor_dof_info);
- face_worker(face_dof_info, neighbor_dof_info, info.face, info.neighbor);
- assembler.assemble(face_dof_info, neighbor_dof_info);
+ dof_info.interior_face_available[face_no] = true;
+ dof_info.exterior_face_available[face_no] = true;
+ dof_info.interior[face_no].reinit(cell, face, face_no);
+ info.face.reinit(dof_info.interior[face_no]);
+ dof_info.exterior[face_no].reinit(
+ neighbor, neighbor->face(neighbor_face_no), neighbor_face_no);
+ info.neighbor.reinit(dof_info.exterior[face_no]);
+
+ face_worker(dof_info.interior[face_no], dof_info.exterior[face_no],
+ info.face, info.neighbor);
+ assembler.assemble(dof_info.interior[face_no], dof_info.exterior[face_no]);
}
}
} // faces
// have to be handled first
if (integrate_cell && !cells_first)
{
- cell_dof_info.reinit(cell);
- info.cell.reinit(cell_dof_info);
- cell_worker(cell_dof_info, info.cell);
- assembler.assemble (cell_dof_info);
+ dof_info.cell.reinit(cell);
+ info.cell.reinit(dof_info.cell);
+ cell_worker(dof_info.cell, info.cell);
+ assembler.assemble(dof_info.cell);
}
}
}
MeshWorker::IntegrationInfoBox<dim> info_box;
MeshWorker::DoFInfo<dim> dof_info(dof_handler);
info_box.initialize(integration_worker, fe, mapping);
- MeshWorker::loop<MeshWorker::DoFInfo<dim>, MeshWorker::IntegrationInfoBox<dim> >(
+ MeshWorker::loop<MeshWorker::IntegrationInfoBox<dim>, dim, dim>(
dof_handler.begin_active(), dof_handler.end(),
dof_info, info_box,
&MatrixIntegrator<dim>::cell,
MeshWorker::IntegrationInfoBox<dim> info_box;
MeshWorker::DoFInfo<dim> dof_info(mg_dof_handler);
info_box.initialize(integration_worker, fe, mapping);
- MeshWorker::loop<MeshWorker::DoFInfo<dim>, MeshWorker::IntegrationInfoBox<dim> > (
+ MeshWorker::loop<MeshWorker::IntegrationInfoBox<dim>, dim, dim> (
mg_dof_handler.begin(), mg_dof_handler.end(),
dof_info, info_box,
&MatrixIntegrator<dim>::cell,
MeshWorker::DoFInfo<dim> dof_info(dof_handler);
info_box.initialize(integration_worker, fe, mapping);
- MeshWorker::loop<MeshWorker::DoFInfo<dim>, MeshWorker::IntegrationInfoBox<dim> >(
+ MeshWorker::loop<MeshWorker::IntegrationInfoBox<dim>, dim, dim>(
dof_handler.begin_active(), dof_handler.end(),
dof_info, info_box,
&RHSIntegrator<dim>::cell,
coarse_matrix.copy_from (mg_matrix[0]);
MGCoarseGridHouseholder<double, Vector<double> > mg_coarse;
mg_coarse.initialize(coarse_matrix);
- typedef PreconditionSOR<SparseMatrix<double> > RELAXATION;
+ typedef PreconditionSSOR<SparseMatrix<double> > RELAXATION;
MGSmootherRelaxation<SparseMatrix<double>, RELAXATION, Vector<double> >
mg_smoother(mem);
RELAXATION::AdditionalData smoother_data(1.);
MeshWorker::IntegrationInfoBox<dim> info_box;
MeshWorker::DoFInfo<dim> dof_info(dof_handler);
info_box.initialize(integration_worker, fe, mapping, solution_data);
- MeshWorker::loop<MeshWorker::DoFInfo<dim>, MeshWorker::IntegrationInfoBox<dim> > (
+ MeshWorker::loop<MeshWorker::IntegrationInfoBox<dim>, dim, dim> (
dof_handler.begin_active(), dof_handler.end(),
dof_info, info_box,
&Estimator<dim>::cell,