/*
 This Source Code Form is subject to the terms of the Mozilla Public
 License, v. 2.0. If a copy of the MPL was not distributed with this
 file, You can obtain one at http://mozilla.org/MPL/2.0/.
*/


namespace Vitelotte {

/*!
  \page vitelotte_user_manual_fem_solver_page The FEM Solver

  \section vitelotte_user_manual_fem_solver_intro_sec Introduction

  FemSolver interpolates values (colors, normals, ...) smoothly through diffusion directly on a VGMesh using the finite element method (FEM). This is a simple way to get smooth color gradients from a sparse set of constraints.

  Vitelotte currently supports two kind of diffusions, with their advantages and drawbacks:

  - **Harmonic** diffusion satisfies the equation

    \f[ \Delta f = 0, \quad f(x_i) = v_i \f]

    in such a way that it interpolates constraints (nodes with a value set). This equation can be solved using very simple elements. Vitelotte currently supports *linear* and *quadratic* triangular elements for this, which uses the same interpolation functions than rendering.

    However, harmonic interpolation has some drawbacks. It does produce have \f$C^1\f$ transitions across constraints, which limits its usefulness in many cases. For instance, to represent smooth shading, biharmonic diffusion produces better results.

  - **Bi-harmonic** diffusion satisfies a similar equation

    \f[ \Delta^2 f = 0, \quad f(x_i) = v_i, \quad \nabla f(x_i) = g_i \f]

    Bi-harmonic diffusion is a bit more costly, but also offers smoother gradients and more control. For instance, it is possible to constrain derivative nodes to add control over gradients. A drawback of bi-harmonic diffusion is that it extrapolates: it can produce values higher or lower than constraints. When diffusing colors, it means that you can have color components below 0 and above the maximum, which may in turn lead to artifacts.

    Theoretically, solving this equation using FEM requires \f$C^1\f$ continuity across the elements, which in turn requires using high order elements. As this is unpractical, we use _non-conforming_ elements that do not meet this requirement but ensure convergence nonetheless: Morley elements and Fraeijs de Veubeke's elements (FV). They are similar to linear and quadratic elements respectively, but add a derivative node over each edge.

    These elements are not even \f$C^0\f$ on their boundaries, they are not suited for rendering. Instead, we interpolate value nodes per-element (thus ignoring gradient nodes), which lead to the linear and quadratic interpolation functions used for harmonic elements.



  \section vitelotte_user_manual_fem_solver_usage_sec Using the solver

  Using the solver is quite easy. You first need to include the header `Patate/vitelotte.h` to get access to the required classes.

  FemSolver takes two template parameters: a VGMesh type and an ElementBuilder. Currently, only triangular elements are supported. Here is a short description of available element builders:

  - **LinearElementBuilder** solves harmonic diffusion over linear elements. Input mesh must have the TO_VERTEX_VALUE attribute enabled.
  - **QuadraticElementBuilder** solves harmonic diffusion over quadratic elements. Input mesh must have the TO_VERTEX_VALUE and EDGE_VALUE attributes enabled.
  - **MorleyElementBuilder** solves biharmonic diffusion using Morley elements. Input mesh must have TO_VERTEX_VALUE and EDGE_GRADIENT attributes enabled.
  - **FVElementBuilder** solves biharmonic diffusion using Fraeijs de Veubeke's (FV) elements. Input mesh must have TO_VERTEX_VALUE, EDGE_VALUE and EDGE_GRADIENT attributes enabled.

  These element builders do not support singularities (in practice, the attachment point FROM_VERTEX_VALUE is ignored). To enable it, you must wrap your element builder into a SingularElementDecorator. This only supports up to one singularity by element, so it may be necessary to refine the input mesh in case 2 singularities are directly connected.

  Here is a typical example:

  \code{.cpp}
  typedef VGMesh<float, [...]> Mesh;

  typedef FVElementBuilder<Mesh, double> FVBuilder;
  typedef SingularElementDecorator<FVBuilder> ElementBuilder;
  typedef FemSolver<Mesh, ElementBuilder> Solver;
  \endcode

  Note that the element builders are parametrized by a scalar type. It will be used internally by the solver instead of the mesh Scalar type. This may be used to avoid numerical problems. In our example, we always use floats for meshes and double in the solver, as above.

  All is left to do is to create a solver and solve the problem:

  \code
  Mesh mesh;
  // [ Fill the mesh... ]

  Solver solver();
  solver.build(mesh);
  solver.solve(mesh);
  \endcode

  FemSolver::build() builds the internal matrices and compute factorizations. It does not depend on constraint values, so if you are using the solver in an interactive application (like our `mvg_editor` example), it is not always necessary to call it. FemSolver::solve() effectively solves the diffusion and sets a value for each unknown node of the mesh.

  The solver may fail, either because input is invalid or because of numerical instability. We try to detect errors and to provide meaningful messages in case of problems:

  \code
  SolverError::Status status = solver.error().status();
  if(status != SolverError::STATUS_OK) {
    cout << ((status == SolverError::STATUS_ERROR)? "Error: ": "Warning: ")
         << solver.error().message() << endl;
  }
  \endcode


 */

 }
