v0.8.23
Public Member Functions | Public Attributes | List of all members
HcurlFaceBase Struct Reference
Collaboration diagram for HcurlFaceBase:
[legend]

Public Member Functions

MoFEMErrorCode operator() (int shift, int p, int nb_integration_pts, int n0f0_idx, int n1f0_idx, int n2f0_idx, double n[], FTensor::Tensor1< double, 3 > t_grad_n[], FTensor::Tensor1< double *, 3 > &t_phi, FTensor::Tensor2< double *, 3, 3 > &t_diff_phi)
 
MoFEMErrorCode operator() (int shift, int p, int nb_integration_pts, int n0f0_idx, int n1f0_idx, int n2f0_idx, int n0f1_idx, int n1f1_idx, int n2f1_idx, double n[], FTensor::Tensor1< double, 3 > t_grad_n[], FTensor::Tensor1< double *, 3 > &t_phi, FTensor::Tensor2< double *, 3, 3 > &t_diff_phi)
 
MoFEMErrorCode operator() (int shift, int p, int nb_integration_pts, int n0f0_idx, int n1f0_idx, int n2f0_idx, int n0f1_idx, int n1f1_idx, int n2f1_idx, double n[], FTensor::Tensor1< double, 3 > t_grad_n[], FTensor::Tensor1< double *, 3 > &t_phi, FTensor::Tensor2< FTensor::PackPtr< double *, 6 >, 3, 2 > &t_diff_phi)
 

Public Attributes

HcurlEdgeBase hCurlBaseOnEdge
 
VectorDouble f0PhiII
 
VectorDouble diffF0PhiII
 
VectorDouble f1PhiII
 
VectorDouble diffF1PhiII
 
VectorDouble iFiF0
 
VectorDouble diffIFiF0
 
VectorDouble iFiF1
 
VectorDouble diffIFiF1
 
FTensor::Index< 'i', 3 > i
 

Detailed Description

Definition at line 2312 of file Hcurl.cpp.

Member Function Documentation

◆ operator()() [1/3]

MoFEMErrorCode HcurlFaceBase::operator() ( int  shift,
int  p,
int  nb_integration_pts,
int  n0f0_idx,
int  n1f0_idx,
int  n2f0_idx,
double  n[],
FTensor::Tensor1< double, 3 >  t_grad_n[],
FTensor::Tensor1< double *, 3 > &  t_phi,
FTensor::Tensor2< double *, 3, 3 > &  t_diff_phi 
)

Definition at line 2321 of file Hcurl.cpp.

