From 12efa07f841a364a96a58774379258f3a256442c Mon Sep 17 00:00:00 2001 From: guido Date: Thu, 1 Aug 2002 08:20:10 +0000 Subject: [PATCH] prepared for lp-norms git-svn-id: https://svn.dealii.org/trunk@6308 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/deal.II/include/numerics/vectors.h | 50 +- deal.II/deal.II/source/numerics/vectors.cc | 513 +++++++++++---------- 2 files changed, 291 insertions(+), 272 deletions(-) diff --git a/deal.II/deal.II/include/numerics/vectors.h b/deal.II/deal.II/include/numerics/vectors.h index 62b20af51a..160619dc76 100644 --- a/deal.II/deal.II/include/numerics/vectors.h +++ b/deal.II/deal.II/include/numerics/vectors.h @@ -305,10 +305,17 @@ class VectorTools * absolute value of the * function is integrated. * @item @p{L2_norm}: the square - * of the function is - * integrated on each cell; - * afterwards the root is - * taken of this value. + * of the function is integrated + * and the the square root of the + * result is computed on each + * cell. + * @item @p{Lp_norm}: the + * absolute value to the pth + * power is integrated and the + * pth root is computed on each + * cell. The exponent @p{p} is + * the last parameter of the + * function. * @item @p{Linfty_norm}: the * maximum absolute value of the * function. @@ -316,8 +323,11 @@ class VectorTools * square of the function * gradient is integrated on * each cell; afterwards the - * root is taken of this * + * root is taken of this * value. + * @item @p{W1p_seminorm}: this + * is the @p{Lp_norm} of the + * gradient. * @item @p{H1_norm}: the square * of the function plus the * square of the function @@ -330,14 +340,21 @@ class VectorTools * square of the * @p{H1_seminorm}. * @end{itemize} - */ + * @item @p{W1p_norm}: like + * @p{H1_norm}, but for + * @p{Lp_norm} instead of + * @p{L2_norm} + */ enum NormType { mean, L1_norm, L2_norm, Linfty_norm, H1_seminorm, - H1_norm + H1_norm, + Lp_norm, + W1p_seminorm, + W1p_norm }; /** @@ -736,18 +753,13 @@ class VectorTools * Compute the error of the finite element solution. * Integrate the difference between * a finite element function and - * the reference function, which + * a reference function, which * is given as a continuous function * object. * - * Note that this function returns - * its results in a vector of @p{float}s, - * rather than in a vector of @p{double}s, - * since accuracy is not that important - * here and to save memory. During - * computation of the results, the full - * accuracy of the @p{double} data type is - * used. + * The value of @p{exponent} is + * used for computing $L^p$-norms + * and $W^{1,p}$-norms. * * The additional argument * @p{weight} allows to evaluate @@ -790,7 +802,8 @@ class VectorTools OutVector &difference, const Quadrature &q, const NormType &norm, - const Function *weight=0); + const Function *weight=0, + const double exponent = 2.); /** * Calls the @p{integrate_difference} @@ -804,7 +817,8 @@ class VectorTools OutVector &difference, const Quadrature &q, const NormType &norm, - const Function *weight=0); + const Function *weight=0, + const double exponent = 2.); /** * Mean-value filter for Stokes. diff --git a/deal.II/deal.II/source/numerics/vectors.cc b/deal.II/deal.II/source/numerics/vectors.cc index 53efe8b20d..821fbdc10a 100644 --- a/deal.II/deal.II/source/numerics/vectors.cc +++ b/deal.II/deal.II/source/numerics/vectors.cc @@ -1058,7 +1058,8 @@ VectorTools::integrate_difference (const Mapping &mapping, OutVector &difference, const Quadrature &q, const NormType &norm, - const Function *weight) + const Function *weight, + double exponent) { const unsigned int n_q_points = q.n_quadrature_points; const FiniteElement &fe = dof.get_fe(); @@ -1070,16 +1071,39 @@ VectorTools::integrate_difference (const Mapping &mapping, Assert ((weight->n_components==1) || (weight->n_components==n_components), ExcDimensionMismatch(weight->n_components, n_components)); } - + difference.reinit (dof.get_tria().n_active_cells()); + switch (norm) + { + case L2_norm: + case H1_seminorm: + case H1_norm: + exponent = 2.; + break; + case L1_norm: + exponent = 1.; + break; + default: + break; + } + UpdateFlags update_flags = UpdateFlags (update_q_points | update_JxW_values); - if (norm != H1_seminorm) - update_flags = UpdateFlags(update_flags | update_values); - - if ((norm==H1_seminorm) || (norm==H1_norm)) - update_flags = UpdateFlags (update_flags | update_gradients); + switch (norm) + { + case H1_seminorm: + case W1p_seminorm: + update_flags |= UpdateFlags (update_gradients); + break; + case H1_norm: + case W1p_norm: + update_flags |= UpdateFlags (update_gradients); + // no break! + default: + update_flags |= UpdateFlags (update_values); + break; + } FEValues fe_values(mapping, fe, q, update_flags); @@ -1115,267 +1139,247 @@ VectorTools::integrate_difference (const Mapping &mapping, { if (weight->n_components>1) weight->vector_value_list (fe_values.get_quadrature_points(), - weight_vectors); + weight_vectors); else - weight->value_list (fe_values.get_quadrature_points(), - weight_values); + { + weight->value_list (fe_values.get_quadrature_points(), + weight_values); +// for (unsigned int k=0;kn_components > 1) + for (unsigned int q=0; qn_components == 1 + for (unsigned int q=0; qn_components > 1) - for (unsigned int q=0; qn_components == 1 - for (unsigned int q=0; qn_components > 1) - for (unsigned int q=0; qn_components == 1 - for (unsigned int q=0; qn_components > 1) - for (unsigned int q=0; qn_components == 1 - for (unsigned int q=0; qn_components > 1) + for (unsigned int q=0; qn_components == 1 + for (unsigned int q=0; qn_components > 1) + for (unsigned int q=0; qn_components == 1 + for (unsigned int q=0; qn_components > 1) + for (unsigned int q=0; qn_components > 1) - for (unsigned int q=0; qn_components == 1 - for (unsigned int q=0; qn_components == 1 + for (unsigned int q=0; qn_components > 1) + for (unsigned int q=0; qn_components == 1 + for (unsigned int q=0; qn_components > 1) - for (unsigned int q=0; qn_components == 1 - for (unsigned int q=0; q @@ -1386,12 +1390,13 @@ VectorTools::integrate_difference (const DoFHandler &dof, OutVector &difference, const Quadrature &q, const NormType &norm, - const Function *weight) + const Function *weight, + double exponent) { Assert (DEAL_II_COMPAT_MAPPING, ExcCompatibility("mapping")); static const MappingQ1 mapping; integrate_difference(mapping, dof, fe_function, exact_solution, - difference, q, norm, weight); + difference, q, norm, weight, exponent); } -- 2.39.5