Flow123d  master-1fea4ce
quadrature_lib.cc
Go to the documentation of this file.
1 /*!
2  *
3  * Copyright (C) 2015 Technical University of Liberec. All rights reserved.
4  *
5  * This program is free software; you can redistribute it and/or modify it under
6  * the terms of the GNU General Public License version 3 as published by the
7  * Free Software Foundation. (http://www.gnu.org/licenses/gpl-3.0.en.html)
8  *
9  * This program is distributed in the hope that it will be useful, but WITHOUT
10  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
11  * FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.
12  *
13  *
14  * @file quadrature_lib.cc
15  * @brief Definitions of particular quadrature rules on simplices.
16  * @author Jan Stebel
17  */
18 
19 #include "system/global_defs.h"
20 #include "system/system.hh"
22 
23 
24 #define FLOAT double
25 #define SHORT int
26 #define _F(n) n
27 #include "quad.c"
28 #undef FLOAT
29 #undef SHORT
30 #undef _F
31 
34  {},
35  { QUAD_1D_P1,
43  },
44  { QUAD_2D_P1,
51  },
52  { QUAD_3D_P1,
57  }
58 };
59 
60 double __unit_cell_volume[] = { 1, 1, 0.5, 1./6 };
61 
62 template<uint dimension>
63 void QGauss::init(uint order) {
64  DimQuadList & quads = __gauss_quadratures[dimension];
65  ASSERT_LT(order, quads.size()).error("Quadrature of given order is not implemented.");
66  auto &point_list = quads[order];
67 
68  this->quadrature_points.reinit(point_list->npoints);
69  this->weights.resize(0);
70  for (int i=0; i<point_list->npoints; i++)
71  {
72  Armor::ArmaVec<double, dimension> p(& point_list->points[i*(dimension + 1)]);
73  this->quadrature_points.append(p);
74  this->weights.push_back(point_list->weights[i] * __unit_cell_volume[dimension]);
75  }
76 }
77 
78 
79 QGauss::QGauss(unsigned int dim, unsigned int order)
80 : Quadrature(dim)
81 {
82 
83  switch (dim)
84  {
85  case 0:
86  // Quadrature on 0-dim element have single quadrature point
87  // with 0 local coordinates.
88  // No way to append 0-dim arma vec, we just resize the Array.
89  this->quadrature_points.reinit(1);
90  this->quadrature_points.resize(1);
91  this->weights.push_back(1);
92  ASSERT_EQ(size(), 1);
93  return;
94  case 1:
95  init<1>(order);
96  break;
97  case 2:
98  init<2>(order);
99  break;
100  case 3:
101  init<3>(order);
102  break;
103  }
104 
105 
106 }
107 
108 
109 
110 
#define ASSERT_LT(a, b)
Definition of comparative assert macro (Less Than) only for debug mode.
Definition: asserts.hh:301
#define ASSERT_EQ(a, b)
Definition of comparative assert macro (EQual) only for debug mode.
Definition: asserts.hh:333
void resize(uint size)
Definition: armor.hh:758
void reinit(uint size)
Definition: armor.hh:746
void append(const ArmaMat< Type, nr, nc > &item)
Definition: armor.hh:784
QGauss(unsigned int dim, const unsigned int order)
Create a formula of given order.
void init(uint order)
Base class for quadrature rules on simplices in arbitrary dimensions.
Definition: quadrature.hh:48
unsigned int size() const
Returns number of quadrature points.
Definition: quadrature.hh:86
unsigned int dim() const
Definition: quadrature.hh:72
std::vector< double > weights
List of weights to the quadrature points.
Definition: quadrature.hh:143
Armor::Array< double > quadrature_points
List of quadrature points.
Definition: quadrature.hh:136
Global macros to enhance readability and debugging, general constants.
unsigned int uint
typename arma::Col< Type >::template fixed< nr > ArmaVec
Definition: armor.hh:505
#define QUAD_3D_P12
Definition: quad.h:143
#define QUAD_2D_P19
Definition: quad.h:113
#define QUAD_2D_P15
Definition: quad.h:105
#define QUAD_1D_P2
Definition: quad.h:45
#define QUAD_1D_P20
Definition: quad.h:72
#define QUAD_1D_P14
Definition: quad.h:63
#define QUAD_2D_P4
Definition: quad.h:83
#define QUAD_3D_P4
Definition: quad.h:127
#define QUAD_2D_P7
Definition: quad.h:89
#define QUAD_2D_P12
Definition: quad.h:99
#define QUAD_1D_P3
Definition: quad.h:46
#define QUAD_1D_P16
Definition: quad.h:66
#define QUAD_2D_P5
Definition: quad.h:85
#define QUAD_2D_P11
Definition: quad.h:97
#define QUAD_3D_P14
Definition: quad.h:147
#define QUAD_2D_P16
Definition: quad.h:107
#define QUAD_2D_P18
Definition: quad.h:111
#define QUAD_2D_P10
Definition: quad.h:95
#define QUAD_2D_P20
Definition: quad.h:115
#define QUAD_2D_P6
Definition: quad.h:87
#define QUAD_1D_P5
Definition: quad.h:49
#define QUAD_3D_P3
Definition: quad.h:125
#define QUAD_2D_P21
Definition: quad.h:117
#define QUAD_3D_P8
Definition: quad.h:135
#define QUAD_2D_P2
Definition: quad.h:79
#define QUAD_1D_P18
Definition: quad.h:69
#define QUAD_1D_P21
Definition: quad.h:73
#define QUAD_3D_P6
Definition: quad.h:131
#define QUAD_3D_P2
Definition: quad.h:123
#define QUAD_1D_P4
Definition: quad.h:48
#define QUAD_2D_P9
Definition: quad.h:93
#define QUAD_1D_P12
Definition: quad.h:60
#define QUAD_3D_P5
Definition: quad.h:129
#define QUAD_1D_P13
Definition: quad.h:61
#define QUAD_1D_P10
Definition: quad.h:57
#define QUAD_1D_P19
Definition: quad.h:70
#define QUAD_1D_P7
Definition: quad.h:52
#define QUAD_1D_P8
Definition: quad.h:54
#define QUAD_1D_P17
Definition: quad.h:67
#define QUAD_2D_P13
Definition: quad.h:101
#define QUAD_1D_P15
Definition: quad.h:64
#define QUAD_1D_P9
Definition: quad.h:55
#define QUAD_3D_P11
Definition: quad.h:141
#define QUAD_2D_P8
Definition: quad.h:91
#define QUAD_2D_P1
Definition: quad.h:77
#define QUAD_3D_P9
Definition: quad.h:137
#define QUAD_2D_P3
Definition: quad.h:81
#define QUAD_2D_P17
Definition: quad.h:109
#define QUAD_1D_P6
Definition: quad.h:51
#define QUAD_2D_P14
Definition: quad.h:103
#define QUAD_1D_P1
Definition: quad.h:43
#define QUAD_3D_P1
Definition: quad.h:121
#define QUAD_3D_P13
Definition: quad.h:145
#define QUAD_1D_P11
Definition: quad.h:58
#define QUAD_3D_P10
Definition: quad.h:139
#define QUAD_3D_P7
Definition: quad.h:133
double __unit_cell_volume[]
std::vector< DimQuadList > __gauss_quadratures
std::vector< QUAD * > DimQuadList
Definitions of particular quadrature rules on simplices.