2325  {
2326 
2328 
2330  f0PhiII.resize(3 * NBEDGE_DEMKOWICZ_HCURL(p) * nb_integration_pts,false);
2331  diffF0PhiII.resize(9 * NBEDGE_DEMKOWICZ_HCURL(p) * nb_integration_pts,false);
2332 
2333  // edge base for family I
2334  double *f0_phi_ii = &*f0PhiII.data().begin();
2335  double *diff_f0_phi_ii = &*diffF0PhiII.data().begin();
2336  FTensor::Tensor1<double *, 3> t_f0_phi_ii(
2337  &f0_phi_ii[HVEC0], &f0_phi_ii[HVEC1], &f0_phi_ii[HVEC2], 3);
2338  FTensor::Tensor2<double *, 3, 3> t_diff_f0_phi_ii(
2339  &diff_f0_phi_ii[HVEC0_0], &diff_f0_phi_ii[HVEC0_1],
2340  &diff_f0_phi_ii[HVEC0_2], &diff_f0_phi_ii[HVEC1_0],
2341  &diff_f0_phi_ii[HVEC1_1], &diff_f0_phi_ii[HVEC1_2],
2342  &diff_f0_phi_ii[HVEC2_0], &diff_f0_phi_ii[HVEC2_1],
2343  &diff_f0_phi_ii[HVEC2_2], 9);
2344  CHKERR hCurlBaseOnEdge(4, p-1, nb_integration_pts, n0f0_idx, n1f0_idx, n,
2345  t_grad_n, t_f0_phi_ii, t_diff_f0_phi_ii);
2346 
2347  FTensor::Tensor1<double, 3> &t_grad_n0f0 = t_grad_n[n0f0_idx];
2348  FTensor::Tensor1<double, 3> &t_grad_n1f0 = t_grad_n[n1f0_idx];
2349  FTensor::Tensor1<double, 3> &t_grad_n2f0 = t_grad_n[n2f0_idx];
2350  FTensor::Tensor1<double, 3> t_grad_n0f0_p_n1f0;
2351  t_grad_n0f0_p_n1f0(i) = t_grad_n0f0(i) + t_grad_n1f0(i) + t_grad_n2f0(i);
2352 
2353  iFiF0.resize(p + 1, false);
2354  diffIFiF0.resize(3 * p + 3, false);
2355  double *ifif0 = &*iFiF0.data().begin();
2356  double *diff_ifif0 = &*diffIFiF0.data().begin();
2357 
2358  for (int gg = 0; gg != nb_integration_pts; ++gg) {
2359 
2360  const int shift_n = shift * gg;
2361  const double n0f0 = n[shift_n + n0f0_idx];
2362  const double n1f0 = n[shift_n + n1f0_idx];
2363  const double n2f0 = n[shift_n + n2f0_idx];
2364 
2365  int phi_shift = 3 * NBEDGE_DEMKOWICZ_HCURL(p - 1) * gg;
2366  int diff_phi_shift = 9 * NBEDGE_DEMKOWICZ_HCURL(p - 1) * gg;
2367 
2368  for (int oo = 2; oo <= p; ++oo) {
2369 
2370  FTensor::Tensor1<double *, 3> t_f0_phi_ii(
2371  &f0_phi_ii[phi_shift + HVEC0], &f0_phi_ii[phi_shift + HVEC1],
2372  &f0_phi_ii[phi_shift + HVEC2], 3);
2373  FTensor::Tensor2<double *, 3, 3> t_diff_f0_phi_ii(
2374  &diff_f0_phi_ii[diff_phi_shift + HVEC0_0],
2375  &diff_f0_phi_ii[diff_phi_shift + HVEC0_1],
2376  &diff_f0_phi_ii[diff_phi_shift + HVEC0_2],
2377  &diff_f0_phi_ii[diff_phi_shift + HVEC1_0],
2378  &diff_f0_phi_ii[diff_phi_shift + HVEC1_1],
2379  &diff_f0_phi_ii[diff_phi_shift + HVEC1_2],
2380  &diff_f0_phi_ii[diff_phi_shift + HVEC2_0],
2381  &diff_f0_phi_ii[diff_phi_shift + HVEC2_1],
2382  &diff_f0_phi_ii[diff_phi_shift + HVEC2_2], 9);
2383 
2384  for (int ii = 0; ii <= oo - 2; ii++) {
2385 
2386  int jj = oo - 2 - ii;
2387 
2388  // family I
2390  jj + 1, 2 * ii + 1, n2f0, n0f0 + n1f0 + n2f0, &t_grad_n2f0(0),
2391  &t_grad_n0f0_p_n1f0(0), ifif0, diff_ifif0, 3);
2392  FTensor::Tensor1<double, 3> t_diff_ifif0(
2393  diff_ifif0[0 + jj], diff_ifif0[(jj + 1) + jj],
2394  diff_ifif0[2 * (jj + 1) + jj]);
2395  t_phi(i) = ifif0[jj] * t_f0_phi_ii(i);
2396  t_diff_phi(i, j) = ifif0[jj] * t_diff_f0_phi_ii(i, j) +
2397  t_diff_ifif0(j) * t_f0_phi_ii(i);
2398  ++t_phi;
2399  ++t_diff_phi;
2400  ++t_f0_phi_ii;
2401  ++t_diff_f0_phi_ii;
2402 
2403  }
2404  }
2405  }
2406 
2408  }
VectorDouble diffF0PhiII
Definition: Hcurl.cpp:2314
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
HcurlEdgeBase hCurlBaseOnEdge
Definition: Hcurl.cpp:2313
VectorDouble diffIFiF0
Definition: Hcurl.cpp:2316
FTensor::Index< 'i', 3 > i
Definition: Hcurl.cpp:2319
VectorDouble f0PhiII
Definition: Hcurl.cpp:2314
#define CHKERR
Inline error check.
Definition: definitions.h:595
#define NBEDGE_DEMKOWICZ_HCURL(P)
VectorDouble iFiF0
Definition: Hcurl.cpp:2316
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406
PetscErrorCode IntegratedJacobi_polynomials(int p, double alpha, double x, double t, double *diff_x, double *diff_t, double *L, double *diffL, const int dim)
Calculate integrated Jacobi approximation basis.

