Logo Search packages:      
Sourcecode: ufc version File versions  Download package

ufc.h

// This is UFC (Unified Form-assembly Code) v. 1.2.
// This code is released into the public domain.
//
// The FEniCS Project (http://www.fenics.org/) 2006-2009.

#ifndef __UFC_H
#define __UFC_H

#define UFC_VERSION_MAJOR 1
#define UFC_VERSION_MINOR 2
#define UFC_VERSION_MAINTENANCE 0

#include <stdexcept>

const char UFC_VERSION[] = "1.2";

00017 namespace ufc
{

  /// Valid cell shapes
00021   enum shape {interval, triangle, quadrilateral, tetrahedron, hexahedron};

  /// This class defines the data structure for a finite element mesh.

00025   class mesh
  {
  public:

    /// Constructor
00030     mesh(): topological_dimension(0), geometric_dimension(0), num_entities(0) {}

    /// Destructor
00033     virtual ~mesh() {}

    /// Topological dimension of the mesh
00036     unsigned int topological_dimension;

    /// Geometric dimension of the mesh
00039     unsigned int geometric_dimension;

    /// Array of the global number of entities of each topological dimension
00042     unsigned int* num_entities;

  };

  /// This class defines the data structure for a cell in a mesh.

00048   class cell
  {
  public:

    /// Constructor
00053     cell(): cell_shape(interval),
            topological_dimension(0), geometric_dimension(0),
            entity_indices(0), coordinates(0) {}

    /// Destructor
00058     virtual ~cell() {}

    /// Shape of the cell
00061     shape cell_shape;

    /// Topological dimension of the mesh
00064     unsigned int topological_dimension;

    /// Geometric dimension of the mesh
00067     unsigned int geometric_dimension;

    /// Array of global indices for the mesh entities of the cell
00070     unsigned int** entity_indices;

    /// Array of coordinates for the vertices of the cell
00073     double** coordinates;

  };

  /// This class defines the interface for a general tensor-valued function.

00079   class function
  {
  public:

    /// Constructor
00084     function() {}

    /// Destructor
00087     virtual ~function() {}

    /// Evaluate function at given point in cell
    virtual void evaluate(double* values,
                          const double* coordinates,
                          const cell& c) const = 0;

  };

  /// This class defines the interface for a finite element.

00098   class finite_element
  {
  public:

    /// Constructor
00103     finite_element() {}

    /// Destructor
00106     virtual ~finite_element() {}

    /// Return a string identifying the finite element
    virtual const char* signature() const = 0;

    /// Return the cell shape
    virtual shape cell_shape() const = 0;

    /// Return the dimension of the finite element function space
    virtual unsigned int space_dimension() const = 0;

    /// Return the rank of the value space
    virtual unsigned int value_rank() const = 0;

    /// Return the dimension of the value space for axis i
    virtual unsigned int value_dimension(unsigned int i) const = 0;

    /// Evaluate basis function i at given point in cell
    virtual void evaluate_basis(unsigned int i,
                                double* values,
                                const double* coordinates,
                                const cell& c) const = 0;

    /// Evaluate all basis functions at given point in cell
00130     virtual void evaluate_basis_all(double* values,
                                    const double* coordinates,
                                    const cell& c) const
    { throw std::runtime_error("Not implemented (introduced in UFC v1.1)."); }

    /// Evaluate order n derivatives of basis function i at given point in cell
    virtual void evaluate_basis_derivatives(unsigned int i,
                                            unsigned int n,
                                            double* values,
                                            const double* coordinates,
                                            const cell& c) const = 0;

    /// Evaluate order n derivatives of all basis functions at given point in cell
00143     virtual void evaluate_basis_derivatives_all(unsigned int n,
                                                double* values,
                                                const double* coordinates,
                                                const cell& c) const
    { throw std::runtime_error("Not implemented (introduced in UFC v1.1)."); }

    /// Evaluate linear functional for dof i on the function f
    virtual double evaluate_dof(unsigned int i,
                                const function& f,
                                const cell& c) const = 0;

    /// Evaluate linear functionals for all dofs on the function f
00155     virtual void evaluate_dofs(double* values,
                               const function& f,
                               const cell& c) const
    { throw std::runtime_error("Not implemented (introduced in UFC v1.1)."); }

    /// Interpolate vertex values from dof values
    virtual void interpolate_vertex_values(double* vertex_values,
                                           const double* dof_values,
                                           const cell& c) const = 0;

    /// Return the number of sub elements (for a mixed element)
    virtual unsigned int num_sub_elements() const = 0;

    /// Create a new finite element for sub element i (for a mixed element)
    virtual finite_element* create_sub_element(unsigned int i) const = 0;

  };

  /// This class defines the interface for a local-to-global mapping of
  /// degrees of freedom (dofs).

00176   class dof_map
  {
  public:

    /// Constructor
00181     dof_map() {}

    /// Destructor
00184     virtual ~dof_map() {}

    /// Return a string identifying the dof map
    virtual const char* signature() const = 0;

    /// Return true iff mesh entities of topological dimension d are needed
    virtual bool needs_mesh_entities(unsigned int d) const = 0;

    /// Initialize dof map for mesh (return true iff init_cell() is needed)
    virtual bool init_mesh(const mesh& mesh) = 0;

    /// Initialize dof map for given cell
    virtual void init_cell(const mesh& m,
                           const cell& c) = 0;

    /// Finish initialization of dof map for cells
    virtual void init_cell_finalize() = 0;

    /// Return the dimension of the global finite element function space
    virtual unsigned int global_dimension() const = 0;

    /// Return the dimension of the local finite element function space for a cell
    virtual unsigned int local_dimension(const cell& c) const = 0;

    /// Return the maximum dimension of the local finite element function space
    virtual unsigned int max_local_dimension() const = 0;

    // Return the geometric dimension of the coordinates this dof map provides
    virtual unsigned int geometric_dimension() const
    { throw std::runtime_error("Not implemented (introduced in UFC v1.1)."); }

    /// Return the number of dofs on each cell facet
    virtual unsigned int num_facet_dofs() const = 0;

    /// Return the number of dofs associated with each cell entity of dimension d
00219     virtual unsigned int num_entity_dofs(unsigned int d) const
    { throw std::runtime_error("Not implemented (introduced in UFC v1.1)."); }

    /// Tabulate the local-to-global mapping of dofs on a cell
    virtual void tabulate_dofs(unsigned int* dofs,
                               const mesh& m,
                               const cell& c) const = 0;

    /// Tabulate the local-to-local mapping from facet dofs to cell dofs
    virtual void tabulate_facet_dofs(unsigned int* dofs,
                                     unsigned int facet) const = 0;

    /// Tabulate the local-to-local mapping of dofs on entity (d, i)
00232     virtual void tabulate_entity_dofs(unsigned int* dofs,
                                      unsigned int d, unsigned int i) const
    { throw std::runtime_error("Not implemented (introduced in UFC v1.1)."); }

    /// Tabulate the coordinates of all dofs on a cell
    virtual void tabulate_coordinates(double** coordinates,
                                      const cell& c) const = 0;

    /// Return the number of sub dof maps (for a mixed element)
    virtual unsigned int num_sub_dof_maps() const = 0;

    /// Create a new dof_map for sub dof map i (for a mixed element)
    virtual dof_map* create_sub_dof_map(unsigned int i) const = 0;

  };

  /// This class defines the interface for the tabulation of the cell
  /// tensor corresponding to the local contribution to a form from
  /// the integral over a cell.

00252   class cell_integral
  {
  public:

    /// Constructor
00257     cell_integral() {}

    /// Destructor
00260     virtual ~cell_integral() {}

    /// Tabulate the tensor for the contribution from a local cell
    virtual void tabulate_tensor(double* A,
                                 const double * const * w,
                                 const cell& c) const = 0;

  };

  /// This class defines the interface for the tabulation of the
  /// exterior facet tensor corresponding to the local contribution to
  /// a form from the integral over an exterior facet.

00273   class exterior_facet_integral
  {
  public:

    /// Constructor
00278     exterior_facet_integral() {}

    /// Destructor
00281     virtual ~exterior_facet_integral() {}

    /// Tabulate the tensor for the contribution from a local exterior facet
    virtual void tabulate_tensor(double* A,
                                 const double * const * w,
                                 const cell& c,
                                 unsigned int facet) const = 0;

  };

  /// This class defines the interface for the tabulation of the
  /// interior facet tensor corresponding to the local contribution to
  /// a form from the integral over an interior facet.

00295   class interior_facet_integral
  {
  public:

    /// Constructor
00300     interior_facet_integral() {}

    /// Destructor
00303     virtual ~interior_facet_integral() {}

    /// Tabulate the tensor for the contribution from a local interior facet
    virtual void tabulate_tensor(double* A,
                                 const double * const * w,
                                 const cell& c0,
                                 const cell& c1,
                                 unsigned int facet0,
                                 unsigned int facet1) const = 0;

  };

  /// This class defines the interface for the assembly of the global
  /// tensor corresponding to a form with r + n arguments, that is, a
  /// mapping
  ///
  ///     a : V1 x V2 x ... Vr x W1 x W2 x ... x Wn -> R
  ///
  /// with arguments v1, v2, ..., vr, w1, w2, ..., wn. The rank r
  /// global tensor A is defined by
  ///
  ///     A = a(V1, V2, ..., Vr, w1, w2, ..., wn),
  ///
  /// where each argument Vj represents the application to the
  /// sequence of basis functions of Vj and w1, w2, ..., wn are given
  /// fixed functions (coefficients).

00330   class form
  {
  public:

    /// Constructor
00335     form() {}

    /// Destructor
00338     virtual ~form() {}

    /// Return a string identifying the form
    virtual const char* signature() const = 0;

    /// Return the rank of the global tensor (r)
    virtual unsigned int rank() const = 0;

    /// Return the number of coefficients (n)
    virtual unsigned int num_coefficients() const = 0;

    /// Return the number of cell integrals
    virtual unsigned int num_cell_integrals() const = 0;

    /// Return the number of exterior facet integrals
    virtual unsigned int num_exterior_facet_integrals() const = 0;

    /// Return the number of interior facet integrals
    virtual unsigned int num_interior_facet_integrals() const = 0;

    /// Create a new finite element for argument function i
    virtual finite_element* create_finite_element(unsigned int i) const = 0;

    /// Create a new dof map for argument function i
    virtual dof_map* create_dof_map(unsigned int i) const = 0;

    /// Create a new cell integral on sub domain i
    virtual cell_integral* create_cell_integral(unsigned int i) const = 0;

    /// Create a new exterior facet integral on sub domain i
    virtual exterior_facet_integral*
    create_exterior_facet_integral(unsigned int i) const = 0;

    /// Create a new interior facet integral on sub domain i
    virtual interior_facet_integral*
    create_interior_facet_integral(unsigned int i) const = 0;

  };

}

#endif

Generated by  Doxygen 1.6.0   Back to index