// This file was automatically generated by FFC, the FEniCS Form Compiler.
// Licensed under the GNU GPL Version 2.

#ifndef __POISSON_BILINEAR_H
#define __POISSON_BILINEAR_H

#include <dolfin/NewFiniteElement.h>
#include <dolfin/BilinearForm.h>

using namespace dolfin;

/// This is the finite element for which the form is generated,
/// providing the information neccessary to do assembly.

class PoissonFiniteElement : public NewFiniteElement
{
public:

  PoissonFiniteElement() : NewFiniteElement() {}

  unsigned int spacedim() const
  {
    return 4;
  }

  unsigned int shapedim() const
  {
    return 3;
  }

  unsigned int tensordim(unsigned int i) const
  {
    unsigned int tensordims[] = {};
    return tensordims[i];
  }

  unsigned int rank() const
  {
    return 0;
  }

  // FIXME: Only works for nodal basis
  unsigned int dof(unsigned int i, const Cell& cell) const
  {
    return cell.id() * spacedim() + i;
  }

  // FIXME: Only works for nodal basis
  const Point& coord(unsigned int i, const Cell& cell) const
  {
    return cell.node(i).coord();
  }

};

/// This class contains the form to be evaluated, including
/// contributions from the interior and boundary of the domain.

class PoissonBilinearForm : public BilinearForm
{
public:

  PoissonBilinearForm(const NewFiniteElement& element) : BilinearForm(element) {}

  bool interior(real** A) const
  {
    // Compute geometry tensors
    real G0_00 = det*(g00*g00 + g01*g01 + g02*g02);
    real G0_01 = det*(g00*g10 + g01*g11 + g02*g12);
    real G0_02 = det*(g00*g20 + g01*g21 + g02*g22);
    real G0_10 = det*(g10*g00 + g11*g01 + g12*g02);
    real G0_11 = det*(g10*g10 + g11*g11 + g12*g12);
    real G0_12 = det*(g10*g20 + g11*g21 + g12*g22);
    real G0_20 = det*(g20*g00 + g21*g01 + g22*g02);
    real G0_21 = det*(g20*g10 + g21*g11 + g22*g12);
    real G0_22 = det*(g20*g20 + g21*g21 + g22*g22);

    // Compute element tensor
    A[0][0] = 0.166666666667*G0_00 + 0.166666666667*G0_01 + 0.166666666667*G0_02 + 0.166666666667*G0_10 + 0.166666666667*G0_11 + 0.166666666667*G0_12 + 0.166666666667*G0_20 + 0.166666666667*G0_21 + 0.166666666667*G0_22;
    A[0][1] = -0.166666666667*G0_00 - 0.166666666667*G0_10 - 0.166666666667*G0_20;
    A[0][2] = -0.166666666667*G0_01 - 0.166666666667*G0_11 - 0.166666666667*G0_21;
    A[0][3] = -0.166666666667*G0_02 - 0.166666666667*G0_12 - 0.166666666667*G0_22;
    A[1][0] = -0.166666666667*G0_00 - 0.166666666667*G0_01 - 0.166666666667*G0_02;
    A[1][1] = 0.166666666667*G0_00;
    A[1][2] = 0.166666666667*G0_01;
    A[1][3] = 0.166666666667*G0_02;
    A[2][0] = -0.166666666667*G0_10 - 0.166666666667*G0_11 - 0.166666666667*G0_12;
    A[2][1] = 0.166666666667*G0_10;
    A[2][2] = 0.166666666667*G0_11;
    A[2][3] = 0.166666666667*G0_12;
    A[3][0] = -0.166666666667*G0_20 - 0.166666666667*G0_21 - 0.166666666667*G0_22;
    A[3][1] = 0.166666666667*G0_20;
    A[3][2] = 0.166666666667*G0_21;
    A[3][3] = 0.166666666667*G0_22;

    return true;
  }

};

#endif