◆ operator()() [2/3]

MoFEMErrorCode HcurlFaceBase::operator() ( int  shift,
int  p,
int  nb_integration_pts,
int  n0f0_idx,
int  n1f0_idx,
int  n2f0_idx,
int  n0f1_idx,
int  n1f1_idx,
int  n2f1_idx,
double  n[],
FTensor::Tensor1< double, 3 >  t_grad_n[],
FTensor::Tensor1< double *, 3 > &  t_phi,
FTensor::Tensor2< double *, 3, 3 > &  t_diff_phi 
)

Definition at line 2410 of file Hcurl.cpp.

2415  {
2416 
2418 
2420 
2421  f0PhiII.resize(3 * NBEDGE_DEMKOWICZ_HCURL(p) * nb_integration_pts,false);
2422  diffF0PhiII.resize(9 * NBEDGE_DEMKOWICZ_HCURL(p) * nb_integration_pts,false);
2423  f1PhiII.resize(3 * NBEDGE_DEMKOWICZ_HCURL(p) * nb_integration_pts,false);
2424  diffF1PhiII.resize(9 * NBEDGE_DEMKOWICZ_HCURL(p) * nb_integration_pts,false);
2425 
2426  // edge base for family I
2427  double *f0_phi_ii = &*f0PhiII.data().begin();
2428  double *diff_f0_phi_ii = &*diffF0PhiII.data().begin();
2429  FTensor::Tensor1<double *, 3> t_f0_phi_ii(
2430  &f0_phi_ii[HVEC0], &f0_phi_ii[HVEC1], &f0_phi_ii[HVEC2], 3);
2431  FTensor::Tensor2<double *, 3, 3> t_diff_f0_phi_ii(
2432  &diff_f0_phi_ii[HVEC0_0], &diff_f0_phi_ii[HVEC0_1],
2433  &diff_f0_phi_ii[HVEC0_2], &diff_f0_phi_ii[HVEC1_0],
2434  &diff_f0_phi_ii[HVEC1_1], &diff_f0_phi_ii[HVEC1_2],
2435  &diff_f0_phi_ii[HVEC2_0], &diff_f0_phi_ii[HVEC2_1],
2436  &diff_f0_phi_ii[HVEC2_2], 9);
2437  CHKERR hCurlBaseOnEdge(4, p-1, nb_integration_pts, n0f0_idx, n1f0_idx, n,
2438  t_grad_n, t_f0_phi_ii, t_diff_f0_phi_ii);
2439 
2440  // edge base for family II
2441  double *f1_phi_ii = &*f1PhiII.data().begin();
2442  double *diff_f1_phi_ii = &*diffF1PhiII.data().begin();
2443  FTensor::Tensor1<double *, 3> t_f1_phi_ii(
2444  &f1_phi_ii[HVEC0], &f1_phi_ii[HVEC1], &f1_phi_ii[HVEC2], 3);
2445  FTensor::Tensor2<double *, 3, 3> t_diff_f1_phi_ii(
2446  &diff_f1_phi_ii[HVEC0_0], &diff_f1_phi_ii[HVEC0_1],
2447  &diff_f1_phi_ii[HVEC0_2], &diff_f1_phi_ii[HVEC1_0],
2448  &diff_f1_phi_ii[HVEC1_1], &diff_f1_phi_ii[HVEC1_2],
2449  &diff_f1_phi_ii[HVEC2_0], &diff_f1_phi_ii[HVEC2_1],
2450  &diff_f1_phi_ii[HVEC2_2], 9);
2451  CHKERR hCurlBaseOnEdge(4, p-1, nb_integration_pts, n0f1_idx, n1f1_idx, n,
2452  t_grad_n, t_f1_phi_ii, t_diff_f1_phi_ii);
2453 
2454  FTensor::Tensor1<double, 3> &t_grad_n0f0 = t_grad_n[n0f0_idx];
2455  FTensor::Tensor1<double, 3> &t_grad_n1f0 = t_grad_n[n1f0_idx];
2456  FTensor::Tensor1<double, 3> &t_grad_n2f0 = t_grad_n[n2f0_idx];
2457  FTensor::Tensor1<double, 3> t_grad_n0f0_p_n1f0;
2458  t_grad_n0f0_p_n1f0(i) = t_grad_n0f0(i) + t_grad_n1f0(i);
2459 
2460  FTensor::Tensor1<double, 3> &t_grad_n0f1 = t_grad_n[n0f1_idx];
2461  FTensor::Tensor1<double, 3> &t_grad_n1f1 = t_grad_n[n1f1_idx];
2462  FTensor::Tensor1<double, 3> &t_grad_n2f1 = t_grad_n[n2f1_idx];
2463  FTensor::Tensor1<double, 3> t_grad_n0f1_p_n1f1;
2464  t_grad_n0f1_p_n1f1(i) = t_grad_n0f1(i) + t_grad_n1f1(i);
2465 
2466  iFiF0.resize(p + 1, false);
2467  diffIFiF0.resize(3 * p + 3, false);
2468  double *ifif0 = &*iFiF0.data().begin();
2469  double *diff_ifif0 = &*diffIFiF0.data().begin();
2470  iFiF1.resize(p + 1, false);
2471  diffIFiF1.resize(3 * p + 3, false);
2472  double *ifif1 = &*iFiF1.data().begin();
2473  double *diff_ifif1 = &*diffIFiF1.data().begin();
2474 
2475  for (int gg = 0; gg != nb_integration_pts; ++gg) {
2476 
2477  const int shift_n = shift * gg;
2478  const double n0f0 = n[shift_n + n0f0_idx];
2479  const double n1f0 = n[shift_n + n1f0_idx];
2480  const double n2f0 = n[shift_n + n2f0_idx];
2481  const double n0f1 = n[shift_n + n0f1_idx];
2482  const double n1f1 = n[shift_n + n1f1_idx];
2483  const double n2f1 = n[shift_n + n2f1_idx];
2484 
2485  int phi_shift = 3 * NBEDGE_DEMKOWICZ_HCURL(p - 1) * gg;
2486  int diff_phi_shift = 9 * NBEDGE_DEMKOWICZ_HCURL(p - 1) * gg;
2487 
2488  int kk = 0;
2489  for (int oo = 2; oo <= p; ++oo) {
2490 
2491  FTensor::Tensor1<double *, 3> t_f0_phi_ii(
2492  &f0_phi_ii[phi_shift + HVEC0], &f0_phi_ii[phi_shift + HVEC1],
2493  &f0_phi_ii[phi_shift + HVEC2], 3);
2494  FTensor::Tensor2<double *, 3, 3> t_diff_f0_phi_ii(
2495  &diff_f0_phi_ii[diff_phi_shift + HVEC0_0],
2496  &diff_f0_phi_ii[diff_phi_shift + HVEC0_1],
2497  &diff_f0_phi_ii[diff_phi_shift + HVEC0_2],
2498  &diff_f0_phi_ii[diff_phi_shift + HVEC1_0],
2499  &diff_f0_phi_ii[diff_phi_shift + HVEC1_1],
2500  &diff_f0_phi_ii[diff_phi_shift + HVEC1_2],
2501  &diff_f0_phi_ii[diff_phi_shift + HVEC2_0],
2502  &diff_f0_phi_ii[diff_phi_shift + HVEC2_1],
2503  &diff_f0_phi_ii[diff_phi_shift + HVEC2_2], 9);
2504 
2505  FTensor::Tensor1<double *, 3> t_f1_phi_ii(
2506  &f1_phi_ii[phi_shift + HVEC0],
2507  &f1_phi_ii[phi_shift + HVEC1],
2508  &f1_phi_ii[phi_shift + HVEC2], 3);
2509  FTensor::Tensor2<double *, 3, 3> t_diff_f1_phi_ii(
2510  &diff_f1_phi_ii[diff_phi_shift + HVEC0_0],
2511  &diff_f1_phi_ii[diff_phi_shift + HVEC0_1],
2512  &diff_f1_phi_ii[diff_phi_shift + HVEC0_2],
2513  &diff_f1_phi_ii[diff_phi_shift + HVEC1_0],
2514  &diff_f1_phi_ii[diff_phi_shift + HVEC1_1],
2515  &diff_f1_phi_ii[diff_phi_shift + HVEC1_2],
2516  &diff_f1_phi_ii[diff_phi_shift + HVEC2_0],
2517  &diff_f1_phi_ii[diff_phi_shift + HVEC2_1],
2518  &diff_f1_phi_ii[diff_phi_shift + HVEC2_2], 9);
2519 
2520  for (int ii = 0; ii <= oo - 2; ii++) {
2521 
2522  int jj = oo - 2 - ii;
2523 
2524  // family I
2526  jj + 1, 2 * ii + 1, n2f0, n0f0 + n1f0, &t_grad_n2f0(0),
2527  &t_grad_n0f0_p_n1f0(0), ifif0, diff_ifif0, 3);
2528  FTensor::Tensor1<double, 3> t_diff_ifif0(
2529  diff_ifif0[0 + jj], diff_ifif0[(jj + 1) + jj],
2530  diff_ifif0[2 * (jj + 1) + jj]);
2531 
2532  t_phi(i) = ifif0[jj] * t_f0_phi_ii(i);
2533  t_diff_phi(i, j) = ifif0[jj] * t_diff_f0_phi_ii(i, j) +
2534  t_diff_ifif0(j) * t_f0_phi_ii(i);
2535 
2536  ++t_phi;
2537  ++t_diff_phi;
2538  ++t_f0_phi_ii;
2539  ++t_diff_f0_phi_ii;
2540  ++kk;
2541 
2542  // family II
2544  jj + 1, 2 * ii + 1, n2f1, n0f1 + n1f1, &t_grad_n2f1(0),
2545  &t_grad_n0f1_p_n1f1(0), ifif1, diff_ifif1, 3);
2546  FTensor::Tensor1<double, 3> t_diff_ifif1(
2547  diff_ifif1[0 + jj], diff_ifif1[(jj + 1) + jj],
2548  diff_ifif1[2 * (jj + 1) + jj]);
2549  t_phi(i) = ifif1[jj] * t_f1_phi_ii(i);
2550  t_diff_phi(i, j) = ifif1[jj] * t_diff_f1_phi_ii(i, j) +
2551  t_diff_ifif1(j) * t_f1_phi_ii(i);
2552  ++t_phi;
2553  ++t_diff_phi;
2554  ++t_f1_phi_ii;
2555  ++t_diff_f1_phi_ii;
2556  ++kk;
2557  }
2558  }
2559  if(kk != NBFACETRI_DEMKOWICZ_HCURL(p)) {
2560  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
2561  "Wrong number of base functions");
2562  }
2563 
2564  }
2566  }
VectorDouble diffF0PhiII
Definition: Hcurl.cpp:2314
VectorDouble f1PhiII
Definition: Hcurl.cpp:2315
VectorDouble diffF1PhiII
Definition: Hcurl.cpp:2315
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
HcurlEdgeBase hCurlBaseOnEdge
Definition: Hcurl.cpp:2313
VectorDouble diffIFiF0
Definition: Hcurl.cpp:2316
FTensor::Index< 'i', 3 > i
Definition: Hcurl.cpp:2319
VectorDouble f0PhiII
Definition: Hcurl.cpp:2314
#define CHKERR
Inline error check.
Definition: definitions.h:595
#define NBFACETRI_DEMKOWICZ_HCURL(P)
#define NBEDGE_DEMKOWICZ_HCURL(P)
VectorDouble iFiF1
Definition: Hcurl.cpp:2317
VectorDouble iFiF0
Definition: Hcurl.cpp:2316
VectorDouble diffIFiF1
Definition: Hcurl.cpp:2317
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406
PetscErrorCode IntegratedJacobi_polynomials(int p, double alpha, double x, double t, double *diff_x, double *diff_t, double *L, double *diffL, const int dim)
Calculate integrated Jacobi approximation basis.

