Function< dim > Class Template Reference
[Functions]

Inheritance diagram for Function< dim >:
Inheritance graph
[legend]

List of all members.

Classes

class  ExcNumberOfComponents

Public Member Functions

 Function (const unsigned int n_components=1, const double initial_time=0.0)
virtual ~Function ()
Functionoperator= (const Function &f)
virtual double value (const Point< dim > &p, const unsigned int component=0) const
virtual void vector_value (const Point< dim > &p, Vector< double > &values) const
virtual void value_list (const std::vector< Point< dim > > &points, std::vector< double > &values, const unsigned int component=0) const
virtual void vector_value_list (const std::vector< Point< dim > > &points, std::vector< Vector< double > > &values) const
virtual Tensor< 1, dim > gradient (const Point< dim > &p, const unsigned int component=0) const
virtual void vector_gradient (const Point< dim > &p, std::vector< Tensor< 1, dim > > &gradients) const
virtual void gradient_list (const std::vector< Point< dim > > &points, std::vector< Tensor< 1, dim > > &gradients, const unsigned int component=0) const
virtual void vector_gradient_list (const std::vector< Point< dim > > &points, std::vector< std::vector< Tensor< 1, dim > > > &gradients) const
virtual double laplacian (const Point< dim > &p, const unsigned int component=0) const
virtual void vector_laplacian (const Point< dim > &p, Vector< double > &values) const
virtual void laplacian_list (const std::vector< Point< dim > > &points, std::vector< double > &values, const unsigned int component=0) const
virtual void vector_laplacian_list (const std::vector< Point< dim > > &points, std::vector< Vector< double > > &values) const
unsigned int memory_consumption () const

Public Attributes

const unsigned int n_components

Static Public Attributes

static const unsigned int dimension = dim

Detailed Description

template<int dim>
class Function< dim >

This class is a model for a general function. It serves the purpose of representing scalar and vector valued functions. To this end, we consider scalar functions as a special case of vector valued functions, in the former case only having a single component return vector. Since handling with vectors is comparatively expensive, functions are provided which only ask for a single component of the function, which is what you will usually need in case you know that your function is scalar-valued.

Access to function objects therefore is through the following methods:

                 // access to one component at one point
   double value        (const Point<dim>   &p,
                        const unsigned int  component = 0) const;

                 // return all components at one point
   void   vector_value (const Point<dim>   &p,
                        Vector<double>     &value) const;

For more efficiency, there are other functions returning one or all components at a list of points at once:

                 // access to one component at several points
   void   value_list (const std::vector<Point<dim> >  &point_list,
                      std::vector<double>             &value_list,
                      const unsigned int  component = 0) const;

                 // return all components at several points
   void   vector_value_list (const std::vector<Point<dim> >    &point_list,
                             std::vector<Vector<double> >      &value_list) const;

Furthermore, there are functions returning the gradient of the function at one or several points.

You will usually only overload those functions you need; the functions returning several values at a time (value_list(), vector_value_list(), and gradient analoga) will call those returning only one value (value(), vector_value(), and gradient analoga), while those ones will throw an exception when called but not overloaded.

Note however, that the functions returning all components of the function at one or several points (i.e. vector_value(), vector_value_list()), will not call the function returning one component at one point repeatedly, once for each point and component. The reason is efficiency: this would amount to too many virtual function calls. If you have vector-valued functions, you should therefore also provide overloads of the virtual functions for all components at a time.

Also note, that unless only called a very small number of times, you should overload all sets of functions (returning only one value, as well as those returning a whole array), since the cost of evaluation of a point value is often less than the virtual function call itself.

Support for time dependant functions can be found in the base class FunctionTime.

Note:
if the functions you are dealing with have sizes which are a priori known (for example, dim elements), you might consider using the TensorFunction class instead.
Author:
Wolfgang Bangerth, 1998, 1999

Constructor & Destructor Documentation

template<int dim>
Function< dim >::Function ( const unsigned int  n_components = 1,
const double  initial_time = 0.0 
)

Constructor. May take an initial value for the number of components (which defaults to one, i.e. a scalar function), and the time variable, which defaults to zero.

template<int dim>
virtual Function< dim >::~Function (  )  [virtual]

Virtual destructor; absolutely necessary in this case.

This destructor is declared pure virtual, such that objects of this class cannot be created. Since all the other virtual functions have a pseudo-implementation to avoid overhead in derived classes, this is the best place to do this.

Nevertheless, since derived classes want to call the destructor of a base class, the destructor is implemented (despite it being pure virtual).

Note: Compaq's cxx compiler does not allow defining a function that was declared pure. It simply refuses to instantiate the function when it later sees it, but also does not generate a respective function itself, which then leads to linker errors claiming that the destructor of this class is missing. We therefore only make the function abstract if the compiler can handle this.


Member Function Documentation

template<int dim>
Function& Function< dim >::operator= ( const Function< dim > &  f  ) 

Assignment operator. This is here only so that you can have objects of derived classes in containers, or assign them otherwise. It will raise an exception if the object from which you assign has a different number of components than the one being assigned to.

Reimplemented from Subscriptor.

template<int dim>
virtual double Function< dim >::value ( const Point< dim > &  p,
const unsigned int  component = 0 
) const [virtual]
template<int dim>
virtual void Function< dim >::vector_value ( const Point< dim > &  p,
Vector< double > &  values 
) const [virtual]
template<int dim>
virtual void Function< dim >::value_list ( const std::vector< Point< dim > > &  points,
std::vector< double > &  values,
const unsigned int  component = 0 
) const [virtual]

Set values to the point values of the specified component of the function at the points. It is assumed that values already has the right size, i.e. the same size as the points array.

