]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Example functions for singularities
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 3 Dec 1999 02:49:56 +0000 (02:49 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Fri, 3 Dec 1999 02:49:56 +0000 (02:49 +0000)
git-svn-id: https://svn.dealii.org/trunk@1971 0785d39b-7218-0410-832d-ea1e28bc413d

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

index b2762dd54f8acd89d632bfd7a751bc7399be8638..c2012f76d86989dcbb762a223c0306827ad43c2d 100644 (file)
@@ -9,9 +9,9 @@
 #include <base/function.h>
 
 /**
- * n-quadratic pillow on the unit square.
+ * d-quadratic pillow on the unit hypercube.
  *
- * This is a function for testing the implementation. It has zero
+ * This is a function for testing the implementation. It has zero Dirichlet
  * boundary values on the domain $(-1,1)^d$. In the inside, it is the
  * product of $x_i^2-1$.
  *
@@ -20,7 +20,7 @@
  * @author: Guido Kanschat, 1999
  */
 template<int dim>
-class PillowFunction : Function<dim>
+class PillowFunction : public Function<dim>
 {
   public:
                                     /**
@@ -63,6 +63,100 @@ class PillowFunction : Function<dim>
                        const unsigned int         component = 0) const;
 };
 
+/**
+ * Singularity on the L-shaped domain in 2D.
+ * @author Guido Kanschat, 1999
+ */
+class LSingularityFunction : public Function<2>
+{
+  public:
+                                    /**
+                                     * The value at a single point.
+                                     */
+    virtual double value (const Point<2>   &p,
+                         const unsigned int  component = 0) const;
+
+                                    /**
+                                     * Values at multiple points.
+                                     */
+    virtual void value_list (const vector<Point<2> > &points,
+                            vector<double>            &values,
+                            const unsigned int         component = 0) const;
+
+                                    /**
+                                     * Gradient at a single point.
+                                     */
+    virtual Tensor<1,2> gradient (const Point<2>   &p,
+                                   const unsigned int  component = 0) const;
+
+                                    /**
+                                     * Gradients at multiple points.
+                                     */
+    virtual void gradient_list (const vector<Point<2> > &points,
+                               vector<Tensor<1,2> >    &gradients,
+                               const unsigned int         component = 0) const;
+
+                                    /**
+                                     * Laplacian at a single point.
+                                     */
+    double laplacian(const Point<2>   &p,
+                    const unsigned int  component = 0) const;
+
+                                    /**
+                                     * Laplacian at multiple points.
+                                     */
+    void laplacian_list(const vector<Point<2> > &points,
+                       vector<double>            &values,
+                       const unsigned int         component = 0) const;
+};
+
+/**
+ * Singularity on the slit domain in 2D.
+ * @author Guido Kanschat, 1999
+ */
+class SlitSingularityFunction : public Function<2>
+{
+  public:
+                                    /**
+                                     * The value at a single point.
+                                     */
+    virtual double value (const Point<2>   &p,
+                         const unsigned int  component = 0) const;
+
+                                    /**
+                                     * Values at multiple points.
+                                     */
+    virtual void value_list (const vector<Point<2> > &points,
+                            vector<double>            &values,
+                            const unsigned int         component = 0) const;
+
+                                    /**
+                                     * Gradient at a single point.
+                                     */
+    virtual Tensor<1,2> gradient (const Point<2>   &p,
+                                   const unsigned int  component = 0) const;
+
+                                    /**
+                                     * Gradients at multiple points.
+                                     */
+    virtual void gradient_list (const vector<Point<2> > &points,
+                               vector<Tensor<1,2> >    &gradients,
+                               const unsigned int         component = 0) const;
+
+                                    /**
+                                     * Laplacian at a single point.
+                                     */
+    double laplacian(const Point<2>   &p,
+                    const unsigned int  component = 0) const;
+
+                                    /**
+                                     * Laplacian at multiple points.
+                                     */
+    void laplacian_list(const vector<Point<2> > &points,
+                       vector<double>            &values,
+                       const unsigned int         component = 0) const;
+};
+
 
 
 /*----------------------------   function_lib.h     ---------------------------*/
index 47cfdfa0333ec7f628919969075a87c7641b1cf0..575730f7679a4483197004394fb15264352529bf 100644 (file)
@@ -4,6 +4,8 @@
 #include <base/point.h>
 #include <base/function_lib.h>
 
+#include <cmath>
+
 //TODO: Derivatives in 3d wrong (GK!)
 
 template<int dim>
@@ -162,7 +164,176 @@ PillowFunction<dim>::gradient_list (const vector<Point<dim> > &points,
     }
 }
 
+//////////////////////////////////////////////////////////////////////
+
+
+double
+LSingularityFunction::value (const Point<2>   &p,
+                           const unsigned int) const
+{
+  double x = p(0);
+  double y = p(1);
+
+  if ((x>=0) && (y>=0))
+    return 0.;
+  
+  double phi = atan2(y,-x)+M_PI;
+  double r2 = x*x+y*y;
+
+  return pow(r2,1./3.) * sin(2./3.*phi);
+}
+
+
+void
+LSingularityFunction::value_list (const vector<Point<2> > &points,
+                                vector<double>            &values,
+                                const unsigned int) const
+{
+  Assert (values.size() == points.size(),
+         ExcVectorHasWrongSize(values.size(), points.size()));
+
+  for (unsigned int i=0;i<points.size();++i)
+    {
+      double x = points[i](0);
+      double y = points[i](1);
+
+      if ((x>=0) && (y>=0))
+       values[i] = 0.;
+      else
+       {
+         double phi = atan2(y,-x)+M_PI;
+         double r2 = x*x+y*y;
+
+         values[i] = pow(r2,1./3.) * sin(2./3.*phi);
+       }
+    }
+}
+
+
+double
+LSingularityFunction::laplacian (const Point<2>   &,
+                           const unsigned int) const
+{
+  return 0.;
+}
+
+
+void
+LSingularityFunction::laplacian_list (const vector<Point<2> > &points,
+                                    vector<double>            &values,
+                                    const unsigned int) const
+{
+  Assert (values.size() == points.size(),
+         ExcVectorHasWrongSize(values.size(), points.size()));
+
+  for (unsigned int i=0;i<points.size();++i)
+    values[i] = 0.;
+}
+
+
+//TODO: Implement derivatives
+
+Tensor<1,2>
+LSingularityFunction::gradient (const Point<2>   &/*p*/,
+                              const unsigned int) const
+{
+  Assert(false, ExcNotImplemented());
+  return Tensor<1,2>();
+}
+
+
+void
+LSingularityFunction::gradient_list (const vector<Point<2> > &points,
+                                   vector<Tensor<1,2> >    &gradients,
+                                   const unsigned int) const
+{
+  Assert (gradients.size() == points.size(),
+         ExcVectorHasWrongSize(gradients.size(), points.size()));
+  Assert(false, ExcNotImplemented());
+}
+
+//////////////////////////////////////////////////////////////////////
+
+
+double
+SlitSingularityFunction::value (const Point<2>   &p,
+                           const unsigned int) const
+{
+  double x = p(0);
+  double y = p(1);
+
+  double phi = atan2(x,-y)+M_PI;
+  double r2 = x*x+y*y;
+
+  return pow(r2,.25) * sin(.5*phi);
+}
+
+
+void
+SlitSingularityFunction::value_list (const vector<Point<2> > &points,
+                                vector<double>            &values,
+                                const unsigned int) const
+{
+  Assert (values.size() == points.size(),
+         ExcVectorHasWrongSize(values.size(), points.size()));
+
+  for (unsigned int i=0;i<points.size();++i)
+    {
+      double x = points[i](0);
+      double y = points[i](1);
+
+      double phi = atan2(x,-y)+M_PI;
+      double r2 = x*x+y*y;
+
+      values[i] = pow(r2,.25) * sin(.5*phi);
+    }
+}
+
+
+double
+SlitSingularityFunction::laplacian (const Point<2>   &,
+                           const unsigned int) const
+{
+  return 0.;
+}
+
+
+void
+SlitSingularityFunction::laplacian_list (const vector<Point<2> > &points,
+                                    vector<double>            &values,
+                                    const unsigned int) const
+{
+  Assert (values.size() == points.size(),
+         ExcVectorHasWrongSize(values.size(), points.size()));
+
+  for (unsigned int i=0;i<points.size();++i)
+    values[i] = 0.;
+}
+
+
+//TODO: Implement derivatives
+
+Tensor<1,2>
+SlitSingularityFunction::gradient (const Point<2>   &/*p*/,
+                              const unsigned int) const
+{
+  Assert(false, ExcNotImplemented());
+  return Tensor<1,2>();
+}
+
+
+void
+SlitSingularityFunction::gradient_list (const vector<Point<2> > &points,
+                                   vector<Tensor<1,2> >    &gradients,
+                                   const unsigned int) const
+{
+  Assert (gradients.size() == points.size(),
+         ExcVectorHasWrongSize(gradients.size(), points.size()));
+  Assert(false, ExcNotImplemented());
+}
+
+//////////////////////////////////////////////////////////////////////
+
 template PillowFunction<1>;
 template PillowFunction<2>;
 template PillowFunction<3>;
-

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.