From: kanschat Date: Fri, 3 Dec 1999 02:49:15 +0000 (+0000) Subject: make_flux_sparsity_pattern in DofTools X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=f3d66220d952b655bc565c7e58b3796c67801906;p=dealii-svn.git make_flux_sparsity_pattern in DofTools git-svn-id: https://svn.dealii.org/trunk@1970 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/dofs/dof_handler.h b/deal.II/deal.II/include/dofs/dof_handler.h index 4eea08f89c..6bf15cb273 100644 --- a/deal.II/deal.II/include/dofs/dof_handler.h +++ b/deal.II/deal.II/include/dofs/dof_handler.h @@ -449,7 +449,7 @@ class DoFHandler : public Subscriptor, * have to compress the matrix yourself, * using #SparseMatrixStruct::compress()#. */ - void make_sparsity_pattern (SparseMatrixStruct &) const; + void make_sparsity_pattern (SparseMatrixStruct &) const; /** * This function does mistly the same as @@ -494,7 +494,7 @@ class DoFHandler : public Subscriptor, */ void make_sparsity_pattern (const vector > &mask, SparseMatrixStruct &) const; - + /** * Write the sparsity structure of the * matrix composed of the basis functions diff --git a/deal.II/deal.II/include/dofs/dof_tools.h b/deal.II/deal.II/include/dofs/dof_tools.h index f024890241..5a0e3afb56 100644 --- a/deal.II/deal.II/include/dofs/dof_tools.h +++ b/deal.II/deal.II/include/dofs/dof_tools.h @@ -24,6 +24,22 @@ class DoFTools { public: + /** + * Generate sparsity pattern for fluxes. + * This is a replacement of the + * function + * #make_sparsity_pattern# for + * discontinuous methods. Since + * the fluxes include couplings + * between neighboring elements, + * the normal couplings and these + * extra matrix entries are + * considered. + */ + template + static void make_flux_sparsity_pattern (const DoFHandler&, + SparseMatrixStruct &); + /** * Extract the indices of the degrees * of freedom belonging to certain diff --git a/deal.II/deal.II/source/dofs/dof_handler.cc b/deal.II/deal.II/source/dofs/dof_handler.cc index 7b480a0e4b..b7b9589a09 100644 --- a/deal.II/deal.II/source/dofs/dof_handler.cc +++ b/deal.II/deal.II/source/dofs/dof_handler.cc @@ -1974,7 +1974,8 @@ void DoFHandler<3>::make_hanging_node_constraints (ConstraintMatrix &constraints template -void DoFHandler::make_sparsity_pattern (SparseMatrixStruct &sparsity) const { +void DoFHandler::make_sparsity_pattern (SparseMatrixStruct &sparsity) const +{ Assert (selected_fe != 0, ExcNoFESelected()); Assert (sparsity.n_rows() == n_dofs(), ExcDifferentDimensions (sparsity.n_rows(), n_dofs())); @@ -1993,10 +1994,8 @@ void DoFHandler::make_sparsity_pattern (SparseMatrixStruct &sparsity) const for (unsigned int j=0; j void diff --git a/deal.II/deal.II/source/dofs/dof_tools.cc b/deal.II/deal.II/source/dofs/dof_tools.cc index 7847191faa..8a30654c59 100644 --- a/deal.II/deal.II/source/dofs/dof_tools.cc +++ b/deal.II/deal.II/source/dofs/dof_tools.cc @@ -9,8 +9,76 @@ #include #include #include +#include +#if deal_II_dimension == 1 +//TODO: Implement make_flux_sparsity_pattern in 1d + +template <> +void +DoFTools::make_flux_sparsity_pattern (const DoFHandler<1>&, + SparseMatrixStruct &) +{ + Assert(false, ExcNotImplemented()); +} + +#else + +//TODO: Check this function for potential of optimization. (G) + +template +void +DoFTools::make_flux_sparsity_pattern (const DoFHandler& dof, + SparseMatrixStruct &sparsity) +{ + const unsigned int n_dofs = dof.n_dofs(); + + Assert (sparsity.n_rows() == n_dofs, + ExcDimensionMismatch (sparsity.n_rows(), n_dofs)); + Assert (sparsity.n_cols() == n_dofs, + ExcDimensionMismatch (sparsity.n_cols(), n_dofs)); + + const unsigned int total_dofs = dof.get_fe().dofs_per_cell; + vector dofs_on_this_cell(total_dofs); + vector dofs_on_other_cell(total_dofs); + DoFHandler::active_cell_iterator cell = dof.begin_active(), + endc = dof.end(); + for (; cell!=endc; ++cell) + { + cell->get_dof_indices (dofs_on_this_cell); + // make sparsity pattern for this cell + for (unsigned int i=0; i::faces_per_cell; + ++face) + { + if (cell->face(face)->boundary_indicator() == 255) + { + DoFHandler::active_cell_iterator neighbor = cell->neighbor(face); + neighbor->get_dof_indices (dofs_on_other_cell); + // only add one direction + // The other is taken care of + // by neighbor. + for (unsigned int i=0; i void @@ -79,6 +147,10 @@ DoFTools::extract_level_dofs(const unsigned int level, // explicit instantiations +#if deal_II_dimension > 1 +template void DoFTools::make_flux_sparsity_pattern (const DoFHandler& dof, + SparseMatrixStruct &sparsity); +#endif template void DoFTools::extract_dofs(const DoFHandler& dof, const vector& local_select, vector& selected_dofs);