p5tri.h

finite element library for mathematic majored research

// This code conforms with the UFC specification version 1.0
// and was automatically generated by FFC version 0.4.3.

#ifndef __P5TRI_H
#define __P5TRI_H

#include <cmath>
#include <stdexcept>
#include <ufc.h>

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

class P5tri_finite_element_0: public ufc::finite_element
{
public:

  /// Constructor
  P5tri_finite_element_0() : ufc::finite_element()
  {
    // Do nothing
  }

  /// Destructor
  virtual ~P5tri_finite_element_0()
  {
    // Do nothing
  }

  /// Return a string identifying the finite element
  virtual const char* signature() const
  {
    return "Lagrange finite element of degree 5 on a triangle";
  }

  /// Return the cell shape
  virtual ufc::shape cell_shape() const
  {
    return ufc::triangle;
  }

  /// Return the dimension of the finite element function space
  virtual unsigned int space_dimension() const
  {
    return 21;
  }

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

  /// Return the dimension of the value space for axis i
  virtual unsigned int value_dimension(unsigned int i) const
  {
    return 1;
  }

  /// Evaluate basis function i at given point in cell
  virtual void evaluate_basis(unsigned int i,
                              double* values,
                              const double* coordinates,
                              const ufc::cell& c) const
  {
    // Extract vertex coordinates
    const double * const * element_coordinates = c.coordinates;
    
    // Compute Jacobian of affine map from reference cell
    const double J_00 = element_coordinates[1][0] - element_coordinates[0][0];
    const double J_01 = element_coordinates[2][0] - element_coordinates[0][0];
    const double J_10 = element_coordinates[1][1] - element_coordinates[0][1];
    const double J_11 = element_coordinates[2][1] - element_coordinates[0][1];
      
    // Compute determinant of Jacobian
    const double detJ = J_00*J_11 - J_01*J_10;
    
    // Get coordinates and map to the reference (UFC) element
    double x = (element_coordinates[0][1]*element_coordinates[2][0] -\
                element_coordinates[0][0]*element_coordinates[2][1] +\
                J_11*coordinates[0] - J_01*coordinates[1]) / detJ;
    double y = (element_coordinates[1][1]*element_coordinates[0][0] -\
                element_coordinates[1][0]*element_coordinates[0][1] -\
                J_10*coordinates[0] + J_00*coordinates[1]) / detJ;
    
    // Map coordinates to the reference square
    if (std::abs(y - 1.0) < 1e-09)
      x = -1.0;
    else
      x = 2.0 *x/(1.0 - y) - 1.0;
    y = 2.0*y - 1.0;
    
    // Reset values
    *values = 0;
    
    // Map degree of freedom to element degree of freedom
    const unsigned int dof = i;
    
    // Generate scalings
    const double scalings_y_0 = 1;
    const double scalings_y_1 = scalings_y_0*(0.5 - 0.5*y);
    const double scalings_y_2 = scalings_y_1*(0.5 - 0.5*y);
    const double scalings_y_3 = scalings_y_2*(0.5 - 0.5*y);
    const double scalings_y_4 = scalings_y_3*(0.5 - 0.5*y);
    const double scalings_y_5 = scalings_y_4*(0.5 - 0.5*y);
    
    // Compute psitilde_a
    const double psitilde_a_0 = 1;
    const double psitilde_a_1 = x;
    const double psitilde_a_2 = 1.5*x*psitilde_a_1 - 0.5*psitilde_a_0;
    const double psitilde_a_3 = 1.666666667*x*psitilde_a_2 - 0.6666666667*psitilde_a_1;
    const double psitilde_a_4 = 1.75*x*psitilde_a_3 - 0.75*psitilde_a_2;
    const double psitilde_a_5 = 1.8*x*psitilde_a_4 - 0.8*psitilde_a_3;
    
    // Compute psitilde_bs
    const double psitilde_bs_0_0 = 1;
    const double psitilde_bs_0_1 = 1.5*y + 0.5;
    const double psitilde_bs_0_2 = 0.1111111111*psitilde_bs_0_1 + 1.666666667*y*psitilde_bs_0_1 - 0.5555555556*psitilde_bs_0_0;
    const double psitilde_bs_0_3 = 0.05*psitilde_bs_0_2 + 1.75*y*psitilde_bs_0_2 - 0.7*psitilde_bs_0_1;
    const double psitilde_bs_0_4 = 0.02857142857*psitilde_bs_0_3 + 1.8*y*psitilde_bs_0_3 - 0.7714285714*psitilde_bs_0_2;
    const double psitilde_bs_0_5 = 0.01851851852*psitilde_bs_0_4 + 1.833333333*y*psitilde_bs_0_4 - 0.8148148148*psitilde_bs_0_3;
    const double psitilde_bs_1_0 = 1;
    const double psitilde_bs_1_1 = 2.5*y + 1.5;
    const double psitilde_bs_1_2 = 0.54*psitilde_bs_1_1 + 2.1*y*psitilde_bs_1_1 - 0.56*psitilde_bs_1_0;
    const double psitilde_bs_1_3 = 0.2857142857*psitilde_bs_1_2 + 2*y*psitilde_bs_1_2 - 0.7142857143*psitilde_bs_1_1;
    const double psitilde_bs_1_4 = 0.1785714286*psitilde_bs_1_3 + 1.964285714*y*psitilde_bs_1_3 - 0.7857142857*psitilde_bs_1_2;
    const double psitilde_bs_2_0 = 1;
    const double psitilde_bs_2_1 = 3.5*y + 2.5;
    const double psitilde_bs_2_2 = 1.020408163*psitilde_bs_2_1 + 2.571428571*y*psitilde_bs_2_1 - 0.5510204082*psitilde_bs_2_0;
    const double psitilde_bs_2_3 = 0.5787037037*psitilde_bs_2_2 + 2.291666667*y*psitilde_bs_2_2 - 0.712962963*psitilde_bs_2_1;
    const double psitilde_bs_3_0 = 1;
    const double psitilde_bs_3_1 = 4.5*y + 3.5;
    const double psitilde_bs_3_2 = 1.512345679*psitilde_bs_3_1 + 3.055555556*y*psitilde_bs_3_1 - 0.5432098765*psitilde_bs_3_0;
    const double psitilde_bs_4_0 = 1;
    const double psitilde_bs_4_1 = 5.5*y + 4.5;
    const double psitilde_bs_5_0 = 1;
    
    // Compute basisvalues
    const double basisvalue0 = 0.7071067812*psitilde_a_0*scalings_y_0*psitilde_bs_0_0;
    const double basisvalue1 = 1.732050808*psitilde_a_1*scalings_y_1*psitilde_bs_1_0;
    const double basisvalue2 = psitilde_a_0*scalings_y_0*psitilde_bs_0_1;
    const double basisvalue3 = 2.738612788*psitilde_a_2*scalings_y_2*psitilde_bs_2_0;
    const double basisvalue4 = 2.121320344*psitilde_a_1*scalings_y_1*psitilde_bs_1_1;
    const double basisvalue5 = 1.224744871*psitilde_a_0*scalings_y_0*psitilde_bs_0_2;
    const double basisvalue6 = 3.741657387*psitilde_a_3*scalings_y_3*psitilde_bs_3_0;
    const double basisvalue7 = 3.16227766*psitilde_a_2*scalings_y_2*psitilde_bs_2_1;
    const double basisvalue8 = 2.449489743*psitilde_a_1*scalings_y_1*psitilde_bs_1_2;
    const double basisvalue9 = 1.414213562*psitilde_a_0*scalings_y_0*psitilde_bs_0_3;
    const double basisvalue10 = 4.74341649*psitilde_a_4*scalings_y_4*psitilde_bs_4_0;
    const double basisvalue11 = 4.183300133*psitilde_a_3*scalings_y_3*psitilde_bs_3_1;
    const double basisvalue12 = 3.535533906*psitilde_a_2*scalings_y_2*psitilde_bs_2_2;
    const double basisvalue13 = 2.738612788*psitilde_a_1*scalings_y_1*psitilde_bs_1_3;
    const double basisvalue14 = 1.58113883*psitilde_a_0*scalings_y_0*psitilde_bs_0_4;