From 6556e696f6805b109b9d217aba7b18ae098664d6 Mon Sep 17 00:00:00 2001 From: Guido Kanschat Date: Thu, 10 Jan 2002 12:37:53 +0000 Subject: [PATCH] vector valued cut-off functions git-svn-id: https://svn.dealii.org/trunk@5376 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/include/base/function_lib.h | 98 ++++++++++- deal.II/base/source/function_lib_cutoff.cc | 193 ++++++++++++++++----- deal.II/doc/news/2001/c-3-2.html | 13 +- 3 files changed, 253 insertions(+), 51 deletions(-) diff --git a/deal.II/base/include/base/function_lib.h b/deal.II/base/include/base/function_lib.h index def9951f84..30374dbe67 100644 --- a/deal.II/base/include/base/function_lib.h +++ b/deal.II/base/include/base/function_lib.h @@ -788,9 +788,10 @@ namespace Functions /** * Cut-off function in L-infinity for an arbitrary ball. This * function is the characteristic function of a ball around @p{center} - * with a specified @p{radius}. + * with a specified @p{radius}. If vector valued, it can be restricted + * to a single component. * - * @author Guido Kanschat, 2001 + * @author Guido Kanschat, 2001, 2002 */ template class CutOffFunctionLinfty : public Function @@ -800,9 +801,17 @@ namespace Functions * Constructor. Arguments are the * center of the ball and its * radius. + * + * If an argument @p{select} is + * given and not -1, the + * cut-off function will be + * non-zero for this component + * only. */ CutOffFunctionLinfty (const double radius = 1., - Point = Point()); + Point = Point(), + const unsigned int n_components = 1, + const unsigned int select = no_component); /** * Function value at one point. @@ -817,6 +826,12 @@ namespace Functions std::vector &values, const unsigned int component = 0) const; + /** + * Function values at multiple points. + */ + virtual void vector_value_list (const typename std::vector > &points, + std::vector > &values) const; + private: /** * Center of the integration ball. @@ -827,15 +842,27 @@ namespace Functions * Radius of the ball. */ const double radius; + + /** + * Component selected. If + * @p{no_component}, the function is + * the same in all components. + */ + const unsigned int selected; + /** + * Value for no selected component. + */ + static const unsigned int no_component = static_cast(-1); }; /** * Cut-off function for an arbitrary ball. This function is a cone * with support in a ball of certain @p{radius} around @p{center}. The - * maximum value is 1. + * maximum value is 1. If vector valued, it can be restricted + * to a single component. * - * @author Guido Kanschat, 2001 + * @author Guido Kanschat, 2001, 2002 */ template class CutOffFunctionW1 : public Function @@ -845,9 +872,17 @@ namespace Functions * Constructor. Arguments are the * center of the ball and its * radius. + * radius. + * + * If an argument @p{select} is + * given, the cut-off function + * will be non-zero for this + * component only. */ CutOffFunctionW1 (const double radius = 1., - Point = Point()); + Point = Point(), + const unsigned int n_components = 1, + const unsigned int select = no_component); /** * Function value at one point. @@ -862,6 +897,12 @@ namespace Functions std::vector &values, const unsigned int component = 0) const; + /** + * Function values at multiple points. + */ + virtual void vector_value_list (const typename std::vector > &points, + std::vector > &values) const; + private: /** * Center of the integration ball. @@ -872,6 +913,18 @@ namespace Functions * Radius of the ball. */ const double radius; + + /** + * Component selected. If + * @p{no_component}, the + * function is the same in all + * components. + */ + const unsigned int selected; + /** + * Value for no selected component. + */ + static const unsigned int no_component = static_cast(-1); }; @@ -879,9 +932,10 @@ namespace Functions * Cut-off function for an arbitrary ball. This is the traditional * cut-off function in C-infinity for a ball of certain @p{radius} * around @p{center}, $f(r)=exp(1-1/(1-r**2/s**2))$, where $r$ is the - * distance to the center, and $s$ is the radius of the sphere. + * distance to the center, and $s$ is the radius of the sphere. If + * vector valued, it can be restricted to a single component. * - * @author Guido Kanschat, 2001 + * @author Guido Kanschat, 2001, 2002 */ template class CutOffFunctionCinfty : public Function @@ -891,9 +945,17 @@ namespace Functions * Constructor. Arguments are the * center of the ball and its * radius. + * radius. + * + * If an argument @p{select} is + * given, the cut-off function + * will be non-zero for this + * component only. */ CutOffFunctionCinfty (const double radius = 1., - Point = Point()); + Point = Point(), + const unsigned int n_components = 1, + const unsigned int select = no_component); /** * Function value at one point. @@ -908,6 +970,12 @@ namespace Functions std::vector &values, const unsigned int component = 0) const; + /** + * Function values at multiple points. + */ + virtual void vector_value_list (const typename std::vector > &points, + std::vector > &values) const; + /** * Function gradient at one point. */ @@ -924,6 +992,18 @@ namespace Functions * Radius of the ball. */ const double radius; + + /** + * Component selected. If + * @p{no_component}, the function is + * the same in all components. + */ + const unsigned int selected; + + /** + * Value for no selected component. + */ + static const unsigned int no_component = static_cast(-1); }; }; diff --git a/deal.II/base/source/function_lib_cutoff.cc b/deal.II/base/source/function_lib_cutoff.cc index 1aaaad1712..776c566318 100644 --- a/deal.II/base/source/function_lib_cutoff.cc +++ b/deal.II/base/source/function_lib_cutoff.cc @@ -15,6 +15,7 @@ #include #include #include +#include #include @@ -41,20 +42,25 @@ namespace Functions template CutOffFunctionLinfty::CutOffFunctionLinfty (const double r, - const Point p) + const Point p, + const unsigned int n_components, + unsigned int select) : - Function (1), - center(p), - radius(r) + Function (n_components), + center(p), + radius(r), + selected(select) {} template double CutOffFunctionLinfty::value (const Point &p, - const unsigned int) const + const unsigned int component) const { - return ((center.distance(p)::value_list (const typename std::vector > &points, std::vector &values, - const unsigned int) const + const unsigned int component) const { Assert (values.size() == points.size(), ExcDimensionMismatch(values.size(), points.size())); + Assert (component < n_components, + ExcIndexRange(component,0,n_components)); + - for (unsigned int i=0;i + void + CutOffFunctionLinfty::vector_value_list ( + const typename std::vector > &points, + std::vector > &values) const + { + Assert (values.size() == points.size(), + ExcDimensionMismatch(values.size(), points.size())); + for (unsigned int k=0;k CutOffFunctionW1::CutOffFunctionW1 (const double r, - const Point p) + const Point p, + const unsigned int n_components, + unsigned int select) : - Function (1), + Function (n_components), center(p), - radius(r) + radius(r), + selected(select) {} template double CutOffFunctionW1::value (const Point &p, - const unsigned int) const + const unsigned int component) const { - const double d = center.distance(p); - return ((d::value_list (const typename std::vector > &points, std::vector &values, - const unsigned int) const + const unsigned int component) const { Assert (values.size() == points.size(), ExcDimensionMismatch(values.size(), points.size())); - for (unsigned int i=0;i + void + CutOffFunctionW1::vector_value_list ( + const typename std::vector > &points, + std::vector > &values) const + { + Assert (values.size() == points.size(), + ExcDimensionMismatch(values.size(), points.size())); + + for (unsigned int k=0;k CutOffFunctionCinfty::CutOffFunctionCinfty (const double r, - const Point p) + const Point p, + const unsigned int n_components, + unsigned int select) : - Function (1), + Function (n_components), center(p), - radius(r) + radius(r), + selected(select) {} template double CutOffFunctionCinfty::value (const Point &p, - const unsigned int) const + const unsigned int component) const { - const double d = center.distance(p); - const double r = radius; - if (d>=r) - return 0.; - const double e = -r*r/(r*r-d*d); - return ((e<-50) ? 0. : M_E * exp(e)); + if (selected==no_component || component==selected) + { + const double d = center.distance(p); + const double r = radius; + if (d>=r) + return 0.; + const double e = -r*r/(r*r-d*d); + return ((e<-50) ? 0. : M_E * exp(e)); + } + return 0.; } - + template void CutOffFunctionCinfty::value_list (const typename std::vector > &points, std::vector &values, - const unsigned int) const + const unsigned int component) const { Assert (values.size() == points.size(), ExcDimensionMismatch(values.size(), points.size())); const double r = radius; - for (unsigned int i=0;i=r) + { + values[i] = 0.; + } else { + const double e = -r*r/(r*r-d*d); + values[i] = (e<-50) ? 0. : M_E * exp(e); + } + } + else + fill (values.begin(), values.end(), 0.); + } + + + template + void + CutOffFunctionCinfty::vector_value_list ( + const typename std::vector > &points, + std::vector > &values) const + { + Assert (values.size() == points.size(), + ExcDimensionMismatch(values.size(), points.size())); + + for (unsigned int k=0;k=r) + const double d = center.distance(points[k]); + const double r = radius; + double e = 0.; + double val = 0.; + if (d-50) + val = M_E * exp(e); + } + + if (selected==no_component) + values[k] = val; + else { - values[i] = 0.; - } else { - const double e = -r*r/(r*r-d*d); - values[i] = (e<-50) ? 0. : M_E * exp(e); + values[k] = 0.; + values[k](selected) = val; } } } - + template diff --git a/deal.II/doc/news/2001/c-3-2.html b/deal.II/doc/news/2001/c-3-2.html index b585ed0fb4..11fd597f4a 100644 --- a/deal.II/doc/news/2001/c-3-2.html +++ b/deal.II/doc/news/2001/c-3-2.html @@ -107,11 +107,22 @@ documentation, etc.

base

    +
  1. + Improved: The cut-off functions Functios::CutOffFunctionLinfty, Functios::CutOffFunctionW1, and Functios::CutOffFunctionCinfty can be + vector-valued now and optionally only a single componente can + be selected. +
    + (GK 2002/01/10) +

    +
  2. New: the deal_II_exceptions::set_additional_assert_output function allows to set additional output to be printed upon - triggering an Assert() call. This + triggering an Assert() call. This is helpful for parallel applications where you only see the text of the message but do not know from which cluster node it stems. -- 2.39.5