]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
new SlitHyperSingularityFunction, yet to be tested
authorguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 7 Nov 2002 18:36:05 +0000 (18:36 +0000)
committerguido <guido@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 7 Nov 2002 18:36:05 +0000 (18:36 +0000)
git-svn-id: https://svn.dealii.org/trunk@6739 0785d39b-7218-0410-832d-ea1e28bc413d

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

index 2a8dcc3c5c8f9dafbbe2caf8da4c03115fe55575..a43c781bbb8050f177aa5ad4b491ca1989c86ee9 100644 (file)
@@ -324,8 +324,6 @@ namespace Functions
 /**
  * Singularity on the L-shaped domain in 2D.
  *
- * Caveat: derivatives of this function are not implemented!
- *
  * @author Guido Kanschat, 1999
  */
   class LSingularityFunction : public Function<2>
@@ -376,8 +374,6 @@ namespace Functions
 /**
  * Singularity on the slit domain in 2D.
  *
- * Caveat: derivatives of this function are not implemented!
- *
  * @author Guido Kanschat, 1999
  */
   class SlitSingularityFunction : public Function<2>
@@ -424,6 +420,55 @@ namespace Functions
   };
   
   
+/**
+ * Singularity on the slit domain with one Neumann boundary in 2D.
+ *
+ * @author Guido Kanschat, 2002
+ */
+  class SlitHyperSingularityFunction : 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 std::vector<Point<2> > &points,
+                              std::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 std::vector<Point<2> > &points,
+                                 std::vector<Tensor<1,2> >    &gradients,
+                                 const unsigned int            component = 0) const;
+      
+                                      /**
+                                       * Laplacian at a single point.
+                                       */
+      virtual double laplacian (const Point<2>   &p,
+                               const unsigned int  component = 0) const;
+      
+                                      /**
+                                       * Laplacian at multiple points.
+                                       */
+      virtual void laplacian_list (const std::vector<Point<2> > &points,
+                                  std::vector<double>          &values,
+                                  const unsigned int            component = 0) const;
+  };
+  
+  
 /**
  * A jump in x-direction transported into some direction.
  *
index 6ba48bc36d48e644150ce4089a0deb521942b7fd..fcdba001258fbadac2c218ea7ab618648fda9cbd 100644 (file)
@@ -985,6 +985,108 @@ namespace Functions
       }
   }
   
+//////////////////////////////////////////////////////////////////////
+  
+  
+  double
+  SlitHyperSingularityFunction::value (const Point<2>   &p,
+                                      const unsigned int) const
+  {
+    double x = p(0);
+    double y = p(1);
+    
+    double phi = std::atan2(x,y)+M_PI;
+    double r2 = x*x+y*y;
+    
+    return std::pow(r2,.125) * std::sin(.25*phi);
+  }
+  
+  
+  void
+  SlitHyperSingularityFunction::value_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)
+      {
+       double x = points[i](0);
+       double y = points[i](1);
+       
+       double phi = std::atan2(x,y)+M_PI;
+       double r2 = x*x+y*y;
+       
+       values[i] = std::pow(r2,.125) * std::sin(.25*phi);
+      }
+  }
+  
+  
+  double
+  SlitHyperSingularityFunction::laplacian (
+    const Point<2>   &,
+    const unsigned int) const
+  {
+    return 0.;
+  }
+  
+  
+  void
+  SlitHyperSingularityFunction::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>
+  SlitHyperSingularityFunction::gradient (
+    const Point<2>   &p,
+    const unsigned int) const
+  {
+    double x = p(0);
+    double y = p(1);
+    double phi = std::atan2(x,y)+M_PI;
+    double r78 = std::pow(x*x+y*y,7./8.);
+    
+    
+    Tensor<1,2> result;
+    result[0] = 1./4.*(std::sin(1./4.*phi)*x + std::cos(1./4.*phi)*y)/r78;
+    result[1] = 1./4.*(std::sin(1./4.*phi)*y - std::cos(1./4.*phi)*x)/r78;
+    return result;
+  }
+  
+  
+  void
+  SlitHyperSingularityFunction::gradient_list (
+    const std::vector<Point<2> > &points,
+    std::vector<Tensor<1,2> >    &gradients,
+    const unsigned int) const
+  {
+    Assert (gradients.size() == points.size(),
+           ExcDimensionMismatch(gradients.size(), points.size()));
+    
+    for (unsigned int i=0;i<points.size();++i)
+      {
+       const Point<2>& p = points[i];
+       double x = p(0);
+       double y = p(1);
+       double phi = std::atan2(x,y)+M_PI;
+       double r78 = std::pow(x*x+y*y,7./8.);
+       
+       gradients[i][0] = 1./4.*(std::sin(1./4.*phi)*x + std::cos(1./4.*phi)*y)/r78;
+       gradients[i][1] = 1./4.*(std::sin(1./4.*phi)*y - std::cos(1./4.*phi)*x)/r78;
+      }
+  }
+  
 //////////////////////////////////////////////////////////////////////
   
   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.