Line data Source code
1 : /*
2 : * Copyright (c) 2013: G-CSC, Goethe University Frankfurt
3 : * Author: Lisa Grau
4 : *
5 : * This file is part of UG4.
6 : *
7 : * UG4 is free software: you can redistribute it and/or modify it under the
8 : * terms of the GNU Lesser General Public License version 3 (as published by the
9 : * Free Software Foundation) with the following additional attribution
10 : * requirements (according to LGPL/GPL v3 §7):
11 : *
12 : * (1) The following notice must be displayed in the Appropriate Legal Notices
13 : * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
14 : *
15 : * (2) The following notice must be displayed at a prominent place in the
16 : * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
17 : *
18 : * (3) The following bibliography is recommended for citation and must be
19 : * preserved in all covered files:
20 : * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
21 : * parallel geometric multigrid solver on hierarchically distributed grids.
22 : * Computing and visualization in science 16, 4 (2013), 151-164"
23 : * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
24 : * flexible software system for simulating pde based models on high performance
25 : * computers. Computing and visualization in science 16, 4 (2013), 165-179"
26 : *
27 : * This program is distributed in the hope that it will be useful,
28 : * but WITHOUT ANY WARRANTY; without even the implied warranty of
29 : * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
30 : * GNU Lesser General Public License for more details.
31 : */
32 :
33 : #include "newton_cotes.h"
34 : #include "../quadrature.h"
35 : #include "common/util/provider.h"
36 :
37 : namespace ug
38 : {
39 :
40 : //constructor
41 0 : NewtonCotes::NewtonCotes(size_t order)
42 : {
43 :
44 : if (order == 0) order = 1;
45 0 : m_order = order;
46 0 : m_numPoints = order+1;
47 0 : position_type* pvPoint= new position_type[m_numPoints];
48 0 : weight_type* pvWeight= new weight_type[m_numPoints];
49 0 : m_pvPoint = pvPoint;
50 0 : m_pvWeight = pvWeight;
51 :
52 0 : switch(m_numPoints)
53 : {
54 : case 0:
55 : case 1:
56 : case 2:
57 0 : pvPoint[0][0] = 0;
58 0 : pvPoint[1][0] = 1.;
59 :
60 0 : pvWeight[0] = 0.5;
61 0 : pvWeight[1] = 0.5;
62 0 : break;
63 :
64 : case 3:
65 0 : pvPoint[0][0] = 0;
66 0 : pvPoint[1][0] = 0.5;
67 0 : pvPoint[2][0] = 1.;
68 :
69 0 : pvWeight[0] = 0.1666666666666666666666666666666666666667;
70 0 : pvWeight[1] = 0.6666666666666666666666666666666666666667;
71 0 : pvWeight[2] = 0.1666666666666666666666666666666666666667;
72 0 : break;
73 :
74 : case 4:
75 0 : pvPoint[0][0] = 0;
76 0 : pvPoint[1][0] = 0.3333333333333333333333333333333333333333;
77 0 : pvPoint[2][0] = 0.6666666666666666666666666666666666666667;
78 0 : pvPoint[3][0] = 1.;
79 :
80 0 : pvWeight[0] = 0.125;
81 0 : pvWeight[1] = 0.375;
82 0 : pvWeight[2] = 0.375;
83 0 : pvWeight[3] = 0.125;
84 0 : break;
85 :
86 : case 5:
87 0 : pvPoint[0][0] = 0;
88 0 : pvPoint[1][0] = 0.25;
89 0 : pvPoint[2][0] = 0.5;
90 0 : pvPoint[3][0] = 0.75;
91 0 : pvPoint[4][0] = 1.;
92 :
93 0 : pvWeight[0] = 0.07777777777777777777777777777777777777778;
94 0 : pvWeight[1] = 0.3555555555555555555555555555555555555556;
95 0 : pvWeight[2] = 0.1333333333333333333333333333333333333333;
96 0 : pvWeight[3] = 0.3555555555555555555555555555555555555556;
97 0 : pvWeight[4] = 0.07777777777777777777777777777777777777778;
98 0 : break;
99 :
100 : case 6:
101 0 : pvPoint[0][0] = 0;
102 0 : pvPoint[1][0] = 0.2;
103 0 : pvPoint[2][0] = 0.4;
104 0 : pvPoint[3][0] = 0.6;
105 0 : pvPoint[4][0] = 0.8;
106 0 : pvPoint[5][0] = 1.;
107 :
108 0 : pvWeight[0] = 0.06597222222222222222222222222222222222222;
109 0 : pvWeight[1] = 0.2604166666666666666666666666666666666667;
110 0 : pvWeight[2] = 0.1736111111111111111111111111111111111111;
111 0 : pvWeight[3] = 0.1736111111111111111111111111111111111111;
112 0 : pvWeight[4] = 0.2604166666666666666666666666666666666667;
113 0 : pvWeight[5] = 0.06597222222222222222222222222222222222222;
114 0 : break;
115 :
116 : case 7:
117 0 : pvPoint[0][0] = 0;
118 0 : pvPoint[1][0] = 0.1666666666666666666666666666666666666667;
119 0 : pvPoint[2][0] = 0.3333333333333333333333333333333333333333;
120 0 : pvPoint[3][0] = 0.5;
121 0 : pvPoint[4][0] = 0.6666666666666666666666666666666666666667;
122 0 : pvPoint[5][0] = 0.8333333333333333333333333333333333333333;
123 0 : pvPoint[6][0] = 1.;
124 :
125 0 : pvWeight[0] = 0.04880952380952380952380952380952380952381;
126 0 : pvWeight[1] = 0.2571428571428571428571428571428571428571;
127 0 : pvWeight[2] = 0.03214285714285714285714285714285714285714;
128 0 : pvWeight[3] = 0.3238095238095238095238095238095238095238;
129 0 : pvWeight[4] = 0.03214285714285714285714285714285714285714;
130 0 : pvWeight[5] = 0.2571428571428571428571428571428571428571;
131 0 : pvWeight[6] = 0.04880952380952380952380952380952380952381;
132 0 : break;
133 :
134 : case 8:
135 0 : pvPoint[0][0] = 0;
136 0 : pvPoint[1][0] = 0.1428571428571428571428571428571428571429;
137 0 : pvPoint[2][0] = 0.2857142857142857142857142857142857142857;
138 0 : pvPoint[3][0] = 0.4285714285714285714285714285714285714286;
139 0 : pvPoint[4][0] = 0.5714285714285714285714285714285714285714;
140 0 : pvPoint[5][0] = 0.7142857142857142857142857142857142857143;
141 0 : pvPoint[6][0] = 0.8571428571428571428571428571428571428571;
142 0 : pvPoint[7][0] = 1.;
143 :
144 0 : pvWeight[0] = 0.04346064814814814814814814814814814814815;
145 0 : pvWeight[1] = 0.2070023148148148148148148148148148148148;
146 0 : pvWeight[2] = 0.0765625;
147 0 : pvWeight[3] = 0.172974537037037037037037037037037037037;
148 0 : pvWeight[4] = 0.172974537037037037037037037037037037037;
149 0 : pvWeight[5] = 0.0765625;
150 0 : pvWeight[6] = 0.2070023148148148148148148148148148148148;
151 0 : pvWeight[7] = 0.04346064814814814814814814814814814814815;
152 0 : break;
153 :
154 : case 9:
155 0 : pvPoint[0][0] = 0;
156 0 : pvPoint[1][0] = 0.125;
157 0 : pvPoint[2][0] = 0.25;
158 0 : pvPoint[3][0] = 0.375;
159 0 : pvPoint[4][0] = 0.5;
160 0 : pvPoint[5][0] = 0.625;
161 0 : pvPoint[6][0] = 0.75;
162 0 : pvPoint[7][0] = 0.875;
163 0 : pvPoint[8][0] = 1.;
164 :
165 0 : pvWeight[0] = 0.03488536155202821869488536155202821869489;
166 0 : pvWeight[1] = 0.2076895943562610229276895943562610229277;
167 0 : pvWeight[2] = -0.03273368606701940035273368606701940035273;
168 0 : pvWeight[3] = 0.3702292768959435626102292768959435626102;
169 0 : pvWeight[4] = -0.1601410934744268077601410934744268077601;
170 0 : pvWeight[5] = 0.3702292768959435626102292768959435626102;
171 0 : pvWeight[6] = -0.03273368606701940035273368606701940035273;
172 0 : pvWeight[7] = 0.2076895943562610229276895943562610229277;
173 0 : pvWeight[8] = 0.03488536155202821869488536155202821869489;
174 0 : break;
175 :
176 : case 10:
177 0 : pvPoint[0][0] = 0;
178 0 : pvPoint[1][0] = 0.1111111111111111111111111111111111111111;
179 0 : pvPoint[2][0] = 0.2222222222222222222222222222222222222222;
180 0 : pvPoint[3][0] = 0.3333333333333333333333333333333333333333;
181 0 : pvPoint[4][0] = 0.4444444444444444444444444444444444444444;
182 0 : pvPoint[5][0] = 0.5555555555555555555555555555555555555556;
183 0 : pvPoint[6][0] = 0.6666666666666666666666666666666666666667;
184 0 : pvPoint[7][0] = 0.7777777777777777777777777777777777777778;
185 0 : pvPoint[8][0] = 0.8888888888888888888888888888888888888889;
186 0 : pvPoint[9][0] = 1.;
187 :
188 0 : pvWeight[0] = 0.03188616071428571428571428571428571428571;
189 0 : pvWeight[1] = 0.1756808035714285714285714285714285714286;
190 0 : pvWeight[2] = 0.01205357142857142857142857142857142857143;
191 0 : pvWeight[3] = 0.2158928571428571428571428571428571428571;
192 0 : pvWeight[4] = 0.06448660714285714285714285714285714285714;
193 0 : pvWeight[5] = 0.06448660714285714285714285714285714285714;
194 0 : pvWeight[6] = 0.2158928571428571428571428571428571428571;
195 0 : pvWeight[7] = 0.01205357142857142857142857142857142857143;
196 0 : pvWeight[8] = 0.1756808035714285714285714285714285714286;
197 0 : pvWeight[9] = 0.03188616071428571428571428571428571428571;
198 0 : break;
199 :
200 : case 11:
201 0 : pvPoint[0][0] = 0;
202 0 : pvPoint[1][0] = 0.1;
203 0 : pvPoint[2][0] = 0.2;
204 0 : pvPoint[3][0] = 0.3;
205 0 : pvPoint[4][0] = 0.4;
206 0 : pvPoint[5][0] = 0.5;
207 0 : pvPoint[6][0] = 0.6;
208 0 : pvPoint[7][0] = 0.7;
209 0 : pvPoint[8][0] = 0.8;
210 0 : pvPoint[9][0] = 0.9;
211 0 : pvPoint[10][0] = 1.;
212 :
213 0 : pvWeight[0] = 0.02683414836192613970391748169525947303725;
214 0 : pvWeight[1] = 0.1775359414248303137192026080914969803859;
215 0 : pvWeight[2] = -0.08104357062690396023729357062690396023729;
216 0 : pvWeight[3] = 0.4549462882796216129549462882796216129549;
217 0 : pvWeight[4] = -0.4351551226551226551226551226551226551227;
218 0 : pvWeight[5] = 0.7137646304312970979637646304312970979638;
219 0 : pvWeight[6] = -0.4351551226551226551226551226551226551227;
220 0 : pvWeight[7] = 0.4549462882796216129549462882796216129549;
221 0 : pvWeight[8] = -0.08104357062690396023729357062690396023729;
222 0 : pvWeight[9] = 0.1775359414248303137192026080914969803859;
223 0 : pvWeight[10] = 0.02683414836192613970391748169525947303725;
224 0 : break;
225 :
226 0 : default:
227 0 : UG_THROW("NewtonCotes: Rule for order " << order << " not supported.");
228 : }
229 0 : }
230 :
231 0 : NewtonCotes::~NewtonCotes()
232 : {
233 0 : delete[] m_pvPoint;
234 0 : delete[] m_pvWeight;
235 0 : };
236 : }
|