discrete_problem_thread_assembler.h

#ifndef __H2D_DISCRETE_PROBLEM_ASSEMBLY_DATA_H
#define __H2D_DISCRETE_PROBLEM_ASSEMBLY_DATA_H

#include "../weakform/weakform.h"
#include "../shapeset/precalc.h"
#include "../function/solution.h"
#include "discrete_problem_helpers.h"
#include "discrete_problem_integration_order_calculator.h"
#include "discrete_problem_selective_assembler.h"

namespace Hermes
{
  namespace Hermes2D
  {
    template<typename Scalar>
    class HERMES_API DiscreteProblemThreadAssembler :
      public Hermes::Hermes2D::Mixins::DiscreteProblemWeakForm<Scalar>,
      public Hermes::Hermes2D::Mixins::DiscreteProblemRungeKutta<Scalar>,
      public Hermes::Hermes2D::Mixins::DiscreteProblemMatrixVector < Scalar >
    {
    public:
      void free();

    private:
      DiscreteProblemThreadAssembler(DiscreteProblemSelectiveAssembler<Scalar>* selectiveAssembler, bool nonlinear);
      ~DiscreteProblemThreadAssembler();

      void init_spaces(const std::vector<SpaceSharedPtr<Scalar> > spaces);
      void set_weak_formulation(WeakFormSharedPtr<Scalar> wf);
      void init_u_ext(const std::vector<SpaceSharedPtr<Scalar> > spaces, Solution<Scalar>** u_ext_sln);

      void init_assembling(Solution<Scalar>** u_ext_sln, const std::vector<SpaceSharedPtr<Scalar> >& spaces, bool add_dirichlet_lift);

      void init_funcs_wf();
      void init_funcs_space();
      pj_pool_t *FuncMemoryPool;
      void init_funcs_memory_pool();

      void deinit_funcs();

      void deinit_funcs_space();
      bool funcs_space_initialized;
      void deinit_funcs_wf();
      bool funcs_wf_initialized;
      void init_u_ext_values(int order);
      template<typename Geom>
      void init_ext_values(Func<Scalar>** target_array, std::vector<MeshFunctionSharedPtr<Scalar> >& ext, std::vector<UExtFunctionSharedPtr<Scalar> >& u_ext_fns, int order, Func<Scalar>** u_ext_func, Geom* geometry);

      void init_assembling_one_state(const std::vector<SpaceSharedPtr<Scalar> >& spaces, Traverse::State* current_state);
      void assemble_one_state();
      template<typename MatrixFormType, typename Geom>
      void assemble_matrix_form(MatrixFormType* form, int order, Func<double>** base_fns, Func<double>** test_fns,
        AsmList<Scalar>* current_als_i, AsmList<Scalar>* current_als_j, int n_quadrature_points, Geom* geometry, double* jacobian_x_weights);
      template<typename VectorFormType, typename Geom>
      void assemble_vector_form(VectorFormType* form, int order, Func<double>** test_fns, AsmList<Scalar>* current_als,
        int n_quadrature_points, Geom* geometry, double* jacobian_x_weights);
      void deinit_assembling_one_state();

      void deinit_assembling();

      void free_spaces();
      void free_weak_formulation();
      void free_u_ext();

      Vector<Scalar>* dirichlet_lift_rhs;

      PrecalcShapesetAssembling** pss;
      RefMap** refmaps;
      RefMap* rep_refmap;
      Solution<Scalar>** u_ext;
      std::vector<Transformable *> fns;

      DiscreteProblemSelectiveAssembler<Scalar>* selectiveAssembler;

      Traverse::State* current_state;
      Scalar local_stiffness_matrix[H2D_MAX_LOCAL_BASIS_SIZE * H2D_MAX_LOCAL_BASIS_SIZE * 4];

      DiscreteProblemIntegrationOrderCalculator<Scalar> integrationOrderCalculator;
      int order;
      int orderSurface[H2D_MAX_NUMBER_EDGES];

      void init_calculation_variables();
      void deinit_calculation_variables();
      Func<double>* funcs[H2D_MAX_COMPONENTS][H2D_MAX_LOCAL_BASIS_SIZE];
      Func<double>* funcsSurface[H2D_MAX_NUMBER_EDGES][H2D_MAX_COMPONENTS][H2D_MAX_LOCAL_BASIS_SIZE];
      GeomVol<double> geometry;
      GeomSurf<double> geometrySurface[H2D_MAX_NUMBER_EDGES];
      double jacobian_x_weights[H2D_MAX_INTEGRATION_POINTS_COUNT];
      double jacobian_x_weightsSurface[H2D_MAX_NUMBER_EDGES][H2D_MAX_INTEGRATION_POINTS_COUNT];
      int n_quadrature_points;
      int n_quadrature_pointsSurface[H2D_MAX_NUMBER_EDGES];

      Func<Scalar>* u_ext_funcs[H2D_MAX_COMPONENTS];
      Func<Scalar>** ext_funcs;
      int ext_funcs_allocated_size;
      Func<Scalar>** ext_funcs_local;
      int ext_funcs_local_allocated_size;

      AsmList<Scalar> als[H2D_MAX_COMPONENTS];
      AsmList<Scalar> alsSurface[H2D_MAX_NUMBER_EDGES][H2D_MAX_COMPONENTS];
      unsigned short spaces_size;
      bool nonlinear, add_dirichlet_lift;

      friend class DiscreteProblem < Scalar > ;
      friend class DiscreteProblemDGAssembler < Scalar > ;

      bool** reusable_DOFs;
      bool** reusable_Dirichlet;
    };
  }
}
#endif