]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
ComponentSelectFunction added, CosineFunction fixed
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 25 May 2000 21:24:21 +0000 (21:24 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 25 May 2000 21:24:21 +0000 (21:24 +0000)
git-svn-id: https://svn.dealii.org/trunk@2946 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/base/include/base/function.h
deal.II/base/source/function.cc
deal.II/base/source/function_lib.cc

index 5b44c49f2d9648e3ee834f7f9ca5bfebc7236002..2be2230aee39bdbc9ac8ece09c3fe44751218110 100644 (file)
@@ -462,5 +462,53 @@ class ConstantFunction : public ZeroFunction<dim>
     const double function_value;
 };
 
+/**
+ * This is a constant vector-valued function, that is different from
+ * zero only in one component.
+ *
+ * It is especially useful as a weight function for
+ * @p{VectorTools::integrate_difference}, where it allows to integrate
+ * only one component.
+ * @author Guido Kanschat, 2000
+ */
+template <int dim>
+class ComponentSelectFunction : public ConstantFunction<dim>
+{
+  public:
+                                    /**
+                                     * Constructor. Provide the
+                                     * component selected, the value
+                                     * for that component and the
+                                     * number of components.
+                                     */
+    ComponentSelectFunction (const unsigned int selected,
+                            const double value,
+                            const unsigned int n_components);
+
+                                    /**
+                                     * Return the value of the function
+                                     * at the given point for all
+                                     * components.
+                                     */
+    virtual void   vector_value (const Point<dim> &p,
+                                Vector<double>   &return_value) const;
+
+                                    /**
+                                     * Set #values# to the point values
+                                     * of the function at the #points#,
+                                     * for all components.
+                                     * It is assumed that #values#
+                                     * already has the right size, i.e.
+                                     * the same size as the #points#
+                                     * array.
+                                     */
+    virtual void vector_value_list (const vector<Point<dim> > &points,
+                                   vector<Vector<double> >   &values) const;
+
+  protected:
+    const unsigned int selected;
+};
+
+
 
 #endif
index 7ce1f9fa73f07a6139bf279582af1265a6ea0639..aa19b5436eb1ee528f926490936640d6b504d591 100644 (file)
@@ -123,17 +123,69 @@ void Function<dim>::vector_gradient_list (const vector<Point<dim> >       &point
       for (unsigned int component=0; component<n_components; ++component)
        gradients[i][component] = gradient(points[i], component);
     };
-};
+}
+
+
+
+template <int dim>
+double Function<dim>::laplacian (const Point<dim> &,
+                            const unsigned int) const
+{
+  Assert (false, ExcPureFunctionCalled());
+  return 0;
+}
+
+
+template <int dim>
+void Function<dim>::vector_laplacian (const Point<dim> &,
+                                 Vector<double>   &) const
+{
+  Assert (false, ExcPureFunctionCalled());
+}
+
+
+
+template <int dim>
+void Function<dim>::laplacian_list (const vector<Point<dim> > &points,
+                               vector<double>            &laplacians,
+                               const unsigned int         component) const
+{
+                                  // check whether component is in
+                                  // the valid range is up to the
+                                  // derived class
+  Assert (laplacians.size() == points.size(),
+         ExcDimensionMismatch(laplacians.size(), points.size()));
+
+  for (unsigned int i=0; i<points.size(); ++i)
+    laplacians[i]  = this->laplacian (points[i], component);
+}
+
+
+template <int dim>
+void Function<dim>::vector_laplacian_list (const vector<Point<dim> > &points,
+                                      vector<Vector<double> >   &laplacians) const
+{
+                                  // check whether component is in
+                                  // the valid range is up to the
+                                  // derived class
+  Assert (laplacians.size() == points.size(),
+         ExcDimensionMismatch(laplacians.size(), points.size()));
 
+  for (unsigned int i=0; i<points.size(); ++i)
+    this->vector_laplacian (points[i], laplacians[i]);
+}
+
+//------------------------------------------------------------//
 
 template <int dim>
 ZeroFunction<dim>::ZeroFunction (const unsigned int n_components) :
                Function<dim> (n_components)
-{};
+{}
 
 
 template <int dim>
-ZeroFunction<dim>::~ZeroFunction () {};
+ZeroFunction<dim>::~ZeroFunction ()
+{}
 
 
 template <int dim>
@@ -230,6 +282,8 @@ void ZeroFunction<dim>::vector_gradient_list (const vector<Point<dim> >       &p
     };
 };
 