◆ operator()() [3/3]

MoFEMErrorCode HcurlFaceBase::operator() ( int  shift,
int  p,
int  nb_integration_pts,
int  n0f0_idx,
int  n1f0_idx,
int  n2f0_idx,
int  n0f1_idx,
int  n1f1_idx,
int  n2f1_idx,
double  n[],
FTensor::Tensor1< double, 3 >  t_grad_n[],
FTensor::Tensor1< double *, 3 > &  t_phi,
FTensor::Tensor2< FTensor::PackPtr< double *, 6 >, 3, 2 > &  t_diff_phi 
)

Definition at line 2568 of file Hcurl.cpp.

2573  {
2575 
2577 
2578  f0PhiII.resize(3 * NBEDGE_DEMKOWICZ_HCURL(p) * nb_integration_pts, false);
2579  diffF0PhiII.resize(6 * NBEDGE_DEMKOWICZ_HCURL(p) * nb_integration_pts,
2580  false);
2581  f1PhiII.resize(3 * NBEDGE_DEMKOWICZ_HCURL(p) * nb_integration_pts, false);
2582  diffF1PhiII.resize(6 * NBEDGE_DEMKOWICZ_HCURL(p) * nb_integration_pts,
2583  false);
2584 
2585  // edge base for family I
2586  double *f0_phi_ii = &*f0PhiII.data().begin();
2587  double *diff_f0_phi_ii = &*diffF0PhiII.data().begin();
2588  FTensor::Tensor1<double *, 3> t_f0_phi_ii(
2589  &f0_phi_ii[HVEC0], &f0_phi_ii[HVEC1], &f0_phi_ii[HVEC2], 3);
2590  FTensor::Tensor2<FTensor::PackPtr<double *, 6>, 3, 2> t_diff_f0_phi_ii(
2591  &diff_f0_phi_ii[HVEC0_0], &diff_f0_phi_ii[HVEC0_1],
2592  &diff_f0_phi_ii[HVEC1_0], &diff_f0_phi_ii[HVEC1_1],
2593  &diff_f0_phi_ii[HVEC2_0], &diff_f0_phi_ii[HVEC2_1]);
2594  CHKERR hCurlBaseOnEdge(3, p - 1, nb_integration_pts, n0f0_idx, n1f0_idx, n,
2595  t_grad_n, t_f0_phi_ii, t_diff_f0_phi_ii);
2596 
2597  // edge base for family II
2598  double *f1_phi_ii = &*f1PhiII.data().begin();
2599  double *diff_f1_phi_ii = &*diffF1PhiII.data().begin();
2600  FTensor::Tensor1<double *, 3> t_f1_phi_ii(
2601  &f1_phi_ii[HVEC0], &f1_phi_ii[HVEC1], &f1_phi_ii[HVEC2], 3);
2602  FTensor::Tensor2<FTensor::PackPtr<double *, 6>, 3, 2> t_diff_f1_phi_ii(
2603  &diff_f1_phi_ii[HVEC0_0], &diff_f1_phi_ii[HVEC0_1],
2604  &diff_f1_phi_ii[HVEC1_0], &diff_f1_phi_ii[HVEC1_1],
2605  &diff_f1_phi_ii[HVEC2_0], &diff_f1_phi_ii[HVEC2_1]);
2606  CHKERR hCurlBaseOnEdge(3, p - 1, nb_integration_pts, n0f1_idx, n1f1_idx, n,
2607  t_grad_n, t_f1_phi_ii, t_diff_f1_phi_ii);
2608 
2609  FTensor::Tensor1<double, 3> &t_grad_n0f0 = t_grad_n[n0f0_idx];
2610  FTensor::Tensor1<double, 3> &t_grad_n1f0 = t_grad_n[n1f0_idx];
2611  FTensor::Tensor1<double, 3> &t_grad_n2f0 = t_grad_n[n2f0_idx];
2612  FTensor::Tensor1<double, 3> t_grad_n0f0_p_n1f0;
2613  t_grad_n0f0_p_n1f0(i) = t_grad_n0f0(i) + t_grad_n1f0(i);
2614 
2615  FTensor::Tensor1<double, 3> &t_grad_n0f1 = t_grad_n[n0f1_idx];
2616  FTensor::Tensor1<double, 3> &t_grad_n1f1 = t_grad_n[n1f1_idx];
2617  FTensor::Tensor1<double, 3> &t_grad_n2f1 = t_grad_n[n2f1_idx];
2618  FTensor::Tensor1<double, 3> t_grad_n0f1_p_n1f1;
2619  t_grad_n0f1_p_n1f1(i) = t_grad_n0f1(i) + t_grad_n1f1(i);
2620 
2621  iFiF0.resize((p+1), false);
2622  diffIFiF0.resize(2 * p + 2, false);
2623  double *ifif0 = &*iFiF0.data().begin();
2624  double *diff_ifif0 = &*diffIFiF0.data().begin();
2625  iFiF1.resize(p + 1, false);
2626  diffIFiF1.resize(2 * p + 2, false);
2627  double *ifif1 = &*iFiF1.data().begin();
2628  double *diff_ifif1 = &*diffIFiF1.data().begin();
2629 
2630  FTensor::Tensor1<double, 2> t_grad_n2f0_j(t_grad_n2f0(0), t_grad_n2f0(1));
2631 
2632  for (int gg = 0; gg != nb_integration_pts; ++gg) {
2633 
2634  const int shift_n = shift * gg;
2635  const double n0f0 = n[shift_n + n0f0_idx];
2636  const double n1f0 = n[shift_n + n1f0_idx];
2637  const double n2f0 = n[shift_n + n2f0_idx];
2638  const double n0f1 = n[shift_n + n0f1_idx];
2639  const double n1f1 = n[shift_n + n1f1_idx];
2640  const double n2f1 = n[shift_n + n2f1_idx];
2641 
2642  int phi_shift = 3 * NBEDGE_DEMKOWICZ_HCURL(p - 1) * gg;
2643  int diff_phi_shift = 6 * NBEDGE_DEMKOWICZ_HCURL(p - 1) * gg;
2644 
2645  int kk = 0;
2646  for (int oo = 2; oo <= p; ++oo) {
2647 
2648  FTensor::Tensor1<double *, 3> t_f0_phi_ii(
2649  &f0_phi_ii[phi_shift + HVEC0], &f0_phi_ii[phi_shift + HVEC1],
2650  &f0_phi_ii[phi_shift + HVEC2], 3);
2651  FTensor::Tensor2<FTensor::PackPtr<double *, 6>, 3, 2> t_diff_f0_phi_ii(
2652  &diff_f0_phi_ii[diff_phi_shift + HVEC0_0],
2653  &diff_f0_phi_ii[diff_phi_shift + HVEC0_1],
2654  &diff_f0_phi_ii[diff_phi_shift + HVEC1_0],
2655  &diff_f0_phi_ii[diff_phi_shift + HVEC1_1],
2656  &diff_f0_phi_ii[diff_phi_shift + HVEC2_0],
2657  &diff_f0_phi_ii[diff_phi_shift + HVEC2_1]);
2658 
2659  FTensor::Tensor1<double *, 3> t_f1_phi_ii(
2660  &f1_phi_ii[phi_shift + HVEC0], &f1_phi_ii[phi_shift + HVEC1],
2661  &f1_phi_ii[phi_shift + HVEC2], 3);
2662  FTensor::Tensor2<FTensor::PackPtr<double *, 6>, 3, 2> t_diff_f1_phi_ii(
2663  &diff_f1_phi_ii[diff_phi_shift + HVEC0_0],
2664  &diff_f1_phi_ii[diff_phi_shift + HVEC0_1],
2665  &diff_f1_phi_ii[diff_phi_shift + HVEC1_0],
2666  &diff_f1_phi_ii[diff_phi_shift + HVEC1_1],
2667  &diff_f1_phi_ii[diff_phi_shift + HVEC2_0],
2668  &diff_f1_phi_ii[diff_phi_shift + HVEC2_1]);
2669 
2670  for (int ii = 0; ii <= oo - 2; ii++) {
2671 
2672  int jj = oo - 2 - ii;
2673 
2674  // family I
2676  jj + 1, 2 * ii + 1, n2f0, n0f0 + n1f0, &t_grad_n2f0(0),
2677  &t_grad_n0f0_p_n1f0(0), ifif0, diff_ifif0, 2);
2678  FTensor::Tensor1<double, 2> t_diff_ifif0(
2679  diff_ifif0[0 + jj], diff_ifif0[(jj + 1) + jj]);
2680  t_phi(i) = ifif0[jj] * t_f0_phi_ii(i);
2681  t_diff_phi(i, j) = ifif0[jj] * t_diff_f0_phi_ii(i, j) +
2682  t_f0_phi_ii(i) * t_diff_ifif0(j);
2683 
2684  ++t_phi;
2685  ++t_diff_phi;
2686  ++t_f0_phi_ii;
2687  ++t_diff_f0_phi_ii;
2688  ++kk;
2689 
2690  // family II
2691  CHKERR IntegratedJacobi_polynomials(jj + 1, 2 * ii + 1, n2f1, n0f1+n1f1,
2692  &t_grad_n2f1(0), &t_grad_n0f1_p_n1f1(0),
2693  ifif1, diff_ifif1, 2);
2694  FTensor::Tensor1<double, 2> t_diff_ifif1(diff_ifif1[0 + jj],
2695  diff_ifif1[(jj + 1) + jj]);
2696  t_phi(i) = ifif1[jj] * t_f1_phi_ii(i);
2697  t_diff_phi(i, j) = ifif1[jj] * t_diff_f1_phi_ii(i, j) +
2698  t_f1_phi_ii(i) * t_diff_ifif1(j);
2699  ++t_phi;
2700  ++t_diff_phi;
2701  ++t_f1_phi_ii;
2702  ++t_diff_f1_phi_ii;
2703  ++kk;
2704  }
2705  }
2706  if(kk != NBFACETRI_DEMKOWICZ_HCURL(p)) {
2707  SETERRQ(PETSC_COMM_SELF, MOFEM_DATA_INCONSISTENCY,
2708  "Wrong number of base functions");
2709  }
2710  }
2712  }
VectorDouble diffF0PhiII
Definition: Hcurl.cpp:2314
VectorDouble f1PhiII
Definition: Hcurl.cpp:2315
VectorDouble diffF1PhiII
Definition: Hcurl.cpp:2315
#define MoFEMFunctionReturn(a)
Last executable line of each PETSc function used for error handling. Replaces return()
Definition: definitions.h:476
HcurlEdgeBase hCurlBaseOnEdge
Definition: Hcurl.cpp:2313
VectorDouble diffIFiF0
Definition: Hcurl.cpp:2316
FTensor::Index< 'i', 3 > i
Definition: Hcurl.cpp:2319
VectorDouble f0PhiII
Definition: Hcurl.cpp:2314
#define CHKERR
Inline error check.
Definition: definitions.h:595
#define NBFACETRI_DEMKOWICZ_HCURL(P)
#define NBEDGE_DEMKOWICZ_HCURL(P)
VectorDouble iFiF1
Definition: Hcurl.cpp:2317
VectorDouble iFiF0
Definition: Hcurl.cpp:2316
VectorDouble diffIFiF1
Definition: Hcurl.cpp:2317
#define MoFEMFunctionBegin
First executable line of each MoFEM function, used for error handling. Final line of MoFEM functions ...
Definition: definitions.h:406
PetscErrorCode IntegratedJacobi_polynomials(int p, double alpha, double x, double t, double *diff_x, double *diff_t, double *L, double *diffL, const int dim)
Calculate integrated Jacobi approximation basis.

