field(2rheolef) rheolef-6.1 field(2rheolef)
NAME
field - piecewise polynomial finite element field
DESCRIPTION
Store degrees of freedom associated to a mesh and a piecewise polynomial approximation, with respect to the numbering defined by the under-
lying space(2).
This class contains two vectors, namely unknown and blocked degrees of freedoms, and the associated finite element space. Blocked and
unknown degrees of freedom can be set by using domain name indexation:
geo omega ("circle");
space Vh (omega, "P1");
Vh.block ("boundary");
field uh (Vh);
uh ["boundary"] = 0;
INTERPOLATION
Interpolation of a function u in a field uh with respect to the interpolation writes:
Float u (const point& x) { return x[0]*x[1]; }
...
field uh = interpolate (Vh, u);
LINEAR ALGEBRA
Linear algebra, such as uh+vh, uh-vh, lambda*uh + mu*vh, ...are supported. The Euclidian dot product writes simply dot(uh,vh). The appli-
cation of a bilinear form (see form(2)) writes m(uh,vh) and is equivalent to dot(m*uh,vh).
NON-LINEAR ALGEBRA
Non-linear operations, such as sqrt(uh) or 1/uh are also available. Notice that non-linear operations do not returns in general picewise
polynomials: the value returned by sqrt(uh) is the Lagrange interpolant of the square root of uh. These operators applies directly to the
set of degrees-of-freedom (dofs, see below). All standard unary and binary math functions abs, cos, sin... are available on fields. Also
sqr(uh), the square of a field, and min(uh,vh), max(uh,vh) are provided. Binary functions can be used also with a scalar, as in
field wh = max (abs(uh), 0);
field zh = pow (abs(uh), 1./3);
For applying a user-provided function to a field, use:
field vh = compose(f, uh);
field wh = compose(f, uh, vh);
The composition supports also general unary and binary class-functions.
For convenience, uh.max(), uh.min() and uh.max_abs() retuns respectively the maximum, minimum and maximum of the absolute value.
ACCESS BY DOMAIN
The restriction of a field to a geometric domain, says "boundary" writes uh["boundary"]: it represents the trace of the field on the bound-
ary:
space Vh (omega, "P1");
uh["boundary"] = 0;
Extraction of the trace as a field is also possible:
field wh = uh["boundary"];
The space associated to the trace writes wh.get_space() and is equivalent to space(omega["boundary"], "P1"). See see space(2).
VECTOR VALUED FIELD
A vector-valued field contains several components, as:
space Vh (omega, "P2", "vector");
field uh (Vh);
field vh = uh[0] - uh[1];
field nh = norm (uh);
The norm function returns the euclidian norm of the vector-valuated field at each degree of freedom: its result is a scalar field.
TENSOR VALUED FIELD
A tensor-valued field can be constructed and used as:
space Th (omega, "P1d", "tensor");
field sigma_h (Vh);
field trace_h = sigma(0,0) + sigma_h(1,1);
field nh = norm (sigma_h);
The norm function returns the euclidian norm of the tensor-valuated field at each degree of freedom: its result is a scalar field. Notice
that, as tensor-valued fields are symmetric, extra-diagonals are counted twice.
GENERAL MULTI-COMPONENT INTERFACE
A general multi-component field writes:
space Th (omega, "P1d", "tensor");
space Vh (omega, "P2", "vector");
space Qh (omega, "P1");
space Xh = Th*Vh*Qh;
field xh (Xh);
field tau_h = xh[0]; // tensor-valued
field uh = xh[1]; // vector-valued
field qh = xh[2]; // scalar
Remark the hierarchical multi-component field structure: the first-component is tensor-valued and the second-one is vector-valued. There
is no limitation upon the hierarchical number of levels in use.
For any field xh, the string xh.valued() returns "scalar" for a scalar field and "vector" for a vector-valued one. Other possible valued
are "tensor" and "other". The xh.size() returns the number of field components. When the field is scalar, it returns zero by convention,
and xh[0] is undefined. A vector-valued field has d components, where d=omega.dimension(). A tensor-valued field has d*(d+1)/2 compo-
nents, where d=omega.dimension().
BLOCKED AND UNBLOCKED ARRAYS
The field class contains two arrays of degrees-of-freedom (dof) associated respectively to blocked and unknown dofs. Blocked dofs corre-
sponds to Dirichlet boundary conditions, as specified by space (See see space(2)). For simpliity, direct public access to these array is
allowed, as uh.b and uh.u: see see vec(2).
LOW-LEVEL DEGREE-OF-FREEDOM ACCESS
The field class provides a STL-like container interface for accessing the degrees-of-freedom (dof) of a finite element field uh. The number
of dofs is uh.ndof() and any dof can be accessed via uh.dof(idof). A non-local dof at the partition interface can be obtain via
uh.dis_dof(dis_idof) where dis_idof is the (global) distribued index assoiated to the distribution uh.ownership().
For performances, a STL-like iterator interface is available, with uh.begin_dof() and uh.end_dof() returns iterators to the arrays of dofs
on the current processor. See see array(2) for more about distributed arrays.
FILE FORMAT
TODO
IMPLEMENTATION NOTE
The field expression use the expression template technics in order to avoid temporaries when evaluating complex expressions.
IMPLEMENTATION
template <class T, class M = rheo_default_memory_model>
class field_basic : public std::unary_function<point_basic<typename scalar_traits<T>::type>,T> {
public :
// typedefs:
typedef typename std::size_t size_type;
typedef T value_type;
typedef typename scalar_traits<T>::type float_type;
typedef geo_basic <float_type,M> geo_type;
typedef space_basic<float_type,M> space_type;
typedef space_constant::valued_type valued_type;
class iterator;
class const_iterator;
// allocator/deallocator:
field_basic();
explicit field_basic (
const space_type& V,
const T& init_value = std::numeric_limits<T>::max());
void resize (
const space_type& V,
const T& init_value = std::numeric_limits<T>::max());
field_basic (const field_indirect<T,M>&);
field_basic<T, M>& operator= (const field_indirect<T,M>&);
field_basic (const field_indirect_const<T,M>&);
field_basic<T, M>& operator= (const field_indirect_const<T,M>&);
field_basic (const field_component<T,M>&);
field_basic<T, M>& operator= (const field_component<T,M>&);
field_basic (const field_component_const<T,M>&);
field_basic<T, M>& operator= (const field_component_const<T,M>&);
template <class Expr> field_basic (const field_expr<Expr>&);
template <class Expr> field_basic<T, M>& operator= (const field_expr<Expr>&);
field_basic<T, M>& operator= (const T&);
// initializer list (c++ 2011):
#ifdef _RHEOLEF_HAVE_STD_INITIALIZER_LIST
field_basic (const std::initializer_list<field_concat_value<T,M> >& init_list);
field_basic<T,M>& operator= (const std::initializer_list<field_concat_value<T,M> >& init_list);
#endif // _RHEOLEF_HAVE_STD_INITIALIZER_LIST
// accessors:
const space_type& get_space() const { return _V; }
const geo_type& get_geo() const { return _V.get_geo(); }
std::string stamp() const { return _V.stamp(); }
std::string get_approx() const { return _V.get_approx(); }
valued_type valued_tag() const { return _V.valued_tag(); }
const std::string& valued() const { return _V.valued(); }
// accessors & modifiers to unknown & blocked parts:
const vec<T,M>& u() const { dis_dof_update_needed(); return _u; }
const vec<T,M>& b() const { dis_dof_update_needed(); return _b; }
vec<T,M>& set_u() { return _u; }
vec<T,M>& set_b() { return _b; }
// accessors to extremas:
T min() const;
T max() const;
T max_abs() const;
T min_abs() const;
// accessors by domains:
field_indirect<T,M> operator[] (const geo_basic<T,M>& dom);
field_indirect_const<T,M> operator[] (const geo_basic<T,M>& dom) const;
field_indirect<T,M> operator[] (std::string dom_name);
field_indirect_const<T,M> operator[] (std::string dom_name) const;
// accessors by components:
size_type size() const { return _V.size(); }
field_component<T,M> operator[] (size_type i_comp);
field_component_const<T,M> operator[] (size_type i_comp) const;
field_component<T,M> operator() (size_type i_comp, size_type j_comp);
field_component_const<T,M> operator() (size_type i_comp, size_type j_comp) const;
// accessors by degrees-of-freedom (dof):
const distributor& ownership() const { return get_space().ownership(); }
const communicator& comm() const { return ownership().comm(); }
size_type ndof() const { return ownership().size(); }
size_type dis_ndof() const { return ownership().dis_size(); }
T& dof (size_type idof);
const T& dof (size_type idof) const;
const T& dis_dof (size_type dis_idof) const;
iterator begin_dof();
iterator end_dof();
const_iterator begin_dof() const;
const_iterator end_dof() const;
// input/output:
idiststream& get (idiststream& ips);
odiststream& put (odiststream& ops) const;
odiststream& put_field (odiststream& ops) const;
// evaluate uh(x) where x is given locally as hat_x in K:
T dis_evaluate (const point_basic<T>& x, size_type i_comp = 0) const;
T operator() (const point_basic<T>& x) const { return dis_evaluate (x,0); }
// internals:
public:
// evaluate uh(x) where x is given locally as hat_x in K:
// requiers to call field::dis_dof_upgrade() before.
T evaluate (const geo_element& K, const point_basic<T>& hat_xq, size_type i_comp = 0) const;
// propagate changed values shared at partition boundaries to others procs
void dis_dof_update() const;
protected:
void dis_dof_update_internal() const;
void dis_dof_update_needed() const;
// data:
space_type _V;
vec<T,M> _u;
vec<T,M> _b;
mutable bool _dis_dof_update_needed;
};
template <class T, class M>
idiststream& operator >> (odiststream& ips, field_basic<T,M>& u);
template <class T, class M>
odiststream& operator << (odiststream& ops, const field_basic<T,M>& uh);
template <class T, class M>
field_basic<T,M> norm (const field_basic<T,M>& uh);
template <class T, class M>
field_basic<T,M> norm2 (const field_basic<T,M>& uh);
typedef field_basic<Float> field;
typedef field_basic<Float,sequential> field_sequential;
SEE ALSO
space(2), form(2), space(2), space(2), vec(2), array(2)
rheolef-6.1 rheolef-6.1 field(2rheolef)