+//------------------------------------------------------------//
+
 
 template <int dim>
 ConstantFunction<dim>::ConstantFunction (const double value,
@@ -287,52 +341,45 @@ void ConstantFunction<dim>::vector_value_list (const vector<Point<dim> > &points
     };
 };
 
+//------------------------------------------------------------//
 
 template <int dim>
-double Function<dim>::laplacian (const Point<dim> &,
-                            const unsigned int) const
-{
-  Assert (false, ExcPureFunctionCalled());
-  return 0;
-}
-
-
-template <int dim>
-void Function<dim>::vector_laplacian (const Point<dim> &,
-                                 Vector<double>   &) const
-{
-  Assert (false, ExcPureFunctionCalled());
-}
+ComponentSelectFunction<dim>::ComponentSelectFunction (const unsigned int selected,
+                                                      const double value,
+                                                      const unsigned int n_components)
+               :
+               ConstantFunction<dim> (value, n_components),
+               selected(selected)
+{}
 
 
 template <int dim>
-void Function<dim>::laplacian_list (const vector<Point<dim> > &points,
-                               vector<double>            &laplacians,
-                               const unsigned int         component) const
+void ComponentSelectFunction<dim>::vector_value (const Point<dim> &,
+                                                Vector<double>   &return_value) const
 {
-                                  // check whether component is in
-                                  // the valid range is up to the
-                                  // derived class
-  Assert (laplacians.size() == points.size(),
-         ExcDimensionMismatch(laplacians.size(), points.size()));
+  Assert (return_value.size() == n_components,
+         ExcDimensionMismatch (return_value.size(), n_components));
 
-  for (unsigned int i=0; i<points.size(); ++i)
-    laplacians[i]  = this->laplacian (points[i], component);
+  fill_n (return_value.begin(), n_components, 0.);
+  return_value(selected) = function_value;
 }
 
 
+
 template <int dim>
-void Function<dim>::vector_laplacian_list (const vector<Point<dim> > &points,
-                                      vector<Vector<double> >   &laplacians) const
+void ComponentSelectFunction<dim>::vector_value_list (const vector<Point<dim> > &points,
+                                              vector<Vector<double> >   &values) const
 {
-                                  // check whether component is in
-                                  // the valid range is up to the
-                                  // derived class
-  Assert (laplacians.size() == points.size(),
-         ExcDimensionMismatch(laplacians.size(), points.size()));
+  Assert (values.size() == points.size(),
+         ExcDimensionMismatch(values.size(), points.size()));
 
   for (unsigned int i=0; i<points.size(); ++i)
-    this->vector_laplacian (points[i], laplacians[i]);
+    {
+      Assert (values[i].size() == n_components,
+             ExcDimensionMismatch(values[i].size(), n_components));
+      fill_n (values[i].begin(), n_components, 0.);
+      values[i](selected) = function_value;
+    }
 }
 
 
@@ -341,9 +388,12 @@ void Function<dim>::vector_laplacian_list (const vector<Point<dim> > &points,
 template class Function<1>;
 template class ZeroFunction<1>;
 template class ConstantFunction<1>;
+template class ComponentSelectFunction<1>;
 template class Function<2>;
 template class ZeroFunction<2>;
 template class ConstantFunction<2>;
+template class ComponentSelectFunction<2>;
 template class Function<3>;
 template class ZeroFunction<3>;
 template class ConstantFunction<3>;
+template class ComponentSelectFunction<3>;
index 13d8225c682c1b882dd371bbd392e66e2b026802..001bc23227b50e0c8b509bc9203c6a62e8ebcd4d 100644 (file)
@@ -319,7 +319,7 @@ CosineFunction<dim>::gradient_list (const vector<Point<dim> > &points,
          case 2:
                gradients[i][0] = -M_PI_2* sin(M_PI_2*p(0)) * cos(M_PI_2*p(1));
                gradients[i][1] = -M_PI_2* cos(M_PI_2*p(0)) * sin(M_PI_2*p(1));
-               return;
+               break;
          case 3:
                gradients[i][0] = -M_PI_2* sin(M_PI_2*p(0)) * cos(M_PI_2*p(1)) * cos(M_PI_2*p(2));
                gradients[i][1] = -M_PI_2* cos(M_PI_2*p(0)) * sin(M_PI_2*p(1)) * cos(M_PI_2*p(2));

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.