16 #ifndef dealii_mapping_q_internal_h
17 #define dealii_mapping_q_internal_h
63 template <
int spacedim>
77 template <
int spacedim>
89 const long double x = p(0);
90 const long double y = p(1);
92 const long double x0 =
vertices[0](0);
93 const long double x1 =
vertices[1](0);
94 const long double x2 =
vertices[2](0);
95 const long double x3 =
vertices[3](0);
97 const long double y0 =
vertices[0](1);
98 const long double y1 =
vertices[1](1);
99 const long double y2 =
vertices[2](1);
100 const long double y3 =
vertices[3](1);
102 const long double a = (x1 - x3) * (y0 - y2) - (x0 - x2) * (y1 - y3);
103 const long double b = -(x0 - x1 - x2 + x3) * y + (x - 2 * x1 + x3) * y0 -
104 (x - 2 * x0 + x2) * y1 - (x - x1) * y2 +
106 const long double c = (x0 - x1) * y - (x - x1) * y0 + (x - x0) * y1;
108 const long double discriminant =
b *
b - 4 * a * c;
117 const long double sqrt_discriminant = std::sqrt(discriminant);
120 if (
b != 0.0 && std::abs(
b) == sqrt_discriminant)
127 else if (std::abs(a) < 1
e-8 * std::abs(
b))
131 eta1 = 2 * c / (-
b - sqrt_discriminant);
132 eta2 = 2 * c / (-
b + sqrt_discriminant);
137 eta1 = (-
b - sqrt_discriminant) / (2 * a);
138 eta2 = (-
b + sqrt_discriminant) / (2 * a);
141 const long double eta =
142 (std::abs(eta1 - 0.5) < std::abs(eta2 - 0.5)) ? eta1 : eta2;
148 const long double subexpr0 = -eta * x2 + x0 * (eta - 1);
149 const long double xi_denominator0 = eta * x3 - x1 * (eta - 1) + subexpr0;
151 std::max(std::abs(x2), std::abs(x3)));
153 if (std::abs(xi_denominator0) > 1
e-10 * max_x)
155 const double xi = (x + subexpr0) / xi_denominator0;
156 return {xi,
static_cast<double>(eta)};
160 const long double max_y =
162 std::max(std::abs(y2), std::abs(y3)));
163 const long double subexpr1 = -eta * y2 + y0 * (eta - 1);
164 const long double xi_denominator1 =
165 eta * y3 - y1 * (eta - 1) + subexpr1;
166 if (std::abs(xi_denominator1) > 1
e-10 * max_y)
168 const double xi = (subexpr1 + y) / xi_denominator1;
169 return {xi,
static_cast<double>(eta)};
176 spacedim>::ExcTransformationFailed()));
182 return {std::numeric_limits<double>::quiet_NaN(),
183 std::numeric_limits<double>::quiet_NaN()};
188 template <
int spacedim>
197 return {std::numeric_limits<double>::quiet_NaN(),
198 std::numeric_limits<double>::quiet_NaN(),
199 std::numeric_limits<double>::quiet_NaN()};
210 namespace MappingQImplementation
217 std::vector<Point<dim>>
219 const std::vector<unsigned int> &renumbering)
223 std::vector<Point<dim>> points(renumbering.size());
224 const unsigned int n1 = line_support_points.size();
225 for (
unsigned int q2 = 0, q = 0; q2 < (dim > 2 ? n1 : 1); ++q2)
226 for (
unsigned int q1 = 0; q1 < (dim > 1 ? n1 : 1); ++q1)
227 for (
unsigned int q0 = 0; q0 < n1; ++q0, ++q)
229 points[renumbering[q]][0] = line_support_points[q0][0];
231 points[renumbering[q]][1] = line_support_points[q1][0];
233 points[renumbering[q]][2] = line_support_points[q2][0];
247 inline ::Table<2, double>
254 if (polynomial_degree == 1)
257 const unsigned int M = polynomial_degree - 1;
258 const unsigned int n_inner_2d = M * M;
259 const unsigned int n_outer_2d = 4 + 4 * M;
262 loqvs.reinit(n_inner_2d, n_outer_2d);
264 for (
unsigned int i = 0; i < M; ++i)
265 for (
unsigned int j = 0; j < M; ++j)
268 gl.
point((i + 1) * (polynomial_degree + 1) + (j + 1));
269 const unsigned int index_table = i * M + j;
270 for (
unsigned int v = 0; v < 4; ++v)
271 loqvs(index_table, v) =
273 loqvs(index_table, 4 + i) = 1. - p[0];
274 loqvs(index_table, 4 + i + M) = p[0];
275 loqvs(index_table, 4 + j + 2 * M) = 1. - p[1];
276 loqvs(index_table, 4 + j + 3 * M) = p[1];
281 for (
unsigned int unit_point = 0; unit_point < n_inner_2d; ++unit_point)
283 loqvs[unit_point].
end(),
285 1) < 1
e-13 * polynomial_degree,
299 inline ::Table<2, double>
306 if (polynomial_degree == 1)
309 const unsigned int M = polynomial_degree - 1;
311 const unsigned int n_inner = Utilities::fixed_power<3>(M);
312 const unsigned int n_outer = 8 + 12 * M + 6 * M * M;
315 lohvs.reinit(n_inner, n_outer);
317 for (
unsigned int i = 0; i < M; ++i)
318 for (
unsigned int j = 0; j < M; ++j)
319 for (
unsigned int k = 0; k < M; ++k)
322 (j + 1) * (M + 2) + (k + 1));
323 const unsigned int index_table = i * M * M + j * M + k;
326 for (
unsigned int v = 0; v < 8; ++v)
327 lohvs(index_table, v) =
332 constexpr std::array<unsigned int, 4> line_coordinates_y(
335 for (
unsigned int l = 0;
l < 4; ++
l)
336 lohvs(index_table, 8 + line_coordinates_y[
l] * M + j) =
341 constexpr std::array<unsigned int, 4> line_coordinates_x(
344 for (
unsigned int l = 0;
l < 4; ++
l)
345 lohvs(index_table, 8 + line_coordinates_x[
l] * M + k) =
350 constexpr std::array<unsigned int, 4> line_coordinates_z(
353 for (
unsigned int l = 0;
l < 4; ++
l)
354 lohvs(index_table, 8 + line_coordinates_z[
l] * M + i) =
359 lohvs(index_table, 8 + 12 * M + 0 * M * M + i * M + j) =
361 lohvs(index_table, 8 + 12 * M + 1 * M * M + i * M + j) = p[0];
362 lohvs(index_table, 8 + 12 * M + 2 * M * M + k * M + i) =
364 lohvs(index_table, 8 + 12 * M + 3 * M * M + k * M + i) = p[1];
365 lohvs(index_table, 8 + 12 * M + 4 * M * M + j * M + k) =
367 lohvs(index_table, 8 + 12 * M + 5 * M * M + j * M + k) = p[2];
372 for (
unsigned int unit_point = 0; unit_point < n_inner; ++unit_point)
374 lohvs[unit_point].
end(),
376 1) < 1
e-13 * polynomial_degree,
388 inline std::vector<::Table<2, double>>
390 const unsigned int polynomial_degree,
391 const unsigned int dim)
394 std::vector<::Table<2, double>> output(dim);
395 if (polynomial_degree <= 1)
400 output[0].reinit(polynomial_degree - 1,
402 for (
unsigned int q = 0; q < polynomial_degree - 1; ++q)
423 inline ::Table<2, double>
427 if (polynomial_degree <= 1)
428 return ::Table<2, double>();
431 const std::vector<unsigned int> h2l =
432 FETools::hierarchic_to_lexicographic_numbering<dim>(polynomial_degree);
437 for (
unsigned int q = 0; q < output.size(0); ++q)
454 template <
int dim,
int spacedim>
457 const typename ::MappingQ<dim, spacedim>::InternalData &data)
460 data.mapping_support_points.size());
464 for (
unsigned int i = 0; i < data.mapping_support_points.size(); ++i)
465 p_real += data.mapping_support_points[i] * data.shape(0, i);
476 template <
int dim,
int spacedim,
typename Number>
483 const std::vector<unsigned int> & renumber,
484 const bool print_iterations_to_deallog =
false)
486 if (print_iterations_to_deallog)
487 deallog <<
"Start MappingQ::do_transform_real_to_unit_cell for real "
488 <<
"point [ " << p <<
" ] " << std::endl;
505 polynomials_1d, points, p_unit, polynomials_1d.size() == 2, renumber);
514 f.
norm_square() - 1
e-24 * p_real.second[0].norm_square()) ==
552 const double eps = 1.e-11;
553 const unsigned int newton_iteration_limit = 20;
556 invalid_point[0] = std::numeric_limits<double>::infinity();
557 bool tried_project_to_unit_cell =
false;
559 unsigned int newton_iteration = 0;
560 Number f_weighted_norm_square = 1.;
561 Number last_f_weighted_norm_square = 1.;
565 if (print_iterations_to_deallog)
566 deallog <<
"Newton iteration " << newton_iteration
567 <<
" for unit point guess " << p_unit << std::endl;
571 for (
unsigned int d = 0;
d < spacedim; ++
d)
572 for (
unsigned int e = 0;
e < dim; ++
e)
573 df[
d][
e] = p_real.second[
e][
d];
586 if (tried_project_to_unit_cell ==
false)
593 polynomials_1d.size() == 2,
595 f = p_real.first - p;
596 f_weighted_norm_square = 1.;
597 last_f_weighted_norm_square = 1;
598 tried_project_to_unit_cell =
true;
602 return invalid_point;
610 if (print_iterations_to_deallog)
611 deallog <<
" delta=" << delta << std::endl;
614 double step_length = 1.0;
622 for (
unsigned int i = 0; i < dim; ++i)
623 p_unit_trial[i] -= step_length * delta[i];
626 const auto p_real_trial =
631 polynomials_1d.size() == 2,
634 p_real_trial.first - p;
635 f_weighted_norm_square = (df_inverse * f_trial).norm_square();
637 if (print_iterations_to_deallog)
639 deallog <<
" step_length=" << step_length << std::endl;
640 if (step_length == 1.0)
641 deallog <<
" ||f || =" << f.norm() << std::endl;
642 deallog <<
" ||f*|| =" << f_trial.
norm() << std::endl
644 << std::sqrt(f_weighted_norm_square) << std::endl;
664 if (
std::max(f_weighted_norm_square - 1
e-6 * 1
e-6, Number(0.)) *
669 p_real = p_real_trial;
670 p_unit = p_unit_trial;
674 else if (step_length > 0.05)
685 if (step_length <= 0.05 && tried_project_to_unit_cell ==
false)
692 polynomials_1d.size() == 2,
694 f = p_real.first - p;
695 f_weighted_norm_square = 1.;
696 last_f_weighted_norm_square = 1;
697 tried_project_to_unit_cell =
true;
700 else if (step_length <= 0.05)
701 return invalid_point;
704 if (newton_iteration > newton_iteration_limit)
705 return invalid_point;
714 std::max(last_f_weighted_norm_square -
715 std::min(f_weighted_norm_square, Number(1
e-6 * 1
e-6)) *
720 if (print_iterations_to_deallog)
721 deallog <<
"Iteration converged for p_unit = [ " << p_unit
722 <<
" ] and iteration error "
723 << std::sqrt(f_weighted_norm_square) << std::endl;
736 const typename ::Triangulation<dim, dim + 1>::cell_iterator &cell,
739 typename ::MappingQ<dim, dim + 1>::InternalData &mdata)
741 const unsigned int spacedim = dim + 1;
743 const unsigned int n_shapes = mdata.shape_values.size();
747 Assert(mdata.shape_second_derivatives.size() == n_shapes,
750 std::vector<Point<spacedim>> &points = mdata.mapping_support_points;
763 mdata.compute_shape_function_values(std::vector<
Point<dim>>(1, p_unit));
765 for (
unsigned int k = 0; k < mdata.n_shape_functions; ++k)
771 for (
unsigned int j = 0; j < dim; ++j)
773 DF[j] += grad_phi_k[j] * point_k;
774 for (
unsigned int l = 0;
l < dim; ++
l)
775 D2F[j][
l] += hessian_k[j][
l] * point_k;
780 p_minus_F -= compute_mapped_location_of_point<dim, spacedim>(mdata);
783 for (
unsigned int j = 0; j < dim; ++j)
784 f[j] = DF[j] * p_minus_F;
786 for (
unsigned int j = 0; j < dim; ++j)
788 f[j] = DF[j] * p_minus_F;
789 for (
unsigned int l = 0;
l < dim; ++
l)
790 df[j][
l] = -DF[j] * DF[
l] + D2F[j][
l] * p_minus_F;
794 const double eps = 1.e-12 * cell->diameter();
795 const unsigned int loop_limit = 10;
797 unsigned int loop = 0;
806 for (
unsigned int j = 0; j < dim; ++j)
809 for (
unsigned int l = 0;
l < dim; ++
l)
813 mdata.compute_shape_function_values(
816 for (
unsigned int k = 0; k < mdata.n_shape_functions; ++k)
822 for (
unsigned int j = 0; j < dim; ++j)
824 DF[j] += grad_phi_k[j] * point_k;
825 for (
unsigned int l = 0;
l < dim; ++
l)
826 D2F[j][
l] += hessian_k[j][
l] * point_k;
833 p_minus_F -= compute_mapped_location_of_point<dim, spacedim>(mdata);
835 for (
unsigned int j = 0; j < dim; ++j)
837 f[j] = DF[j] * p_minus_F;
838 for (
unsigned int l = 0;
l < dim; ++
l)
839 df[j][
l] = -DF[j] * DF[
l] + D2F[j][
l] * p_minus_F;
873 template <
int dim,
int spacedim>
881 (spacedim == 1 ? 3 : (spacedim == 2 ? 6 : 10));
899 1. / real_support_points[0].distance(real_support_points[1]))
912 Assert(dim == spacedim || real_support_points.size() ==
917 const auto affine = GridTools::affine_cell_approximation<dim>(
920 affine.first.covariant_form().transpose();
927 for (
unsigned int d = 0;
d < spacedim; ++
d)
928 for (
unsigned int e = 0;
e < dim; ++
e)
936 std::array<double, n_functions> shape_values;
942 shape_values[0] = 1.;
946 for (
unsigned int d = 0;
d < spacedim; ++
d)
947 shape_values[1 +
d] = p_scaled[
d];
948 for (
unsigned int d = 0, c = 0;
d < spacedim; ++
d)
949 for (
unsigned int e = 0;
e <=
d; ++
e, ++c)
950 shape_values[1 + spacedim + c] = p_scaled[
d] * p_scaled[
e];
959 matrix[i][j] += shape_values[i] * shape_values[j];
972 for (
unsigned int j = 0; j < i; ++j)
974 double Lik_Ljk_sum = 0;
975 for (
unsigned int k = 0; k < j; ++k)
981 ExcMessage(
"Matrix of normal equations not positive "
994 for (
unsigned int j = 0; j < i; ++j)
1004 for (
unsigned int j = i + 1; j <
n_functions; ++j)
1012 for (
unsigned int i = dim + 1; i <
n_functions; ++i)
1029 template <
typename Number>
1034 for (
unsigned int d = 0;
d < dim; ++
d)
1042 for (
unsigned int d = 0;
d < spacedim; ++
d)
1045 for (
unsigned int d = 0;
d < spacedim; ++
d)
1051 for (
unsigned int d = 0, c = 0;
d < spacedim; ++
d)
1052 for (
unsigned int e = 0;
e <=
d; ++
e, ++c)
1065 const Number affine_distance_to_unit_cell =
1068 for (
unsigned int d = 0;
d < dim; ++
d)
1069 result[
d] = compare_and_apply_mask<SIMDComparison::greater_than>(
1070 distance_to_unit_cell,
1071 affine_distance_to_unit_cell + 0.5,
1113 template <
int dim,
int spacedim>
1117 const typename ::MappingQ<dim, spacedim>::InternalData &data,
1121 const UpdateFlags update_flags = data.update_each;
1123 using VectorizedArrayType =
1124 typename ::MappingQ<dim,
1125 spacedim>::InternalData::VectorizedArrayType;
1126 const unsigned int n_shape_values = data.n_shape_functions;
1127 const unsigned int n_q_points = data.shape_info.n_q_points;
1128 constexpr
unsigned int n_lanes = VectorizedArrayType::size();
1129 constexpr
unsigned int n_comp = 1 + (spacedim - 1) / n_lanes;
1130 constexpr
unsigned int n_hessians = (dim * (dim + 1)) / 2;
1150 data.n_shape_functions > 0,
1153 n_q_points == data.contravariant.size(),
1156 n_q_points == jacobian_grads.size(),
1162 data.shape_info.element_type ==
1165 for (
unsigned int q = 0; q < n_q_points; ++q)
1167 data.mapping_support_points[data.shape_info
1168 .lexicographic_numbering[q]];
1181 for (
unsigned int i = 0; i < n_shape_values * n_comp; ++i)
1184 const std::vector<unsigned int> &renumber_to_lexicographic =
1185 data.shape_info.lexicographic_numbering;
1186 for (
unsigned int i = 0; i < n_shape_values; ++i)
1187 for (
unsigned int d = 0;
d < spacedim; ++
d)
1189 const unsigned int in_comp =
d % n_lanes;
1190 const unsigned int out_comp =
d / n_lanes;
1193 data.mapping_support_points[renumber_to_lexicographic[i]][
d];
1204 for (
unsigned int out_comp = 0; out_comp < n_comp; ++out_comp)
1205 for (
unsigned int i = 0; i < n_q_points; ++i)
1206 for (
unsigned int in_comp = 0;
1207 in_comp < n_lanes && in_comp < spacedim - out_comp * n_lanes;
1210 eval.
begin_values()[out_comp * n_q_points + i][in_comp];
1215 std::fill(data.contravariant.begin(),
1216 data.contravariant.end(),
1219 for (
unsigned int out_comp = 0; out_comp < n_comp; ++out_comp)
1221 for (
unsigned int j = 0; j < dim; ++j)
1222 for (
unsigned int in_comp = 0;
1223 in_comp < n_lanes &&
1224 in_comp < spacedim - out_comp * n_lanes;
1227 const unsigned int total_number =
point * dim + j;
1228 const unsigned int new_comp = total_number / n_q_points;
1229 const unsigned int new_point = total_number % n_q_points;
1230 data.contravariant[new_point][out_comp * n_lanes + in_comp]
1240 data.covariant[
point] =
1241 (data.contravariant[
point]).covariant_form();
1246 data.volume_elements[
point] =
1247 data.contravariant[
point].determinant();
1251 constexpr
int desymmetrize_3d[6][2] = {
1252 {0, 0}, {1, 1}, {2, 2}, {0, 1}, {0, 2}, {1, 2}};
1253 constexpr
int desymmetrize_2d[3][2] = {{0, 0}, {1, 1}, {0, 1}};
1256 for (
unsigned int out_comp = 0; out_comp < n_comp; ++out_comp)
1258 for (
unsigned int j = 0; j < n_hessians; ++j)
1259 for (
unsigned int in_comp = 0;
1260 in_comp < n_lanes &&
1261 in_comp < spacedim - out_comp * n_lanes;
1264 const unsigned int total_number =
point * n_hessians + j;
1265 const unsigned int new_point = total_number % n_q_points;
1266 const unsigned int new_hessian_comp =
1267 total_number / n_q_points;
1268 const unsigned int new_hessian_comp_i =
1269 dim == 2 ? desymmetrize_2d[new_hessian_comp][0] :
1270 desymmetrize_3d[new_hessian_comp][0];
1271 const unsigned int new_hessian_comp_j =
1272 dim == 2 ? desymmetrize_2d[new_hessian_comp][1] :
1273 desymmetrize_3d[new_hessian_comp][1];
1274 const double value =
1278 jacobian_grads[new_point][out_comp * n_lanes + in_comp]
1279 [new_hessian_comp_i][new_hessian_comp_j] =
1281 jacobian_grads[new_point][out_comp * n_lanes + in_comp]
1282 [new_hessian_comp_j][new_hessian_comp_i] =
1295 template <
int dim,
int spacedim>
1299 const typename ::MappingQ<dim, spacedim>::InternalData &data,
1302 const UpdateFlags update_flags = data.update_each;
1307 const double * shape = &data.shape(
point + data_set, 0);
1309 (shape[0] * data.mapping_support_points[0]);
1310 for (
unsigned int k = 1; k < data.n_shape_functions; ++k)
1311 for (
unsigned int i = 0; i < spacedim; ++i)
1312 result[i] += shape[k] * data.mapping_support_points[k][i];
1327 template <
int dim,
int spacedim>
1331 const typename ::QProjector<dim>::DataSetDescriptor data_set,
1332 const typename ::MappingQ<dim, spacedim>::InternalData &data)
1334 const UpdateFlags update_flags = data.update_each;
1342 const unsigned int n_q_points = data.contravariant.size();
1344 std::fill(data.contravariant.begin(),
1345 data.contravariant.end(),
1352 double result[spacedim][dim];
1356 for (
unsigned int i = 0; i < spacedim; ++i)
1357 for (
unsigned int j = 0; j < dim; ++j)
1358 result[i][j] = data.derivative(
point + data_set, 0)[j] *
1359 data.mapping_support_points[0][i];
1360 for (
unsigned int k = 1; k < data.n_shape_functions; ++k)
1361 for (
unsigned int i = 0; i < spacedim; ++i)
1362 for (
unsigned int j = 0; j < dim; ++j)
1363 result[i][j] += data.derivative(
point + data_set, k)[j] *
1364 data.mapping_support_points[k][i];
1371 for (
unsigned int i = 0; i < spacedim; ++i)
1372 for (
unsigned int j = 0; j < dim; ++j)
1373 data.contravariant[
point][i][j] = result[i][j];
1380 const unsigned int n_q_points = data.contravariant.size();
1383 data.covariant[
point] =
1384 (data.contravariant[
point]).covariant_form();
1391 const unsigned int n_q_points = data.contravariant.size();
1393 data.volume_elements[
point] =
1394 data.contravariant[
point].determinant();
1406 template <
int dim,
int spacedim>
1411 const typename ::MappingQ<dim, spacedim>::InternalData &data,
1414 const UpdateFlags update_flags = data.update_each;
1417 const unsigned int n_q_points = jacobian_grads.size();
1423 &data.second_derivative(
point + data_set, 0);
1424 double result[spacedim][dim][dim];
1425 for (
unsigned int i = 0; i < spacedim; ++i)
1426 for (
unsigned int j = 0; j < dim; ++j)
1427 for (
unsigned int l = 0;
l < dim; ++
l)
1429 (
second[0][j][
l] * data.mapping_support_points[0][i]);
1430 for (
unsigned int k = 1; k < data.n_shape_functions; ++k)
1431 for (
unsigned int i = 0; i < spacedim; ++i)
1432 for (
unsigned int j = 0; j < dim; ++j)
1433 for (
unsigned int l = 0;
l < dim; ++
l)
1435 (
second[k][j][
l] * data.mapping_support_points[k][i]);
1437 for (
unsigned int i = 0; i < spacedim; ++i)
1438 for (
unsigned int j = 0; j < dim; ++j)
1439 for (
unsigned int l = 0;
l < dim; ++
l)
1440 jacobian_grads[
point][i][j][
l] = result[i][j][
l];
1453 template <
int dim,
int spacedim>
1458 const typename ::MappingQ<dim, spacedim>::InternalData &data,
1461 const UpdateFlags update_flags = data.update_each;
1464 const unsigned int n_q_points = jacobian_pushed_forward_grads.size();
1468 double tmp[spacedim][spacedim][spacedim];
1472 &data.second_derivative(
point + data_set, 0);
1473 double result[spacedim][dim][dim];
1474 for (
unsigned int i = 0; i < spacedim; ++i)
1475 for (
unsigned int j = 0; j < dim; ++j)
1476 for (
unsigned int l = 0;
l < dim; ++
l)
1478 (
second[0][j][
l] * data.mapping_support_points[0][i]);
1479 for (
unsigned int k = 1; k < data.n_shape_functions; ++k)
1480 for (
unsigned int i = 0; i < spacedim; ++i)
1481 for (
unsigned int j = 0; j < dim; ++j)
1482 for (
unsigned int l = 0;
l < dim; ++
l)
1485 data.mapping_support_points[k][i]);
1488 for (
unsigned int i = 0; i < spacedim; ++i)
1489 for (
unsigned int j = 0; j < spacedim; ++j)
1490 for (
unsigned int l = 0;
l < dim; ++
l)
1493 result[i][0][
l] * data.covariant[
point][j][0];
1494 for (
unsigned int jr = 1; jr < dim; ++jr)
1497 result[i][jr][
l] * data.covariant[
point][j][jr];
1502 for (
unsigned int i = 0; i < spacedim; ++i)
1503 for (
unsigned int j = 0; j < spacedim; ++j)
1504 for (
unsigned int l = 0;
l < spacedim; ++
l)
1506 jacobian_pushed_forward_grads[
point][i][j][
l] =
1507 tmp[i][j][0] * data.covariant[
point][
l][0];
1508 for (
unsigned int lr = 1; lr < dim; ++lr)
1510 jacobian_pushed_forward_grads[
point][i][j][
l] +=
1511 tmp[i][j][lr] * data.covariant[
point][
l][lr];
1527 template <
int dim,
int spacedim>
1532 const typename ::MappingQ<dim, spacedim>::InternalData &data,
1535 const UpdateFlags update_flags = data.update_each;
1538 const unsigned int n_q_points = jacobian_2nd_derivatives.size();
1545 &data.third_derivative(
point + data_set, 0);
1546 double result[spacedim][dim][dim][dim];
1547 for (
unsigned int i = 0; i < spacedim; ++i)
1548 for (
unsigned int j = 0; j < dim; ++j)
1549 for (
unsigned int l = 0;
l < dim; ++
l)
1550 for (
unsigned int m = 0; m < dim; ++m)
1551 result[i][j][
l][m] =
1552 (third[0][j][
l][m] *
1553 data.mapping_support_points[0][i]);
1554 for (
unsigned int k = 1; k < data.n_shape_functions; ++k)
1555 for (
unsigned int i = 0; i < spacedim; ++i)
1556 for (
unsigned int j = 0; j < dim; ++j)
1557 for (
unsigned int l = 0;
l < dim; ++
l)
1558 for (
unsigned int m = 0; m < dim; ++m)
1559 result[i][j][
l][m] +=
1560 (third[k][j][
l][m] *
1561 data.mapping_support_points[k][i]);
1563 for (
unsigned int i = 0; i < spacedim; ++i)
1564 for (
unsigned int j = 0; j < dim; ++j)
1565 for (
unsigned int l = 0;
l < dim; ++
l)
1566 for (
unsigned int m = 0; m < dim; ++m)
1567 jacobian_2nd_derivatives[
point][i][j][
l][m] =
1583 template <
int dim,
int spacedim>
1588 const typename ::MappingQ<dim, spacedim>::InternalData &data,
1591 const UpdateFlags update_flags = data.update_each;
1594 const unsigned int n_q_points =
1595 jacobian_pushed_forward_2nd_derivatives.size();
1599 double tmp[spacedim][spacedim][spacedim][spacedim];
1603 &data.third_derivative(
point + data_set, 0);
1604 double result[spacedim][dim][dim][dim];
1605 for (
unsigned int i = 0; i < spacedim; ++i)
1606 for (
unsigned int j = 0; j < dim; ++j)
1607 for (
unsigned int l = 0;
l < dim; ++
l)
1608 for (
unsigned int m = 0; m < dim; ++m)
1609 result[i][j][
l][m] =
1610 (third[0][j][
l][m] *
1611 data.mapping_support_points[0][i]);
1612 for (
unsigned int k = 1; k < data.n_shape_functions; ++k)
1613 for (
unsigned int i = 0; i < spacedim; ++i)
1614 for (
unsigned int j = 0; j < dim; ++j)
1615 for (
unsigned int l = 0;
l < dim; ++
l)
1616 for (
unsigned int m = 0; m < dim; ++m)
1617 result[i][j][
l][m] +=
1618 (third[k][j][
l][m] *
1619 data.mapping_support_points[k][i]);
1622 for (
unsigned int i = 0; i < spacedim; ++i)
1623 for (
unsigned int j = 0; j < spacedim; ++j)
1624 for (
unsigned int l = 0;
l < dim; ++
l)
1625 for (
unsigned int m = 0; m < dim; ++m)
1627 jacobian_pushed_forward_2nd_derivatives
1628 [
point][i][j][
l][m] = result[i][0][
l][m] *
1629 data.covariant[
point][j][0];
1630 for (
unsigned int jr = 1; jr < dim; ++jr)
1631 jacobian_pushed_forward_2nd_derivatives[
point][i]
1634 result[i][jr][
l][m] *
1635 data.covariant[
point][j][jr];
1639 for (
unsigned int i = 0; i < spacedim; ++i)
1640 for (
unsigned int j = 0; j < spacedim; ++j)
1641 for (
unsigned int l = 0;
l < spacedim; ++
l)
1642 for (
unsigned int m = 0; m < dim; ++m)
1645 jacobian_pushed_forward_2nd_derivatives[
point][i]
1647 data.covariant[
point][
l][0];
1648 for (
unsigned int lr = 1; lr < dim; ++lr)
1650 jacobian_pushed_forward_2nd_derivatives[
point]
1653 data.covariant[
point][
l][lr];
1657 for (
unsigned int i = 0; i < spacedim; ++i)
1658 for (
unsigned int j = 0; j < spacedim; ++j)
1659 for (
unsigned int l = 0;
l < spacedim; ++
l)
1660 for (
unsigned int m = 0; m < spacedim; ++m)
1662 jacobian_pushed_forward_2nd_derivatives
1664 tmp[i][j][
l][0] * data.covariant[
point][m][0];
1665 for (
unsigned int mr = 1; mr < dim; ++mr)
1666 jacobian_pushed_forward_2nd_derivatives[
point][i]
1669 tmp[i][j][
l][mr] * data.covariant[
point][m][mr];
1684 template <
int dim,
int spacedim>
1689 const typename ::MappingQ<dim, spacedim>::InternalData &data,
1692 const UpdateFlags update_flags = data.update_each;
1695 const unsigned int n_q_points = jacobian_3rd_derivatives.size();
1702 &data.fourth_derivative(
point + data_set, 0);
1703 double result[spacedim][dim][dim][dim][dim];
1704 for (
unsigned int i = 0; i < spacedim; ++i)
1705 for (
unsigned int j = 0; j < dim; ++j)
1706 for (
unsigned int l = 0;
l < dim; ++
l)
1707 for (
unsigned int m = 0; m < dim; ++m)
1708 for (
unsigned int n = 0; n < dim; ++n)
1709 result[i][j][
l][m][n] =
1710 (fourth[0][j][
l][m][n] *
1711 data.mapping_support_points[0][i]);
1712 for (
unsigned int k = 1; k < data.n_shape_functions; ++k)
1713 for (
unsigned int i = 0; i < spacedim; ++i)
1714 for (
unsigned int j = 0; j < dim; ++j)
1715 for (
unsigned int l = 0;
l < dim; ++
l)
1716 for (
unsigned int m = 0; m < dim; ++m)
1717 for (
unsigned int n = 0; n < dim; ++n)
1718 result[i][j][
l][m][n] +=
1719 (fourth[k][j][
l][m][n] *
1720 data.mapping_support_points[k][i]);
1722 for (
unsigned int i = 0; i < spacedim; ++i)
1723 for (
unsigned int j = 0; j < dim; ++j)
1724 for (
unsigned int l = 0;
l < dim; ++
l)
1725 for (
unsigned int m = 0; m < dim; ++m)
1726 for (
unsigned int n = 0; n < dim; ++n)
1727 jacobian_3rd_derivatives[
point][i][j][
l][m][n] =
1728 result[i][j][
l][m][n];
1743 template <
int dim,
int spacedim>
1748 const typename ::MappingQ<dim, spacedim>::InternalData &data,
1751 const UpdateFlags update_flags = data.update_each;
1754 const unsigned int n_q_points =
1755 jacobian_pushed_forward_3rd_derivatives.size();
1759 double tmp[spacedim][spacedim][spacedim][spacedim][spacedim];
1763 &data.fourth_derivative(
point + data_set, 0);
1764 double result[spacedim][dim][dim][dim][dim];
1765 for (
unsigned int i = 0; i < spacedim; ++i)
1766 for (
unsigned int j = 0; j < dim; ++j)
1767 for (
unsigned int l = 0;
l < dim; ++
l)
1768 for (
unsigned int m = 0; m < dim; ++m)
1769 for (
unsigned int n = 0; n < dim; ++n)
1770 result[i][j][
l][m][n] =
1771 (fourth[0][j][
l][m][n] *
1772 data.mapping_support_points[0][i]);
1773 for (
unsigned int k = 1; k < data.n_shape_functions; ++k)
1774 for (
unsigned int i = 0; i < spacedim; ++i)
1775 for (
unsigned int j = 0; j < dim; ++j)
1776 for (
unsigned int l = 0;
l < dim; ++
l)
1777 for (
unsigned int m = 0; m < dim; ++m)
1778 for (
unsigned int n = 0; n < dim; ++n)
1779 result[i][j][
l][m][n] +=
1780 (fourth[k][j][
l][m][n] *
1781 data.mapping_support_points[k][i]);
1784 for (
unsigned int i = 0; i < spacedim; ++i)
1785 for (
unsigned int j = 0; j < spacedim; ++j)
1786 for (
unsigned int l = 0;
l < dim; ++
l)
1787 for (
unsigned int m = 0; m < dim; ++m)
1788 for (
unsigned int n = 0; n < dim; ++n)
1790 tmp[i][j][
l][m][n] = result[i][0][
l][m][n] *
1791 data.covariant[
point][j][0];
1792 for (
unsigned int jr = 1; jr < dim; ++jr)
1793 tmp[i][j][
l][m][n] +=
1794 result[i][jr][
l][m][n] *
1795 data.covariant[
point][j][jr];
1799 for (
unsigned int i = 0; i < spacedim; ++i)
1800 for (
unsigned int j = 0; j < spacedim; ++j)
1801 for (
unsigned int l = 0;
l < spacedim; ++
l)
1802 for (
unsigned int m = 0; m < dim; ++m)
1803 for (
unsigned int n = 0; n < dim; ++n)
1805 jacobian_pushed_forward_3rd_derivatives
1807 tmp[i][j][0][m][n] *
1808 data.covariant[
point][
l][0];
1809 for (
unsigned int lr = 1; lr < dim; ++lr)
1810 jacobian_pushed_forward_3rd_derivatives[
point]
1813 tmp[i][j][lr][m][n] *
1814 data.covariant[
point][
l][lr];
1818 for (
unsigned int i = 0; i < spacedim; ++i)
1819 for (
unsigned int j = 0; j < spacedim; ++j)
1820 for (
unsigned int l = 0;
l < spacedim; ++
l)
1821 for (
unsigned int m = 0; m < spacedim; ++m)
1822 for (
unsigned int n = 0; n < dim; ++n)
1824 tmp[i][j][
l][m][n] =
1825 jacobian_pushed_forward_3rd_derivatives[
point]
1828 data.covariant[
point][m][0];
1829 for (
unsigned int mr = 1; mr < dim; ++mr)
1830 tmp[i][j][
l][m][n] +=
1831 jacobian_pushed_forward_3rd_derivatives
1832 [
point][i][j][
l][mr][n] *
1833 data.covariant[
point][m][mr];
1837 for (
unsigned int i = 0; i < spacedim; ++i)
1838 for (
unsigned int j = 0; j < spacedim; ++j)
1839 for (
unsigned int l = 0;
l < spacedim; ++
l)
1840 for (
unsigned int m = 0; m < spacedim; ++m)
1841 for (
unsigned int n = 0; n < spacedim; ++n)
1843 jacobian_pushed_forward_3rd_derivatives
1845 tmp[i][j][
l][m][0] *
1846 data.covariant[
point][n][0];
1847 for (
unsigned int nr = 1; nr < dim; ++nr)
1848 jacobian_pushed_forward_3rd_derivatives[
point]
1851 tmp[i][j][
l][m][nr] *
1852 data.covariant[
point][n][nr];
1870 template <
int dim,
int spacedim>
1873 const ::MappingQ<dim, spacedim> &mapping,
1874 const typename ::Triangulation<dim, spacedim>::cell_iterator &cell,
1875 const unsigned int face_no,
1876 const unsigned int subface_no,
1877 const unsigned int n_q_points,
1878 const std::vector<double> & weights,
1879 const typename ::MappingQ<dim, spacedim>::InternalData &data,
1883 const UpdateFlags update_flags = data.update_each;
1904 for (
unsigned int d = 0;
d != dim - 1; ++
d)
1907 data.unit_tangentials.size(),
1910 data.aux[
d].size() <=
1912 .unit_tangentials[face_no +
1919 data.unit_tangentials[face_no +
1930 if (dim == spacedim)
1932 for (
unsigned int i = 0; i < n_q_points; ++i)
1942 (face_no == 0 ? -1 : +1);
1946 cross_product_2d(data.aux[0][i]);
1950 cross_product_3d(data.aux[0][i], data.aux[1][i]);
1972 data.contravariant[
point].transpose()[0];
1974 (face_no == 0 ? -1. : +1.) *
1984 cross_product_3d(DX_t[0], DX_t[1]);
1985 cell_normal /= cell_normal.
norm();
1990 cross_product_3d(data.aux[0][
point], cell_normal);
1997 for (
unsigned int i = 0; i < output_data.
boundary_forms.size(); ++i)
2005 cell->subface_case(face_no), subface_no);
2011 for (
unsigned int i = 0; i < output_data.
normal_vectors.size(); ++i)
2023 data.covariant[
point].transpose();
2034 template <
int dim,
int spacedim>
2037 const ::MappingQ<dim, spacedim> &mapping,
2038 const typename ::Triangulation<dim, spacedim>::cell_iterator &cell,
2039 const unsigned int face_no,
2040 const unsigned int subface_no,
2043 const typename ::MappingQ<dim, spacedim>::InternalData &data,
2047 if (dim > 1 && data.tensor_product_quadrature)
2049 maybe_update_q_points_Jacobians_and_grads_tensor<dim, spacedim>(
2057 maybe_compute_q_points<dim, spacedim>(data_set,
2063 maybe_update_jacobian_grads<dim, spacedim>(
2066 maybe_update_jacobian_pushed_forward_grads<dim, spacedim>(
2071 maybe_update_jacobian_2nd_derivatives<dim, spacedim>(
2076 maybe_update_jacobian_pushed_forward_2nd_derivatives<dim, spacedim>(
2081 maybe_update_jacobian_3rd_derivatives<dim, spacedim>(
2086 maybe_update_jacobian_pushed_forward_3rd_derivatives<dim, spacedim>(
2107 template <
int dim,
int spacedim,
int rank>
2117 const typename ::MappingQ<dim, spacedim>::InternalData *
>(
2118 &mapping_data) !=
nullptr),
2120 const typename ::MappingQ<dim, spacedim>::InternalData &data =
2122 const typename ::MappingQ<dim, spacedim>::InternalData &
>(
2125 switch (mapping_kind)
2131 "update_contravariant_transformation"));
2133 for (
unsigned int i = 0; i < output.size(); ++i)
2144 "update_contravariant_transformation"));
2147 "update_volume_elements"));
2152 for (
unsigned int i = 0; i < output.size(); ++i)
2156 output[i] /= data.volume_elements[i];
2167 "update_covariant_transformation"));
2169 for (
unsigned int i = 0; i < output.size(); ++i)
2185 template <
int dim,
int spacedim,
int rank>
2195 const typename ::MappingQ<dim, spacedim>::InternalData *
>(
2196 &mapping_data) !=
nullptr),
2198 const typename ::MappingQ<dim, spacedim>::InternalData &data =
2200 const typename ::MappingQ<dim, spacedim>::InternalData &
>(
2203 switch (mapping_kind)
2209 "update_covariant_transformation"));
2212 "update_contravariant_transformation"));
2215 for (
unsigned int i = 0; i < output.size(); ++i)
2231 "update_covariant_transformation"));
2234 for (
unsigned int i = 0; i < output.size(); ++i)
2250 "update_covariant_transformation"));
2253 "update_contravariant_transformation"));
2256 "update_volume_elements"));
2259 for (
unsigned int i = 0; i < output.size(); ++i)
2267 output[i] /= data.volume_elements[i];
2283 template <
int dim,
int spacedim>
2293 const typename ::MappingQ<dim, spacedim>::InternalData *
>(
2294 &mapping_data) !=
nullptr),
2296 const typename ::MappingQ<dim, spacedim>::InternalData &data =
2298 const typename ::MappingQ<dim, spacedim>::InternalData &
>(
2301 switch (mapping_kind)
2307 "update_covariant_transformation"));
2310 "update_contravariant_transformation"));
2312 for (
unsigned int q = 0; q < output.size(); ++q)
2313 for (
unsigned int i = 0; i < spacedim; ++i)
2315 double tmp1[dim][dim];
2316 for (
unsigned int J = 0; J < dim; ++J)
2317 for (
unsigned int K = 0;
K < dim; ++
K)
2320 data.contravariant[q][i][0] * input[q][0][J][
K];
2321 for (
unsigned int I = 1; I < dim; ++I)
2323 data.contravariant[q][i][I] * input[q][I][J][
K];
2325 for (
unsigned int j = 0; j < spacedim; ++j)
2328 for (
unsigned int K = 0;
K < dim; ++
K)
2330 tmp2[
K] = data.covariant[q][j][0] * tmp1[0][
K];
2331 for (
unsigned int J = 1; J < dim; ++J)
2332 tmp2[
K] += data.covariant[q][j][J] * tmp1[J][
K];
2334 for (
unsigned int k = 0; k < spacedim; ++k)
2336 output[q][i][j][k] =
2337 data.covariant[q][k][0] * tmp2[0];
2338 for (
unsigned int K = 1;
K < dim; ++
K)
2339 output[q][i][j][k] +=
2340 data.covariant[q][k][
K] * tmp2[
K];
2351 "update_covariant_transformation"));
2353 for (
unsigned int q = 0; q < output.size(); ++q)
2354 for (
unsigned int i = 0; i < spacedim; ++i)
2356 double tmp1[dim][dim];
2357 for (
unsigned int J = 0; J < dim; ++J)
2358 for (
unsigned int K = 0;
K < dim; ++
K)
2361 data.covariant[q][i][0] * input[q][0][J][
K];
2362 for (
unsigned int I = 1; I < dim; ++I)
2364 data.covariant[q][i][I] * input[q][I][J][
K];
2366 for (
unsigned int j = 0; j < spacedim; ++j)
2369 for (
unsigned int K = 0;
K < dim; ++
K)
2371 tmp2[
K] = data.covariant[q][j][0] * tmp1[0][
K];
2372 for (
unsigned int J = 1; J < dim; ++J)
2373 tmp2[
K] += data.covariant[q][j][J] * tmp1[J][
K];
2375 for (
unsigned int k = 0; k < spacedim; ++k)
2377 output[q][i][j][k] =
2378 data.covariant[q][k][0] * tmp2[0];
2379 for (
unsigned int K = 1;
K < dim; ++
K)
2380 output[q][i][j][k] +=
2381 data.covariant[q][k][
K] * tmp2[
K];
2393 "update_covariant_transformation"));
2396 "update_contravariant_transformation"));
2399 "update_volume_elements"));
2401 for (
unsigned int q = 0; q < output.size(); ++q)
2402 for (
unsigned int i = 0; i < spacedim; ++i)
2405 for (
unsigned int I = 0; I < dim; ++I)
2407 data.contravariant[q][i][I] / data.volume_elements[q];
2408 double tmp1[dim][dim];
2409 for (
unsigned int J = 0; J < dim; ++J)
2410 for (
unsigned int K = 0;
K < dim; ++
K)
2412 tmp1[J][
K] = factor[0] * input[q][0][J][
K];
2413 for (
unsigned int I = 1; I < dim; ++I)
2414 tmp1[J][
K] += factor[I] * input[q][I][J][
K];
2416 for (
unsigned int j = 0; j < spacedim; ++j)
2419 for (
unsigned int K = 0;
K < dim; ++
K)
2421 tmp2[
K] = data.covariant[q][j][0] * tmp1[0][
K];
2422 for (
unsigned int J = 1; J < dim; ++J)
2423 tmp2[
K] += data.covariant[q][j][J] * tmp1[J][
K];
2425 for (
unsigned int k = 0; k < spacedim; ++k)
2427 output[q][i][j][k] =
2428 data.covariant[q][k][0] * tmp2[0];
2429 for (
unsigned int K = 1;
K < dim; ++
K)
2430 output[q][i][j][k] +=
2431 data.covariant[q][k][
K] * tmp2[
K];
2450 template <
int dim,
int spacedim,
int rank>
2460 const typename ::MappingQ<dim, spacedim>::InternalData *
>(
2461 &mapping_data) !=
nullptr),
2463 const typename ::MappingQ<dim, spacedim>::InternalData &data =
2465 const typename ::MappingQ<dim, spacedim>::InternalData &
>(
2468 switch (mapping_kind)
2474 "update_covariant_transformation"));
2476 for (
unsigned int i = 0; i < output.size(); ++i)
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
void set_data_pointers(AlignedVector< Number > *scratch_data, const unsigned int n_components)
const Number * begin_values() const
const Number * begin_hessians() const
const Number * begin_gradients() const
const Number * begin_dof_values() const
Abstract base class for mapping classes.
numbers::NumberTraits< Number >::real_type distance_square(const Point< dim, Number > &p) const
const std::vector< double > & get_weights() const
const Point< dim > & point(const unsigned int i) const
unsigned int size() const
constexpr numbers::NumberTraits< Number >::real_type norm_square() const
numbers::NumberTraits< Number >::real_type norm() const
const Point< spacedim > normalization_shift
const double normalization_length
Point< dim, Number > compute(const Point< spacedim, Number > &p) const
InverseQuadraticApproximation(const InverseQuadraticApproximation &)=default
static constexpr unsigned int n_functions
InverseQuadraticApproximation(const std::vector< Point< spacedim >> &real_support_points, const std::vector< Point< dim >> &unit_support_points)
std::array< Point< dim >, n_functions > coefficients
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
static ::ExceptionBase & ExcNotImplemented()
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
void loop(ITERATOR begin, std_cxx20::type_identity_t< ITERATOR > end, DOFINFO &dinfo, INFOBOX &info, const std::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &cell_worker, const std::function< void(DOFINFO &, typename INFOBOX::CellInfo &)> &boundary_worker, const std::function< void(DOFINFO &, DOFINFO &, typename INFOBOX::CellInfo &, typename INFOBOX::CellInfo &)> &face_worker, ASSEMBLER &assembler, const LoopControl &lctrl=LoopControl())
@ update_jacobian_pushed_forward_2nd_derivatives
@ update_volume_elements
Determinant of the Jacobian.
@ update_contravariant_transformation
Contravariant transformation.
@ update_jacobian_pushed_forward_grads
@ update_jacobian_3rd_derivatives
@ update_jacobian_grads
Gradient of volume element.
@ update_normal_vectors
Normal vectors.
@ update_JxW_values
Transformed quadrature weights.
@ update_covariant_transformation
Covariant transformation.
@ update_jacobians
Volume element.
@ update_inverse_jacobians
Volume element.
@ update_quadrature_points
Transformed quadrature points.
@ update_jacobian_pushed_forward_3rd_derivatives
@ update_boundary_forms
Outer normal vector, not normalized.
@ update_jacobian_2nd_derivatives
@ mapping_covariant_gradient
@ mapping_contravariant_hessian
@ mapping_covariant_hessian
@ mapping_contravariant_gradient
@ tensor_symmetric_collocation
CGAL::Exact_predicates_exact_constructions_kernel_with_sqrt K
Expression fabs(const Expression &x)
EvaluationFlags
The EvaluationFlags enum.
@ matrix
Contents is actually a matrix.
Point< spacedim > point(const gp_Pnt &p, const double tolerance=1e-10)
void quadrature_points(const Triangulation< dim, spacedim > &triangulation, const Quadrature< dim > &quadrature, const std::vector< std::vector< BoundingBox< spacedim >>> &global_bounding_boxes, ParticleHandler< dim, spacedim > &particle_handler, const Mapping< dim, spacedim > &mapping=(ReferenceCells::get_hypercube< dim >() .template get_default_linear_mapping< dim, spacedim >()), const std::vector< std::vector< double >> &properties={})
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
VectorType::value_type * begin(VectorType &V)
VectorType::value_type * end(VectorType &V)
T sum(const T &t, const MPI_Comm &mpi_communicator)
constexpr T pow(const T base, const int iexp)
Point< 1 > transform_real_to_unit_cell(const std::array< Point< spacedim >, GeometryInfo< 1 >::vertices_per_cell > &vertices, const Point< spacedim > &p)
void transform_gradients(const ArrayView< const Tensor< rank, dim >> &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< rank, spacedim >> &output)
std::vector< Point< dim > > unit_support_points(const std::vector< Point< 1 >> &line_support_points, const std::vector< unsigned int > &renumbering)
inline ::Table< 2, double > compute_support_point_weights_on_hex(const unsigned int polynomial_degree)
void transform_differential_forms(const ArrayView< const DerivativeForm< rank, dim, spacedim >> &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< rank+1, spacedim >> &output)
void maybe_update_jacobian_pushed_forward_3rd_derivatives(const CellSimilarity::Similarity cell_similarity, const typename QProjector< dim >::DataSetDescriptor data_set, const typename ::MappingQ< dim, spacedim >::InternalData &data, std::vector< Tensor< 5, spacedim >> &jacobian_pushed_forward_3rd_derivatives)
Point< spacedim > compute_mapped_location_of_point(const typename ::MappingQ< dim, spacedim >::InternalData &data)
void maybe_update_q_points_Jacobians_and_grads_tensor(const CellSimilarity::Similarity cell_similarity, const typename ::MappingQ< dim, spacedim >::InternalData &data, std::vector< Point< spacedim >> &quadrature_points, std::vector< DerivativeForm< 2, dim, spacedim >> &jacobian_grads)
Point< dim > do_transform_real_to_unit_cell_internal_codim1(const typename ::Triangulation< dim, dim+1 >::cell_iterator &cell, const Point< dim+1 > &p, const Point< dim > &initial_p_unit, typename ::MappingQ< dim, dim+1 >::InternalData &mdata)
void maybe_update_jacobian_pushed_forward_2nd_derivatives(const CellSimilarity::Similarity cell_similarity, const typename QProjector< dim >::DataSetDescriptor data_set, const typename ::MappingQ< dim, spacedim >::InternalData &data, std::vector< Tensor< 4, spacedim >> &jacobian_pushed_forward_2nd_derivatives)
void maybe_compute_q_points(const typename QProjector< dim >::DataSetDescriptor data_set, const typename ::MappingQ< dim, spacedim >::InternalData &data, std::vector< Point< spacedim >> &quadrature_points)
void maybe_update_Jacobians(const CellSimilarity::Similarity cell_similarity, const typename ::QProjector< dim >::DataSetDescriptor data_set, const typename ::MappingQ< dim, spacedim >::InternalData &data)
void do_fill_fe_face_values(const ::MappingQ< dim, spacedim > &mapping, const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const typename QProjector< dim >::DataSetDescriptor data_set, const Quadrature< dim - 1 > &quadrature, const typename ::MappingQ< dim, spacedim >::InternalData &data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data)
void maybe_update_jacobian_grads(const CellSimilarity::Similarity cell_similarity, const typename QProjector< dim >::DataSetDescriptor data_set, const typename ::MappingQ< dim, spacedim >::InternalData &data, std::vector< DerivativeForm< 2, dim, spacedim >> &jacobian_grads)
void maybe_update_jacobian_3rd_derivatives(const CellSimilarity::Similarity cell_similarity, const typename QProjector< dim >::DataSetDescriptor data_set, const typename ::MappingQ< dim, spacedim >::InternalData &data, std::vector< DerivativeForm< 4, dim, spacedim >> &jacobian_3rd_derivatives)
Point< dim, Number > do_transform_real_to_unit_cell_internal(const Point< spacedim, Number > &p, const Point< dim, Number > &initial_p_unit, const std::vector< Point< spacedim >> &points, const std::vector< Polynomials::Polynomial< double >> &polynomials_1d, const std::vector< unsigned int > &renumber, const bool print_iterations_to_deallog=false)
inline ::Table< 2, double > compute_support_point_weights_cell(const unsigned int polynomial_degree)
std::vector<::Table< 2, double > > compute_support_point_weights_perimeter_to_interior(const unsigned int polynomial_degree, const unsigned int dim)
inline ::Table< 2, double > compute_support_point_weights_on_quad(const unsigned int polynomial_degree)
void maybe_update_jacobian_pushed_forward_grads(const CellSimilarity::Similarity cell_similarity, const typename QProjector< dim >::DataSetDescriptor data_set, const typename ::MappingQ< dim, spacedim >::InternalData &data, std::vector< Tensor< 3, spacedim >> &jacobian_pushed_forward_grads)
void maybe_compute_face_data(const ::MappingQ< dim, spacedim > &mapping, const typename ::Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const unsigned int subface_no, const unsigned int n_q_points, const std::vector< double > &weights, const typename ::MappingQ< dim, spacedim >::InternalData &data, internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &output_data)
void transform_hessians(const ArrayView< const Tensor< 3, dim >> &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< 3, spacedim >> &output)
void transform_fields(const ArrayView< const Tensor< rank, dim >> &input, const MappingKind mapping_kind, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_data, const ArrayView< Tensor< rank, spacedim >> &output)
void maybe_update_jacobian_2nd_derivatives(const CellSimilarity::Similarity cell_similarity, const typename QProjector< dim >::DataSetDescriptor data_set, const typename ::MappingQ< dim, spacedim >::InternalData &data, std::vector< DerivativeForm< 3, dim, spacedim >> &jacobian_2nd_derivatives)
std::pair< typename ProductTypeNoPoint< Number, Number2 >::type, Tensor< 1, dim, typename ProductTypeNoPoint< Number, Number2 >::type > > evaluate_tensor_product_value_and_gradient(const std::vector< Polynomials::Polynomial< double >> &poly, const std::vector< Number > &values, const Point< dim, Number2 > &p, const bool d_linear=false, const std::vector< unsigned int > &renumber={})
static const unsigned int invalid_unsigned_int
static double subface_ratio(const internal::SubfaceCase< dim > &subface_case, const unsigned int subface_no)
static double d_linear_shape_function(const Point< dim > &xi, const unsigned int i)
static Point< dim, Number > project_to_unit_cell(const Point< dim, Number > &p)
static void evaluate(const unsigned int n_components, const EvaluationFlags::EvaluationFlags evaluation_flag, const Number *values_dofs, FEEvaluationData< dim, Number, false > &fe_eval)
constexpr DEAL_II_HOST SymmetricTensor< 2, dim, Number > invert(const SymmetricTensor< 2, dim, Number > &)
constexpr DEAL_II_HOST Number determinant(const SymmetricTensor< 2, dim, Number > &)