// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+// Copyright (C) 1998 - 2005 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
#include <set>
#include <map>
+template<int dim, class T> class Table;
class SparsityPattern;
template <typename number> class Vector;
-template <typename number> class FullMatrix;
template <int dim> class Function;
template <int dim> class Point;
+template <int dim> class FiniteElement;
template <int dim> class DoFHandler;
template <int dim> class MGDoFHandler;
class ConstraintMatrix;
* don't do it here) and not relevant in this context.
*
*
- * @author Wolfgang Bangerth and others, 1998, 1999, 2000, 2001
+ * @author Wolfgang Bangerth, Guido Kanschat and others, 1998 - 2005
*/
class DoFTools
{
public:
+ /**
+ * The flags used in tables by certain
+ * <tt>make_*_pattern</tt>
+ * functions to determine whether
+ * two components of the solution
+ * couple in the interior of mesh
+ * cells or at the boundary.
+ */
+ enum Coupling
+ {
+ /// The two components do not couple
+ none,
+ /// The two components do couple
+ always,
+/// The two components couple only if their shape functions can be nonzero on this face
+ nonzero
+ };
+
+
/**
* Locate non-zero entries of the
* system matrix.
static void
make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
SparsityPattern &sparsity,
- const FullMatrix<double>& int_mask,
- const FullMatrix<double>& flux_mask);
+ const Table<2,Coupling>& int_mask,
+ const Table<2,Coupling>& flux_mask);
/**
* Make up the constraints which
std::vector<std::map<unsigned int, float> > &weights,
const typename DoFHandler<dim>::active_cell_iterator &begin,
const typename DoFHandler<dim>::active_cell_iterator &end);
+
+ /**
+ * Compute coupling of dofs from
+ * coupling of components.
+ */
+ template <int dim>
+ static void compute_dof_couplings (
+ Table<2,Coupling>& dof_couplings,
+ const Table<2,Coupling>& component_couplings,
+ const FiniteElement<dim>& fe);
+
+ friend class MGTools;
};
-//---------------------------- mg_dof_tools.h ---------------------------
+//---------------------------------------------------------------------------
// $Id$
// Version: $Name$
//
// to the file deal.II/doc/license.html for the text and
// further information on this license.
//
-//---------------------------- mg_dof_tools.h ---------------------------
+//---------------------------------------------------------------------------
#ifndef __deal2__mg_dof_tools_h
#define __deal2__mg_dof_tools_h
#include <base/config.h>
+#include <dofs/dof_tools.h>
#include <multigrid/mg_dof_handler.h>
#include <vector>
* for more information.
*
* All member functions are static, so there is no need to create an
- * object of class @p MGTools.
+ * object of class MGTools.
*
- * @author Wolfgang Bangerth, Guido Kanschat, 1999-2003
+ * @author Wolfgang Bangerth, Guido Kanschat, 1999 - 2005
*/
class MGTools
{
make_flux_sparsity_pattern (const MGDoFHandler<dim> &dof,
SparsityPattern &sparsity,
const unsigned int level,
- const FullMatrix<double> &int_mask,
- const FullMatrix<double> &flux_mask);
+ const Table<2,DoFTools::Coupling> &int_mask,
+ const Table<2,DoFTools::Coupling> &flux_mask);
+
+ /**
+ * Create sparsity pattern for
+ * the fluxes at refinement
+ * edges. The matrix maps a
+ * function of the fine level
+ * space @p level to the coarser
+ * space. This is the version
+ * restricting the pattern to the
+ * elements actually needed.
+ *
+ * make_flux_sparsity_pattern()
+ */
+ template <int dim, class SparsityPattern>
+ static void
+ make_flux_sparsity_pattern_edge (const MGDoFHandler<dim> &dof_handler,
+ SparsityPattern &sparsity,
+ const unsigned int level,
+ const Table<2,DoFTools::Coupling> &flux_mask);
/**
* Count the dofs component-wise
// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004, 2005 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
#include <base/multithread_info.h>
#include <base/thread_management.h>
#include <base/quadrature_lib.h>
+#include <base/table.h>
#include <grid/tria.h>
#include <grid/tria_iterator.h>
#include <grid/intergrid_map.h>
template <int dim, class SparsityPattern>
void
-DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim> &dof,
- SparsityPattern &sparsity,
- const FullMatrix<double> &int_mask,
- const FullMatrix<double> &flux_mask)
+DoFTools::make_flux_sparsity_pattern (const DoFHandler<dim>& dof,
+ SparsityPattern& sparsity,
+ const Table<2,Coupling>& int_mask,
+ const Table<2,Coupling>& flux_mask)
{
const unsigned int n_dofs = dof.n_dofs();
const FiniteElement<dim>& fe = dof.get_fe();
ExcDimensionMismatch (sparsity.n_rows(), n_dofs));
Assert (sparsity.n_cols() == n_dofs,
ExcDimensionMismatch (sparsity.n_cols(), n_dofs));
- Assert (int_mask.m() == n_comp,
- ExcDimensionMismatch (int_mask.m(), n_comp));
- Assert (int_mask.n() == n_comp,
- ExcDimensionMismatch (int_mask.n(), n_comp));
- Assert (flux_mask.m() == n_comp,
- ExcDimensionMismatch (flux_mask.m(), n_comp));
- Assert (flux_mask.n() == n_comp,
- ExcDimensionMismatch (flux_mask.n(), n_comp));
+ Assert (int_mask.n_rows() == n_comp,
+ ExcDimensionMismatch (int_mask.n_rows(), n_comp));
+ Assert (int_mask.n_cols() == n_comp,
+ ExcDimensionMismatch (int_mask.n_cols(), n_comp));
+ Assert (flux_mask.n_rows() == n_comp,
+ ExcDimensionMismatch (flux_mask.n_rows(), n_comp));
+ Assert (flux_mask.n_cols() == n_comp,
+ ExcDimensionMismatch (flux_mask.n_cols(), n_comp));
const unsigned int total_dofs = fe.dofs_per_cell;
std::vector<unsigned int> dofs_on_this_cell(total_dofs);
typename DoFHandler<dim>::active_cell_iterator cell = dof.begin_active(),
endc = dof.end();
-
-
-//TODO[GK]: would someone explain what the letters 'f', 'a', and 'n' mean?
-//Even better would, of course, be to change this to an enum
- Table<2,char> int_dof_mask(total_dofs, total_dofs);
- Table<2,char> flux_dof_mask(total_dofs, total_dofs);
-
+
+ Table<2,Coupling> int_dof_mask(total_dofs, total_dofs);
+ Table<2,Coupling> flux_dof_mask(total_dofs, total_dofs);
+
+ compute_dof_couplings(int_dof_mask, int_mask, fe);
+ compute_dof_couplings(flux_dof_mask, flux_mask, fe);
+
for (unsigned int i=0; i<total_dofs; ++i)
- {
- const unsigned int ii
- = (dof.get_fe().is_primitive(i) ?
- dof.get_fe().system_to_component_index(i).first
- :
- (std::find (dof.get_fe().get_nonzero_components(i).begin(),
- dof.get_fe().get_nonzero_components(i).end(),
- true)
- -
- dof.get_fe().get_nonzero_components(i).begin())
- );
- Assert (ii < dof.get_fe().n_components(),
- ExcInternalError());
-
- for (unsigned int j=0; j<total_dofs; ++j)
- {
- const unsigned int jj
- = (dof.get_fe().is_primitive(j) ?
- dof.get_fe().system_to_component_index(j).first
- :
- (std::find (dof.get_fe().get_nonzero_components(j).begin(),
- dof.get_fe().get_nonzero_components(j).end(),
- true)
- -
- dof.get_fe().get_nonzero_components(j).begin())
- );
- Assert (jj < dof.get_fe().n_components(),
- ExcInternalError());
-
-//TODO[GK]: the documentation only says that non-zeroness counts, but here
-//individual values are tested. this should be documented. Even better, rather
-//than passing in a matrix of doubles, we should pass in a Table<2,some-enum>
- if (int_mask(ii,jj) != 0)
- int_dof_mask(i,j) = 'f';
- if (flux_mask(ii,jj) == 1.)
- flux_dof_mask(i,j) = 'a';
- if (flux_mask(ii,jj) == 2.)
- flux_dof_mask(i,j) = 'f';
- }
- for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell;++f)
- support_on_face(i,f) = fe.has_support_on_face(i,f);
- }
+ for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell;++f)
+ support_on_face(i,f) = fe.has_support_on_face(i,f);
// Clear user flags because we will
// need them. But first we save
// make sparsity pattern for this cell
for (unsigned int i=0; i<total_dofs; ++i)
for (unsigned int j=0; j<total_dofs; ++j)
- if (int_dof_mask(i,j))
+ if (int_dof_mask(i,j) != none)
sparsity.add (dofs_on_this_cell[i],
dofs_on_this_cell[j]);
{
const bool j_non_zero_i = support_on_face (j, face);
- if (flux_dof_mask(i,j) == 'f')
+ if (flux_dof_mask(i,j) == always)
sparsity.add (dofs_on_this_cell[i],
dofs_on_this_cell[j]);
- if (flux_dof_mask(i,j) == 'a')
- if (i_non_zero_i && j_non_zero_i)
- sparsity.add (dofs_on_this_cell[i],
- dofs_on_this_cell[j]);
+ if (flux_dof_mask(i,j) == nonzero
+ && i_non_zero_i && j_non_zero_i)
+ sparsity.add (dofs_on_this_cell[i],
+ dofs_on_this_cell[j]);
}
}
}
// by coarser cells
if (neighbor->level() < cell->level())
continue;
-
+
const unsigned int
neighbor_face = cell->neighbor_of_neighbor(face);
-
+
if (neighbor->has_children())
{
for (unsigned int sub_nr = 0;
{
const bool j_non_zero_i = support_on_face (j, face);
const bool j_non_zero_e =support_on_face (j, neighbor_face);
-
- if (flux_dof_mask(i,j)==0)
- flux_dof_mask(i,j) = 'n';
-
- if (flux_dof_mask(i,j) == 'f')
+ if (flux_dof_mask(i,j) == always)
{
sparsity.add (dofs_on_this_cell[i],
dofs_on_other_cell[j]);
sparsity.add (dofs_on_other_cell[i],
dofs_on_other_cell[j]);
}
- if (flux_dof_mask(i,j) == 'a')
+ if (flux_dof_mask(i,j) == nonzero)
{
if (i_non_zero_i && j_non_zero_e)
sparsity.add (dofs_on_this_cell[i],
dofs_on_other_cell[j]);
}
- if (flux_dof_mask(j,i) == 'f')
+ if (flux_dof_mask(j,i) == always)
{
sparsity.add (dofs_on_this_cell[j],
dofs_on_other_cell[i]);
sparsity.add (dofs_on_other_cell[j],
dofs_on_other_cell[i]);
}
- if (flux_dof_mask(j,i) == 'a')
+ if (flux_dof_mask(j,i) == nonzero)
{
if (j_non_zero_i && i_non_zero_e)
sparsity.add (dofs_on_this_cell[j],
{
const bool j_non_zero_i = support_on_face (j, face);
const bool j_non_zero_e = support_on_face (j, neighbor_face);
- if (flux_dof_mask(i,j)==0)
- flux_dof_mask(i,j) = 'n';
-
- if (flux_dof_mask(i,j) == 'f')
+ if (flux_dof_mask(i,j) == always)
{
sparsity.add (dofs_on_this_cell[i],
dofs_on_other_cell[j]);
sparsity.add (dofs_on_other_cell[i],
dofs_on_other_cell[j]);
}
- if (flux_dof_mask(i,j) == 'a')
+ if (flux_dof_mask(i,j) == nonzero)
{
if (i_non_zero_i && j_non_zero_e)
sparsity.add (dofs_on_this_cell[i],
dofs_on_other_cell[j]);
}
- if (flux_dof_mask(j,i) == 'f')
+ if (flux_dof_mask(j,i) == always)
{
sparsity.add (dofs_on_this_cell[j],
dofs_on_other_cell[i]);
sparsity.add (dofs_on_other_cell[j],
dofs_on_other_cell[i]);
}
- if (flux_dof_mask(j,i) == 'a')
+ if (flux_dof_mask(j,i) == nonzero)
{
if (j_non_zero_i && i_non_zero_e)
sparsity.add (dofs_on_this_cell[j],
}
+template <int dim>
+void
+DoFTools::compute_dof_couplings(
+ Table<2,Coupling>& dof_couplings,
+ const Table<2,Coupling>& component_couplings,
+ const FiniteElement<dim>& fe)
+{
+ Assert(component_couplings.n_rows() == fe.n_components(),
+ ExcDimensionMismatch(component_couplings.n_rows(),
+ fe.n_components()));
+ Assert(component_couplings.n_cols() == fe.n_components(),
+ ExcDimensionMismatch(component_couplings.n_cols(),
+ fe.n_components()));
+
+ const unsigned int n_dofs = fe.dofs_per_cell;
+
+ Assert(dof_couplings.n_rows() == n_dofs,
+ ExcDimensionMismatch(dof_couplings.n_rows(), n_dofs));
+ Assert(dof_couplings.n_cols() == n_dofs,
+ ExcDimensionMismatch(dof_couplings.n_cols(), n_dofs));
+
+ for (unsigned int i=0; i<n_dofs; ++i)
+ {
+ const unsigned int ii
+ = (fe.is_primitive(i) ?
+ fe.system_to_component_index(i).first
+ :
+ (std::find (fe.get_nonzero_components(i).begin(),
+ fe.get_nonzero_components(i).end(),
+ true)
+ -
+ fe.get_nonzero_components(i).begin())
+ );
+ Assert (ii < fe.n_components(),
+ ExcInternalError());
+
+ for (unsigned int j=0; j<n_dofs; ++j)
+ {
+ const unsigned int jj
+ = (fe.is_primitive(j) ?
+ fe.system_to_component_index(j).first
+ :
+ (std::find (fe.get_nonzero_components(j).begin(),
+ fe.get_nonzero_components(j).end(),
+ true)
+ -
+ fe.get_nonzero_components(j).begin())
+ );
+ Assert (jj < fe.n_components(),
+ ExcInternalError());
+
+ dof_couplings(i,j) = component_couplings(ii,jj);
+ }
+ }
+}
+
// explicit instantiations
template void
DoFTools::make_flux_sparsity_pattern<deal_II_dimension,SparsityPattern>
(const DoFHandler<deal_II_dimension>& dof,
SparsityPattern &,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,Coupling>&,
+ const Table<2,Coupling>&);
template void
DoFTools::make_flux_sparsity_pattern<deal_II_dimension,CompressedSparsityPattern>
(const DoFHandler<deal_II_dimension>& dof,
CompressedSparsityPattern &,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,Coupling>&,
+ const Table<2,Coupling>&);
template void
DoFTools::make_flux_sparsity_pattern<deal_II_dimension,BlockSparsityPattern>
(const DoFHandler<deal_II_dimension>& dof,
BlockSparsityPattern &,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,Coupling>&,
+ const Table<2,Coupling>&);
template void
DoFTools::make_flux_sparsity_pattern<deal_II_dimension,CompressedBlockSparsityPattern>
(const DoFHandler<deal_II_dimension>& dof,
CompressedBlockSparsityPattern &,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,Coupling>&,
+ const Table<2,Coupling>&);
#endif
(const Mapping<deal_II_dimension> &mapping,
const DoFHandler<deal_II_dimension> &dof_handler,
std::vector<Point<deal_II_dimension> > &support_points);
+
+template
+void
+DoFTools::compute_dof_couplings<deal_II_dimension>(
+ Table<2,Coupling>& dof_couplings,
+ const Table<2,Coupling>& component_couplings,
+ const FiniteElement<deal_II_dimension>& fe);
-//---------------------------- mg_dof_tools.cc ---------------------------
+//---------------------------------------------------------------------------
// $Id$
// Version: $Name$
//
-// Copyright (C) 1998, 1999, 2000, 2001, 2002, 2003, 2004 by the deal.II authors
+// Copyright (C) 1998 - 2005 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
// to the file deal.II/doc/license.html for the text and
// further information on this license.
//
-//---------------------------- mg_dof_tools.cc ---------------------------
+//---------------------------------------------------------------------------
#include <base/multithread_info.h>
-//TODO[GK]: Replace FullMatrix by Table<2,char>
+
+
+
template<int dim, class SparsityPattern>
void
MGTools::make_flux_sparsity_pattern (
const MGDoFHandler<dim> &dof,
SparsityPattern &sparsity,
const unsigned int level,
- const FullMatrix<double> &int_mask,
- const FullMatrix<double> &flux_mask)
+ const Table<2,DoFTools::Coupling> &int_mask,
+ const Table<2,DoFTools::Coupling> &flux_mask)
{
+ const FiniteElement<dim>& fe = dof.get_fe();
const unsigned int n_dofs = dof.n_dofs(level);
- const unsigned int n_comp = dof.get_fe().n_components();
+ const unsigned int n_comp = fe.n_components();
Assert (sparsity.n_rows() == n_dofs,
ExcDimensionMismatch (sparsity.n_rows(), n_dofs));
Assert (sparsity.n_cols() == n_dofs,
ExcDimensionMismatch (sparsity.n_cols(), n_dofs));
- Assert (int_mask.m() == n_comp,
- ExcDimensionMismatch (int_mask.m(), n_comp));
- Assert (int_mask.n() == n_comp,
- ExcDimensionMismatch (int_mask.n(), n_comp));
- Assert (flux_mask.m() == n_comp,
- ExcDimensionMismatch (flux_mask.m(), n_comp));
- Assert (flux_mask.n() == n_comp,
- ExcDimensionMismatch (flux_mask.n(), n_comp));
+ Assert (int_mask.n_rows() == n_comp,
+ ExcDimensionMismatch (int_mask.n_rows(), n_comp));
+ Assert (int_mask.n_cols() == n_comp,
+ ExcDimensionMismatch (int_mask.n_cols(), n_comp));
+ Assert (flux_mask.n_rows() == n_comp,
+ ExcDimensionMismatch (flux_mask.n_rows(), n_comp));
+ Assert (flux_mask.n_cols() == n_comp,
+ ExcDimensionMismatch (flux_mask.n_cols(), n_comp));
- const unsigned int total_dofs = dof.get_fe().dofs_per_cell;
+ const unsigned int total_dofs = fe.dofs_per_cell;
std::vector<unsigned int> dofs_on_this_cell(total_dofs);
std::vector<unsigned int> dofs_on_other_cell(total_dofs);
+ Table<2,bool> support_on_face(total_dofs, GeometryInfo<dim>::faces_per_cell);
+
typename MGDoFHandler<dim>::cell_iterator cell = dof.begin(level),
endc = dof.end(level);
+ Table<2,DoFTools::Coupling> int_dof_mask(total_dofs, total_dofs);
+ Table<2,DoFTools::Coupling> flux_dof_mask(total_dofs, total_dofs);
- std::vector<std::vector<bool> > int_dof_mask(total_dofs,
- std::vector<bool>(total_dofs, false));
- std::vector<std::vector<bool> > flux_dof_mask(total_dofs,
- std::vector<bool>(total_dofs, false));
+ DoFTools::compute_dof_couplings(int_dof_mask, int_mask, fe);
+ DoFTools::compute_dof_couplings(flux_dof_mask, flux_mask, fe);
+
for (unsigned int i=0; i<total_dofs; ++i)
- for (unsigned int j=0; j<total_dofs; ++j)
- {
- unsigned int ii = dof.get_fe().system_to_component_index(i).first;
- unsigned int jj = dof.get_fe().system_to_component_index(j).first;
-
- if (int_mask(ii,jj) != 0)
- int_dof_mask[i][j] = true;
- if (flux_mask(ii,jj) != 0)
- flux_dof_mask[i][j] = true;
- }
+ for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell;++f)
+ support_on_face(i,f) = fe.has_support_on_face(i,f);
// Clear user flags because we will
// need them. But first we save
// make sparsity pattern for this cell
for (unsigned int i=0; i<total_dofs; ++i)
for (unsigned int j=0; j<total_dofs; ++j)
- if (int_dof_mask[i][j])
+ if (int_dof_mask[i][j] != DoFTools::none)
sparsity.add (dofs_on_this_cell[i],
dofs_on_this_cell[j]);
if (cell_face->user_flag_set ())
continue;
- if (! cell->at_boundary (face) )
+ if (cell->at_boundary (face) )
+ {
+ for (unsigned int i=0; i<total_dofs; ++i)
+ {
+ const bool i_non_zero_i = support_on_face (i, face);
+ for (unsigned int j=0; j<total_dofs; ++j)
+ {
+ const bool j_non_zero_i = support_on_face (j, face);
+
+ if (flux_dof_mask(i,j) == DoFTools::always)
+ sparsity.add (dofs_on_this_cell[i],
+ dofs_on_this_cell[j]);
+ if (flux_dof_mask(i,j) == DoFTools::nonzero
+ && i_non_zero_i && j_non_zero_i)
+ sparsity.add (dofs_on_this_cell[i],
+ dofs_on_this_cell[j]);
+ }
+ }
+ }
+ else
{
typename MGDoFHandler<dim>::cell_iterator
neighbor = cell->neighbor(face);
-
+
if (neighbor->level() < cell->level())
continue;
-
+
unsigned int neighbor_face = cell->neighbor_of_neighbor(face);
-
+
neighbor->get_mg_dof_indices (dofs_on_other_cell);
for (unsigned int i=0; i<total_dofs; ++i)
{
+ const bool i_non_zero_i = support_on_face (i, face);
+ const bool i_non_zero_e = support_on_face (i, neighbor_face);
for (unsigned int j=0; j<total_dofs; ++j)
- if (flux_dof_mask[i][j])
- {
- sparsity.add (dofs_on_this_cell[i],
- dofs_on_other_cell[j]);
- sparsity.add (dofs_on_other_cell[i],
- dofs_on_this_cell[j]);
- }
+ {
+ const bool j_non_zero_i = support_on_face (j, face);
+ const bool j_non_zero_e = support_on_face (j, neighbor_face);
+ if (flux_dof_mask(i,j) == DoFTools::always)
+ {
+ sparsity.add (dofs_on_this_cell[i],
+ dofs_on_other_cell[j]);
+ sparsity.add (dofs_on_other_cell[i],
+ dofs_on_this_cell[j]);
+ sparsity.add (dofs_on_this_cell[i],
+ dofs_on_this_cell[j]);
+ sparsity.add (dofs_on_other_cell[i],
+ dofs_on_other_cell[j]);
+ }
+ if (flux_dof_mask(i,j) == DoFTools::nonzero)
+ {
+ if (i_non_zero_i && j_non_zero_e)
+ sparsity.add (dofs_on_this_cell[i],
+ dofs_on_other_cell[j]);
+ if (i_non_zero_e && j_non_zero_i)
+ sparsity.add (dofs_on_other_cell[i],
+ dofs_on_this_cell[j]);
+ if (i_non_zero_i && j_non_zero_i)
+ sparsity.add (dofs_on_this_cell[i],
+ dofs_on_this_cell[j]);
+ if (i_non_zero_e && j_non_zero_e)
+ sparsity.add (dofs_on_other_cell[i],
+ dofs_on_other_cell[j]);
+ }
+
+ if (flux_dof_mask(j,i) == DoFTools::always)
+ {
+ sparsity.add (dofs_on_this_cell[j],
+ dofs_on_other_cell[i]);
+ sparsity.add (dofs_on_other_cell[j],
+ dofs_on_this_cell[i]);
+ sparsity.add (dofs_on_this_cell[j],
+ dofs_on_this_cell[i]);
+ sparsity.add (dofs_on_other_cell[j],
+ dofs_on_other_cell[i]);
+ }
+ if (flux_dof_mask(j,i) == DoFTools::nonzero)
+ {
+ if (j_non_zero_i && i_non_zero_e)
+ sparsity.add (dofs_on_this_cell[j],
+ dofs_on_other_cell[i]);
+ if (j_non_zero_e && i_non_zero_i)
+ sparsity.add (dofs_on_other_cell[j],
+ dofs_on_this_cell[i]);
+ if (j_non_zero_i && i_non_zero_i)
+ sparsity.add (dofs_on_this_cell[j],
+ dofs_on_this_cell[i]);
+ if (j_non_zero_e && i_non_zero_e)
+ sparsity.add (dofs_on_other_cell[j],
+ dofs_on_other_cell[i]);
+ }
+ }
}
neighbor->face(neighbor_face)->set_user_flag ();
}
}
}
-
+
// finally restore the user flags
const_cast<Triangulation<dim> &>(dof.get_tria()).load_user_flags(user_flags);
}
+template<int dim, class SparsityPattern>
+void
+MGTools::make_flux_sparsity_pattern_edge (
+ const MGDoFHandler<dim> &dof,
+ SparsityPattern &sparsity,
+ const unsigned int level,
+ const Table<2,DoFTools::Coupling> &flux_mask)
+{
+ const FiniteElement<dim>& fe = dof.get_fe();
+ const unsigned int n_comp = fe.n_components();
+
+ Assert ((level>=1) && (level<dof.get_tria().n_levels()),
+ ExcIndexRange(level, 1, dof.get_tria().n_levels()));
+
+ const unsigned int fine_dofs = dof.n_dofs(level);
+ const unsigned int coarse_dofs = dof.n_dofs(level-1);
+
+ // Matrix maps from fine level to coarse level
+
+ Assert (sparsity.n_rows() == coarse_dofs,
+ ExcDimensionMismatch (sparsity.n_rows(), coarse_dofs));
+ Assert (sparsity.n_cols() == fine_dofs,
+ ExcDimensionMismatch (sparsity.n_cols(), fine_dofs));
+ Assert (flux_mask.n_rows() == n_comp,
+ ExcDimensionMismatch (flux_mask.n_rows(), n_comp));
+ Assert (flux_mask.n_cols() == n_comp,
+ ExcDimensionMismatch (flux_mask.n_cols(), n_comp));
+
+ const unsigned int dofs_per_cell = dof.get_fe().dofs_per_cell;
+ std::vector<unsigned int> dofs_on_this_cell(dofs_per_cell);
+ std::vector<unsigned int> dofs_on_other_cell(dofs_per_cell);
+ Table<2,bool> support_on_face(total_dofs, GeometryInfo<dim>::faces_per_cell);
+
+ typename MGDoFHandler<dim>::cell_iterator cell = dof.begin(level),
+ endc = dof.end(level);
+
+ Table<2,DoFTools::Coupling> flux_dof_mask(total_dofs, total_dofs);
+
+ DoFTools::compute_dof_couplings(flux_dof_mask, flux_mask, fe);
+
+ for (unsigned int i=0; i<total_dofs; ++i)
+ for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell;++f)
+ support_on_face(i,f) = fe.has_support_on_face(i,f);
+
+ for (; cell!=endc; ++cell)
+ {
+ cell->get_mg_dof_indices (dofs_on_this_cell);
+ // Loop over all interior neighbors
+ for (unsigned int face = 0;
+ face < GeometryInfo<dim>::faces_per_cell;
+ ++face)
+ {
+ // Neighbor is coarser
+
+ if ( (! cell->at_boundary(face)) &&
+ (static_cast<unsigned int>(cell->neighbor_level(face)) != level) )
+ {
+ typename MGDoFHandler<dim>::cell_iterator
+ neighbor = cell->neighbor(face);
+ neighbor->get_mg_dof_indices (dofs_on_other_cell);
+
+ for (unsigned int i=0; i<dofs_per_cell; ++i)
+ {
+ for (unsigned int j=0; j<dofs_per_cell; ++j)
+ {
+ if (flux_dof_mask(i,j) != DoFTools::none)
+ {
+ sparsity.add (dofs_on_other_cell[i],
+ dofs_on_this_cell[j]);
+ sparsity.add (dofs_on_other_cell[j],
+ dofs_on_this_cell[i]);
+ }
+ }
+ }
+ }
+ }
+ }
+}
+
+
+
template <int dim>
void
MGTools::count_dofs_per_component (
MGTools::make_flux_sparsity_pattern<deal_II_dimension> (const MGDoFHandler<deal_II_dimension> &,
SparsityPattern &,
const unsigned int,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,DoFTools::Coupling>&,
+ const Table<2,DoFTools::Coupling>&);
template void
MGTools::make_flux_sparsity_pattern<deal_II_dimension> (const MGDoFHandler<deal_II_dimension> &,
CompressedSparsityPattern &,
const unsigned int,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,DoFTools::Coupling>&,
+ const Table<2,DoFTools::Coupling>&);
template void
MGTools::make_flux_sparsity_pattern<deal_II_dimension> (const MGDoFHandler<deal_II_dimension> &,
BlockSparsityPattern &,
const unsigned int,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,DoFTools::Coupling>&,
+ const Table<2,DoFTools::Coupling>&);
template void
MGTools::make_flux_sparsity_pattern<deal_II_dimension> (const MGDoFHandler<deal_II_dimension> &,
CompressedBlockSparsityPattern &,
const unsigned int,
- const FullMatrix<double>&,
- const FullMatrix<double>&);
+ const Table<2,DoFTools::Coupling>&,
+ const Table<2,DoFTools::Coupling>&);
#endif
// $Id$
// Version: $Name$
//
-// Copyright (C) 2003, 2004 by the deal.II authors
+// Copyright (C) 2003, 2004, 2005 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
// something)
void
make_masks (const unsigned int n,
- FullMatrix<double> &m1,
- FullMatrix<double> &m2)
+ Table<2,DoFTools::Coupling> &m1,
+ Table<2,DoFTools::Coupling> &m2)
{
m1.reinit (n,n);
m2.reinit (n,n);
for (unsigned int i=0; i<n; ++i)
- m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = 1.;
+ m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = DoFTools::nonzero;
}
void
check_this (const DoFHandler<dim> &dof_handler)
{
- FullMatrix<double> mask_int, mask_ext;
+ Table<2,DoFTools::Coupling> mask_int;
+ Table<2,DoFTools::Coupling> mask_ext;
make_masks (dof_handler.get_fe().n_components(),
mask_int, mask_ext);
// $Id$
// Version: $Name$
//
-// Copyright (C) 2003, 2004 by the deal.II authors
+// Copyright (C) 2003, 2004, 2005 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
// something)
void
make_masks (const unsigned int n,
- FullMatrix<double> &m1,
- FullMatrix<double> &m2)
+ Table<2, DoFTools::Coupling> &m1,
+ Table<2, DoFTools::Coupling> &m2)
{
m1.reinit (n,n);
m2.reinit (n,n);
for (unsigned int i=0; i<n; ++i)
- m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = 1.;
+ m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = DoFTools::nonzero;
}
void
check_this (const DoFHandler<dim> &dof_handler)
{
- FullMatrix<double> mask_int, mask_ext;
+ Table<2, DoFTools::Coupling> mask_int;
+ Table<2, DoFTools::Coupling> mask_ext;
make_masks (dof_handler.get_fe().n_components(),
mask_int, mask_ext);
// $Id$
// Version: $Name$
//
-// Copyright (C) 2003, 2004 by the deal.II authors
+// Copyright (C) 2003, 2004, 2005 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
// something)
void
make_masks (const unsigned int n,
- FullMatrix<double> &m1,
- FullMatrix<double> &m2)
+ Table<2, DoFTools::Coupling> &m1,
+ Table<2, DoFTools::Coupling> &m2)
{
m1.reinit (n,n);
m2.reinit (n,n);
for (unsigned int i=0; i<n; ++i)
- m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = 1.;
+ m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = DoFTools::nonzero;
}
if (dof_handler.get_fe().is_primitive() != true)
return;
- FullMatrix<double> mask_int, mask_ext;
+ Table<2, DoFTools::Coupling> mask_int;
+ Table<2, DoFTools::Coupling> mask_ext;
make_masks (dof_handler.get_fe().n_components(),
mask_int, mask_ext);
// $Id$
// Version: $Name$
//
-// Copyright (C) 2003, 2004 by the deal.II authors
+// Copyright (C) 2003, 2004, 2005 by the deal.II authors
//
// This file is subject to QPL and may not be distributed
// without copyright and license information. Please refer
// something)
void
make_masks (const unsigned int n,
- FullMatrix<double> &m1,
- FullMatrix<double> &m2)
+ Table<2, DoFTools::Coupling> &m1,
+ Table<2, DoFTools::Coupling> &m2)
{
m1.reinit (n,n);
m2.reinit (n,n);
for (unsigned int i=0; i<n; ++i)
- m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = 1.;
+ m1(i,0) = m1(0,i) = m2(i,0) = m2(0,i) = m1(i,i) = m2(i,i) = DoFTools::nonzero;
}
if (dof_handler.get_fe().is_primitive() != true)
return;
- FullMatrix<double> mask_int, mask_ext;
+ Table<2, DoFTools::Coupling> mask_int;
+ Table<2, DoFTools::Coupling> mask_ext;
make_masks (dof_handler.get_fe().n_components(),
mask_int, mask_ext);