Hermes2D  3.0
proj_based_selector.h
1 // This file is part of Hermes2D.
2 //
3 // Hermes2D is free software: you can redistribute it and/or modify
4 // it under the terms of the GNU General Public License as published by
5 // the Free Software Foundation, either version 2 of the License, or
6 // (at your option) any later version.
7 //
8 // Hermes2D is distributed in the hope that it will be useful,
9 // but WITHOUT ANY WARRANTY; without even the implied warranty of
10 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
11 // GNU General Public License for more details.
12 //
13 // You should have received a copy of the GNU General Public License
14 // along with Hermes2D. If not, see <http://www.gnu.org/licenses/>.
15 
16 #ifndef __H2D_REFINEMENT_PROJ_BASED_SELECTOR_H
17 #define __H2D_REFINEMENT_PROJ_BASED_SELECTOR_H
18 
19 #include "optimum_selector.h"
20 
22 namespace Hermes
23 {
24  namespace Hermes2D
25  {
26  namespace RefinementSelectors {
28 
31 
33 
47  template<typename Scalar>
48  class HERMES_API ProjBasedSelector : public OptimumSelector < Scalar > {
49  protected:
50  class TrfShapeExp;
51 
52  public: //API
54  virtual ~ProjBasedSelector();
55 
57 
63  void set_error_weights(double weight_h = H2DRS_DEFAULT_ERR_WEIGHT_H, double weight_p = H2DRS_DEFAULT_ERR_WEIGHT_P, double weight_aniso = H2DRS_DEFAULT_ERR_WEIGHT_ANISO);
64 
65  double get_error_weight_h() const;
66  double get_error_weight_p() const;
67  double get_error_weight_aniso() const;
68 
70  typedef std::vector<TrfShapeExp> TrfShape[H2D_TRF_NUM];
71 
77  TrfShape* cached_shape_vals;
78 
79  protected: //evaluated shape basis
81 
90  class TrfShapeExp {
91  public:
93  TrfShapeExp();
94 
96  virtual ~TrfShapeExp();
97 
99  const TrfShapeExp& operator = (const TrfShapeExp& other)
100  {
101  free_with_check(values);
102  values = nullptr;
103 
104  if (other.values == nullptr)
105  throw Exceptions::Exception("Unable to assign a non-empty values. Use references instead.");
106  return *this;
107  }
108  private:
110  int num_gip;
112  int num_expansion;
114  double** values;
115 
117 
119  void allocate(int num_expansion, int num_gip);
120 
122 
124  inline double* operator[](int inx_expansion);
125 
127 
128  inline bool empty() const;
129 
130  template<typename T> friend class ProjBasedSelector;
131  template<typename T> friend class L2ProjBasedSelector;
132  template<typename T> friend class H1ProjBasedSelector;
133  template<typename T> friend class HcurlProjBasedSelector;
134  template<typename T> friend class Adapt;
135  };
136 
138 
152  virtual void precalc_shapes(const double3* gip_points, const int num_gip_points, const Trf* trfs, const int num_noni_trfs, const std::vector<typename OptimumSelector<Scalar>::ShapeInx>& shapes, const int max_shape_inx, TrfShape& svals, ElementMode2D mode) {};
153 
155 
172  virtual void precalc_ortho_shapes(const double3* gip_points, const int num_gip_points, const Trf* trfs, const int num_noni_trfs, const std::vector<typename OptimumSelector<Scalar>::ShapeInx>& shapes, const int max_shape_inx, TrfShape& ortho_svals, ElementMode2D mode) {};
173 
174  protected:
176 
182  ProjBasedSelector(CandList cand_list, int max_order, Shapeset* shapeset, const Range& vertex_order, const Range& edge_bubble_order);
183 
184  protected: //internal logic
186 
191 
192  protected: //error evaluation
194  static const int H2DRS_VALCACHE_INVALID = 0;
196  static const int H2DRS_VALCACHE_VALID = 1;
198  static const int H2DRS_VALCACHE_USER = 2;
199 
201  template<typename T>
202  struct ValueCacheItem {
204 
205  bool is_valid() const { return state != H2DRS_VALCACHE_INVALID; };
206 
208 
209  void mark(int new_state = H2DRS_VALCACHE_VALID) { state = new_state; };
210 
212 
213  void set(T new_value) { value = new_value; };
214 
216 
217  T get() const { return value; };
218 
220 
223  ValueCacheItem(const T& value = 0, const int state = H2DRS_VALCACHE_INVALID) : value(value), state(state) {};
225  private:
227  T value;
229  int state;
230  };
232 
233  typedef double** ProjMatrixCache[H2DRS_MAX_ORDER + 2][H2DRS_MAX_ORDER + 2];
234 
236 
241  ProjMatrixCache proj_matrix_cache[H2D_NUM_MODES];
242 
249 
251 
252  virtual void evaluate_cands_error(std::vector<Cand>& candidates, Element* e, MeshFunction<Scalar>* rsln);
253 
255 
268  virtual void calc_projection_errors(Element* e, const typename OptimumSelector<Scalar>::CandsInfo& info_h, const typename OptimumSelector<Scalar>::CandsInfo& info_p, const typename OptimumSelector<Scalar>::CandsInfo& info_aniso, MeshFunction<Scalar>* rsln, CandElemProjError herr[H2D_MAX_ELEMENT_SONS], CandElemProjError perr, CandElemProjError anisoerr[H2D_MAX_ELEMENT_SONS]);
269 
271 
285  void calc_error_cand_element(const ElementMode2D mode, double3* gip_points, int num_gip_points, const int num_sub, Element** sub_domains, Trf** sub_trfs, int* sons, std::vector<TrfShapeExp>** sub_nonortho_svals, std::vector<TrfShapeExp>** sub_ortho_svals, const typename OptimumSelector<Scalar>::CandsInfo& info, CandElemProjError errors_squared, Scalar* rval[H2D_MAX_ELEMENT_SONS][MAX_NUMBER_FUNCTION_VALUES_FOR_SELECTORS]);
286 
287  protected: //projection
289  struct ElemProj {
295  std::vector<TrfShapeExp>& svals;
297  Scalar* shape_coeffs;
300  };
301 
303 
304  struct ElemGIP {
306  double3* gip_points;
309  };
310 
312  struct ElemSubTrf {
316  double coef_mx;
318  double coef_my;
319  };
320 
324  int inx;
327  };
328 
330 
342  virtual void precalc_ref_solution(int inx_son, MeshFunction<Scalar>* rsln, Element* element, int intr_gip_order, Scalar* rval[H2D_MAX_ELEMENT_SONS][MAX_NUMBER_FUNCTION_VALUES_FOR_SELECTORS]) = 0;
343 
345  virtual void free_ref_solution_data(int inx_son, Scalar* rval[H2D_MAX_ELEMENT_SONS][MAX_NUMBER_FUNCTION_VALUES_FOR_SELECTORS]) = 0;
346 
348 
354  virtual double** build_projection_matrix(double3* gip_points, int num_gip_points, const int* shape_inx, const int num_shapes, ElementMode2D) = 0;
355 
357 
363  virtual Scalar evaluate_rhs_subdomain(Element* sub_elem, const ElemGIP& sub_gip, int son, const ElemSubTrf& sub_trf, const ElemSubShapeFunc& sub_shape, Scalar* rval[H2D_MAX_ELEMENT_SONS][MAX_NUMBER_FUNCTION_VALUES_FOR_SELECTORS]) = 0;
364 
366 
372  virtual double evaluate_error_squared_subdomain(Element* sub_elem, const ElemGIP& sub_gip, int son, const ElemSubTrf& sub_trf, const ElemProj& elem_proj, Scalar* rval[H2D_MAX_ELEMENT_SONS][MAX_NUMBER_FUNCTION_VALUES_FOR_SELECTORS]) = 0;
373  };
374  }
375  }
376 }
377 #endif
double coef_my
A coefficient that scales df/dy for each subdomain. A coefficient represents effects of a transformat...
Definition: adapt.h:24
virtual void precalc_shapes(const double3 *gip_points, const int num_gip_points, const Trf *trfs, const int num_noni_trfs, const std::vector< typename OptimumSelector< Scalar >::ShapeInx > &shapes, const int max_shape_inx, TrfShape &svals, ElementMode2D mode)
Calculates values of shape function at GIP for all transformations.
A projection-based selector for Hcurl space.
Definition: function.h:36
Stores one element of a mesh.
Definition: element.h:107
#define H2DRS_DEFAULT_ERR_WEIGHT_P
A default multiplicative coefficient of an error of a P-candidate.
Definition: global.h:74
#define H2DRS_MAX_ORDER
A maximum order suported by refinement selectors.
Definition: global.h:104
double error_weight_h
A coefficient that multiplies error of H-candidate. The default value is H2DRS_DEFAULT_ERR_WEIGHT_H.
bool is_valid() const
Returns true if value is mared as valid.
Represents a function defined on a mesh.
Definition: mesh_function.h:56
Integration points in the reference domain of an element of a candidate.
A projection-based selector for H1 space.
Definition: function.h:35
double coef_mx
A coefficient that scales df/dx for each subdomain. A coefficient represents effects of a transformat...
double error_weight_aniso
A coefficient that multiplies error of ANISO-candidate. The default value is H2DRS_DEFAULT_ERR_WEIGHT...
bool warn_uniform_orders
True if the selector has already warned about possible inefficiency.
double3 * gip_points
Integration points and weights. The first index is an index of an integration point, the second index is defined through the enum GIP2DIndices.
A general projection-based selector.
Definition: function.h:33
#define H2D_NUM_MODES
Internal.
Definition: global.h:35
TrfShape * cached_shape_ortho_vals
Precalculated valus of orthogonalized shape functions.
CandList
Predefined list of candidates.
Definition: candidates.h:46
double CandElemProjError[H2DRS_MAX_ORDER+2][H2DRS_MAX_ORDER+2]
Error of an element of a candidate for various permutations of orders.
bool * cached_shape_vals_valid
True if shape values were already initialized.
Scalar * shape_coeffs
Coefficients of shape indices of a projection.
2D transformation.
Definition: transformable.h:29
std::vector< TrfShapeExp > & svals
A precalculated shape-function values. Empty is not defined.
TrfShapeExp & svals
Evaluate values of a shape function. If TrfShapeExp::empty(), no precalculated values are available...
Should be exactly the same as is the count of enum ShapesetType.
Definition: shapeset.h:95
TrfShape * cached_shape_vals
Precalculate values of shape functions.
#define H2DRS_DEFAULT_ERR_WEIGHT_H
A default multiplicative coefficient of an error of a H-candidate.
Definition: global.h:73
double error_weight_p
A coefficient that multiplies error of P-candidate. The default value is H2DRS_DEFAULT_ERR_WEIGHT_P.
int max_quad_order
An encoded maximum order of the projection. If triangle, the vertical order is equal to the horizonta...
A transformation from a reference domain of a subdomain to a reference domain of an element of a cand...
A selector that chooses an optimal candidates based on a score.
Definition: function.h:32
#define H2D_MAX_ELEMENT_SONS
Macros.
Definition: global.h:30
#define H2DRS_DEFAULT_ERR_WEIGHT_ANISO
A default multiplicative coefficient of an error of a ANISO-candidate.
Definition: global.h:75
virtual void precalc_ortho_shapes(const double3 *gip_points, const int num_gip_points, const Trf *trfs, const int num_noni_trfs, const std::vector< typename OptimumSelector< Scalar >::ShapeInx > &shapes, const int max_shape_inx, TrfShape &ortho_svals, ElementMode2D mode)
Calculates values of orthogonalized shape function at GIP for all transformations.
A projection-based selector for L2 space.
Definition: function.h:34
void mark(int new_state=H2DRS_VALCACHE_VALID)
Marks a value.