Flow123d  JS_before_hm-1626-gde32303
quadrature.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.hh
15  * @brief Basic definitions of numerical quadrature rules.
16  * @author Jan Stebel
17  */
18 
19 #include "quadrature/quadrature.hh"
20 
21 
22 
24 {
25  ASSERT_DBG( dim_ == q.dim_ );
27  weights = q.weights;
28  return *this;
29 }
30 
31 
32 Quadrature::Quadrature(unsigned int dimension, unsigned int n_q)
33 : dim_(dimension),
34  quadrature_points(dimension, 1, n_q),
35  weights(n_q, 0)
36 {}
37 
38 
40  dim_(q.dim_),
42  weights(q.weights)
43 {}
44 
45 
46 template<unsigned int bulk_dim>
47 Quadrature Quadrature::make_from_side(unsigned int sid, unsigned int pid) const
48 {
49  ASSERT_DBG( bulk_dim == dim_ + 1 );
50 
51  // Below we permute point coordinates according to permutation
52  // of nodes on side. We just check that these numbers equal.
54 
55  Quadrature q(dim_+1, size());
56  Armor::ArmaVec<double, bulk_dim+1> el_bar_coords, final_bar;
57 
58  for (unsigned int k=0; k<size(); k++)
59  {
60  //compute barycentric coordinates on element
63 
64  //permute
65  for (unsigned int i=0; i<RefElement<bulk_dim>::n_nodes_per_side; i++) {
66  pp(RefElement<bulk_dim>::side_permutations[pid][i]) = p(i);
67  }
68 
69  el_bar_coords = RefElement<bulk_dim>::template interpolate<bulk_dim-1>(pp, sid);
70 
71  //get local coordinates and set
73  q.weights[k] = weight(k);
74  }
75 
76  return q;
77 }
78 
79 // Specialized subquadrature consructor for dim=1.
80 template<> Quadrature Quadrature::make_from_side<1>(unsigned int sid, unsigned int) const
81 {
82  ASSERT_EQ_DBG(size(), 1);
83  Quadrature q(1, 1);
84  q.quadrature_points.set(0) = Armor::ArmaVec<double, 1>({ (double)sid });
85  q.weight(0) = 1;
86 
87  return q;
88 }
89 
90 template Quadrature Quadrature::make_from_side<2>(unsigned int sid, unsigned int pid) const;
91 template Quadrature Quadrature::make_from_side<3>(unsigned int sid, unsigned int pid) const;
92 
93 
Armor::Array< double > quadrature_points
List of quadrature points.
Definition: quadrature.hh:141
#define ASSERT_EQ_DBG(a, b)
Definition of comparative assert macro (EQual) only for debug mode.
Definition: asserts.hh:332
Quadrature(const Quadrature &q)
Copy constructor.
Definition: quadrature.cc:39
Quadrature make_from_side(unsigned int sid, unsigned int pid) const
Definition: quadrature.cc:47
double weight(unsigned int i) const
Returns the ith weight.
Definition: quadrature.hh:108
Base class for quadrature rules on simplices in arbitrary dimensions.
Definition: quadrature.hh:48
static LocalPoint bary_to_local(const BaryPoint &bp)
Converts from barycentric to local coordinates.
Definition: ref_element.cc:319
std::vector< double > weights
List of weights to the quadrature points.
Definition: quadrature.hh:148
Basic definitions of numerical quadrature rules.
ArrayMatSet set(uint index)
Definition: armor.hh:838
typename arma::Col< Type >::template fixed< nr > ArmaVec
Definition: armor.hh:505
const unsigned int dim_
Dimension of quadrature points.
Definition: quadrature.hh:134
#define ASSERT_DBG(expr)
static BaryPoint local_to_bary(const LocalPoint &lp)
Converts from local to barycentric coordinates.
Definition: ref_element.cc:305
unsigned int size() const
Returns number of quadrature points.
Definition: quadrature.hh:87
Quadrature & operator=(const Quadrature &q)
Definition: quadrature.cc:23