base_interface.h

Go to the documentation of this file.

#ifndef _pidomus_base_interface_h_
#define _pidomus_base_interface_h_

#include <deal.II/lac/linear_operator.h>
/* #include <deal.II/lac/linear_operator.h> */
/* #include <deal.II/lac/block_linear_operator.h> */
/* #include <deal.II/lac/packaged_operation.h> */
/* #include <deal.II/lac/solver_cg.h> */
#include <deal.II/base/sacado_product_type.h>


#include <deal2lkit/parsed_finite_element.h>
#include <deal2lkit/parsed_function.h>
#include <deal2lkit/parsed_mapped_functions.h>
#include <deal2lkit/parsed_dirichlet_bcs.h>

#include "copy_data.h"
#include "lac/lac_type.h"

using namespace pidomus;
template <int dim,int spacedim=dim, typename LAC=LATrilinos>
class BaseInterface : public ParameterAcceptor
{

public:

  virtual ~BaseInterface() {};

  BaseInterface(const std::string &name="",
                const unsigned int &n_comp=0,
                const unsigned int &n_matrices=0,
                const std::string &default_fe="FE_Q(1)",
                const std::string &default_component_names="u",
                const std::string &default_differential_components="");

  void init ();


  virtual void declare_parameters (ParameterHandler &prm);

  const std::vector<unsigned int> get_differential_blocks() const;

  virtual void postprocess_newly_created_triangulation(Triangulation<dim, spacedim> &tria) const;

  virtual void initialize_data(const typename LAC::VectorType &solution,
                               const typename LAC::VectorType &solution_dot,
                               const typename LAC::VectorType &explicit_solution,
                               const double t,
                               const double alpha) const;


  virtual void assemble_energies_and_residuals(const typename DoFHandler<dim,spacedim>::active_cell_iterator &,
                                               FEValuesCache<dim,spacedim> &,
                                               std::vector<Sdouble> &energies,
                                               std::vector<std::vector<double> > &local_residuals,
                                               bool compute_only_system_terms) const;

  virtual void assemble_energies_and_residuals(const typename DoFHandler<dim,spacedim>::active_cell_iterator &,
                                               FEValuesCache<dim,spacedim> &,
                                               std::vector<SSdouble> &energies,
                                               std::vector<std::vector<Sdouble> > &local_residuals,
                                               bool compute_only_system_terms) const;

  virtual void assemble_local_matrices (const typename DoFHandler<dim,spacedim>::active_cell_iterator &cell,
                                        FEValuesCache<dim,spacedim> &scratch,
                                        CopyData &data) const;

  virtual void assemble_local_system_residual (const typename DoFHandler<dim,spacedim>::active_cell_iterator &cell,
                                               FEValuesCache<dim,spacedim> &scratch,
                                               CopyData &data) const;

  virtual void compute_system_operators(const DoFHandler<dim,spacedim> &,
                                        const std::vector<shared_ptr<typename LATrilinos::BlockMatrix> >,
                                        LinearOperator<typename LATrilinos::VectorType> &,
                                        LinearOperator<typename LATrilinos::VectorType> &) const;

  void compute_system_operators(const DoFHandler<dim,spacedim> &,
                                const std::vector<shared_ptr<typename LADealII::BlockMatrix> >,
                                LinearOperator<typename LADealII::VectorType> &,
                                LinearOperator<typename LADealII::VectorType> &) const
  {};


  virtual const Mapping<dim,spacedim> &get_mapping() const;

  virtual UpdateFlags get_face_update_flags() const;

  virtual UpdateFlags get_cell_update_flags() const;

  void fix_solution_dot_derivative(FEValuesCache<dim,spacedim> &, double) const;

  void fix_solution_dot_derivative(FEValuesCache<dim,spacedim> &fe_cache, Sdouble alpha) const;


  void fix_solution_dot_derivative(FEValuesCache<dim,spacedim> &fe_cache, SSdouble alpha) const;

  template<typename Number>
  void reinit(const Number &alpha,
              const typename DoFHandler<dim,spacedim>::active_cell_iterator &cell,
              FEValuesCache<dim,spacedim> &fe_cache) const;


  template<typename Number>
  void reinit(const Number &alpha,
              const typename DoFHandler<dim,spacedim>::active_cell_iterator &cell,
              const unsigned int face_no,
              FEValuesCache<dim,spacedim> &fe_cache) const;



  const Table<2,DoFTools::Coupling> &get_matrix_coupling(const unsigned int &i) const;

  const unsigned int n_components;
  const unsigned int n_matrices;

  ParsedFiniteElement<dim,spacedim> pfe;

  std::string get_component_names() const;


protected:

  void build_couplings();

  virtual void set_matrix_couplings(std::vector<std::string> &couplings) const;

  void convert_string_to_int(const std::string &str_coupling,
                             std::vector<std::vector<unsigned int> > &int_coupling) const;

  Table<2, DoFTools::Coupling> to_coupling(const std::vector<std::vector<unsigned int> > &table) const;

  std::string str_diff_comp;

  std::vector<unsigned int> _diff_comp;

  mutable const typename LAC::VectorType *solution;

  mutable const typename LAC::VectorType *explicit_solution;

  mutable const typename LAC::VectorType *solution_dot;

  mutable double t;
  mutable double alpha;
  unsigned int dofs_per_cell;
  unsigned int n_q_points;
  unsigned int n_face_q_points;

  std::vector<Table<2,DoFTools::Coupling> > matrix_couplings;


};

template <int dim, int spacedim, typename LAC>
template<typename Number>
void
BaseInterface<dim,spacedim,LAC>::reinit(const Number &,
                                        const typename DoFHandler<dim,spacedim>::active_cell_iterator &cell,
                                        FEValuesCache<dim,spacedim> &fe_cache) const
{
  double dummy=0;
  Number alpha = this->alpha;
  fe_cache.reinit(cell);
  fe_cache.cache_local_solution_vector("explicit_solution", *this->explicit_solution, dummy);
  fe_cache.cache_local_solution_vector("solution", *this->solution, alpha);
  fe_cache.cache_local_solution_vector("solution_dot", *this->solution_dot, alpha);
  this->fix_solution_dot_derivative(fe_cache, alpha);
}


template <int dim, int spacedim, typename LAC>
template<typename Number>
void
BaseInterface<dim,spacedim,LAC>::reinit(const Number &,
                                        const typename DoFHandler<dim,spacedim>::active_cell_iterator &cell,
                                        const unsigned int face_no,
                                        FEValuesCache<dim,spacedim> &fe_cache) const
{
  double dummy=0;
  Number alpha = this->alpha;
  fe_cache.reinit(cell, face_no);
  fe_cache.cache_local_solution_vector("explicit_solution", *this->explicit_solution, dummy);
  fe_cache.cache_local_solution_vector("solution", *this->solution, alpha);
  fe_cache.cache_local_solution_vector("solution_dot", *this->solution_dot, alpha);
  this->fix_solution_dot_derivative(fe_cache, alpha);
}

#endif