Be default, this function repeatedly calls value() for each point separately, to fill the output array.

Reimplemented in ZeroFunction< dim >, ConstantFunction< dim >, FunctionDerivative< dim >, Functions::SquareFunction< dim >, Functions::Q1WedgeFunction< dim >, Functions::PillowFunction< dim >, Functions::CosineFunction< dim >, Functions::ExpFunction< dim >, Functions::SlitSingularityFunction< dim >, Functions::JumpFunction< dim >, Functions::CutOffFunctionLinfty< dim >, Functions::CutOffFunctionW1< dim >, Functions::CutOffFunctionCinfty< dim >, Functions::Monomial< dim >, and Functions::FEFieldFunction< dim, DH, VECTOR >.

template<int dim>
virtual void Function< dim >::vector_value_list ( const std::vector< Point< dim > > &  points,
std::vector< Vector< double > > &  values 
) const [virtual]

Set values to the point values of the function at the points. It is assumed that values already has the right size, i.e. the same size as the points array, and that all elements be vectors with the same number of components as this function has.

Be default, this function repeatedly calls vector_value() for each point separately, to fill the output array.

Reimplemented in Functions::FlowFunction< dim >, ZeroFunction< dim >, ConstantFunction< dim >, ComponentSelectFunction< dim >, Functions::Q1WedgeFunction< dim >, Functions::SlitSingularityFunction< dim >, Functions::CutOffFunctionLinfty< dim >, Functions::CutOffFunctionW1< dim >, Functions::CutOffFunctionCinfty< dim >, Functions::FEFieldFunction< dim, DH, VECTOR >, and Functions::FlowFunction< 2 >.

template<int dim>
virtual Tensor<1,dim> Function< dim >::gradient ( const Point< dim > &  p,
const unsigned int  component = 0 
) const [virtual]
template<int dim>
virtual void Function< dim >::vector_gradient ( const Point< dim > &  p,
std::vector< Tensor< 1, dim > > &  gradients 
) const [virtual]

Return the gradient of all components of the function at the given point.

Reimplemented in AutoDerivativeFunction< dim >, ZeroFunction< dim >, and Functions::FEFieldFunction< dim, DH, VECTOR >.

template<int dim>
virtual void Function< dim >::gradient_list ( const std::vector< Point< dim > > &  points,
std::vector< Tensor< 1, dim > > &  gradients,
const unsigned int  component = 0 
) const [virtual]

Set gradients to the gradients of the specified component of the function at the points. It is assumed that gradients already has the right size, i.e. the same size as the points array.

Reimplemented in AutoDerivativeFunction< dim >, ZeroFunction< dim >, Functions::SquareFunction< dim >, Functions::Q1WedgeFunction< dim >, Functions::PillowFunction< dim >, Functions::CosineFunction< dim >, Functions::ExpFunction< dim >, Functions::SlitSingularityFunction< dim >, Functions::JumpFunction< dim >, and Functions::FEFieldFunction< dim, DH, VECTOR >.

template<int dim>
virtual void Function< dim >::vector_gradient_list ( const std::vector< Point< dim > > &  points,
std::vector< std::vector< Tensor< 1, dim > > > &  gradients 
) const [virtual]

Set gradients to the gradients of the function at the points, for all components. It is assumed that gradients already has the right size, i.e. the same size as the points array.

The outer loop over gradients is over the points in the list, the inner loop over the different components of the function.

Reimplemented in AutoDerivativeFunction< dim >, Functions::FlowFunction< dim >, ZeroFunction< dim >, Functions::Q1WedgeFunction< dim >, Functions::SlitSingularityFunction< dim >, Functions::FEFieldFunction< dim, DH, VECTOR >, and Functions::FlowFunction< 2 >.

template<int dim>
virtual double Function< dim >::laplacian ( const Point< dim > &  p,
const unsigned int  component = 0 
) const [virtual]
template<int dim>
virtual void Function< dim >::vector_laplacian ( const Point< dim > &  p,
Vector< double > &  values 
) const [virtual]

Compute the Laplacian of all components at point p and store them in values.

template<int dim>
virtual void Function< dim >::laplacian_list ( const std::vector< Point< dim > > &  points,
std::vector< double > &  values,
const unsigned int  component = 0 
) const [virtual]
template<int dim>
virtual void Function< dim >::vector_laplacian_list ( const std::vector< Point< dim > > &  points,
std::vector< Vector< double > > &  values 
) const [virtual]

Compute the Laplacians of all components at a set of points.

Reimplemented in Functions::FlowFunction< dim >, and Functions::FlowFunction< 2 >.

template<int dim>
unsigned int Function< dim >::memory_consumption (  )  const

Determine an estimate for the memory consumption (in bytes) of this object. Since sometimes the size of objects can not be determined exactly (for example: what is the memory consumption of an STL std::map type with a certain number of elements?), this is only an estimate. however often quite close to the true value.

Reimplemented in Functions::FlowFunction< dim >, ConstantFunction< dim >, ComponentSelectFunction< dim >, FunctionDerivative< dim >, Functions::JumpFunction< dim >, and Functions::FlowFunction< 2 >.


Member Data Documentation

template<int dim>
const unsigned int Function< dim >::dimension = dim [static]

Export the value of the template parameter as a static member constant. Sometimes useful for some expression template programming.

template<int dim>
const unsigned int Function< dim >::n_components

Number of vector components.

Reimplemented in Functions::FEFieldFunction< dim, DH, VECTOR >.


The documentation for this class was generated from the following file:

deal.II documentation generated on Mon Nov 23 22:57:47 2009 by doxygen 1.6.1