auto_derivative_function.h

// ---------------------------------------------------------------------
//
// Copyright (C) 2001 - 2017 by the deal.II authors
//
// This file is part of the deal.II library.
//
// The deal.II library is free software; you can use it, redistribute
// it, and/or modify it under the terms of the GNU Lesser General
// Public License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// The full text of the license can be found in the file LICENSE.md at
// the top level directory of deal.II.
//
// ---------------------------------------------------------------------

#ifndef dealii_auto_derivative_function_h
#define dealii_auto_derivative_function_h


#include <deal.II/base/exceptions.h>
#include <deal.II/base/function.h>

DEAL_II_NAMESPACE_OPEN

template <int dim>
class AutoDerivativeFunction : public Function<dim>
{
public:
  enum DifferenceFormula
  {
    Euler,
    UpwindEuler,
    FourthOrder
  };

  AutoDerivativeFunction(const double       h,
                         const unsigned int n_components = 1,
                         const double       initial_time = 0.0);

  virtual ~AutoDerivativeFunction() override = default;

  void
  set_formula(const DifferenceFormula formula = Euler);

  void
  set_h(const double h);

  virtual Tensor<1, dim>
  gradient(const Point<dim> & p,
           const unsigned int component = 0) const override;

  virtual void
  vector_gradient(const Point<dim> &           p,
                  std::vector<Tensor<1, dim>> &gradients) const override;

  virtual void
  gradient_list(const std::vector<Point<dim>> &points,
                std::vector<Tensor<1, dim>> &  gradients,
                const unsigned int             component = 0) const override;

  virtual void
  vector_gradient_list(
    const std::vector<Point<dim>> &           points,
    std::vector<std::vector<Tensor<1, dim>>> &gradients) const override;

  static DifferenceFormula
  get_formula_of_order(const unsigned int ord);


private:
  double h;

  std::vector<Tensor<1, dim>> ht;

  DifferenceFormula formula;
};


DEAL_II_NAMESPACE_CLOSE

#endif