]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
add gradient of LSingularityFunction
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 22 Sep 2010 21:46:57 +0000 (21:46 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 22 Sep 2010 21:46:57 +0000 (21:46 +0000)
git-svn-id: https://svn.dealii.org/trunk@22124 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 9e1e2430f65e39dcb543da4ae812bf54009a2b9c..e873f1d4dd3285f557ab2f53150776de2752e5a6 100644 (file)
@@ -1,8 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
-//    Version: $Name$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -323,10 +322,11 @@ namespace Functions
   
   
 /**
- * Singularity on the L-shaped domain in 2D.
+ * Harmonic singularity on the L-shaped domain in 2D.
  *
  * @ingroup functions
- * @author Guido Kanschat, 1999
+ * @author Guido Kanschat
+ * @date 1999
  */
   class LSingularityFunction : public Function<2>
   {
@@ -361,6 +361,53 @@ namespace Functions
   
   
   
+/**
+ * Gradient of the harmonic singularity on the L-shaped domain in 2D.
+ *
+ * The gradient of LSingularityFunction, which is a vector valued
+ * function with vanishing curl and divergence.
+ *
+ * @ingroup functions
+ * @author Guido Kanschat, 2010
+ */
+  class LSingularityGradFunction : public Function<2>
+  {
+    public:
+                                      /**
+                                       * Default constructor setting
+                                       * the dimension to 2.
+                                       */
+      LSingularityGradFunction ();
+      virtual double value (const Point<2>   &p,
+                           const unsigned int  component) const;
+      
+      virtual void value_list (const std::vector<Point<2> > &points,
+                              std::vector<double>            &values,
+                              const unsigned int              component) const;
+
+      virtual void vector_value_list (const std::vector<Point<2> >& points,
+                                     std::vector<Vector<double> >& values) const;
+      
+      virtual Tensor<1,2> gradient (const Point<2>     &p,
+                                   const unsigned int  component) const;
+      
+      virtual void gradient_list (const std::vector<Point<2> > &points,
+                                 std::vector<Tensor<1,2> >    &gradients,
+                                 const unsigned int            component) const;
+      
+      virtual void vector_gradient_list (const std::vector<Point<2> > &,
+                                        std::vector<std::vector<Tensor<1,2> > >&) const;
+      
+      virtual double laplacian (const Point<2>   &p,
+                               const unsigned int  component) const;
+      
+      virtual void laplacian_list (const std::vector<Point<2> > &points,
+                                  std::vector<double>          &values,
+                                  const unsigned int            component) const;
+  };
+  
+  
+  
 /**
  * Singularity on the slit domain in 2D and 3D.
  *
index f26398e5a4897bdac76060360a6d63990eab927c..67fe4af5ecf6ded2f8ba5b0634f50ab1a1a6b908 100644 (file)
@@ -1,8 +1,7 @@
 //---------------------------------------------------------------------------
 //    $Id$
-//    Version: $Name$
 //
-//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009 by the deal.II authors
+//    Copyright (C) 1999, 2000, 2001, 2002, 2003, 2004, 2005, 2006, 2007, 2009, 2010 by the deal.II authors
 //
 //    This file is subject to QPL and may not be  distributed
 //    without copyright and license information. Please refer
@@ -916,10 +915,10 @@ namespace Functions
       values[i] = 0.;
   }
   
-  
+
   Tensor<1,2>
   LSingularityFunction::gradient (const Point<2>   &p,
-                                 const unsigned int) const
+                                     const unsigned int d) const
   {
     double x = p(0);
     double y = p(1);
@@ -978,6 +977,132 @@ namespace Functions
       }
   }
   
+//////////////////////////////////////////////////////////////////////
+  
+  LSingularityGradFunction::LSingularityGradFunction ()
+                 :
+                 Function<2> (2)
+  {}
+
+  
+  double
+  LSingularityGradFunction::value (const Point<2>   &p,
+                                  const unsigned int d) const
+  {
+    AssertIndexRange(d,2);
+    
+    const double x = p(0);
+    const double y = p(1);
+    const double phi = std::atan2(y,-x)+M_PI;
+    const double r43 = std::pow(x*x+y*y,2./3.);
+
+    return 2./3.*(std::sin(2./3.*phi)*p(d) +
+                 (d==0
+                  ? (std::cos(2./3.*phi)*y)
+                  : (-std::cos(2./3.*phi)*x)))
+                 /r43;
+  }
+  
+  
+  void
+  LSingularityGradFunction::value_list (
+    const std::vector<Point<2> >& points,
+    std::vector<double>& values,
+    const unsigned int d) const
+  {
+    AssertIndexRange(d, 2);
+    AssertDimension(values.size(), points.size());
+    
+    for (unsigned int i=0;i<points.size();++i)
+      {
+       const Point<2>& p = points[i];
+       const double x = p(0);
+       const double y = p(1);
+       const double phi = std::atan2(y,-x)+M_PI;
+       const double r43 = std::pow(x*x+y*y,2./3.);
+       
+       values[i] = 2./3.*(std::sin(2./3.*phi)*p(d) +
+                          (d==0
+                           ? (std::cos(2./3.*phi)*y)
+                           : (-std::cos(2./3.*phi)*x)))
+                          /r43;
+      }
+  }
+  
+  
+  void
+  LSingularityGradFunction::vector_value_list (
+    const std::vector<Point<2> >& points,
+    std::vector<Vector<double> >& values) const
+  {
+    Assert (values.size() == points.size(),
+           ExcDimensionMismatch(values.size(), points.size()));
+    
+    for (unsigned int i=0;i<points.size();++i)
+      {
+       AssertDimension(values[i].size(), 2);
+       const Point<2>& p = points[i];
+       const double x = p(0);
+       const double y = p(1);
+       const double phi = std::atan2(y,-x)+M_PI;
+       const double r43 = std::pow(x*x+y*y,2./3.);
+
+       values[i](0) = 2./3.*(std::sin(2./3.*phi)*x + std::cos(2./3.*phi)*y)/r43;
+       values[i](1) = 2./3.*(std::sin(2./3.*phi)*y - std::cos(2./3.*phi)*x)/r43;
+      }
+  }
+  
+  
+  double
+  LSingularityGradFunction::laplacian (const Point<2>   &,
+                                  const unsigned int) const
+  {
+    return 0.;
+  }
+  
+  
+  void
+  LSingularityGradFunction::laplacian_list (const std::vector<Point<2> > &points,
+                                       std::vector<double>            &values,
+                                       const unsigned int) const
+  {
+    Assert (values.size() == points.size(),
+           ExcDimensionMismatch(values.size(), points.size()));
+    
+    for (unsigned int i=0;i<points.size();++i)
+      values[i] = 0.;
+  }
+  
+  
+  
+  Tensor<1,2>
+  LSingularityGradFunction::gradient (
+    const Point<2>   &/*p*/,
+    const unsigned int /*component*/) const
+  {
+    Assert(false, ExcNotImplemented());
+    return Tensor<1,2>();
+  }
+  
+  
+  void
+  LSingularityGradFunction::gradient_list (
+    const std::vector<Point<2> >& /*points*/,
+    std::vector<Tensor<1,2> >& /*gradients*/,
+    const unsigned int /*component*/) const
+  {
+    Assert(false, ExcNotImplemented());
+  }
+  
+  
+  void
+  LSingularityGradFunction::vector_gradient_list (
+    const std::vector<Point<2> >& /*points*/,
+    std::vector<std::vector<Tensor<1,2> > >& /*gradients*/) const
+  {
+    Assert(false, ExcNotImplemented());
+  }
+  
 //////////////////////////////////////////////////////////////////////
   
   template <int dim>

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.