Member Data Documentation

◆ diffF0PhiII

VectorDouble HcurlFaceBase::diffF0PhiII

Definition at line 2314 of file Hcurl.cpp.

◆ diffF1PhiII

VectorDouble HcurlFaceBase::diffF1PhiII

Definition at line 2315 of file Hcurl.cpp.

◆ diffIFiF0

VectorDouble HcurlFaceBase::diffIFiF0

Definition at line 2316 of file Hcurl.cpp.

◆ diffIFiF1

VectorDouble HcurlFaceBase::diffIFiF1

Definition at line 2317 of file Hcurl.cpp.

◆ f0PhiII

VectorDouble HcurlFaceBase::f0PhiII

Definition at line 2314 of file Hcurl.cpp.

◆ f1PhiII

VectorDouble HcurlFaceBase::f1PhiII

Definition at line 2315 of file Hcurl.cpp.

◆ hCurlBaseOnEdge

HcurlEdgeBase HcurlFaceBase::hCurlBaseOnEdge

Definition at line 2313 of file Hcurl.cpp.

◆ i

FTensor::Index<'i', 3> HcurlFaceBase::i

Definition at line 2319 of file Hcurl.cpp.

◆ iFiF0

VectorDouble HcurlFaceBase::iFiF0

Definition at line 2316 of file Hcurl.cpp.

◆ iFiF1

VectorDouble HcurlFaceBase::iFiF1

Definition at line 2317 of file Hcurl.cpp.


The documentation for this struct was generated from the following file: