From: Wolfgang Bangerth Date: Mon, 8 Mar 2004 17:02:30 +0000 (+0000) Subject: Work around a nasty bug in gcc 2.95 that doesn't allow us to specialize an outer... X-Git-Tag: v8.0.0~15667 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c027664038e9e0ff36c9e661dcfb519069b058d6;p=dealii.git Work around a nasty bug in gcc 2.95 that doesn't allow us to specialize an outer class template while retaining a member function template unspecialized. git-svn-id: https://svn.dealii.org/trunk@8679 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/deal.II/include/numerics/error_estimator.h b/deal.II/deal.II/include/numerics/error_estimator.h index d2ff000c1a..34fdf6293f 100644 --- a/deal.II/deal.II/include/numerics/error_estimator.h +++ b/deal.II/deal.II/include/numerics/error_estimator.h @@ -622,4 +622,184 @@ class KellyErrorEstimator +/** + * This is a specialization of the general template for 1d. The implementation + * is sufficiently different for 1d to justify this specialization. The basic + * difference between 1d and all other space dimensions is that in 1d, there + * are no faces of cells, just the vertices between line segments, so we have + * to compute the jump terms differently. However, this class offers exactly + * the same public functions as the general template, so that a user will not + * see any difference. + * + * @author Wolfgang Bangerth, 1998, 2004. + */ +template <> +class KellyErrorEstimator<1> +{ + public: + /** + * Implementation of the error + * estimator described above. You + * may give a coefficient, but + * there is a default value which + * denotes the constant + * coefficient with value + * one. The coefficient function + * may either be a scalar one, in + * which case it is used for all + * components of the finite + * element, or a vector-valued + * one with as many components as + * there are in the finite + * element; in the latter case, + * each component is weighted by + * the respective component in + * the coefficient. + * + * You might give a list of + * components you want to + * evaluate, in case the finite + * element used by the + * @ref{DoFHandler} object is + * vector-valued. You then have + * to set those entries to true + * in the bit-vector + * @p{component_mask} for which the + * respective component is to be + * used in the error + * estimator. The default is to + * use all components, which is + * done by either providing a + * bit-vector with all-set + * entries, or an empty + * bit-vector. + * + * The estimator supports + * multithreading and splits the + * cells to + * @p{multithread_info.n_default_threads} + * (default) threads. The number + * of threads to be used in + * multithreaded mode can be set + * with the last parameter of the + * error estimator. + * Multithreading is only + * implemented in two and three + * dimensions. + */ + template + static void estimate (const Mapping<1> &mapping, + const DoFHandler<1> &dof, + const Quadrature<0> &quadrature, + const FunctionMap<1>::type &neumann_bc, + const InputVector &solution, + Vector &error, + const std::vector &component_mask = std::vector(), + const Function<1> *coefficients = 0, + const unsigned int n_threads = multithread_info.n_default_threads); + + /** + * Calls the @p{estimate} + * function, see above, with + * @p{mapping=MappingQ1<1>()}. + */ + template + static void estimate (const DoFHandler<1> &dof, + const Quadrature<0> &quadrature, + const FunctionMap<1>::type &neumann_bc, + const InputVector &solution, + Vector &error, + const std::vector &component_mask = std::vector(), + const Function<1> *coefficients = 0, + const unsigned int n_threads = multithread_info.n_default_threads); + + /** + * Same function as above, but + * accepts more than one solution + * vectors and returns one error + * vector for each solution + * vector. For the reason of + * existence of this function, + * see the general documentation + * of this class. + * + * Since we do not want to force + * the user of this function to + * copy around their solution + * vectors, the vector of + * solution vectors takes + * pointers to the solutions, + * rather than being a vector of + * vectors. This makes it simpler + * to have the solution vectors + * somewhere in memory, rather + * than to have them collected + * somewhere special. (Note that + * it is not possible to + * construct of vector of + * references, so we had to use a + * vector of pointers.) + */ + template + static void estimate (const Mapping<1> &mapping, + const DoFHandler<1> &dof, + const Quadrature<0> &quadrature, + const FunctionMap<1>::type &neumann_bc, + const std::vector &solutions, + std::vector*> &errors, + const std::vector &component_mask = std::vector(), + const Function<1> *coefficients = 0, + const unsigned int n_threads = multithread_info.n_default_threads); + + /** + * Calls the @p{estimate} + * function, see above, with + * @p{mapping=MappingQ1<1>()}. + */ + template + static void estimate (const DoFHandler<1> &dof, + const Quadrature<0> &quadrature, + const FunctionMap<1>::type &neumann_bc, + const std::vector &solutions, + std::vector*> &errors, + const std::vector &component_mask = std::vector(), + const Function<1> *coefficients = 0, + const unsigned int n_threads = multithread_info.n_default_threads); + + + /** + * Exception + */ + DeclException0 (ExcInvalidBoundaryIndicator); + /** + * Exception + */ + DeclException0 (ExcInvalidComponentMask); + /** + * Exception + */ + DeclException0 (ExcInvalidCoefficient); + /** + * Exception + */ + DeclException0 (ExcInvalidBoundaryFunction); + /** + * Exception + */ + DeclException2 (ExcIncompatibleNumberOfElements, + int, int, + << "The number of elements " << arg1 << " and " << arg2 + << " of the vectors do not match!"); + /** + * Exception + */ + DeclException0 (ExcInvalidSolutionVector); + /** + * Exception + */ + DeclException0 (ExcNoSolutions); +}; + + + #endif diff --git a/deal.II/deal.II/source/numerics/error_estimator.cc b/deal.II/deal.II/source/numerics/error_estimator.cc index 5d218ad2a0..52983956a1 100644 --- a/deal.II/deal.II/source/numerics/error_estimator.cc +++ b/deal.II/deal.II/source/numerics/error_estimator.cc @@ -43,59 +43,19 @@ double sqr (const double x) } - -template -KellyErrorEstimator::PerThreadData:: -PerThreadData (const unsigned int n_solution_vectors, - const unsigned int n_components, - const unsigned int n_q_points) -{ - // Init the size of a lot of vectors - // needed in the calculations once - // per thread. - normal_vectors.resize(n_q_points); - coefficient_values1.resize(n_q_points); - coefficient_values.resize(n_q_points); - JxW_values.resize(n_q_points); - - phi.resize(n_solution_vectors); - psi.resize(n_solution_vectors); - neighbor_psi.resize(n_solution_vectors); - - for (unsigned int i=0; i template -void KellyErrorEstimator::estimate (const Mapping &mapping, - const DoFHandler &dof_handler, - const Quadrature &quadrature, - const typename FunctionMap::type &neumann_bc, - const InputVector &solution, - Vector &error, - const std::vector &component_mask, - const Function *coefficients, - const unsigned int n_threads) +void KellyErrorEstimator<1>::estimate (const Mapping<1> &mapping, + const DoFHandler<1> &dof_handler, + const Quadrature<0> &quadrature, + const FunctionMap<1>::type &neumann_bc, + const InputVector &solution, + Vector &error, + const std::vector &component_mask, + const Function<1> *coefficients, + const unsigned int n_threads) { // just pass on to the other function const std::vector solutions (1, &solution); @@ -105,49 +65,42 @@ void KellyErrorEstimator::estimate (const Mapping &mapping, } -template template -void KellyErrorEstimator::estimate (const DoFHandler &dof_handler, - const Quadrature &quadrature, - const typename FunctionMap::type &neumann_bc, - const InputVector &solution, - Vector &error, - const std::vector &component_mask, - const Function *coefficients, - const unsigned int n_threads) +void KellyErrorEstimator<1>::estimate (const DoFHandler<1> &dof_handler, + const Quadrature<0> &quadrature, + const FunctionMap<1>::type &neumann_bc, + const InputVector &solution, + Vector &error, + const std::vector &component_mask, + const Function<1> *coefficients, + const unsigned int n_threads) { Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping")); - static const MappingQ1 mapping; + static const MappingQ1<1> mapping; estimate(mapping, dof_handler, quadrature, neumann_bc, solution, error, component_mask, coefficients, n_threads); } -#if deal_II_dimension == 1 - -template <> template -void KellyErrorEstimator<1>:: -estimate_some (const Mapping<1> &, - const DoFHandler<1> &, - const Quadrature<0> &, - const FunctionMap<1>::type &, - const std::vector &, - const std::vector &, - const Function<1> *, - const std::pair , - FaceIntegrals &, - PerThreadData &) +void KellyErrorEstimator<1>::estimate (const DoFHandler<1> &dof_handler, + const Quadrature<0> &quadrature, + const FunctionMap<1>::type &neumann_bc, + const std::vector &solutions, + std::vector*> &errors, + const std::vector &component_mask, + const Function<1> *coefficients, + const unsigned int n_threads) { - // in 1d, the @p{estimate} function - // does all the work - Assert (false, ExcInternalError() ); + Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping")); + static const MappingQ1<1> mapping; + estimate(mapping, dof_handler, quadrature, neumann_bc, solutions, + errors, component_mask, coefficients, n_threads); } -template <> template void KellyErrorEstimator<1>:: estimate (const Mapping<1> &mapping, @@ -246,7 +199,7 @@ estimate (const Mapping<1> &mapping, // vertices of each cell. QTrapez<1> quadrature; FEValues<1> fe_values (mapping, dof_handler.get_fe(), quadrature, update_gradients); - active_cell_iterator cell = dof_handler.begin_active(); + DoFHandler<1>::active_cell_iterator cell = dof_handler.begin_active(); for (unsigned int cell_index=0; cell != dof_handler.end(); ++cell, ++cell_index) { for (unsigned int n=0; n &mapping, #else // #if deal_II_dimension !=1 +template +KellyErrorEstimator::PerThreadData:: +PerThreadData (const unsigned int n_solution_vectors, + const unsigned int n_components, + const unsigned int n_q_points) +{ + // Init the size of a lot of vectors + // needed in the calculations once + // per thread. + normal_vectors.resize(n_q_points); + coefficient_values1.resize(n_q_points); + coefficient_values.resize(n_q_points); + JxW_values.resize(n_q_points); + + phi.resize(n_solution_vectors); + psi.resize(n_solution_vectors); + neighbor_psi.resize(n_solution_vectors); + + for (unsigned int i=0; i +template +void KellyErrorEstimator::estimate (const Mapping &mapping, + const DoFHandler &dof_handler, + const Quadrature &quadrature, + const typename FunctionMap::type &neumann_bc, + const InputVector &solution, + Vector &error, + const std::vector &component_mask, + const Function *coefficients, + const unsigned int n_threads) +{ + // just pass on to the other function + const std::vector solutions (1, &solution); + std::vector*> errors (1, &error); + estimate (mapping, dof_handler, quadrature, neumann_bc, solutions, errors, + component_mask, coefficients, n_threads); +} + + +template +template +void KellyErrorEstimator::estimate (const DoFHandler &dof_handler, + const Quadrature &quadrature, + const typename FunctionMap::type &neumann_bc, + const InputVector &solution, + Vector &error, + const std::vector &component_mask, + const Function *coefficients, + const unsigned int n_threads) +{ + Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping")); + static const MappingQ1 mapping; + estimate(mapping, dof_handler, quadrature, neumann_bc, solution, + error, component_mask, coefficients, n_threads); +} + + + + template template void KellyErrorEstimator:: @@ -722,7 +756,6 @@ KellyErrorEstimator::estimate (const Mapping &mapping }; } -#endif template template @@ -743,49 +776,6 @@ void KellyErrorEstimator::estimate (const DoFHandler &do -#if deal_II_dimension == 1 - -template <> -template -void KellyErrorEstimator<1>:: -integrate_over_regular_face (const DoFHandler<1> &, - const Quadrature<0> &, - const FunctionMap<1>::type &, - const std::vector &, - const std::vector &, - const Function<1> *, - FaceIntegrals &, - PerThreadData &, - const active_cell_iterator &, - const unsigned int , - FEFaceValues<1> &, - FEFaceValues<1> &) -{ - Assert (false, ExcInternalError()); -} - - -template <> -template -void KellyErrorEstimator<1>:: -integrate_over_irregular_face (const DoFHandler<1> &, - const Quadrature<0> &, - const std::vector &, - const std::vector &, - const Function<1> *, - FaceIntegrals &, - PerThreadData &, - const active_cell_iterator &, - const unsigned int , - FEFaceValues<1> &, - FESubfaceValues<1> &) -{ - Assert (false, ExcInternalError()); -} - -#endif - - template template void KellyErrorEstimator:: @@ -1165,9 +1155,14 @@ integrate_over_irregular_face (const DoFHandler &dof_handler, face_integrals[face] = sum; } +#endif + + // explicit instantiations +#if deal_II_dimension != 1 template class KellyErrorEstimator; +#endif // instantiate the externally visible functions. define a list of functions // for vectors, where the vector/matrix can be replaced using a preprocessor