field(2rheolef)

πŸ“¦
version 7.2 Β· 4 versions
, in debian trixie, ubuntu noble, ubuntu questing

✨ NAME ✨

field β€” piecewise polynomial finite element function. this page lives in the rheolef-7.2 manual, section 2rheolef, and she's a section-2 girlie (library calls, not commands, not kernel β€” iykyk πŸ’…)


πŸŽ€ DESCRIPTION πŸŽ€

okay so bestie, listen. this class is a piecewise polynomial finite element function. periodt.

and the way it pulls that off is genuinely so elegant: since the function spans onto a finite dimensional basis, it doesn't need to store anything fancy β€” it just stores its coefficients on that basis. those coefficients? they have a name, they're the degrees of freedom (go peek space(2) for the full lore 🫢).


πŸ¦‹ INTERPOLATION πŸ¦‹

for any function u, its interpolation on the finite element space Xh β€” as a field uh in Xh β€” goes down simply via the interpolate(3) function:

Float u (const point& x) { return exp(x[0]*x[1]); }
...
field uh = interpolate (Xh, u);

like. that's it. that's the whole move. iconic ✨


πŸ”₯ LINEAR ALGEBRA πŸ”₯

linear algebra is fully supported, no cap:

  • uh+vh β€” slay
  • uh-vh β€” also slay
  • lambda*uh + mu*vh β€” where lambda and mu are of type Float β€” slay, but scaled

the duality product between two fields lh and uh writes simply dual(lh,uh). and here's the cute part: since we're considering finite dimensional spaces, this duality product coincides with the usual Euclidean dot product in IR^dim(Xh).

and the application of a bilinear form(2) writes a(uh,vh), which is equivalent to dual(m*uh,vh).


πŸ’– COMMON ACCESSORS πŸ’–

for convenience (and she is convenient):

  • uh.max() β€” the maximum of the degrees of freedom
  • uh.min() β€” the minimum of the degrees of freedom
  • uh.max_abs() β€” the maximum of the absolute value of the degrees of freedom

three little accessors, all eating, no notes 🌸


🌈 NON-LINEAR ALGEBRA 🌈

non-linear operations β€” think sqrt(uh) or 1/uh β€” are also available. but bestie, hear me:

note that non-linear operations do not return in general piecewise polynomials. 😭

