bezier.hh

Go to the documentation of this file.

/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
/*                                                                           */
/*  This file is part of the library KASKADE 7                               */
/*  https://www.zib.de/research/projects/kaskade7-finite-element-toolbox     */
/*                                                                           */
/*  Copyright (C) 2019-2020 Zuse Institute Berlin                            */
/*                                                                           */
/*  KASKADE 7 is distributed under the terms of the ZIB Academic License.    */
/*    see $KASKADE/academic.txt                                              */
/*                                                                           */
/* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
 
#ifndef BEZIER_HH
#define BEZIER_HH
 
#include "dune/common/fvector.hh"
 
#include "fem/barycentric.hh"
#include "utilities/combinatorics.hh"
#include "utilities/power.hh"
 
 
namespace Kaskade
{
  template <int d, class Real>
  Real bezier(int p, int i, Dune::FieldVector<Real,d> const& x)
  {
    auto b = barycentric2(x);
    if (d==1)
      return binomial(p,i)*power(b[0],p-i)*power(b[1],i);
    else if (d==2)
    {
      // First get the i-th multiindex of integers such that k[0]+...+k[d-1] <= p.
      auto ks = multiIndex<d>(p,i);
      // Then compute k[d] such that k[0]+...+k[d] == p
      std::array<size_t,d+1> Ks;
      std::copy(begin(ks),end(ks),begin(Ks));
      Ks[d] = p - std::accumulate(begin(ks),end(ks),0);
 
      // Now the Bezier function is multinomial(p,Ks) * lambda_0^Ks[0] * ... * lambda_d^Ks[d].
 
      Real r = multinomial(Ks);
      for (int j=0; j<=d; ++j)
        r *= power(b[j],Ks[j]);
      return r;
    }
 
    static_assert(d<=2,"not yet implemented");
    return -1; // never get here
  }
 
 
 
}
 
#endif