gmsh-TingyuanDoc  0.1
An Open-Source Timing-driven Analytical Mixed-size FPGA Placer
HierarchicalBasisH1Quad.cpp
Go to the documentation of this file.
1 // Gmsh - Copyright (C) 1997-2022 C. Geuzaine, J.-F. Remacle
2 //
3 // See the LICENSE.txt file in the Gmsh root directory for license information.
4 // Please report all issues on https://gitlab.onelab.info/gmsh/gmsh/issues.
5 //
6 // Contributed by Ismail Badia.
7 
8 #include <stdexcept>
9 #include "HierarchicalBasisH1Quad.h"
10 
11 HierarchicalBasisH1Quad::HierarchicalBasisH1Quad(int pf1, int pf2, int pe0,
12  int pe1, int pe2, int pe3)
13 
14 {
15  _nvertex = 4;
16  _nedge = 4;
17  _nfaceQuad = 1;
18  _nfaceTri = 0;
19  _nVertexFunction = 4;
20  _nEdgeFunction = pe0 + pe1 + pe2 + pe3 - 4;
21  _nQuadFaceFunction = (pf1 - 1) * (pf2 - 1);
22  _nTriFaceFunction = 0;
23  _nBubbleFunction = 0;
24  _pf1 = pf1;
25  _pf2 = pf2;
26  if(pe1 > pf2 || pe3 > pf2) {
27  throw std::runtime_error("pe1 and pe3 must be <=pf2");
28  }
29  if(pe0 > pf1 || pe2 > pf1) {
30  throw std::runtime_error("pe0 and pe2 must be <=pf1");
31  }
32  _pOrderEdge[0] = pe0;
33  _pOrderEdge[1] = pe1;
34  _pOrderEdge[2] = pe2;
35  _pOrderEdge[3] = pe3;
36 }
37 HierarchicalBasisH1Quad::HierarchicalBasisH1Quad(int order)
38 
39 {
40  _nvertex = 4;
41  _nedge = 4;
42  _nfaceQuad = 1;
43  _nfaceTri = 0;
44  _nVertexFunction = 4;
45  _nEdgeFunction = 4 * order - 4;
46  _nQuadFaceFunction = (order - 1) * (order - 1);
47  _nTriFaceFunction = 0;
48  _nBubbleFunction = 0;
49  _pf1 = order;
50  _pf2 = order;
51  _pOrderEdge[0] = order;
52  _pOrderEdge[1] = order;
53  _pOrderEdge[2] = order;
54  _pOrderEdge[3] = order;
55 }
56 
57 HierarchicalBasisH1Quad::~HierarchicalBasisH1Quad() {}
58 
59 unsigned int HierarchicalBasisH1Quad::getNumberOfOrientations() const
60 {
61  return 24; // factorial 4
62 }
63 
64 double HierarchicalBasisH1Quad::_affineCoordinate(int const &j, double const &u,
65  double const &v)
66 {
67  switch(j) {
68  case(1): return 0.5 * (1 + u);
69  case(2): return 0.5 * (1 - u);
70  case(3): return 0.5 * (1 + v);
71  case(4): return 0.5 * (1 - v);
72  default: throw std::runtime_error("j must be : 1<=j<=4");
73  }
74 }
75 
76 void HierarchicalBasisH1Quad::generateBasis(double const &u, double const &v,
77  double const &w,
78  std::vector<double> &vertexBasis,
79  std::vector<double> &edgeBasis,
80  std::vector<double> &faceBasis,
81  std::vector<double> &bubbleBasis)
82 {
83  double lambda1 = _affineCoordinate(1, u, v);
84  double lambda2 = _affineCoordinate(2, u, v);
85  double lambda3 = _affineCoordinate(3, u, v);
86  double lambda4 = _affineCoordinate(4, u, v);
87  // vertex shape functions:
88  vertexBasis[0] = lambda2 * lambda4;
89  vertexBasis[1] = lambda1 * lambda4;
90  vertexBasis[2] = lambda1 * lambda3;
91  vertexBasis[3] = lambda2 * lambda3;
92  // edge 0 & 2 shape functions and a part of face functions (all lk(u)):
93  int minpe0pe2 = std::min(_pOrderEdge[0], _pOrderEdge[2]);
94  int const1 = _pOrderEdge[0] + _pOrderEdge[1] - 4;
95  for(int k = 2; k <= minpe0pe2; k++) {
96  double lk = OrthogonalPoly::EvalLobatto(k, u);
97  double phie0 = lambda4 * lk;
98  double phie2 = lambda3 * lk;
99  int const2 = k - 2;
100  edgeBasis[const2] = phie0;
101  edgeBasis[k + const1] = phie2;
102  int const3 = _pf2 - 1;
103  int iterator = 0;
104  while(iterator < const3) {
105  int nbr = const2 * const3 + iterator;
106  faceBasis[nbr] = lk;
107  iterator++;
108  }
109  }
110  if(_pOrderEdge[0] > minpe0pe2) {
111  for(int k = minpe0pe2 + 1; k <= _pOrderEdge[0]; k++) {
112  double lk = OrthogonalPoly::EvalLobatto(k, u);
113  double phie0 = lambda4 * lk;
114  int const2 = k - 2;
115  edgeBasis[const2] = phie0;
116  int const3 = _pf2 - 1;
117  int iterator = 0;
118  while(iterator < const3) {
119  int nbr = const2 * const3 + iterator;
120  faceBasis[nbr] = lk;
121  iterator++;
122  }
123  }
124  }
125  else {
126  for(int k = minpe0pe2 + 1; k <= _pOrderEdge[2]; k++) {
127  double lk = OrthogonalPoly::EvalLobatto(k, u);
128  double phie2 = lambda3 * lk;
129  edgeBasis[k + const1] = phie2;
130  int iterator = 0;
131  int const3 = _pf2 - 1;
132  while(iterator < const3) {
133  int nbr = (k - 2) * const3 + iterator;
134  faceBasis[nbr] = lk;
135  iterator++;
136  }
137  }
138  }
139  int maxpe0pe2 = std::max(_pOrderEdge[0], _pOrderEdge[2]);
140  for(int k = maxpe0pe2 + 1; k <= _pf1; k++) {
141  double lk = OrthogonalPoly::EvalLobatto(k, u);
142  int iterator = 0;
143  int const3 = _pf2 - 1;
144  while(iterator < const3) {
145  int nbr = (k - 2) * const3 + iterator;
146  faceBasis[nbr] = lk;
147  iterator++;
148  }
149  }
150  // edge 1 & 3 shape functions and a part of face functions (all lk(v)) :
151  int minpe1pe3 = std::min(_pOrderEdge[1], _pOrderEdge[3]);
152  int const3 = _pOrderEdge[0] - 3;
153  int const4 = _pOrderEdge[0] + _pOrderEdge[1] + _pOrderEdge[2] - 5;
154  for(int k = 2; k <= minpe1pe3; k++) {
155  double lk = OrthogonalPoly::EvalLobatto(k, v);
156  double phie3 = lambda2 * lk;
157  double phie1 = lambda1 * lk;
158  edgeBasis[k + const3] = phie1;
159  edgeBasis[k + const4] = phie3;
160  int s = k - 2;
161  int iterator = 1;
162  while(iterator <= _pf1 - 1) {
163  faceBasis[s] = faceBasis[s] * lk;
164  s = s + _pf2 - 1;
165  iterator++;
166  }
167  }
168  if(_pOrderEdge[1] > minpe1pe3) {
169  for(int k = minpe1pe3 + 1; k <= _pOrderEdge[1]; k++) {
170  double lk = OrthogonalPoly::EvalLobatto(k, v);
171  double phie1 = lambda1 * lk;
172  edgeBasis[k + const3] = phie1;
173  int s = k - 2;
174  int iterator = 1;
175  while(iterator <= _pf1 - 1) {
176  faceBasis[s] = faceBasis[s] * lk;
177  s = s + _pf2 - 1;
178  iterator++;
179  }
180  }
181  }
182  else {
183  for(int k = minpe1pe3 + 1; k <= _pOrderEdge[3]; k++) {
184  double lk = OrthogonalPoly::EvalLobatto(k, v);
185  double phie3 = lambda2 * lk;
186  edgeBasis[k + const4] = phie3;
187  int s = k - 2;
188  int iterator = 1;
189  while(iterator <= _pf1 - 1) {
190  faceBasis[s] = faceBasis[s] * lk;
191  s = s + _pf2 - 1;
192  iterator++;
193  }
194  }
195  }
196  int maxpe1pe3 = std::max(_pOrderEdge[1], _pOrderEdge[3]);
197  for(int k = maxpe1pe3 + 1; k <= _pf2; k++) {
198  double lk = OrthogonalPoly::EvalLobatto(k, v);
199  int s = k - 2;
200  int iterator = 1;
201  while(iterator <= _pf1 - 1) {
202  faceBasis[s] = faceBasis[s] * lk;
203  s = s + _pf2 - 1;
204  iterator++;
205  }
206  }
207 }
208 
209 void HierarchicalBasisH1Quad::generateGradientBasis(
210  double const &u, double const &v, double const &w,
211  std::vector<std::vector<double> > &gradientVertex,
212  std::vector<std::vector<double> > &gradientEdge,
213  std::vector<std::vector<double> > &gradientFace,
214  std::vector<std::vector<double> > &gradientBubble)
215 {
216  double dlambda1 = 0.5;
217  double dlambda2 = -0.5;
218  double dlambda3 = 0.5;
219  double dlambda4 = -0.5;
220  double lambda1 = _affineCoordinate(1, u, v);
221  double lambda2 = _affineCoordinate(2, u, v);
222  double lambda3 = _affineCoordinate(3, u, v);
223  double lambda4 = _affineCoordinate(4, u, v);
224  // vertex gradient functions:
225  gradientVertex[0][0] = dlambda2 * lambda4;
226  gradientVertex[0][1] = lambda2 * dlambda4;
227  gradientVertex[1][0] = dlambda1 * lambda4;
228  gradientVertex[1][1] = lambda1 * dlambda4;
229  gradientVertex[2][0] = dlambda1 * lambda3;
230  gradientVertex[2][1] = lambda1 * dlambda3;
231  gradientVertex[3][0] = dlambda2 * lambda3;
232  gradientVertex[3][1] = lambda2 * dlambda3;
233  // edge 0 & 2 gradient and a part of face gradient ( in the direction u
234  // ):
235  int const1 = _pOrderEdge[0] + _pOrderEdge[1] - 4;
236  int minpe0pe2 = std::min(_pOrderEdge[0], _pOrderEdge[2]);
237  for(int k = 2; k <= minpe0pe2; k++) {
238  double lk = OrthogonalPoly::EvalLobatto(k, u);
239  double dlk = OrthogonalPoly::EvalDLobatto(k, u);
240  double dphie0U = lambda4 * dlk;
241  double dphie0V = dlambda4 * lk;
242  double dphie2U = lambda3 * dlk;
243  double dphie2V = dlambda3 * lk;
244  int const2 = k - 2;
245  int const3 = _pf2 - 1;
246  gradientEdge[const2][0] = dphie0U;
247  gradientEdge[const2][1] = dphie0V;
248  gradientEdge[k + const1][0] = dphie2U;
249  gradientEdge[k + const1][1] = dphie2V;
250  int iterator = 0;
251  while(iterator < const3) {
252  int nbr = const2 * const3 + iterator;
253  gradientFace[nbr][0] = dlk;
254  gradientFace[nbr][1] = lk;
255  iterator++;
256  }
257  }
258  if(_pOrderEdge[0] > minpe0pe2) {
259  for(int k = minpe0pe2 + 1; k <= _pOrderEdge[0]; k++) {
260  double lk = OrthogonalPoly::EvalLobatto(k, u);
261  double dlk = OrthogonalPoly::EvalDLobatto(k, u);
262  double dphie0U = lambda4 * dlk;
263  double dphie0V = dlambda4 * lk;
264  int const2 = k - 2;
265  int const3 = _pf2 - 1;
266  gradientEdge[const2][0] = dphie0U;
267  gradientEdge[const2][1] = dphie0V;
268  int iterator = 0;
269  while(iterator < const3) {
270  int nbr = const2 * const3 + iterator;
271  gradientFace[nbr][0] = dlk;
272  gradientFace[nbr][1] = lk;
273  iterator++;
274  }
275  }
276  }
277  else {
278  for(int k = minpe0pe2 + 1; k <= _pOrderEdge[2]; k++) {
279  double lk = OrthogonalPoly::EvalLobatto(k, u);
280  double dlk = OrthogonalPoly::EvalDLobatto(k, u);
281  double dphie2U = lambda3 * dlk;
282  double dphie2V = dlambda3 * lk;
283  gradientEdge[k + const1][0] = dphie2U;
284  gradientEdge[k + const1][1] = dphie2V;
285  int const3 = _pf2 - 1;
286  int iterator = 0;
287  while(iterator < const3) {
288  int nbr = (k - 2) * (const3) + iterator;
289  gradientFace[nbr][0] = dlk;
290  gradientFace[nbr][1] = lk;
291  iterator++;
292  }
293  }
294  }
295  int maxpe0pe2 = std::max(_pOrderEdge[0], _pOrderEdge[2]);
296  for(int k = maxpe0pe2 + 1; k <= _pf1; k++) {
297  double lk = OrthogonalPoly::EvalLobatto(k, u);
298  double dlk = OrthogonalPoly::EvalDLobatto(k, u);
299  int iterator = 1;
300  while(iterator < _pf2) {
301  int nbr = (k - 2) * (_pf2 - 1) + iterator;
302  gradientFace[nbr][0] = dlk;
303  gradientFace[nbr][1] = lk;
304  iterator++;
305  }
306  }
307 
308  // edge 1 & 3 shape functions and a part of face functions (all lk(v)) :
309  int minpe1pe3 = std::min(_pOrderEdge[1], _pOrderEdge[3]);
310  int const3 = _pOrderEdge[0] - 3;
311  int const4 = _pOrderEdge[0] + _pOrderEdge[1] + _pOrderEdge[2] - 5;
312  for(int k = 2; k <= minpe1pe3; k++) {
313  double lk = OrthogonalPoly::EvalLobatto(k, v);
314  double dlk = OrthogonalPoly::EvalDLobatto(k, v);
315  double dphie3U = dlambda2 * lk;
316  double dphie3V = lambda2 * dlk;
317  double dphie1U = dlambda1 * lk;
318  double dphie1V = lambda1 * dlk;
319  gradientEdge[k + const3][0] = dphie1U;
320  gradientEdge[k + const3][1] = dphie1V;
321  gradientEdge[k + const4][0] = dphie3U;
322  gradientEdge[k + const4][1] = dphie3V;
323  int s = k - 2;
324  int iterator = 1;
325  while(iterator <= _pf1 - 1) {
326  gradientFace[s][0] = gradientFace[s][0] * lk;
327  gradientFace[s][1] = gradientFace[s][1] * dlk;
328  s = s + _pf2 - 1;
329  iterator++;
330  }
331  }
332  if(_pOrderEdge[1] > minpe1pe3) {
333  for(int k = minpe1pe3 + 1; k <= _pOrderEdge[1]; k++) {
334  double lk = OrthogonalPoly::EvalLobatto(k, v);
335  double dlk = OrthogonalPoly::EvalDLobatto(k, v);
336  double dphie1U = dlambda1 * lk;
337  double dphie1V = lambda1 * dlk;
338  gradientEdge[k + const3][0] = dphie1U;
339  gradientEdge[k + const3][1] = dphie1V;
340  int s = k - 2;
341  int iterator = 1;
342  while(iterator <= _pf1 - 1) {
343  gradientFace[s][0] = gradientFace[s][0] * lk;
344  gradientFace[s][1] = gradientFace[s][1] * dlk;
345  s = s + _pf2 - 1;
346  iterator++;
347  }
348  }
349  }
350  else {
351  for(int k = minpe1pe3 + 1; k <= _pOrderEdge[3]; k++) {
352  double lk = OrthogonalPoly::EvalLobatto(k, v);
353  double dlk = OrthogonalPoly::EvalDLobatto(k, v);
354  double dphie3U = dlambda2 * lk;
355  double dphie3V = lambda2 * dlk;
356  gradientEdge[k + const4][0] = dphie3U;
357  gradientEdge[k + const4][1] = dphie3V;
358  int s = k - 2;
359  int iterator = 1;
360  while(iterator <= _pf1 - 1) {
361  gradientFace[s][0] = gradientFace[s][0] * lk;
362  gradientFace[s][1] = gradientFace[s][1] * dlk;
363  s = s + _pf2 - 1;
364  iterator++;
365  }
366  }
367  }
368  int maxpe1pe3 = std::max(_pOrderEdge[1], _pOrderEdge[3]);
369  for(int k = maxpe1pe3 + 1; k <= _pf2; k++) {
370  double lk = OrthogonalPoly::EvalLobatto(k, v);
371  double dlk = OrthogonalPoly::EvalDLobatto(k, v);
372  int s = k - 2;
373  int iterator = 1;
374  while(iterator <= _pf1 - 1) {
375  gradientFace[s][0] = gradientFace[s][0] * lk;
376  gradientFace[s][1] = gradientFace[s][1] * dlk;
377  s = s + _pf2 - 1;
378  iterator++;
379  }
380  }
381 }
382 
383 void HierarchicalBasisH1Quad::orientEdge(
384  int const &flagOrientation, int const &edgeNumber,
385  std::vector<double> &edgeFunctions,
386  const std::vector<double> &eTablePositiveFlag,
387  const std::vector<double> &eTableNegativeFlag)
388 {
389  if(flagOrientation == -1) {
390  int constant1 = 0;
391  int constant2 = 0;
392  for(int i = 0; i <= edgeNumber; i++) { constant2 += _pOrderEdge[i] - 1; }
393  constant2 = constant2 - 1;
394  constant1 = constant2 - _pOrderEdge[edgeNumber] + 2;
395  for(int k = constant1; k <= constant2; k++) {
396  edgeFunctions[k] = eTableNegativeFlag[k];
397  }
398  }
399  else {
400  int constant1 = 0;
401  int constant2 = 0;
402  for(int i = 0; i <= edgeNumber; i++) { constant2 += _pOrderEdge[i] - 1; }
403  constant2 = constant2 - 1;
404  constant1 = constant2 - _pOrderEdge[edgeNumber] + 2;
405  for(int k = constant1; k <= constant2; k++) {
406  edgeFunctions[k] = eTablePositiveFlag[k];
407  }
408  }
409 }
410 
411 void HierarchicalBasisH1Quad::orientEdge(
412  int const &flagOrientation, int const &edgeNumber,
413  std::vector<std::vector<double> > &edgeFunctions,
414  const std::vector<std::vector<double> > &eTablePositiveFlag,
415  const std::vector<std::vector<double> > &eTableNegativeFlag)
416 {
417  if(flagOrientation == -1) {
418  int constant1 = 0;
419  int constant2 = 0;
420  for(int i = 0; i <= edgeNumber; i++) { constant2 += _pOrderEdge[i] - 1; }
421  constant2 = constant2 - 1;
422  constant1 = constant2 - _pOrderEdge[edgeNumber] + 2;
423  for(int k = constant1; k <= constant2; k++) {
424  edgeFunctions[k][0] = eTableNegativeFlag[k][0];
425  edgeFunctions[k][1] = eTableNegativeFlag[k][1];
426  }
427  }
428  else {
429  int constant1 = 0;
430  int constant2 = 0;
431  for(int i = 0; i <= edgeNumber; i++) { constant2 += _pOrderEdge[i] - 1; }
432  constant2 = constant2 - 1;
433  constant1 = constant2 - _pOrderEdge[edgeNumber] + 2;
434  for(int k = constant1; k <= constant2; k++) {
435  edgeFunctions[k][0] = eTablePositiveFlag[k][0];
436  edgeFunctions[k][1] = eTablePositiveFlag[k][1];
437  }
438  }
439 }
440 
441 void HierarchicalBasisH1Quad::orientEdgeFunctionsForNegativeFlag(
442  std::vector<double> &edgeFunctions)
443 {
444  int constant1 = 0;
445  int constant2 = 0;
446  for(int edgeNumber = 0; edgeNumber < _nedge; edgeNumber++) {
447  constant2 = 0;
448  constant2 = 0;
449  for(int i = 0; i <= edgeNumber; i++) { constant2 += _pOrderEdge[i] - 1; }
450  constant2 = constant2 - 1;
451  constant1 = constant2 - _pOrderEdge[edgeNumber] + 2;
452  for(int k = constant1; k <= constant2; k++) {
453  if((k - constant1) % 2 != 0) {
454  edgeFunctions[k] = edgeFunctions[k] * (-1);
455  }
456  }
457  }
458 }
459 
460 void HierarchicalBasisH1Quad::orientEdgeFunctionsForNegativeFlag(
461  std::vector<std::vector<double> > &edgeFunctions)
462 {
463  int constant1 = 0;
464  int constant2 = 0;
465  for(int edgeNumber = 0; edgeNumber < _nedge; edgeNumber++) {
466  constant2 = 0;
467  constant2 = 0;
468  for(int i = 0; i <= edgeNumber; i++) { constant2 += _pOrderEdge[i] - 1; }
469  constant2 = constant2 - 1;
470  constant1 = constant2 - _pOrderEdge[edgeNumber] + 2;
471  for(int k = constant1; k <= constant2; k++) {
472  if((k - constant1) % 2 != 0) {
473  edgeFunctions[k][0] = edgeFunctions[k][0] * (-1);
474  edgeFunctions[k][1] = edgeFunctions[k][1] * (-1);
475  edgeFunctions[k][2] = edgeFunctions[k][2] * (-1);
476  }
477  }
478  }
479 }
480 void HierarchicalBasisH1Quad::orientOneFace(double const &u, double const &v,
481  double const &w, int const &flag1,
482  int const &flag2, int const &flag3,
483  int const &faceNumber,
484  std::vector<double> &faceBasis)
485 {
486  if(!(flag1 == 1 && flag2 == 1 && flag3 == 1)) {
487  int iterator = 0;
488  if(flag3 == 1) {
489  for(int it1 = 2; it1 <= _pf1; it1++) {
490  for(int it2 = 2; it2 <= _pf2; it2++) {
491  int impactFlag1 = 1;
492  int impactFlag2 = 1;
493  if(flag1 == -1 && it1 % 2 != 0) { impactFlag1 = -1; }
494  if(flag2 == -1 && it2 % 2 != 0) { impactFlag2 = -1; }
495  faceBasis[iterator] = faceBasis[iterator] * impactFlag1 * impactFlag2;
496  iterator++;
497  }
498  }
499  }
500  else {
501  std::vector<double> lkVector1(_pf2 - 1);
502  std::vector<double> lkVector2(_pf1 - 1);
503  for(int it = 2; it <= _pf1; it++) {
504  lkVector1[it - 2] = OrthogonalPoly::EvalLobatto(it, u);
505  }
506  for(int it = 2; it <= _pf2; it++) {
507  lkVector2[it - 2] = OrthogonalPoly::EvalLobatto(it, v);
508  }
509 
510  for(int it1 = 2; it1 <= _pf2; it1++) {
511  for(int it2 = 2; it2 <= _pf1; it2++) {
512  int impactFlag1 = 1;
513  int impactFlag2 = 1;
514  if(flag2 == -1 && it1 % 2 != 0) { impactFlag1 = -1; }
515  if(flag1 == -1 && it2 % 2 != 0) { impactFlag2 = -1; }
516  faceBasis[iterator] =
517  lkVector1[it2 - 2] * lkVector2[it1 - 2] * impactFlag1 * impactFlag2;
518  iterator++;
519  }
520  }
521  }
522  }
523 }
524 
525 void HierarchicalBasisH1Quad::orientOneFace(
526  double const &u, double const &v, double const &w, int const &flag1,
527  int const &flag2, int const &flag3, int const &faceNumber,
528  std::vector<std::vector<double> > &gradientFace, std::string typeFunction)
529 {
530  if(!(flag1 == 1 && flag2 == 1 && flag3 == 1)) {
531  int iterator = 0;
532  if(flag3 == 1) {
533  for(int it1 = 2; it1 <= _pf1; it1++) {
534  for(int it2 = 2; it2 <= _pf2; it2++) {
535  int impactFlag1 = 1;
536  int impactFlag2 = 1;
537  if(flag1 == -1 && it1 % 2 != 0) { impactFlag1 = -1; }
538  if(flag2 == -1 && it2 % 2 != 0) { impactFlag2 = -1; }
539  gradientFace[iterator][0] =
540  gradientFace[iterator][0] * impactFlag1 * impactFlag2;
541  gradientFace[iterator][1] =
542  gradientFace[iterator][1] * impactFlag1 * impactFlag2;
543 
544  iterator++;
545  }
546  }
547  }
548  else {
549  std::vector<double> lkVector1(_pf1 - 1);
550  std::vector<double> lkVector2(_pf2 - 1);
551  std::vector<double> dlkVector1(_pf1 - 1);
552  std::vector<double> dlkVector2(_pf2 - 1);
553  for(int it = 2; it <= _pf1; it++) {
554  lkVector1[it - 2] = OrthogonalPoly::EvalLobatto(it, u);
555  dlkVector1[it - 2] = OrthogonalPoly::EvalDLobatto(it, u);
556  }
557  for(int it = 2; it <= _pf2; it++) {
558  lkVector2[it - 2] = OrthogonalPoly::EvalLobatto(it, v);
559  dlkVector2[it - 2] = OrthogonalPoly::EvalDLobatto(it, v);
560  }
561  for(int it1 = 2; it1 <= _pf2; it1++) {
562  for(int it2 = 2; it2 <= _pf1; it2++) {
563  int impactFlag1 = 1;
564  int impactFlag2 = 1;
565  if(flag2 == -1 && it1 % 2 != 0) { impactFlag1 = -1; }
566  if(flag1 == -1 && it2 % 2 != 0) { impactFlag2 = -1; }
567 
568  gradientFace[iterator][0] = dlkVector1[it2 - 2] * lkVector2[it1 - 2] *
569  impactFlag1 * impactFlag2;
570  gradientFace[iterator][1] = lkVector1[it2 - 2] * dlkVector2[it1 - 2] *
571  impactFlag1 * impactFlag2;
572  iterator++;
573  }
574  }
575  }
576  }
577 }
578 void HierarchicalBasisH1Quad::orientFace(
579  int const &flag1, int const &flag2, int const &flag3, int const &faceNumber,
580  const std::vector<double> &quadFaceFunctionsAllOrientation,
581  const std::vector<double> &triFaceFunctionsAllOrientation,
582  std::vector<double> &fTableCopy)
583 {
584  int iOrientation = numberOrientationQuadFace(flag1, flag2, flag3);
585  int offset = iOrientation * _nQuadFaceFunction;
586  for(int i = 0; i < _nQuadFaceFunction; i++) {
587  fTableCopy[i] = quadFaceFunctionsAllOrientation[i + offset];
588  }
589 }
590 
591 void HierarchicalBasisH1Quad::orientFace(
592  int const &flag1, int const &flag2, int const &flag3, int const &faceNumber,
593  const std::vector<std::vector<double> > &quadFaceFunctionsAllOrientation,
594  const std::vector<std::vector<double> > &triFaceFunctionsAllOrientation,
595  std::vector<std::vector<double> > &fTableCopy)
596 {
597  int iOrientation = numberOrientationQuadFace(flag1, flag2, flag3);
598  int offset = iOrientation * _nQuadFaceFunction;
599  for(int i = 0; i < _nQuadFaceFunction; i++) {
600  fTableCopy[i][0] = quadFaceFunctionsAllOrientation[i + offset][0];
601  fTableCopy[i][1] = quadFaceFunctionsAllOrientation[i + offset][1];
602  fTableCopy[i][2] = quadFaceFunctionsAllOrientation[i + offset][2];
603  }
604 }
605 
606 void HierarchicalBasisH1Quad::getKeysInfo(std::vector<int> &functionTypeInfo,
607  std::vector<int> &orderInfo)
608 {
609  functionTypeInfo[0] = 0;
610  functionTypeInfo[1] = 0;
611  functionTypeInfo[2] = 0;
612  functionTypeInfo[3] = 0;
613  orderInfo[0] = 1;
614  orderInfo[1] = 1;
615  orderInfo[2] = 1;
616  orderInfo[3] = 1;
617  int it = 4;
618  for(int numEdge = 0; numEdge < 4; numEdge++) {
619  for(int i = 2; i <= _pOrderEdge[numEdge]; i++) {
620  functionTypeInfo[it] = 1;
621  orderInfo[it] = i;
622  it++;
623  }
624  }
625  for(int n1 = 2; n1 <= _pf1; n1++) {
626  for(int n2 = 2; n2 <= _pf2; n2++) {
627  functionTypeInfo[it] = 2;
628  orderInfo[it] = std::max(n1, n2);
629  it++;
630  }
631  }
632 }
HierarchicalBasisH1Quad::_pf2
int _pf2
Definition: HierarchicalBasisH1Quad.h:92
HierarchicalBasisH1Quad::generateGradientBasis
void generateGradientBasis(double const &u, double const &v, double const &w, std::vector< std::vector< double > > &gradientVertex, std::vector< std::vector< double > > &gradientEdge, std::vector< std::vector< double > > &gradientFace, std::vector< std::vector< double > > &gradientBubble)
Definition: HierarchicalBasisH1Quad.cpp:209
HierarchicalBasis::_nTriFaceFunction
int _nTriFaceFunction
Definition: HierarchicalBasis.h:27
HierarchicalBasisH1Quad::orientOneFace
virtual void orientOneFace(double const &u, double const &v, double const &w, int const &flag1, int const &flag2, int const &flag3, int const &faceNumber, std::vector< double > &faceBasis)
Definition: HierarchicalBasisH1Quad.cpp:480
HierarchicalBasisH1Quad::_pOrderEdge
int _pOrderEdge[4]
Definition: HierarchicalBasisH1Quad.h:93
HierarchicalBasisH1Quad::orientEdgeFunctionsForNegativeFlag
virtual void orientEdgeFunctionsForNegativeFlag(std::vector< double > &edgeFunctions)
Definition: HierarchicalBasisH1Quad.cpp:441
HierarchicalBasisH1Quad::HierarchicalBasisH1Quad
HierarchicalBasisH1Quad(int pf1, int pf2, int pe0, int pe1, int pe2, int pe3)
Definition: HierarchicalBasisH1Quad.cpp:11
HierarchicalBasis::_nvertex
int _nvertex
Definition: HierarchicalBasis.h:20
HierarchicalBasis::_nBubbleFunction
int _nBubbleFunction
Definition: HierarchicalBasis.h:28
OrthogonalPoly::EvalLobatto
double EvalLobatto(int order, double x)
Definition: OrthogonalPoly.cpp:9
OrthogonalPoly::EvalDLobatto
double EvalDLobatto(int order, double x)
Definition: OrthogonalPoly.cpp:127
pe2
const double pe2
Definition: GaussQuadratureQuad.cpp:76
HierarchicalBasisH1Quad::getNumberOfOrientations
virtual unsigned int getNumberOfOrientations() const
Definition: HierarchicalBasisH1Quad.cpp:59
HierarchicalBasisH1Quad::~HierarchicalBasisH1Quad
virtual ~HierarchicalBasisH1Quad()
Definition: HierarchicalBasisH1Quad.cpp:57
HierarchicalBasisH1Quad::generateBasis
virtual void generateBasis(double const &u, double const &v, double const &w, std::vector< double > &vertexBasis, std::vector< double > &edgeBasis, std::vector< double > &faceBasis, std::vector< double > &bubbleBasis)
Definition: HierarchicalBasisH1Quad.cpp:76
HierarchicalBasisH1Quad::_affineCoordinate
static double _affineCoordinate(int const &j, double const &u, double const &v)
Definition: HierarchicalBasisH1Quad.cpp:64
pf2
const double pf2
Definition: GaussQuadratureQuad.cpp:77
HierarchicalBasisH1Quad::getKeysInfo
virtual void getKeysInfo(std::vector< int > &functionTypeInfo, std::vector< int > &orderInfo)
Definition: HierarchicalBasisH1Quad.cpp:606
HierarchicalBasisH1Quad.h
HierarchicalBasis::_nQuadFaceFunction
int _nQuadFaceFunction
Definition: HierarchicalBasis.h:26
HierarchicalBasis::_nfaceTri
int _nfaceTri
Definition: HierarchicalBasis.h:23
HierarchicalBasis::_nfaceQuad
int _nfaceQuad
Definition: HierarchicalBasis.h:22
HierarchicalBasis::_nVertexFunction
int _nVertexFunction
Definition: HierarchicalBasis.h:24
HierarchicalBasisH1Quad::orientFace
virtual void orientFace(int const &flag1, int const &flag2, int const &flag3, int const &faceNumber, const std::vector< double > &quadFaceFunctionsAllOrientation, const std::vector< double > &triFaceFunctionsAllOrientation, std::vector< double > &fTableCopy)
Definition: HierarchicalBasisH1Quad.cpp:578
HierarchicalBasisH1Quad::_pf1
int _pf1
Definition: HierarchicalBasisH1Quad.h:91
HierarchicalBasisH1Quad::orientEdge
virtual void orientEdge(int const &flagOrientation, int const &edgeNumber, std::vector< double > &edgeFunctions, const std::vector< double > &eTablePositiveFlag, const std::vector< double > &eTableNegativeFlag)
Definition: HierarchicalBasisH1Quad.cpp:383
HierarchicalBasis::_nEdgeFunction
int _nEdgeFunction
Definition: HierarchicalBasis.h:25
HierarchicalBasis::_nedge
int _nedge
Definition: HierarchicalBasis.h:21
HierarchicalBasis::numberOrientationQuadFace
int numberOrientationQuadFace(int const &flag1, int const &flag2, int const &flag3)
Definition: HierarchicalBasis.h:105