the value returned by sqrt(uh) may be filtered by the interpolate(3) function`:

field vh = interpolate (Xh, sqrt(uh));

also, the multiplication uh*vh and the division uh/vh return a result that is not in the same discrete finite element space β€” so that result also may be filtered by interpolate(3):

field wh = interpolate(Xh, uh*vh);

all standard unary and binary math functions β€” abs, cos, sin… β€” are extended to scalar fields. plus you get sqr(uh) (the square of a field) and min(uh,vh), max(uh,vh). and binary functions can be used with a scalar too, as in:

field vh = interpolate (Xh, max (abs(uh), 0));
field wh = interpolate (Xh, pow (abs(uh), 1./3));

for applying a user-provided function to a field, please see the compose(3) function 🫢


🧚 ACCESS BY DOMAIN 🧚

the restriction of a field to a geometric domain β€” say, 'boundary' β€” writes uh['boundary']. it represents the trace of the field on the boundary:

space Xh (omega, 'P1');
uh['boundary'] = 0;

a possible alternative uses a geo(2) domain as the index:

geo boundary = omega['boundary'];
uh[boundary] = 0;

two paths, same destination. main character energy πŸ‘‘


πŸ“ MULTI-VALUED FIELDS πŸ“

a vector-valued field contains several components, like:

space Xh (omega, 'P2', 'vector');
field uh (Xh);

and conversely, for a tensor-valued field:

space Th (omega, 'P1d', 'tensor');
field sigma_h (Xh);

πŸ’œ GENERAL SPACE PRODUCT 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

now remark the hierarchical multi-component field structure: the first component is tensor-valued and the second one is vector-valued. and there is no limitation upon the hierarchical number of levels in use β€” the hierarchy is not flattened. the nesting stays nested. periodt.

xh.size() returns the number of field components. when the field is scalar, it returns zero by convention, and xh[0] is undefined ⚠️


⭐ DIRECT ACCESS TO DEGREES-OF-FREEDOM ⭐

the field class provides a STL-like container interface for accessing the degrees-of-freedom (dofs) of a finite element field uh.

  • the number of dofs is uh.ndof()
  • any dof can be accessed via uh.dof(idof)
  • in a distributed memory environment, a non-local dof at the partition interface can be obtained via uh.dis_dof(dis_idof) β€” where dis_idof is the (global) distributed index associated to the distribution uh.ownership(). see distributor(4).

and for better performances, a STL-like iterator interface is available, with uh.begin_dof() and uh.end_dof() returning iterators to the dofs on the current processor πŸš€


πŸŽ€ REPRESENTATION πŸŽ€

the degrees of freedom (see space(2)) are split between unknowns and blocked, i.e.

uh=[uh.u,uh.b]

for any field uh in Xh. access to these vectors is allowed via some accessors:

  • a read-only one, as uh.u() and uh.b()
  • a read-and-write one, as uh.set_u() and uh.set_b()

see vec(4).

and note that blocked and unknown degrees of freedom could also be elegantly set by using a domain name indexation (see geo(2)):

geo omega ('circle');
space Xh (omega, 'P1');
Xh.block ('boundary');
field uh (Xh);
uh ['boundary'] = 0;

elegant. she said what she said πŸ’…


🌸 IMPLEMENTATION 🌸

this documentation has been generated from file main/lib/field.h

the field class is simply an alias to the field_basic class:

typedef field_basic<Float> field;

the field_basic class provides an interface to a vector data container:

template <class T, class M = rheo_default_memory_model>
class field_basic: public details::field_wdof_base<field_basic<T,M>> {
public :
// typedefs:

  using wdof_base     = details::field_wdof_base<field_basic<T,M>>;
  using rdof_base     = details::field_rdof_base<field_basic<T,M>>;
  using size_type     = std::size_t;
  using scalar_type   = T;
  using memory_type   = M;
  using float_type    = typename float_traits<T>::type;
  using geo_type      = geo_basic  <float_type,memory_type>;
  using space_type    = space_basic<float_type,memory_type>;
  using dis_reference = typename vec<scalar_type,memory_type>::dis_reference;
  using valued_type   = space_constant::valued_type;
  using value_type    = T; // TODO: run-time dependent, set it to undeterminated<T>
  using result_type   = T; // TODO: run-time dependent, set it to undeterminated<T>
  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());

  // expressions: field_expr or field_lazy are accepted here
  // TODO: merge with FieldLazy piecewise_poly
  template <class Expr, class Sfinae =
    typename std::enable_if<
      (      details::is_field_expr_affine_homogeneous<Expr>::value
        && ! details::is_field<Expr>::value
      )
      ||   details::is_field_lazy<Expr>::value
    >::type>
  field_basic (const Expr& expr);

  field_basic (const std::initializer_list<details::field_concat_value<T,M> >& init_list);

// assignments:

  field_basic<T,M>& operator= (const field_basic<T,M>&);
  field_basic<T,M>& operator= (const std::initializer_list<details::field_concat_value<T,M> >& init_list);

// accessors:

  const space_type&  get_space()  const { return _V; }
  const geo_type&    get_geo()    const { return _V.get_geo(); }
  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(); }
#ifdef TO_CLEAN
  std::string        name()      const { return _V.name(); }
  bool have_homogeneous_space (space_basic<T,M>& Xh) const { Xh = get_space(); return true; }
#endif // TO_CLEAN

// accessors & modifiers to unknown & blocked parts:

  const vec<T,M>&     u() const { return _u; }
  const vec<T,M>&     b() const { return _b; }
        vec<T,M>& set_u()       { dis_dof_indexes_requires_update(); return _u; }
        vec<T,M>& set_b()       { dis_dof_indexes_requires_update(); return _b; }

// accessors to extremas:

  T min() const;
  T max() const;
  T max_abs() const;
  T min_abs() 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;
  // write access to non-local dis_idof changes to others procs
  dis_reference dis_dof_entry (size_type dis_idof);

  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); }
  point_basic<T> dis_vector_evaluate (const point_basic<T>& x) const;

};

namespace details {

// concepts:
template<class T, class M>
struct is_field<field_basic<T,M>> : std::true_type {};

// field class is not an ordinary function/functor : for compose(field,args...) filtering
template <typename T, typename M>  struct function_traits <field_basic<T,M> > {};

template<class T, class M>
struct field_traits<field_basic<T,M>> {
  using size_type     = std::size_t;
  using scalar_type   = T;
  using memory_type   = M;
};

} // namespace details

look at all those using declarations. she's organized. she's structured. she's that girl πŸ’…βœ¨


🫢 AUTHOR 🫢

Pierre Saramito <Pierre.Saramito@imag.fr>


Copyright (C) 2000-2018 Pierre Saramito <Pierre.Saramito@imag.fr>

GPLv3+: GNU GPL version 3 or later <http://gnu.org/licenses/gpl.html>.

This is free software: you are free to change and redistribute it.

There is NO WARRANTY, to the extent permitted by law. 🫣

and that's the whole page bestie, we made it πŸŽ‰ this was rheolef, Version 7.2, field(2rheolef) β€” go interpolate something cute, ily, bye πŸ’–πŸ¦‹πŸŒΈ

πŸ“¦ where this page comes from