16 #include <deal.II/base/std_cxx14/memory.h> 18 #include <deal.II/fe/fe_nedelec_sz.h> 20 DEAL_II_NAMESPACE_OPEN
23 template <
int dim,
int spacedim>
30 std::vector<bool>(compute_num_dofs(degree), true),
32 std::vector<bool>(dim, true)))
40 for (
unsigned int comp = 0; comp < this->
n_components(); ++comp)
50 template <
int dim,
int spacedim>
62 template <
int dim,
int spacedim>
67 const unsigned int )
const 76 template <
int dim,
int spacedim>
88 template <
int dim,
int spacedim>
93 const unsigned int )
const 101 template <
int dim,
int spacedim>
113 template <
int dim,
int spacedim>
118 const unsigned int )
const 125 template <
int dim,
int spacedim>
126 std::unique_ptr<typename ::FiniteElement<dim, spacedim>::InternalDataBase>
135 auto data = std_cxx14::make_unique<InternalData>();
145 const unsigned int n_line_dofs = this->
dofs_per_line * lines_per_cell;
146 const unsigned int n_face_dofs = this->
dofs_per_quad * faces_per_cell;
149 const unsigned int n_q_points = quadrature.
size();
152 data->sigma_imj_values.resize(
154 std::vector<std::vector<double>>(vertices_per_cell,
155 std::vector<double>(vertices_per_cell)));
157 data->sigma_imj_grads.resize(vertices_per_cell,
158 std::vector<std::vector<double>>(
159 vertices_per_cell, std::vector<double>(dim)));
180 std::vector<Point<dim>> p_list(n_q_points);
190 std::vector<std::vector<double>> sigma(
191 n_q_points, std::vector<double>(lines_per_cell));
192 std::vector<std::vector<double>> lambda(
193 n_q_points, std::vector<double>(lines_per_cell));
195 for (
unsigned int q = 0; q < n_q_points; ++q)
197 sigma[q][0] = (1.0 - p_list[q][0]) + (1.0 - p_list[q][1]);
198 sigma[q][1] = p_list[q][0] + (1.0 - p_list[q][1]);
199 sigma[q][2] = (1.0 - p_list[q][0]) + p_list[q][1];
200 sigma[q][3] = p_list[q][0] + p_list[q][1];
202 lambda[q][0] = (1.0 - p_list[q][0]) * (1.0 - p_list[q][1]);
203 lambda[q][1] = p_list[q][0] * (1.0 - p_list[q][1]);
204 lambda[q][2] = (1.0 - p_list[q][0]) * p_list[q][1];
205 lambda[q][3] = p_list[q][0] * p_list[q][1];
206 for (
unsigned int i = 0; i < vertices_per_cell; ++i)
208 for (
unsigned int j = 0; j < vertices_per_cell; ++j)
210 data->sigma_imj_values[q][i][j] =
211 sigma[q][i] - sigma[q][j];
224 {-1, -1}, {1, -1}, {-1, 1}, {1, 1}};
225 int sigma_imj_sign[vertices_per_cell][vertices_per_cell];
226 unsigned int sigma_imj_component[vertices_per_cell]
229 for (
unsigned int i = 0; i < vertices_per_cell; ++i)
231 for (
unsigned int j = 0; j < vertices_per_cell; ++j)
238 sigma_imj_sign[i][j] = (i < j) ? -1 : 1;
239 sigma_imj_sign[i][j] = (i == j) ? 0 : sigma_imj_sign[i][j];
243 sigma_imj_component[i][j] = 0;
244 for (
unsigned int d = 0; d < dim; ++d)
247 sigma_comp_signs[i][d] - sigma_comp_signs[j][d];
252 sigma_imj_component[i][j] = d;
259 data->sigma_imj_grads[i][j][sigma_imj_component[i][j]] =
260 2.0 * (double)sigma_imj_sign[i][j];
268 data->edge_sigma_values.resize(lines_per_cell);
269 data->edge_sigma_grads.resize(lines_per_cell);
270 for (
unsigned int m = 0; m < lines_per_cell; ++m)
272 data->edge_sigma_values[m].resize(n_q_points);
275 data->edge_sigma_grads[m].resize(dim);
285 data->edge_lambda_values.resize(lines_per_cell,
286 std::vector<double>(n_q_points));
287 data->edge_lambda_grads_2d.resize(lines_per_cell,
288 std::vector<double>(dim));
289 for (
unsigned int m = 0; m < lines_per_cell; ++m)
294 const unsigned int e1(
296 const unsigned int e2(
298 for (
unsigned int q = 0; q < n_q_points; ++q)
300 data->edge_sigma_values[m][q] =
301 data->sigma_imj_values[q][e2][e1];
302 data->edge_lambda_values[m][q] =
303 lambda[q][e1] + lambda[q][e2];
306 data->edge_sigma_grads[m][edge_sigma_direction[m]] = -2.0;
308 data->edge_lambda_grads_2d[0] = {-1.0, 0.0};
309 data->edge_lambda_grads_2d[1] = {1.0, 0.0};
310 data->edge_lambda_grads_2d[2] = {0.0, -1.0};
311 data->edge_lambda_grads_2d[3] = {0.0, 1.0};
379 const unsigned int cell_type1_offset = n_line_dofs;
386 const unsigned int cell_type2_offset =
387 cell_type1_offset + degree *
degree;
395 const unsigned int cell_type3_offset1 =
396 cell_type2_offset + degree *
degree;
397 const unsigned int cell_type3_offset2 = cell_type3_offset1 +
degree;
399 if (flags & (update_values | update_gradients))
402 std::vector<Point<dim>> cell_points(n_q_points);
403 for (
unsigned int q = 0; q < n_q_points; ++q)
405 for (
unsigned int d = 0; d < dim; ++d)
407 cell_points[q][d] = 2.0 * p_list[q][d] - 1.0;
412 for (
unsigned int q = 0; q < n_q_points; ++q)
426 const unsigned int poly_length(
427 (flags & update_gradients) ? 3 : 2);
429 std::vector<std::vector<double>> polyx(
430 degree, std::vector<double>(poly_length));
431 std::vector<std::vector<double>> polyy(
432 degree, std::vector<double>(poly_length));
433 for (
unsigned int i = 0; i <
degree; ++i)
439 cell_points[q][0], polyx[i]);
441 cell_points[q][1], polyy[i]);
444 if (flags & update_values)
446 for (
unsigned int j = 0; j <
degree; ++j)
448 const unsigned int shift_j(j * degree);
449 for (
unsigned int i = 0; i <
degree; ++i)
451 const unsigned int shift_ij(i + shift_j);
454 const unsigned int dof_index1(cell_type1_offset +
456 data->shape_values[dof_index1][q][0] =
457 2.0 * polyx[i][1] * polyy[j][0];
458 data->shape_values[dof_index1][q][1] =
459 2.0 * polyx[i][0] * polyy[j][1];
462 const unsigned int dof_index2(cell_type2_offset +
464 data->shape_values[dof_index2][q][0] =
465 data->shape_values[dof_index1][q][0];
466 data->shape_values[dof_index2][q][1] =
467 -1.0 * data->shape_values[dof_index1][q][1];
470 const unsigned int dof_index3_1(cell_type3_offset1 +
472 data->shape_values[dof_index3_1][q][0] = polyy[j][0];
473 data->shape_values[dof_index3_1][q][1] = 0.0;
475 const unsigned int dof_index3_2(cell_type3_offset2 +
477 data->shape_values[dof_index3_2][q][0] = 0.0;
478 data->shape_values[dof_index3_2][q][1] = polyx[j][0];
481 if (flags & update_gradients)
483 for (
unsigned int j = 0; j <
degree; ++j)
485 const unsigned int shift_j(j * degree);
486 for (
unsigned int i = 0; i <
degree; ++i)
488 const unsigned int shift_ij(i + shift_j);
491 const unsigned int dof_index1(cell_type1_offset +
493 data->shape_grads[dof_index1][q][0][0] =
494 4.0 * polyx[i][2] * polyy[j][0];
495 data->shape_grads[dof_index1][q][0][1] =
496 4.0 * polyx[i][1] * polyy[j][1];
497 data->shape_grads[dof_index1][q][1][0] =
498 data->shape_grads[dof_index1][q][0][1];
499 data->shape_grads[dof_index1][q][1][1] =
500 4.0 * polyx[i][0] * polyy[j][2];
503 const unsigned int dof_index2(cell_type2_offset +
505 data->shape_grads[dof_index2][q][0][0] =
506 data->shape_grads[dof_index1][q][0][0];
507 data->shape_grads[dof_index2][q][0][1] =
508 data->shape_grads[dof_index1][q][0][1];
509 data->shape_grads[dof_index2][q][1][0] =
510 -1.0 * data->shape_grads[dof_index1][q][1][0];
511 data->shape_grads[dof_index2][q][1][1] =
512 -1.0 * data->shape_grads[dof_index1][q][1][1];
515 const unsigned int dof_index3_1(cell_type3_offset1 +
517 data->shape_grads[dof_index3_1][q][0][0] = 0.0;
518 data->shape_grads[dof_index3_1][q][0][1] =
520 data->shape_grads[dof_index3_1][q][1][0] = 0.0;
521 data->shape_grads[dof_index3_1][q][1][1] = 0.0;
523 const unsigned int dof_index3_2(cell_type3_offset2 +
525 data->shape_grads[dof_index3_2][q][0][0] = 0.0;
526 data->shape_grads[dof_index3_2][q][0][1] = 0.0;
527 data->shape_grads[dof_index3_2][q][1][0] =
529 data->shape_grads[dof_index3_2][q][1][1] = 0.0;
540 std::vector<std::vector<double>> sigma(
541 n_q_points, std::vector<double>(lines_per_cell));
542 std::vector<std::vector<double>> lambda(
543 n_q_points, std::vector<double>(lines_per_cell));
544 for (
unsigned int q = 0; q < n_q_points; ++q)
546 sigma[q][0] = (1.0 - p_list[q][0]) + (1.0 - p_list[q][1]) +
549 p_list[q][0] + (1.0 - p_list[q][1]) + (1 - p_list[q][2]);
551 (1.0 - p_list[q][0]) + p_list[q][1] + (1 - p_list[q][2]);
552 sigma[q][3] = p_list[q][0] + p_list[q][1] + (1 - p_list[q][2]);
554 (1.0 - p_list[q][0]) + (1.0 - p_list[q][1]) + p_list[q][2];
555 sigma[q][5] = p_list[q][0] + (1.0 - p_list[q][1]) + p_list[q][2];
556 sigma[q][6] = (1.0 - p_list[q][0]) + p_list[q][1] + p_list[q][2];
557 sigma[q][7] = p_list[q][0] + p_list[q][1] + p_list[q][2];
559 lambda[q][0] = (1.0 - p_list[q][0]) * (1.0 - p_list[q][1]) *
560 (1.0 - p_list[q][2]);
562 p_list[q][0] * (1.0 - p_list[q][1]) * (1.0 - p_list[q][2]);
564 (1.0 - p_list[q][0]) * p_list[q][1] * (1.0 - p_list[q][2]);
565 lambda[q][3] = p_list[q][0] * p_list[q][1] * (1.0 - p_list[q][2]);
567 (1.0 - p_list[q][0]) * (1.0 - p_list[q][1]) * p_list[q][2];
568 lambda[q][5] = p_list[q][0] * (1.0 - p_list[q][1]) * p_list[q][2];
569 lambda[q][6] = (1.0 - p_list[q][0]) * p_list[q][1] * p_list[q][2];
570 lambda[q][7] = p_list[q][0] * p_list[q][1] * p_list[q][2];
574 for (
unsigned int i = 0; i < vertices_per_cell; ++i)
576 for (
unsigned int j = 0; j < vertices_per_cell; ++j)
578 data->sigma_imj_values[q][i][j] =
579 sigma[q][i] - sigma[q][j];
612 int sigma_imj_sign[vertices_per_cell][vertices_per_cell];
613 unsigned int sigma_imj_component[vertices_per_cell]
616 for (
unsigned int i = 0; i < vertices_per_cell; ++i)
618 for (
unsigned int j = 0; j < vertices_per_cell; ++j)
625 sigma_imj_sign[i][j] = (i < j) ? -1 : 1;
626 sigma_imj_sign[i][j] = (i == j) ? 0 : sigma_imj_sign[i][j];
630 sigma_imj_component[i][j] = 0;
631 for (
unsigned int d = 0; d < dim; ++d)
634 sigma_comp_signs[i][d] - sigma_comp_signs[j][d];
639 sigma_imj_component[i][j] = d;
646 data->sigma_imj_grads[i][j][sigma_imj_component[i][j]] =
647 2.0 * (double)sigma_imj_sign[i][j];
654 data->edge_sigma_values.resize(lines_per_cell);
655 data->edge_lambda_values.resize(lines_per_cell);
656 data->edge_sigma_grads.resize(lines_per_cell);
657 data->edge_lambda_grads_3d.resize(lines_per_cell);
658 data->edge_lambda_gradgrads_3d.resize(lines_per_cell);
659 for (
unsigned int m = 0; m < lines_per_cell; ++m)
661 data->edge_sigma_values[m].resize(n_q_points);
662 data->edge_lambda_values[m].resize(n_q_points);
665 data->edge_sigma_grads[m].resize(dim);
667 data->edge_lambda_grads_3d[m].resize(n_q_points);
668 for (
unsigned int q = 0; q < n_q_points; ++q)
670 data->edge_lambda_grads_3d[m][q].resize(dim);
674 data->edge_lambda_gradgrads_3d[m].resize(dim);
675 for (
unsigned int d = 0; d < dim; ++d)
677 data->edge_lambda_gradgrads_3d[m][d].resize(dim);
684 1, 1, 0, 0, 1, 1, 0, 0, 2, 2, 2, 2};
686 for (
unsigned int m = 0; m < lines_per_cell; ++m)
691 const unsigned int e1(
693 const unsigned int e2(
696 for (
unsigned int q = 0; q < n_q_points; ++q)
698 data->edge_sigma_values[m][q] =
699 data->sigma_imj_values[q][e2][e1];
700 data->edge_lambda_values[m][q] =
701 lambda[q][e1] + lambda[q][e2];
704 data->edge_sigma_grads[m][edge_sigma_direction[m]] = -2.0;
707 for (
unsigned int q = 0; q < n_q_points; ++q)
709 double x(p_list[q][0]);
710 double y(p_list[q][1]);
711 double z(p_list[q][2]);
712 data->edge_lambda_grads_3d[0][q] = {z - 1.0, 0.0, x - 1.0};
713 data->edge_lambda_grads_3d[1][q] = {1.0 - z, 0.0, -x};
714 data->edge_lambda_grads_3d[2][q] = {0.0, z - 1.0, y - 1.0};
715 data->edge_lambda_grads_3d[3][q] = {0.0, 1.0 - z, -y};
716 data->edge_lambda_grads_3d[4][q] = {-z, 0.0, 1.0 - x};
717 data->edge_lambda_grads_3d[5][q] = {z, 0.0, x};
718 data->edge_lambda_grads_3d[6][q] = {0.0, -z, 1.0 - y};
719 data->edge_lambda_grads_3d[7][q] = {0.0, z, y};
720 data->edge_lambda_grads_3d[8][q] = {y - 1.0, x - 1.0, 0.0};
721 data->edge_lambda_grads_3d[9][q] = {1.0 - y, -x, 0.0};
722 data->edge_lambda_grads_3d[10][q] = {-y, 1.0 - x, 0.0};
723 data->edge_lambda_grads_3d[11][q] = {y, x, 0.0};
753 for (
unsigned int m = 0; m < lines_per_cell; ++m)
755 data->edge_lambda_gradgrads_3d[m][edge_lambda_directions[m][0]]
756 [edge_lambda_directions[m][1]] =
757 (double)edge_lambda_sign[m] * 1.0;
758 data->edge_lambda_gradgrads_3d[m][edge_lambda_directions[m][1]]
759 [edge_lambda_directions[m][0]] =
760 (double)edge_lambda_sign[m] * 1.0;
767 data->face_lambda_values.resize(faces_per_cell);
768 data->face_lambda_grads.resize(faces_per_cell);
770 for (
unsigned int m = 0; m < faces_per_cell; ++m)
772 data->face_lambda_values[m].resize(n_q_points);
773 data->face_lambda_grads[m].resize(3);
776 for (
unsigned int q = 0; q < n_q_points; ++q)
778 double x(p_list[q][0]);
779 double y(p_list[q][1]);
780 double z(p_list[q][2]);
781 data->face_lambda_values[0][q] = 1.0 - x;
782 data->face_lambda_values[1][q] = x;
783 data->face_lambda_values[2][q] = 1.0 - y;
784 data->face_lambda_values[3][q] = y;
785 data->face_lambda_values[4][q] = 1.0 - z;
786 data->face_lambda_values[5][q] = z;
789 data->face_lambda_grads[0] = {-1.0, 0.0, 0.0};
790 data->face_lambda_grads[1] = {1.0, 0.0, 0.0};
791 data->face_lambda_grads[2] = {0.0, -1.0, 0.0};
792 data->face_lambda_grads[3] = {0.0, 1.0, 0.0};
793 data->face_lambda_grads[4] = {0.0, 0.0, -1.0};
794 data->face_lambda_grads[5] = {0.0, 0.0, 1.0};
798 if (flags & (update_values | update_gradients))
807 const unsigned int cell_type1_offset(n_line_dofs +
820 const unsigned int cell_type2_offset1(
821 cell_type1_offset + degree * degree * degree);
822 const unsigned int cell_type2_offset2(
823 cell_type2_offset1 + degree * degree * degree);
833 const unsigned int cell_type3_offset1(
834 cell_type2_offset2 + degree * degree * degree);
835 const unsigned int cell_type3_offset2(cell_type3_offset1 +
837 const unsigned int cell_type3_offset3(cell_type3_offset2 +
841 std::vector<Point<dim>> cell_points(n_q_points);
842 for (
unsigned int q = 0; q < n_q_points; ++q)
844 for (
unsigned int d = 0; d < dim; ++d)
846 cell_points[q][d] = 2.0 * p_list[q][d] - 1.0;
853 const unsigned int poly_length(
854 (flags & update_gradients) ? 3 : 2);
856 for (
unsigned int q = 0; q < n_q_points; ++q)
864 std::vector<std::vector<double>> polyx(
865 degree, std::vector<double>(poly_length));
866 std::vector<std::vector<double>> polyy(
867 degree, std::vector<double>(poly_length));
868 std::vector<std::vector<double>> polyz(
869 degree, std::vector<double>(poly_length));
870 for (
unsigned int i = 0; i <
degree; ++i)
874 cell_points[q][0], polyx[i]);
876 cell_points[q][1], polyy[i]);
878 cell_points[q][2], polyz[i]);
881 if (flags & update_values)
883 for (
unsigned int k = 0; k <
degree; ++k)
885 const unsigned int shift_k(k * degree * degree);
886 const unsigned int shift_j(
889 for (
unsigned int j = 0; j <
degree; ++j)
891 const unsigned int shift_jk(j * degree +
893 for (
unsigned int i = 0; i <
degree; ++i)
895 const unsigned int shift_ijk(shift_jk +
899 const unsigned int dof_index1(
900 cell_type1_offset + shift_ijk);
902 data->shape_values[dof_index1][q][0] =
903 2.0 * polyx[i][1] * polyy[j][0] *
905 data->shape_values[dof_index1][q][1] =
906 2.0 * polyx[i][0] * polyy[j][1] *
908 data->shape_values[dof_index1][q][2] =
909 2.0 * polyx[i][0] * polyy[j][0] *
913 const unsigned int dof_index2_1(
914 cell_type2_offset1 + shift_ijk);
915 const unsigned int dof_index2_2(
916 cell_type2_offset2 + shift_ijk);
918 data->shape_values[dof_index2_1][q][0] =
919 data->shape_values[dof_index1][q][0];
920 data->shape_values[dof_index2_1][q][1] =
922 data->shape_values[dof_index1][q][1];
923 data->shape_values[dof_index2_1][q][2] =
924 data->shape_values[dof_index1][q][2];
926 data->shape_values[dof_index2_2][q][0] =
927 data->shape_values[dof_index1][q][0];
928 data->shape_values[dof_index2_2][q][1] =
930 data->shape_values[dof_index1][q][1];
931 data->shape_values[dof_index2_2][q][2] =
933 data->shape_values[dof_index1][q][2];
937 const unsigned int shift_ij(
940 const unsigned int dof_index3_1(
941 cell_type3_offset1 + shift_ij);
942 const unsigned int dof_index3_2(
943 cell_type3_offset2 + shift_ij);
944 const unsigned int dof_index3_3(
945 cell_type3_offset3 + shift_ij);
947 data->shape_values[dof_index3_1][q][0] =
948 polyy[j][0] * polyz[k][0];
949 data->shape_values[dof_index3_1][q][1] = 0.0;
950 data->shape_values[dof_index3_1][q][2] = 0.0;
952 data->shape_values[dof_index3_2][q][0] = 0.0;
953 data->shape_values[dof_index3_2][q][1] =
954 polyx[j][0] * polyz[k][0];
955 data->shape_values[dof_index3_2][q][2] = 0.0;
957 data->shape_values[dof_index3_3][q][0] = 0.0;
958 data->shape_values[dof_index3_3][q][1] = 0.0;
959 data->shape_values[dof_index3_3][q][2] =
960 polyx[j][0] * polyy[k][0];
964 if (flags & update_gradients)
966 for (
unsigned int k = 0; k <
degree; ++k)
968 const unsigned int shift_k(k * degree * degree);
969 const unsigned int shift_j(
972 for (
unsigned int j = 0; j <
degree; ++j)
974 const unsigned int shift_jk(j * degree +
976 for (
unsigned int i = 0; i <
degree; ++i)
978 const unsigned int shift_ijk(shift_jk +
982 const unsigned int dof_index1(
983 cell_type1_offset + shift_ijk);
985 data->shape_grads[dof_index1][q][0][0] =
986 4.0 * polyx[i][2] * polyy[j][0] *
988 data->shape_grads[dof_index1][q][0][1] =
989 4.0 * polyx[i][1] * polyy[j][1] *
991 data->shape_grads[dof_index1][q][0][2] =
992 4.0 * polyx[i][1] * polyy[j][0] *
995 data->shape_grads[dof_index1][q][1][0] =
996 data->shape_grads[dof_index1][q][0][1];
997 data->shape_grads[dof_index1][q][1][1] =
998 4.0 * polyx[i][0] * polyy[j][2] *
1000 data->shape_grads[dof_index1][q][1][2] =
1001 4.0 * polyx[i][0] * polyy[j][1] *
1004 data->shape_grads[dof_index1][q][2][0] =
1005 data->shape_grads[dof_index1][q][0][2];
1006 data->shape_grads[dof_index1][q][2][1] =
1007 data->shape_grads[dof_index1][q][1][2];
1008 data->shape_grads[dof_index1][q][2][2] =
1009 4.0 * polyx[i][0] * polyy[j][0] *
1013 const unsigned int dof_index2_1(
1014 cell_type2_offset1 + shift_ijk);
1015 const unsigned int dof_index2_2(
1016 cell_type2_offset2 + shift_ijk);
1018 for (
unsigned int d = 0; d < dim; ++d)
1020 data->shape_grads[dof_index2_1][q][0]
1022 data->shape_grads[dof_index1][q][0]
1024 data->shape_grads[dof_index2_1][q][1]
1026 -1.0 * data->shape_grads[dof_index1]
1028 data->shape_grads[dof_index2_1][q][2]
1030 data->shape_grads[dof_index1][q][2]
1033 data->shape_grads[dof_index2_2][q][0]
1035 data->shape_grads[dof_index1][q][0]
1037 data->shape_grads[dof_index2_2][q][1]
1039 -1.0 * data->shape_grads[dof_index1]
1041 data->shape_grads[dof_index2_2][q][2]
1043 -1.0 * data->shape_grads[dof_index1]
1049 const unsigned int shift_ij(
1052 const unsigned int dof_index3_1(
1053 cell_type3_offset1 + shift_ij);
1054 const unsigned int dof_index3_2(
1055 cell_type3_offset2 + shift_ij);
1056 const unsigned int dof_index3_3(
1057 cell_type3_offset3 + shift_ij);
1058 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1060 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1062 data->shape_grads[dof_index3_1][q][d1]
1064 data->shape_grads[dof_index3_2][q][d1]
1066 data->shape_grads[dof_index3_3][q][d1]
1070 data->shape_grads[dof_index3_1][q][0][1] =
1071 2.0 * polyy[j][1] * polyz[k][0];
1072 data->shape_grads[dof_index3_1][q][0][2] =
1073 2.0 * polyy[j][0] * polyz[k][1];
1075 data->shape_grads[dof_index3_2][q][1][0] =
1076 2.0 * polyx[j][1] * polyz[k][0];
1077 data->shape_grads[dof_index3_2][q][1][2] =
1078 2.0 * polyx[j][0] * polyz[k][1];
1080 data->shape_grads[dof_index3_3][q][2][0] =
1081 2.0 * polyx[j][1] * polyy[k][0];
1082 data->shape_grads[dof_index3_3][q][2][1] =
1083 2.0 * polyx[j][0] * polyy[k][1];
1097 return std::move(data);
1100 template <
int dim,
int spacedim>
1121 const unsigned int n_q_points = quadrature.
size();
1127 Assert(!(flags & update_values) ||
1132 const unsigned int degree(
1137 const unsigned int vertices_per_line(2);
1141 std::vector<std::vector<unsigned int>> edge_order(
1142 lines_per_cell, std::vector<unsigned int>(vertices_per_line));
1154 std::vector<int> edge_sign(lines_per_cell);
1155 for (
unsigned int m = 0; m < lines_per_cell; ++m)
1157 unsigned int v0_loc =
1159 unsigned int v1_loc =
1161 unsigned int v0_glob = cell->vertex_index(v0_loc);
1162 unsigned int v1_glob = cell->vertex_index(v1_loc);
1164 if (v0_glob > v1_glob)
1167 edge_sign[m] = -1.0;
1211 std::vector<std::vector<double>> edge_sigma_values(
1213 std::vector<std::vector<double>> edge_sigma_grads(
1216 std::vector<std::vector<double>> edge_lambda_values(
1218 std::vector<std::vector<double>> edge_lambda_grads(
1222 for (
unsigned int m = 0; m < lines_per_cell; ++m)
1224 std::transform(edge_sigma_values[m].begin(),
1225 edge_sigma_values[m].end(),
1226 edge_sigma_values[m].begin(),
1227 [&](
const double edge_sigma_value) {
1228 return edge_sign[m] * edge_sigma_value;
1231 std::transform(edge_sigma_grads[m].begin(),
1232 edge_sigma_grads[m].end(),
1233 edge_sigma_grads[m].begin(),
1234 [&](
const double edge_sigma_grad) {
1235 return edge_sign[m] * edge_sigma_grad;
1245 for (
unsigned int m = 0; m < lines_per_cell; ++m)
1248 for (
unsigned int q = 0; q < n_q_points; ++q)
1251 std::vector<std::vector<double>> poly(
1252 degree, std::vector<double>(poly_length));
1253 for (
unsigned int i = 1; i < degree + 1; ++i)
1259 edge_sigma_values[m][q], poly[i - 1]);
1261 if (flags & update_values)
1264 for (
unsigned int d = 0; d < dim; ++d)
1267 0.5 * edge_sigma_grads[m][d] *
1268 edge_lambda_values[m][q];
1271 for (
unsigned int i = 1; i < degree + 1; ++i)
1273 const unsigned int poly_index = i - 1;
1274 const unsigned int dof_index(i + shift_m);
1275 for (
unsigned int d = 0; d < dim; ++d)
1278 edge_sigma_grads[m][d] *
1279 poly[poly_index][1] *
1280 edge_lambda_values[m][q] +
1281 poly[poly_index][0] *
1282 edge_lambda_grads[m][d];
1286 if (flags & update_gradients)
1289 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1291 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1296 0.5 * edge_sigma_grads[m][d1] *
1297 edge_lambda_grads[m][d2];
1301 for (
unsigned int i = 1; i < degree + 1; ++i)
1303 const unsigned int poly_index = i - 1;
1304 const unsigned int dof_index(i + shift_m);
1307 edge_sigma_grads[m][0] *
1308 edge_sigma_grads[m][0] *
1309 edge_lambda_values[m][q] * poly[poly_index][2];
1312 (edge_sigma_grads[m][0] *
1313 edge_lambda_grads[m][1] +
1314 edge_sigma_grads[m][1] *
1315 edge_lambda_grads[m][0]) *
1316 poly[poly_index][1];
1322 edge_sigma_grads[m][1] *
1323 edge_sigma_grads[m][1] *
1324 edge_lambda_values[m][q] * poly[poly_index][2];
1339 std::vector<int> edge_sign(lines_per_cell);
1340 for (
unsigned int m = 0; m < lines_per_cell; ++m)
1342 unsigned int v0_loc =
1344 unsigned int v1_loc =
1346 unsigned int v0_glob = cell->vertex_index(v0_loc);
1347 unsigned int v1_glob = cell->vertex_index(v1_loc);
1349 if (v0_glob > v1_glob)
1352 edge_sign[m] = -1.0;
1398 std::vector<std::vector<double>> edge_sigma_values(
1400 std::vector<std::vector<double>> edge_lambda_values(
1402 std::vector<std::vector<double>> edge_sigma_grads(
1404 std::vector<std::vector<std::vector<double>>> edge_lambda_grads(
1406 std::vector<std::vector<std::vector<double>>>
1410 for (
unsigned int m = 0; m < lines_per_cell; ++m)
1412 std::transform(edge_sigma_values[m].begin(),
1413 edge_sigma_values[m].end(),
1414 edge_sigma_values[m].begin(),
1415 [&](
const double edge_sigma_value) {
1416 return edge_sign[m] * edge_sigma_value;
1418 std::transform(edge_sigma_grads[m].begin(),
1419 edge_sigma_grads[m].end(),
1420 edge_sigma_grads[m].begin(),
1421 [&](
const double edge_sigma_grad) {
1422 return edge_sign[m] * edge_sigma_grad;
1432 std::vector<std::vector<double>> poly(
1433 degree, std::vector<double>(poly_length));
1434 for (
unsigned int m = 0; m < lines_per_cell; ++m)
1437 for (
unsigned int q = 0; q < n_q_points; ++q)
1442 for (
unsigned int i = 0; i <
degree; ++i)
1445 edge_sigma_values[m][q], poly[i]);
1448 if (flags & update_values)
1451 for (
unsigned int d = 0; d < dim; ++d)
1454 0.5 * edge_sigma_grads[m][d] *
1455 edge_lambda_values[m][q];
1459 for (
unsigned int i = 0; i <
degree; ++i)
1461 const unsigned int dof_index(i + 1 + shift_m);
1462 for (
unsigned int d = 0; d < dim; ++d)
1465 edge_sigma_grads[m][d] * poly[i][1] *
1466 edge_lambda_values[m][q] +
1467 poly[i][0] * edge_lambda_grads[m][q][d];
1472 if (flags & update_gradients)
1475 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1477 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1480 0.5 * edge_sigma_grads[m][d1] *
1481 edge_lambda_grads[m][q][d2];
1486 for (
unsigned int i = 0; i <
degree; ++i)
1488 const unsigned int dof_index(i + 1 + shift_m);
1490 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1492 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1496 edge_sigma_grads[m][d1] *
1497 edge_sigma_grads[m][d2] *
1499 edge_lambda_values[m][q] +
1500 (edge_sigma_grads[m][d1] *
1501 edge_lambda_grads[m][q][d2] +
1502 edge_sigma_grads[m][d2] *
1503 edge_lambda_grads[m][q][d1]) *
1505 edge_lambda_gradgrads_3d[m][d1]
1525 template <
int dim,
int spacedim>
1545 const unsigned int degree(
1556 const unsigned int n_q_points = quadrature.
size();
1562 Assert(!(flags & update_values) ||
1568 const unsigned int vertices_per_face(
1573 const unsigned int n_line_dofs =
1577 std::vector<std::vector<unsigned int>> face_orientation(
1578 faces_per_cell, std::vector<unsigned int>(vertices_per_face));
1588 {1, 2}, {0, 3}, {0, 3}, {1, 2}};
1590 for (
unsigned int m = 0; m < faces_per_cell; ++m)
1594 unsigned int current_max = 0;
1595 unsigned int current_glob = cell->vertex_index(
1597 for (
unsigned int v = 1; v < vertices_per_face; ++v)
1604 current_glob = cell->vertex_index(
1608 face_orientation[m][0] =
1613 m, vertex_opposite_on_face[current_max]);
1618 m, vertices_adjacent_on_face[current_max][0])) >
1620 m, vertices_adjacent_on_face[current_max][1])))
1622 face_orientation[m][1] =
1624 m, vertices_adjacent_on_face[current_max][0]);
1625 face_orientation[m][3] =
1627 m, vertices_adjacent_on_face[current_max][1]);
1631 face_orientation[m][1] =
1633 m, vertices_adjacent_on_face[current_max][1]);
1634 face_orientation[m][3] =
1636 m, vertices_adjacent_on_face[current_max][0]);
1642 std::vector<std::vector<double>> face_xi_values(
1643 faces_per_cell, std::vector<double>(n_q_points));
1644 std::vector<std::vector<double>> face_xi_grads(
1645 faces_per_cell, std::vector<double>(dim));
1646 std::vector<std::vector<double>> face_eta_values(
1647 faces_per_cell, std::vector<double>(n_q_points));
1648 std::vector<std::vector<double>> face_eta_grads(
1649 faces_per_cell, std::vector<double>(dim));
1651 std::vector<std::vector<double>> face_lambda_values(
1652 faces_per_cell, std::vector<double>(n_q_points));
1653 std::vector<std::vector<double>> face_lambda_grads(
1654 faces_per_cell, std::vector<double>(dim));
1655 for (
unsigned int m = 0; m < faces_per_cell; ++m)
1657 for (
unsigned int q = 0; q < n_q_points; ++q)
1659 face_xi_values[m][q] =
1661 [face_orientation[m][1]];
1662 face_eta_values[m][q] =
1664 [face_orientation[m][3]];
1667 for (
unsigned int d = 0; d < dim; ++d)
1669 face_xi_grads[m][d] =
1671 [face_orientation[m][1]][d];
1672 face_eta_grads[m][d] =
1674 [face_orientation[m][3]][d];
1681 std::vector<std::vector<double>> polyxi(
1682 degree, std::vector<double>(poly_length));
1683 std::vector<std::vector<double>> polyeta(
1684 degree, std::vector<double>(poly_length));
1687 for (
unsigned int m = 0; m < faces_per_cell; ++m)
1697 const unsigned int face_type1_offset(n_line_dofs + shift_m);
1707 const unsigned int face_type2_offset(face_type1_offset +
1719 const unsigned int face_type3_offset1(face_type2_offset +
1721 const unsigned int face_type3_offset2(face_type3_offset1 +
1725 for (
unsigned int q = 0; q < n_q_points; ++q)
1734 for (
unsigned int i = 0; i <
degree; ++i)
1738 face_xi_values[m][q], polyxi[i]);
1740 face_eta_values[m][q], polyeta[i]);
1743 if (flags & update_values)
1745 for (
unsigned int j = 0; j <
degree; ++j)
1747 const unsigned int shift_j(j * degree);
1748 for (
unsigned int i = 0; i <
degree; ++i)
1750 const unsigned int shift_ij(shift_j + i);
1752 const unsigned int dof_index1(face_type1_offset +
1754 for (
unsigned int d = 0; d < dim; ++d)
1757 (face_xi_grads[m][d] * polyxi[i][1] *
1759 face_eta_grads[m][d] * polyxi[i][0] *
1761 face_lambda_values[m][q] +
1762 face_lambda_grads[m][d] * polyxi[i][0] *
1766 const unsigned int dof_index2(face_type2_offset +
1768 for (
unsigned int d = 0; d < dim; ++d)
1771 (face_xi_grads[m][d] * polyxi[i][1] *
1773 face_eta_grads[m][d] * polyxi[i][0] *
1775 face_lambda_values[m][q];
1779 const unsigned int dof_index3_1(face_type3_offset1 +
1781 const unsigned int dof_index3_2(face_type3_offset2 +
1783 for (
unsigned int d = 0; d < dim; ++d)
1786 face_xi_grads[m][d] * polyeta[j][0] *
1787 face_lambda_values[m][q];
1790 face_eta_grads[m][d] * polyxi[j][0] *
1791 face_lambda_values[m][q];
1795 if (flags & update_gradients)
1797 for (
unsigned int j = 0; j <
degree; ++j)
1799 const unsigned int shift_j(j * degree);
1800 for (
unsigned int i = 0; i <
degree; ++i)
1802 const unsigned int shift_ij(shift_j + i);
1804 const unsigned int dof_index1(face_type1_offset +
1806 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1808 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1812 (face_xi_grads[m][d1] *
1813 face_xi_grads[m][d2] * polyxi[i][2] *
1815 (face_xi_grads[m][d1] *
1816 face_eta_grads[m][d2] +
1817 face_xi_grads[m][d2] *
1818 face_eta_grads[m][d1]) *
1819 polyxi[i][1] * polyeta[j][1] +
1820 face_eta_grads[m][d1] *
1821 face_eta_grads[m][d2] *
1822 polyxi[i][0] * polyeta[j][2]) *
1823 face_lambda_values[m][q] +
1824 (face_xi_grads[m][d2] * polyxi[i][1] *
1826 face_eta_grads[m][d2] * polyxi[i][0] *
1828 face_lambda_grads[m][d1] +
1829 (face_xi_grads[m][d1] * polyxi[i][1] *
1831 face_eta_grads[m][d1] * polyxi[i][0] *
1833 face_lambda_grads[m][d2];
1837 const unsigned int dof_index2(face_type2_offset +
1839 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1841 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1845 (face_xi_grads[m][d1] *
1846 face_xi_grads[m][d2] * polyxi[i][2] *
1848 (face_xi_grads[m][d1] *
1849 face_eta_grads[m][d2] -
1850 face_xi_grads[m][d2] *
1851 face_eta_grads[m][d1]) *
1852 polyxi[i][1] * polyeta[j][1] -
1853 face_eta_grads[m][d1] *
1854 face_eta_grads[m][d2] *
1855 polyxi[i][0] * polyeta[j][2]) *
1856 face_lambda_values[m][q] +
1857 (face_xi_grads[m][d1] * polyxi[i][1] *
1859 face_eta_grads[m][d1] * polyxi[i][0] *
1861 face_lambda_grads[m][d2];
1866 const unsigned int dof_index3_1(face_type3_offset1 +
1868 const unsigned int dof_index3_2(face_type3_offset2 +
1870 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1872 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1875 face_xi_grads[m][d1] *
1876 (face_eta_grads[m][d2] * polyeta[j][1] *
1877 face_lambda_values[m][q] +
1878 face_lambda_grads[m][d2] * polyeta[j][0]);
1881 face_eta_grads[m][d1] *
1882 (face_xi_grads[m][d2] * polyxi[j][1] *
1883 face_lambda_values[m][q] +
1884 face_lambda_grads[m][d2] * polyxi[j][0]);
1899 template <
int dim,
int spacedim>
1907 const ::internal::FEValuesImplementation::MappingRelatedData<dim, dim>
1914 Assert(dynamic_cast<const InternalData *>(&fe_internal) != 0,
1922 if (dim == 3 && this->
degree > 1)
1928 const unsigned int n_q_points = quadrature.
size();
1934 Assert(!(flags & update_values) ||
1938 if (flags & update_values)
1942 std::vector<Tensor<1, dim>> transformed_shape_values(n_q_points);
1943 for (
unsigned int dof = 0; dof < this->
dofs_per_cell; ++dof)
1945 const unsigned int first =
1953 make_array_view(transformed_shape_values));
1954 for (
unsigned int q = 0; q < n_q_points; ++q)
1956 for (
unsigned int d = 0; d < dim; ++d)
1959 transformed_shape_values[q][d];
1968 std::vector<Tensor<2, dim>> input(n_q_points);
1969 std::vector<Tensor<2, dim>> transformed_shape_grads(n_q_points);
1970 for (
unsigned int dof = 0; dof < this->
dofs_per_cell; ++dof)
1972 for (
unsigned int q = 0; q < n_q_points; ++q)
1976 mapping.
transform(make_array_view(input),
1979 make_array_view(transformed_shape_grads));
1981 const unsigned int first =
1986 for (
unsigned int q = 0; q < n_q_points; ++q)
1988 for (
unsigned int d1 = 0; d1 < dim; ++d1)
1990 for (
unsigned int d2 = 0; d2 < dim; ++d2)
1992 transformed_shape_grads[q][d1] -=
1994 mapping_data.jacobian_pushed_forward_grads[q][d2][d1];
1999 for (
unsigned int q = 0; q < n_q_points; ++q)
2001 for (
unsigned int d = 0; d < dim; ++d)
2004 transformed_shape_grads[q][d];
2011 template <
int dim,
int spacedim>
2015 const unsigned int face_no,
2019 const ::internal::FEValuesImplementation::MappingRelatedData<dim, dim>
2038 Assert(dynamic_cast<const InternalData *>(&fe_internal) != 0,
2048 if (dim == 3 && this->
degree > 1)
2056 const unsigned int n_q_points = quadrature.
size();
2059 cell->face_orientation(face_no),
2060 cell->face_flip(face_no),
2061 cell->face_rotation(face_no),
2068 std::vector<Tensor<1, dim>> transformed_shape_values(n_q_points);
2069 for (
unsigned int dof = 0; dof < this->
dofs_per_cell; ++dof)
2076 make_array_view(transformed_shape_values));
2078 const unsigned int first =
2083 for (
unsigned int q = 0; q < n_q_points; ++q)
2085 for (
unsigned int d = 0; d < dim; ++d)
2088 transformed_shape_values[q][d];
2097 std::vector<Tensor<2, dim>> transformed_shape_grads(n_q_points);
2098 for (
unsigned int dof = 0; dof < this->
dofs_per_cell; ++dof)
2105 make_array_view(transformed_shape_grads));
2107 const unsigned int first =
2112 for (
unsigned int q = 0; q < n_q_points; ++q)
2114 for (
unsigned int d1 = 0; d1 < dim; ++d1)
2116 for (
unsigned int d2 = 0; d2 < dim; ++d2)
2118 transformed_shape_grads[q][d1] -=
2120 mapping_data.jacobian_pushed_forward_grads[q][d2][d1];
2125 for (
unsigned int q = 0; q < n_q_points; ++q)
2127 for (
unsigned int d = 0; d < dim; ++d)
2130 transformed_shape_grads[q][d];
2137 template <
int dim,
int spacedim>
2141 const unsigned int ,
2142 const unsigned int ,
2146 const ::internal::FEValuesImplementation::MappingRelatedData<dim, dim>
2155 template <
int dim,
int spacedim>
2163 template <
int dim,
int spacedim>
2171 out |= update_values;
2176 template <
int dim,
int spacedim>
2186 out |= update_gradients | update_values |
2192 out |= update_hessians | update_values | update_gradients |
2200 template <
int dim,
int spacedim>
2211 std::ostringstream namebuf;
2212 namebuf <<
"FE_NedelecSZ<" << dim <<
">(" << this->
degree - 1 <<
")";
2214 return namebuf.str();
2217 template <
int dim,
int spacedim>
2218 std::unique_ptr<FiniteElement<dim, dim>>
2221 return std_cxx14::make_unique<FE_NedelecSZ<dim, spacedim>>(*this);
2224 template <
int dim,
int spacedim>
2225 std::vector<unsigned int>
2234 std::vector<unsigned int> dpo(dim + 1);
2236 dpo[1] = degree + 1;
2237 dpo[2] = 2 * degree * (degree + 1);
2240 dpo[3] = 3 * degree * degree * (degree + 1);
2245 template <
int dim,
int spacedim>
2254 return 2 * (degree + 1) * (degree + 2);
2257 return 3 * (degree + 1) * (degree + 2) * (degree + 2);
2267 template <
int dim,
int spacedim>
2277 #include "fe_nedelec_sz.inst" 2279 DEAL_II_NAMESPACE_CLOSE
unsigned int compute_num_dofs(const unsigned int degree) const
virtual void fill_fe_subface_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const unsigned int face_no, const unsigned int sub_no, const Quadrature< dim-1 > &quadrature, const Mapping< dim, dim > &mapping, const typename Mapping< dim, dim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, dim > &mapping_data, const typename FiniteElement< dim, dim >::InternalDataBase &fedata,::internal::FEValuesImplementation::FiniteElementRelatedData< dim, dim > &data) const override
static std::vector< Polynomials::Polynomial< double > > generate_complete_basis(const unsigned int degree)
std::vector< std::vector< double > > edge_lambda_grads_2d
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
static unsigned int face_to_cell_vertices(const unsigned int face, const unsigned int vertex, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false)
std::vector< std::vector< double > > face_lambda_values
static unsigned int line_to_cell_vertices(const unsigned int line, const unsigned int vertex)
std::vector< std::vector< double > > edge_sigma_values
const unsigned int dofs_per_quad
const unsigned int degree
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
virtual void fill_fe_face_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const unsigned int face_no, const Quadrature< dim-1 > &quadrature, const Mapping< dim, dim > &mapping, const typename Mapping< dim, dim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, dim > &mapping_data, const typename FiniteElement< dim, dim >::InternalDataBase &fedata,::internal::FEValuesImplementation::FiniteElementRelatedData< dim, dim > &data) const override
const unsigned int dofs_per_line
virtual void transform(const ArrayView< const Tensor< 1, dim >> &input, const MappingType type, const typename Mapping< dim, spacedim >::InternalDataBase &internal, const ArrayView< Tensor< 1, spacedim >> &output) const =0
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
virtual double shape_value(const unsigned int i, const Point< dim > &p) const override
void fill_edge_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const Quadrature< dim > &quadrature, const InternalData &fedata) const
static::ExceptionBase & ExcFENotPrimitive()
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const override
virtual std::unique_ptr< typename::FiniteElement< dim, spacedim >::InternalDataBase > get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature,::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const override
virtual std::string get_name() const override
std::vector< std::vector< DerivativeForm< 1, dim, dim > > > shape_grads
const std::vector< Point< dim > > & get_points() const
unsigned int size() const
static::ExceptionBase & ExcImpossibleInDim(int arg1)
std::vector< std::vector< std::vector< double > > > edge_lambda_gradgrads_3d
#define Assert(cond, exc)
static::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
UpdateFlags update_once(const UpdateFlags flags) const
void fill_face_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const Quadrature< dim > &quadrature, const InternalData &fedata) const
Abstract base class for mapping classes.
UpdateFlags update_each(const UpdateFlags flags) const
std::vector< Polynomials::Polynomial< double > > IntegratedLegendrePolynomials
const ComponentMask & get_nonzero_components(const unsigned int i) const
std::vector< std::vector< Tensor< 1, dim > > > shape_values
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
unsigned int n_components() const
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const override
Second derivatives of shape functions.
virtual std::unique_ptr< FiniteElement< dim, dim > > clone() const override
unsigned int first_selected_component(const unsigned int overall_number_of_components=numbers::invalid_unsigned_int) const
const unsigned int dofs_per_cell
virtual UpdateFlags requires_update_flags(const UpdateFlags update_flags) const override
void create_polynomials(const unsigned int degree)
std::vector< std::vector< double > > edge_lambda_values
virtual void fill_fe_values(const typename Triangulation< dim, dim >::cell_iterator &cell, const CellSimilarity::Similarity cell_similarity, const Quadrature< dim > &quadrature, const Mapping< dim, dim > &mapping, const typename Mapping< dim, dim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, dim > &mapping_data, const typename FiniteElement< dim, dim >::InternalDataBase &fedata,::internal::FEValuesImplementation::FiniteElementRelatedData< dim, dim > &data) const override
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
Shape function gradients.
std::vector< std::vector< std::vector< double > > > edge_lambda_grads_3d
FE_NedelecSZ(const unsigned int degree)
std::vector< std::vector< double > > edge_sigma_grads
std::vector< std::vector< std::vector< double > > > sigma_imj_values
static::ExceptionBase & ExcNotImplemented()
std::vector< std::vector< double > > face_lambda_grads
std::vector< std::vector< std::vector< double > > > sigma_imj_grads
Covariant transformation.
static::ExceptionBase & ExcInternalError()
static DataSetDescriptor face(const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points)