From 3a86d1d0e8bda1849eb2f5e12b7298b80ad5ff9c Mon Sep 17 00:00:00 2001 From: kanschat Date: Tue, 2 Nov 1999 16:14:25 +0000 Subject: [PATCH] new function_lib git-svn-id: https://svn.dealii.org/trunk@1818 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/Todo | 8 +- deal.II/base/include/base/function_lib.h | 63 +++++++++ deal.II/base/source/function_lib.cc | 168 +++++++++++++++++++++++ 3 files changed, 237 insertions(+), 2 deletions(-) create mode 100644 deal.II/base/include/base/function_lib.h create mode 100644 deal.II/base/source/function_lib.cc diff --git a/deal.II/base/Todo b/deal.II/base/Todo index 3f13b27f5a..785d4d2507 100644 --- a/deal.II/base/Todo +++ b/deal.II/base/Todo @@ -10,7 +10,10 @@ Move the exceptions in the Tensor package (tensor_base.h and tensor.h) Think about the determinant function in the tensor package. Is it useful, can it be generalized? - +> There are only non-zero alternating real valued n-forms, so a determinant +> exists only for Tensor +> Alternating exterior forms of lower index return tensors. +> Guido Put Tensor<1,dim>::array_size directly into the array type. At present, egcs 1.1.2 chokes on that. @@ -24,8 +27,9 @@ Check the pattern matching in the parameter handler module using Fill in docs for the timer class. Hopefully finally find a way to let it measure times larger than half an hour. - +> Operating system dependence. Grrr. Guido Review the TensorIndex class. Better documentation, remove general template. Move constructors to the back of the file, rather than inline in the classes. Find out whether it is really needed. +> Remove !? Guido diff --git a/deal.II/base/include/base/function_lib.h b/deal.II/base/include/base/function_lib.h new file mode 100644 index 0000000000..011ccd1051 --- /dev/null +++ b/deal.II/base/include/base/function_lib.h @@ -0,0 +1,63 @@ +// $Id$ +// Copyright Guido Kanschat, University of Minnesota, 1999 + +#ifndef __base__function_lib_h +#define __base__function_lib_h + +#include + +/** + * n-quadratic pillow on the unit square. + * + * This is a function for testing the implementation. It has zero + * boundary values on the domain $(-1,1)^d$. In the inside, it is the + * product of $x_i^2-1$. + * + * Together with the function, its derivatives and Laplacian are defined. + * @author: Guido Kanschat, 1999 + */ +template +class PillowFunction : Function +{ + public: + /** + * The value at a single point. + */ + virtual double value (const Point &p, + const unsigned int component = 0) const; + + /** + * Values at multiple points. + */ + virtual void value_list (const vector > &points, + vector &values, + const unsigned int component = 0) const; + + /** + * Gradient at a single point. + */ + virtual Tensor<1,dim> gradient (const Point &p, + const unsigned int component = 0) const; + + /** + * Gradients at multiple points. + */ + virtual void gradient_list (const vector > &points, + vector > &gradients, + const unsigned int component = 0) const; + + /** + * Laplacian at a single point. + */ + double laplacian(const Point &p, + const unsigned int component = 0) const; + + /** + * Laplacian at multiple points. + */ + void laplacian_list(const vector > &points, + vector &values, + const unsigned int component = 0) const; +}; + +#endif diff --git a/deal.II/base/source/function_lib.cc b/deal.II/base/source/function_lib.cc new file mode 100644 index 0000000000..8dd51ed8dc --- /dev/null +++ b/deal.II/base/source/function_lib.cc @@ -0,0 +1,168 @@ +// $Id$ +// (c) Guido Kanschat, 1999 + +#include +#include + +//TODO: Derivatives in 3d wrong + +template +double +PillowFunction::value (const Point &p, + const unsigned int) const +{ + switch(dim) + { + case 1: + return 1.-p(0)*p(0); + case 2: + return (1.-p(0)*p(0))*(1.-p(1)*p(1)); + case 3: + return (1.-p(0)*p(0))*(1.-p(1)*p(1))*(1.-p(2)*p(2)); + default: + Assert(false, ExcNotImplemented()); + } + return 0.; +} + +template +void +PillowFunction::value_list (const vector > &points, + vector &values, + const unsigned int) const +{ + Assert (values.size() == points.size(), + ExcVectorHasWrongSize(values.size(), points.size())); + + for (unsigned int i=0;i& p = points[i]; + switch(dim) + { + case 1: + values[i] = 1.-p(0)*p(0); + break; + case 2: + values[i] = (1.-p(0)*p(0))*(1.-p(1)*p(1)); + break; + case 3: + values[i] = (1.-p(0)*p(0))*(1.-p(1)*p(1))*(1.-p(2)*p(2)); + break; + default: + Assert(false, ExcNotImplemented()); + } + } +} + +template +double +PillowFunction::laplacian (const Point &p, + const unsigned int) const +{ + switch(dim) + { + case 1: + return 2.; + case 2: + return 2.*((1.-p(0)*p(0))+(1.-p(1)*p(1))); + case 3: + return 2.*((1.-p(0)*p(0))*(1.-p(1)*p(1)) + +(1.-p(1)*p(1))*(1.-p(2)*p(2)) + +(1.-p(2)*p(2))*(1.-p(0)*p(0))); + default: + Assert(false, ExcNotImplemented()); + } + return 0.; +} + +template +void +PillowFunction::laplacian_list (const vector > &points, + vector &values, + const unsigned int) const +{ + Assert (values.size() == points.size(), + ExcVectorHasWrongSize(values.size(), points.size())); + + for (unsigned int i=0;i& p = points[i]; + switch(dim) + { + case 1: + values[i] = 2.; + break; + case 2: + values[i] = 2.*((1.-p(0)*p(0))+(1.-p(1)*p(1))); + break; + case 3: + values[i] = 2.*((1.-p(0)*p(0))+(1.-p(1)*p(1))+(1.-p(2)*p(2))); + break; + default: + Assert(false, ExcNotImplemented()); + } + } +} + +template +Tensor<1,dim> +PillowFunction::gradient (const Point &p, + const unsigned int) const +{ + Tensor<1,dim> result; + switch(dim) + { + case 1: + result[0] = 2.*p(0); + break; + case 2: + result[0] = 2.*p(0)*(1.-p(1)*p(1)); + result[1] = 2.*p(1)*(1.-p(0)*p(0)); + break; + case 3: + result[0] = 2.*p(0)*(1.-p(1)*p(1))*(1.-p(2)*p(2)); + result[1] = 2.*p(1)*(1.-p(0)*p(0))*(1.-p(2)*p(2)); + result[2] = 2.*p(2)*(1.-p(0)*p(0))*(1.-p(1)*p(1)); + break; + default: + Assert(false, ExcNotImplemented()); + } + return result; +} + +template +void +PillowFunction::gradient_list (const vector > &points, + vector > &gradients, + const unsigned int) const +{ + Assert (gradients.size() == points.size(), + ExcVectorHasWrongSize(gradients.size(), points.size())); + + for (unsigned int i=0;i& p = points[i]; + switch(dim) + { + case 1: + gradients[i][0] = 2.*p(0); + break; + case 2: + gradients[i][0] = 2.*p(0)*(1.-p(1)*p(1)); + gradients[i][1] = 2.*p(1)*(1.-p(0)*p(0)); + return; + case 3: + gradients[i][0] = 2.*p(0)*(1.-p(1)*p(1))*(1.-p(2)*p(2)); + gradients[i][1] = 2.*p(1)*(1.-p(0)*p(0))*(1.-p(2)*p(2)); + gradients[i][2] = 2.*p(2)*(1.-p(0)*p(0))*(1.-p(1)*p(1)); + break; + default: + Assert(false, ExcNotImplemented()); + } + } +} + +template PillowFunction<1>; +template PillowFunction<2>; +template PillowFunction<3>; + -- 2.39.5