--- /dev/null
+Added: Evaluation and integration of hessians on faces is now available
+ in the matrix free framework.
+<br>
+(Maximilian Bergbauer, 2021/07/28)
buffer.data(),
phi_m.begin_values(),
false,
+ false,
face);
// Check if the face is an internal or a boundary face and
phi_m.begin_values(),
phi.begin_values(),
false,
+ false,
face);
}
Number * values_dofs_actual,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool add_into_values_array);
};
Number * values_dofs_actual,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool add_into_values_array);
};
Number * values_dofs_actual,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool add_into_values_array)
{
+ // TODO: implement hessians
+ (void)hessians_quad;
+ AssertThrow(!(integration_flag & EvaluationFlags::hessians),
+ ExcNotImplemented());
+
const EvaluatorVariant variant =
EvaluatorSelector<type, (fe_degree + n_q_points_1d > 4)>::variant;
using Eval = EvaluatorTensorProduct<variant,
Number * values_dofs_actual,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool add_into_values_array)
{
+ // TODO: implement hessians
+ (void)hessians_quad;
+ AssertThrow(!(integration_flag & EvaluationFlags::hessians),
+ ExcNotImplemented());
(void)scratch_data;
const unsigned int n_dofs = shape_info.dofs_per_component_on_cell;
Number * values_dofs,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool add_into_values_array);
};
Number * values_dofs,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number *,
const bool add_into_values_array)
{
+ // TODO: implement hessians
+ (void)hessians_quad;
+ AssertThrow(!(integration_flag & EvaluationFlags::hessians),
+ ExcNotImplemented());
+
+
AssertDimension(
shape_info.data.front().shape_gradients_collocation_eo.size(),
(fe_degree + 2) / 2 * (fe_degree + 1));
Number * values_dofs,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool add_into_values_array);
};
Number *values_dofs,
Number *values_quad,
Number *gradients_quad,
+ Number *hessians_quad,
Number *,
const bool add_into_values_array)
{
if (integration_flag & EvaluationFlags::gradients)
FEEvaluationImplCollocation<dim, n_q_points_1d - 1, Number>::
integrate(1,
- integration_flag & EvaluationFlags::gradients,
+ integration_flag & (EvaluationFlags::gradients |
+ EvaluationFlags::hessians),
shape_info,
values_quad,
nullptr,
gradients_quad,
+ hessians_quad,
nullptr,
/*add_into_values_array=*/integration_flag &
EvaluationFlags::values);
Number * values_dofs_actual,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool sum_into_values_array)
{
values_dofs_actual,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
sum_into_values_array);
}
values_dofs_actual,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
sum_into_values_array);
}
values_dofs_actual,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
sum_into_values_array);
}
values_dofs_actual,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
sum_into_values_array);
}
values_dofs_actual,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
sum_into_values_array);
}
values_dofs_actual,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
sum_into_values_array);
}
values_dofs_actual,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
sum_into_values_array);
}
Number * values_dofs,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool evaluate_val,
const bool evaluate_grad,
+ const bool evaluate_hessian,
const unsigned int subface_index)
{
- const AlignedVector<Number> &val1 =
+ const AlignedVector<Number> &val0 =
symmetric_evaluate ?
data.data.front().shape_values_eo :
(subface_index >= GeometryInfo<dim>::max_children_per_cell ?
data.data.front().shape_values :
data.data.front().values_within_subface[subface_index % 2]);
- const AlignedVector<Number> &val2 =
+ const AlignedVector<Number> &val1 =
symmetric_evaluate ?
data.data.front().shape_values_eo :
(subface_index >= GeometryInfo<dim>::max_children_per_cell ?
data.data.front().shape_values :
data.data.front().values_within_subface[subface_index / 2]);
- const AlignedVector<Number> &grad1 =
+ const AlignedVector<Number> &grad0 =
symmetric_evaluate ?
data.data.front().shape_gradients_eo :
(subface_index >= GeometryInfo<dim>::max_children_per_cell ?
data.data.front().shape_gradients :
data.data.front().gradients_within_subface[subface_index % 2]);
- const AlignedVector<Number> &grad2 =
+ const AlignedVector<Number> &grad1 =
symmetric_evaluate ?
data.data.front().shape_gradients_eo :
(subface_index >= GeometryInfo<dim>::max_children_per_cell ?
data.data.front().shape_gradients :
data.data.front().gradients_within_subface[subface_index / 2]);
+ const AlignedVector<Number> &hessian0 =
+ symmetric_evaluate ?
+ data.data.front().shape_hessians_eo :
+ (subface_index >= GeometryInfo<dim>::max_children_per_cell ?
+ data.data.front().shape_hessians :
+ data.data.front().hessians_within_subface[subface_index % 2]);
+ const AlignedVector<Number> &hessian1 =
+ symmetric_evaluate ?
+ data.data.front().shape_hessians_eo :
+ (subface_index >= GeometryInfo<dim>::max_children_per_cell ?
+ data.data.front().shape_hessians :
+ data.data.front().hessians_within_subface[subface_index / 2]);
+
using Eval =
internal::EvaluatorTensorProduct<symmetric_evaluate ?
internal::evaluate_evenodd :
fe_degree + 1,
n_q_points_1d,
Number>;
- Eval eval1(val1,
- grad1,
- AlignedVector<Number>(),
+ Eval eval0(val0,
+ grad0,
+ hessian0,
data.data.front().fe_degree + 1,
data.data.front().n_q_points_1d);
- Eval eval2(val2,
- grad2,
- AlignedVector<Number>(),
+ Eval eval1(val1,
+ grad1,
+ hessian1,
data.data.front().fe_degree + 1,
data.data.front().n_q_points_1d);
Utilities::pow(n_q_points_1d, dim - 1) :
data.n_q_points_face;
- if (evaluate_grad == false)
+ // keep a copy of the original pointer for the case of the Hessians
+ Number *values_dofs_ptr = values_dofs;
+
+ if (evaluate_val && !evaluate_grad)
for (unsigned int c = 0; c < n_components; ++c)
{
switch (dim)
{
case 3:
- eval1.template values<0, true, false>(values_dofs,
+ eval0.template values<0, true, false>(values_dofs,
values_quad);
- eval2.template values<1, true, false>(values_quad,
+ eval1.template values<1, true, false>(values_quad,
values_quad);
break;
case 2:
- eval1.template values<0, true, false>(values_dofs,
+ eval0.template values<0, true, false>(values_dofs,
values_quad);
break;
case 1:
default:
Assert(false, ExcNotImplemented());
}
- values_dofs += 2 * size_deg;
+ values_dofs += 3 * size_deg;
values_quad += n_q_points;
}
- else
+ else if (evaluate_grad)
for (unsigned int c = 0; c < n_components; ++c)
{
switch (dim)
case 3:
if (use_collocation)
{
- eval1.template values<0, true, false>(values_dofs,
+ eval0.template values<0, true, false>(values_dofs,
values_quad);
- eval1.template values<1, true, false>(values_quad,
+ eval0.template values<1, true, false>(values_quad,
values_quad);
internal::EvaluatorTensorProduct<
internal::evaluate_evenodd,
}
else
{
- eval1.template gradients<0, true, false>(values_dofs,
+ // grad x
+ eval0.template gradients<0, true, false>(values_dofs,
scratch_data);
- eval2.template values<1, true, false>(scratch_data,
+ eval1.template values<1, true, false>(scratch_data,
gradients_quad);
- eval1.template values<0, true, false>(values_dofs,
+ // grad y
+ eval0.template values<0, true, false>(values_dofs,
scratch_data);
- eval2.template gradients<1, true, false>(scratch_data,
+ eval1.template gradients<1, true, false>(scratch_data,
gradients_quad +
n_q_points);
if (evaluate_val == true)
- eval2.template values<1, true, false>(scratch_data,
+ eval1.template values<1, true, false>(scratch_data,
values_quad);
}
- eval1.template values<0, true, false>(values_dofs + size_deg,
+ // grad z
+ eval0.template values<0, true, false>(values_dofs + size_deg,
scratch_data);
- eval2.template values<1, true, false>(
+ eval1.template values<1, true, false>(
scratch_data, gradients_quad + (dim - 1) * n_q_points);
break;
case 2:
- eval1.template values<0, true, false>(values_dofs + size_deg,
+ eval0.template values<0, true, false>(values_dofs + size_deg,
gradients_quad +
(dim - 1) *
n_q_points);
- eval1.template gradients<0, true, false>(values_dofs,
+ eval0.template gradients<0, true, false>(values_dofs,
gradients_quad);
if (evaluate_val == true)
- eval1.template values<0, true, false>(values_dofs,
+ eval0.template values<0, true, false>(values_dofs,
values_quad);
break;
case 1:
default:
AssertThrow(false, ExcNotImplemented());
}
- values_dofs += 2 * size_deg;
+ values_dofs += 3 * size_deg;
values_quad += n_q_points;
gradients_quad += dim * n_q_points;
}
+
+ if (evaluate_hessian)
+ {
+ values_dofs = values_dofs_ptr;
+ constexpr unsigned int hdim = (dim * (dim + 1)) / 2;
+ for (unsigned int c = 0; c < n_components; ++c)
+ {
+ switch (dim)
+ {
+ case 3:
+ // grad xx
+ eval0.template hessians<0, true, false>(values_dofs,
+ scratch_data);
+ eval1.template values<1, true, false>(scratch_data,
+ hessians_quad);
+
+ // grad yy
+ eval0.template values<0, true, false>(values_dofs,
+ scratch_data);
+ eval1.template hessians<1, true, false>(scratch_data,
+ hessians_quad +
+ n_q_points);
+
+ // grad zz
+ eval0.template values<0, true, false>(values_dofs +
+ 2 * size_deg,
+ scratch_data);
+ eval1.template values<1, true, false>(scratch_data,
+ hessians_quad +
+ 2 * n_q_points);
+
+ // grad xy
+ eval0.template gradients<0, true, false>(values_dofs,
+ scratch_data);
+ eval1.template gradients<1, true, false>(scratch_data,
+ hessians_quad +
+ 3 * n_q_points);
+
+ // grad xz
+ eval0.template gradients<0, true, false>(values_dofs +
+ size_deg,
+ scratch_data);
+ eval1.template values<1, true, false>(scratch_data,
+ hessians_quad +
+ 4 * n_q_points);
+
+ // grad yz
+ eval0.template values<0, true, false>(values_dofs +
+ size_deg,
+ scratch_data);
+ eval1.template gradients<1, true, false>(scratch_data,
+ hessians_quad +
+ 5 * n_q_points);
+
+ break;
+ case 2:
+ // grad xx
+ eval0.template hessians<0, true, false>(values_dofs,
+ hessians_quad);
+ // grad yy
+ eval0.template values<0, true, false>(
+ values_dofs + 2 * size_deg, hessians_quad + n_q_points);
+ // grad xy
+ eval0.template gradients<0, true, false>(
+ values_dofs + size_deg, hessians_quad + 2 * n_q_points);
+ break;
+ case 1:
+ hessians_quad[0] = values_dofs[2];
+ break;
+ default:
+ AssertThrow(false, ExcNotImplemented());
+ }
+ values_dofs += 3 * size_deg;
+ hessians_quad += hdim * n_q_points;
+ }
+ }
}
static void
Number * values_dofs,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool integrate_val,
const bool integrate_grad,
+ const bool integrate_hessian,
const unsigned int subface_index)
{
- const AlignedVector<Number> &val1 =
+ const AlignedVector<Number> &val0 =
symmetric_evaluate ?
data.data.front().shape_values_eo :
(subface_index >= GeometryInfo<dim>::max_children_per_cell ?
data.data.front().shape_values :
data.data.front().values_within_subface[subface_index % 2]);
- const AlignedVector<Number> &val2 =
+ const AlignedVector<Number> &val1 =
symmetric_evaluate ?
data.data.front().shape_values_eo :
(subface_index >= GeometryInfo<dim>::max_children_per_cell ?
data.data.front().shape_values :
data.data.front().values_within_subface[subface_index / 2]);
- const AlignedVector<Number> &grad1 =
+ const AlignedVector<Number> &grad0 =
symmetric_evaluate ?
data.data.front().shape_gradients_eo :
(subface_index >= GeometryInfo<dim>::max_children_per_cell ?
data.data.front().shape_gradients :
data.data.front().gradients_within_subface[subface_index % 2]);
- const AlignedVector<Number> &grad2 =
+ const AlignedVector<Number> &grad1 =
symmetric_evaluate ?
data.data.front().shape_gradients_eo :
(subface_index >= GeometryInfo<dim>::max_children_per_cell ?
data.data.front().shape_gradients :
data.data.front().gradients_within_subface[subface_index / 2]);
+ const AlignedVector<Number> &hessian0 =
+ symmetric_evaluate ?
+ data.data.front().shape_hessians_eo :
+ (subface_index >= GeometryInfo<dim>::max_children_per_cell ?
+ data.data.front().shape_hessians :
+ data.data.front().hessians_within_subface[subface_index % 2]);
+ const AlignedVector<Number> &hessian1 =
+ symmetric_evaluate ?
+ data.data.front().shape_hessians_eo :
+ (subface_index >= GeometryInfo<dim>::max_children_per_cell ?
+ data.data.front().shape_hessians :
+ data.data.front().hessians_within_subface[subface_index / 2]);
+
using Eval =
internal::EvaluatorTensorProduct<symmetric_evaluate ?
internal::evaluate_evenodd :
fe_degree + 1,
n_q_points_1d,
Number>;
- Eval eval1(val1,
- grad1,
- val1,
+ Eval eval0(val0,
+ grad0,
+ hessian0,
data.data.front().fe_degree + 1,
data.data.front().n_q_points_1d);
- Eval eval2(val2,
- grad2,
- val1,
+ Eval eval1(val1,
+ grad1,
+ hessian1,
data.data.front().fe_degree + 1,
data.data.front().n_q_points_1d);
Utilities::pow(n_q_points_1d, dim - 1) :
data.n_q_points_face;
- if (integrate_grad == false)
+ // keep a copy of the original pointer for the case of the Hessians
+ Number *values_dofs_ptr = values_dofs;
+
+ if (integrate_val && !integrate_grad)
for (unsigned int c = 0; c < n_components; ++c)
{
switch (dim)
{
case 3:
- eval2.template values<1, false, false>(values_quad,
+ eval1.template values<1, false, false>(values_quad,
values_quad);
- eval1.template values<0, false, false>(values_quad,
+ eval0.template values<0, false, false>(values_quad,
values_dofs);
break;
case 2:
- eval1.template values<0, false, false>(values_quad,
+ eval0.template values<0, false, false>(values_quad,
values_dofs);
break;
case 1:
default:
Assert(false, ExcNotImplemented());
}
- values_dofs += 2 * size_deg;
+ values_dofs += 3 * size_deg;
values_quad += n_q_points;
}
- else
+ else if (integrate_grad)
for (unsigned int c = 0; c < n_components; ++c)
{
switch (dim)
{
case 3:
- eval2.template values<1, false, false>(gradients_quad +
+ // grad z
+ eval1.template values<1, false, false>(gradients_quad +
2 * n_q_points,
gradients_quad +
2 * n_q_points);
- eval1.template values<0, false, false>(
+ eval0.template values<0, false, false>(
gradients_quad + 2 * n_q_points, values_dofs + size_deg);
if (use_collocation)
{
gradients_quad + n_q_points, values_quad);
eval_grad.template gradients<0, false, true>(
gradients_quad, values_quad);
- eval1.template values<1, false, false>(values_quad,
+ eval0.template values<1, false, false>(values_quad,
values_quad);
- eval1.template values<0, false, false>(values_quad,
+ eval0.template values<0, false, false>(values_quad,
values_dofs);
}
else
{
if (integrate_val)
{
- eval2.template values<1, false, false>(values_quad,
+ eval1.template values<1, false, false>(values_quad,
scratch_data);
- eval2.template gradients<1, false, true>(
+ eval1.template gradients<1, false, true>(
gradients_quad + n_q_points, scratch_data);
}
else
- eval2.template gradients<1, false, false>(
+ eval1.template gradients<1, false, false>(
gradients_quad + n_q_points, scratch_data);
- eval1.template values<0, false, false>(scratch_data,
+ // grad y
+ eval0.template values<0, false, false>(scratch_data,
values_dofs);
- eval2.template values<1, false, false>(gradients_quad,
+
+ // grad x
+ eval1.template values<1, false, false>(gradients_quad,
scratch_data);
- eval1.template gradients<0, false, true>(scratch_data,
+ eval0.template gradients<0, false, true>(scratch_data,
values_dofs);
}
break;
case 2:
- eval1.template values<0, false, false>(
+ eval0.template values<0, false, false>(
gradients_quad + n_q_points, values_dofs + size_deg);
- eval1.template gradients<0, false, false>(gradients_quad,
+ eval0.template gradients<0, false, false>(gradients_quad,
values_dofs);
if (integrate_val == true)
- eval1.template values<0, false, true>(values_quad,
+ eval0.template values<0, false, true>(values_quad,
values_dofs);
break;
case 1:
default:
AssertThrow(false, ExcNotImplemented());
}
- values_dofs += 2 * size_deg;
+ values_dofs += 3 * size_deg;
values_quad += n_q_points;
gradients_quad += dim * n_q_points;
}
+
+ if (integrate_hessian)
+ {
+ values_dofs = values_dofs_ptr;
+ constexpr unsigned int hdim = (dim * (dim + 1)) / 2;
+ for (unsigned int c = 0; c < n_components; ++c)
+ {
+ switch (dim)
+ {
+ case 3:
+ // grad xx
+ eval1.template values<1, false, false>(hessians_quad,
+ scratch_data);
+ if (integrate_val || integrate_grad)
+ eval0.template hessians<0, false, true>(scratch_data,
+ values_dofs);
+ else
+ eval0.template hessians<0, false, false>(scratch_data,
+ values_dofs);
+
+ // grad yy
+ eval1.template hessians<1, false, false>(hessians_quad +
+ n_q_points,
+ scratch_data);
+ eval0.template values<0, false, true>(scratch_data,
+ values_dofs);
+
+ // grad zz
+ eval1.template values<1, false, false>(hessians_quad +
+ 2 * n_q_points,
+ scratch_data);
+ eval0.template values<0, false, false>(scratch_data,
+ values_dofs +
+ 2 * size_deg);
+
+ // grad xy
+ eval1.template gradients<1, false, false>(hessians_quad +
+ 3 * n_q_points,
+ scratch_data);
+ eval0.template gradients<0, false, true>(scratch_data,
+ values_dofs);
+
+ // grad xz
+ eval1.template values<1, false, false>(hessians_quad +
+ 4 * n_q_points,
+ scratch_data);
+ if (integrate_grad)
+ eval0.template gradients<0, false, true>(scratch_data,
+ values_dofs +
+ size_deg);
+ else
+ eval0.template gradients<0, false, false>(scratch_data,
+ values_dofs +
+ size_deg);
+
+ // grad yz
+ eval1.template gradients<1, false, false>(hessians_quad +
+ 5 * n_q_points,
+ scratch_data);
+ eval0.template values<0, false, true>(scratch_data,
+ values_dofs +
+ size_deg);
+
+ break;
+ case 2:
+ // grad xx
+ if (integrate_val || integrate_grad)
+ eval0.template hessians<0, false, true>(hessians_quad,
+ values_dofs);
+ else
+ eval0.template hessians<0, false, false>(hessians_quad,
+ values_dofs);
+
+ // grad yy
+ eval0.template values<0, false, false>(
+ hessians_quad + n_q_points, values_dofs + 2 * size_deg);
+ // grad xy
+ if (integrate_grad)
+ eval0.template gradients<0, false, true>(
+ hessians_quad + 2 * n_q_points, values_dofs + size_deg);
+ else
+ eval0.template gradients<0, false, false>(
+ hessians_quad + 2 * n_q_points, values_dofs + size_deg);
+ break;
+ case 1:
+ values_dofs[2] = hessians_quad[0];
+ if (!integrate_val)
+ values_dofs[0] = 0;
+ if (!integrate_grad)
+ values_dofs[1] = 0;
+ break;
+ default:
+ AssertThrow(false, ExcNotImplemented());
+ }
+ values_dofs += 3 * size_deg;
+ hessians_quad += hdim * n_q_points;
+ }
+ }
}
};
const Number * input,
Number * output,
const bool do_gradients,
+ const bool do_hessians,
const unsigned int face_no)
{
Assert(static_cast<unsigned int>(fe_degree) ==
input,
output,
do_gradients,
+ do_hessians,
face_no,
data.data.front().fe_degree + 1,
data.data.front().shape_data_on_face,
data.dofs_per_component_on_cell,
- 2 * data.dofs_per_component_on_face);
+ 3 * data.dofs_per_component_on_face);
}
/**
const Number * input,
Number * output,
const bool do_gradients,
+ const bool do_hessians,
const unsigned int face_no)
{
Assert(static_cast<unsigned int>(fe_degree + 1) ==
input,
output,
do_gradients,
+ do_hessians,
face_no,
data.data.front().quadrature.size(),
data.data.front().quadrature_data_on_face,
const Number * input,
Number * output,
const bool do_gradients,
+ const bool do_hessians,
const unsigned int face_no,
const unsigned int n_points_1d,
const std::array<AlignedVector<Number>, 2> &shape_data,
for (unsigned int c = 0; c < n_components; c++)
{
- if (do_gradients)
+ if (do_hessians)
+ evalf.template apply_face<face_direction,
+ do_evaluate,
+ add_into_output,
+ 2,
+ lex_faces>(input, output);
+ else if (!do_hessians && do_gradients)
evalf.template apply_face<face_direction,
do_evaluate,
add_into_output,
input,
output,
do_gradients,
+ do_hessians,
face_no,
n_points_1d,
shape_data,
const bool integrate,
const bool values,
const bool gradients,
+ const bool hessians,
const unsigned int n_q_points,
Number * tmp_values,
Number * values_quad,
- Number * gradients_quad)
+ Number * gradients_quad,
+ Number * hessians_quad)
{
Assert(face_orientation, ExcInternalError());
const unsigned int *orientation = &orientation_map[face_orientation][0];
for (unsigned int q = 0; q < n_q_points; ++q)
gradients_quad[(c * dim + d) * n_q_points + q] = tmp_values[q];
}
+ if (hessians == true)
+ {
+ const unsigned int hdim = (dim * (dim + 1)) / 2;
+ for (unsigned int d = 0; d < hdim; ++d)
+ {
+ if (integrate)
+ for (unsigned int q = 0; q < n_q_points; ++q)
+ tmp_values[q] = hessians_quad[(c * hdim + d) * n_q_points +
+ orientation[q]];
+ else
+ for (unsigned int q = 0; q < n_q_points; ++q)
+ tmp_values[orientation[q]] =
+ hessians_quad[(c * hdim + d) * n_q_points + q];
+ for (unsigned int q = 0; q < n_q_points; ++q)
+ hessians_quad[(c * hdim + d) * n_q_points + q] =
+ tmp_values[q];
+ }
+ }
}
}
const VectorizedArrayType * values_array,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const bool evaluate_values,
const bool evaluate_gradients,
+ const bool evaluate_hessians,
const unsigned int face_no,
const unsigned int subface_index,
const unsigned int face_orientation,
}
}
+ if (evaluate_hessians)
+ {
+ AssertThrow(false, ExcNotImplemented());
+ }
return true;
}
VectorizedArrayType *temp1 = scratch_data;
FEFaceNormalEvaluationImpl<dim, fe_degree, VectorizedArrayType>::
- template interpolate<true, false>(
- n_components, data, values_array, temp1, evaluate_gradients, face_no);
+ template interpolate<true, false>(n_components,
+ data,
+ values_array,
+ temp1,
+ evaluate_gradients,
+ evaluate_hessians,
+ face_no);
const unsigned int n_q_points_1d_actual =
fe_degree > -1 ? n_q_points_1d : 0;
temp1,
values_quad,
gradients_quad,
- scratch_data + 2 *
+ hessians_quad,
+ scratch_data + 3 *
n_components *
dofs_per_face,
evaluate_values,
evaluate_gradients,
+ evaluate_hessians,
subface_index);
else
FEFaceEvaluationImpl<
temp1,
values_quad,
gradients_quad,
- scratch_data + 2 *
+ hessians_quad,
+ scratch_data + 3 *
n_components *
dofs_per_face,
evaluate_values,
evaluate_gradients,
+ evaluate_hessians,
subface_index);
if (face_orientation)
false,
evaluate_values,
evaluate_gradients,
+ evaluate_hessians,
data.n_q_points_face,
scratch_data,
values_quad,
- gradients_quad);
+ gradients_quad,
+ hessians_quad);
return false;
}
VectorizedArrayType * values_array,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const bool integrate_values,
const bool integrate_gradients,
+ const bool integrate_hessians,
const unsigned int face_no,
const unsigned int subface_index,
const unsigned int face_orientation,
}
}
+ if (integrate_hessians)
+ {
+ AssertThrow(false, ExcNotImplemented());
+ }
return true;
}
true,
integrate_values,
integrate_gradients,
+ integrate_hessians,
data.n_q_points_face,
scratch_data,
values_quad,
- gradients_quad);
+ gradients_quad,
+ hessians_quad);
constexpr unsigned int static_dofs_per_face =
fe_degree > -1 ? Utilities::pow(fe_degree + 1, dim - 1) :
temp1,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data +
- 2 * n_components *
+ 3 * n_components *
dofs_per_face,
integrate_values,
integrate_gradients,
+ integrate_hessians,
subface_index);
else
FEFaceEvaluationImpl<
temp1,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data +
- 2 * n_components *
+ 3 * n_components *
dofs_per_face,
integrate_values,
integrate_gradients,
+ integrate_hessians,
subface_index);
FEFaceNormalEvaluationImpl<dim, fe_degree, VectorizedArrayType>::
temp1,
values_array,
integrate_gradients,
+ integrate_hessians,
face_no);
return false;
}
auto &dof_info = proc.dof_info;
auto values_quad = proc.values_quad;
auto gradients_quad = proc.gradients_quad;
+ auto hessians_quad = proc.hessians_quad;
auto scratch_data = proc.scratch_data;
auto do_values = proc.do_values;
auto do_gradients = proc.do_gradients;
+ auto do_hessians = proc.do_hessians;
auto active_fe_index = proc.active_fe_index;
auto first_selected_component = proc.first_selected_component;
auto cells = proc.cells;
true,
do_values,
do_gradients,
+ do_hessians,
data.n_q_points_face,
scratch_data,
values_quad,
- gradients_quad);
+ gradients_quad,
+ hessians_quad);
}
// we know that the gradient weights for the Hermite case on the
}
// case 4: contiguous indices without interleaving
- else if (n_face_orientations > 1 ||
- dof_info.index_storage_variants[dof_access_index][cell] ==
- MatrixFreeFunctions::DoFInfo::IndexStorageVariants::
- contiguous)
+ else if ((n_face_orientations > 1 ||
+ dof_info.index_storage_variants[dof_access_index][cell] ==
+ MatrixFreeFunctions::DoFInfo::IndexStorageVariants::
+ contiguous) &&
+ !do_hessians)
{
const unsigned int *indices =
&dof_info.dof_indices_contiguous[dof_access_index]
false,
do_values,
do_gradients,
+ do_hessians,
data.n_q_points_face,
scratch_data,
values_quad,
- gradients_quad);
+ gradients_quad,
+ hessians_quad);
}
return accesses_global_vector;
const MatrixFreeFunctions::DoFInfo & dof_info,
VectorizedArrayType * values_quad,
VectorizedArrayType *gradients_quad,
+ VectorizedArrayType *hessians_quad,
VectorizedArrayType *scratch_data,
const bool evaluate_values,
const bool evaluate_gradients,
+ const bool evaluate_hessians,
const unsigned int active_fe_index,
const unsigned int first_selected_component,
const std::array<unsigned int, VectorizedArrayType::size()> cells,
dof_info,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
evaluate_values,
evaluate_gradients,
+ evaluate_hessians,
active_fe_index,
first_selected_component,
cells,
const MatrixFreeFunctions::DoFInfo & dof_info,
VectorizedArrayType * values_quad,
VectorizedArrayType *gradients_quad,
+ VectorizedArrayType *hessians_quad,
VectorizedArrayType *scratch_data,
const bool do_values,
const bool do_gradients,
+ const bool do_hessians,
const unsigned int active_fe_index,
const unsigned int first_selected_component,
const std::array<unsigned int, VectorizedArrayType::size()> cells,
, dof_info(dof_info)
, values_quad(values_quad)
, gradients_quad(gradients_quad)
+ , hessians_quad(hessians_quad)
, scratch_data(scratch_data)
, do_values(do_values)
, do_gradients(do_gradients)
+ , do_hessians(do_hessians)
, active_fe_index(active_fe_index)
, first_selected_component(first_selected_component)
, cells(cells)
const unsigned int n_q_points =
fe_degree > -1 ? Utilities::pow(n_q_points_1d, dim - 1) :
data.n_q_points_face;
+ constexpr unsigned int hdim = (dim * (dim + 1)) / 2;
if (fe_degree > -1 &&
subface_index >= GeometryInfo<dim>::max_children_per_cell &&
data.element_type <= MatrixFreeFunctions::tensor_symmetric)
temp1,
values_quad + comp * n_q_points,
gradients_quad + comp * dim * n_q_points,
- scratch_data + 2 * dofs_per_face,
+ hessians_quad + comp * hdim * n_q_points,
+ scratch_data + 3 * dofs_per_face,
do_values,
do_gradients,
+ do_hessians,
subface_index);
else
FEFaceEvaluationImpl<false,
temp1,
values_quad + comp * n_q_points,
gradients_quad + comp * dim * n_q_points,
- scratch_data + 2 * dofs_per_face,
+ hessians_quad + comp * hdim * n_q_points,
+ scratch_data + 3 * dofs_per_face,
do_values,
do_gradients,
+ do_hessians,
subface_index);
}
const MatrixFreeFunctions::DoFInfo & dof_info;
VectorizedArrayType * values_quad;
VectorizedArrayType * gradients_quad;
+ VectorizedArrayType * hessians_quad;
VectorizedArrayType * scratch_data;
const bool do_values;
const bool do_gradients;
+ const bool do_hessians;
const unsigned int active_fe_index;
const unsigned int first_selected_component;
const std::array<unsigned int, VectorizedArrayType::size()> cells;
VectorizedArrayType * values_array,
VectorizedArrayType * values_quad,
VectorizedArrayType *gradients_quad,
+ VectorizedArrayType *hessians_quad,
VectorizedArrayType *scratch_data,
const bool integrate_values,
const bool integrate_gradients,
+ const bool integrate_hessians,
const unsigned int active_fe_index,
const unsigned int first_selected_component,
const std::array<unsigned int, VectorizedArrayType::size()> cells,
values_array,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
integrate_values,
integrate_gradients,
+ integrate_hessians,
face_nos[0],
subface_index,
face_orientations[0],
dof_info,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
integrate_values,
integrate_gradients,
+ integrate_hessians,
active_fe_index,
first_selected_component,
cells,
const MatrixFreeFunctions::DoFInfo & dof_info,
VectorizedArrayType * values_quad,
VectorizedArrayType *gradients_quad,
+ VectorizedArrayType *hessians_quad,
VectorizedArrayType *scratch_data,
const bool do_values,
const bool do_gradients,
+ const bool do_hessians,
const unsigned int active_fe_index,
const unsigned int first_selected_component,
const std::array<unsigned int, VectorizedArrayType::size()> cells,
, dof_info(dof_info)
, values_quad(values_quad)
, gradients_quad(gradients_quad)
+ , hessians_quad(hessians_quad)
, scratch_data(scratch_data)
, do_values(do_values)
, do_gradients(do_gradients)
+ , do_hessians(do_hessians)
, active_fe_index(active_fe_index)
, first_selected_component(first_selected_component)
, cells(cells)
temp1,
values_array + comp * data.dofs_per_component_on_cell,
do_gradients,
+ do_hessians,
face_nos[0]);
}
const unsigned int n_q_points =
fe_degree > -1 ? Utilities::pow(n_q_points_1d, dim - 1) :
data.n_q_points_face;
+ constexpr unsigned int hdim = (dim * (dim + 1)) / 2;
if (fe_degree > -1 &&
subface_index >= GeometryInfo<dim>::max_children_per_cell &&
data.element_type <=
temp1,
values_quad + comp * n_q_points,
gradients_quad + dim * comp * n_q_points,
- scratch_data + 2 * dofs_per_face,
+ hessians_quad + hdim * comp * n_q_points,
+ scratch_data + 3 * dofs_per_face,
do_values,
do_gradients,
+ do_hessians,
subface_index);
else
internal::FEFaceEvaluationImpl<false,
temp1,
values_quad + comp * n_q_points,
gradients_quad + dim * comp * n_q_points,
- scratch_data + 2 * dofs_per_face,
+ hessians_quad + hdim * comp * n_q_points,
+ scratch_data + 3 * dofs_per_face,
do_values,
do_gradients,
+ do_hessians,
subface_index);
}
const MatrixFreeFunctions::DoFInfo & dof_info;
VectorizedArrayType * values_quad;
VectorizedArrayType * gradients_quad;
+ VectorizedArrayType * hessians_quad;
VectorizedArrayType * scratch_data;
const bool do_values;
const bool do_gradients;
+ const bool do_hessians;
const unsigned int active_fe_index;
const unsigned int first_selected_component;
const std::array<unsigned int, VectorizedArrayType::size()> cells;
Number * values_dofs_actual,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool sum_into_values_array = false);
};
Number * values_dofs_actual,
Number * values_quad,
Number * gradients_quad,
+ Number * hessians_quad,
Number * scratch_data,
const bool sum_into_values_array)
{
values_dofs_actual,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
sum_into_values_array);
}
VectorizedArrayType *values_dofs_actual,
VectorizedArrayType *values_quad,
VectorizedArrayType *gradients_quad,
+ VectorizedArrayType *hessians_quad,
VectorizedArrayType *scratch_data,
const bool sum_into_values_array);
const VectorizedArrayType * values_array,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const bool evaluate_values,
const bool evaluate_gradients,
+ const bool evaluate_hessians,
const unsigned int face_no,
const unsigned int subface_index,
const unsigned int face_orientation,
VectorizedArrayType * values_array,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const bool integrate_values,
const bool integrate_gradients,
+ const bool integrate_hessians,
const unsigned int face_no,
const unsigned int subface_index,
const unsigned int face_orientation,
const MatrixFreeFunctions::DoFInfo & dof_info,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const bool evaluate_values,
const bool evaluate_gradients,
+ const bool evaluate_hessians,
const unsigned int active_fe_index,
const unsigned int first_selected_component,
const std::array<unsigned int, VectorizedArrayType::size()> cells,
VectorizedArrayType * values_array,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const bool integrate_values,
const bool integrate_gradients,
+ const bool integrate_hessians,
const unsigned int active_fe_index,
const unsigned int first_selected_component,
const std::array<unsigned int, VectorizedArrayType::size()> cells,
VectorizedArrayType *values_dofs_actual,
VectorizedArrayType *values_quad,
VectorizedArrayType *gradients_quad,
+ VectorizedArrayType *hessians_quad,
VectorizedArrayType *scratch_data,
const bool sum_into_values_array)
{
values_dofs_actual,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
sum_into_values_array);
}
const VectorizedArrayType * values_array,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const bool evaluate_values,
const bool evaluate_gradients,
+ const bool evaluate_hessians,
const unsigned int face_no,
const unsigned int subface_index,
const unsigned int face_orientation,
values_array,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
evaluate_values,
evaluate_gradients,
+ evaluate_hessians,
face_no,
subface_index,
face_orientation,
VectorizedArrayType * values_array,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const bool integrate_values,
const bool integrate_gradients,
+ const bool integrate_hessians,
const unsigned int face_no,
const unsigned int subface_index,
const unsigned int face_orientation,
values_array,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
integrate_values,
integrate_gradients,
+ integrate_hessians,
face_no,
subface_index,
face_orientation,
const MatrixFreeFunctions::DoFInfo & dof_info,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const bool evaluate_values,
const bool evaluate_gradients,
+ const bool evaluate_hessians,
const unsigned int active_fe_index,
const unsigned int first_selected_component,
const std::array<unsigned int, VectorizedArrayType::size()> cells,
dof_info,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
evaluate_values,
evaluate_gradients,
+ evaluate_hessians,
active_fe_index,
first_selected_component,
cells,
VectorizedArrayType * values_array,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const bool integrate_values,
const bool integrate_gradients,
+ const bool integrate_hessians,
const unsigned int active_fe_index,
const unsigned int first_selected_component,
const std::array<unsigned int, VectorizedArrayType::size()> cells,
values_array,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
integrate_values,
integrate_gradients,
+ integrate_hessians,
active_fe_index,
first_selected_component,
cells,
using value_type = Tensor<1, n_components_, VectorizedArrayType>;
using gradient_type =
Tensor<1, n_components_, Tensor<1, dim, VectorizedArrayType>>;
+ using hessian_type =
+ Tensor<1, n_components_, Tensor<2, dim, VectorizedArrayType>>;
static constexpr unsigned int dimension = dim;
static constexpr unsigned int n_components = n_components_;
submit_normal_derivative(const value_type grad_in,
const unsigned int q_point);
+ /**
+ * Write a contribution that is tested by the Hessian to the field
+ * containing the values at quadrature points with component @p q_point.
+ * Access to the same field as through get_hessian(). If applied before the
+ * function FEEvaluation::integrate(EvaluationFlags::hessians) is called,
+ * this specifies what is tested by the Hessians of all basis functions on the
+ * current cell and integrated over.
+ *
+ * Note that the derived class FEEvaluationAccess overloads this operation
+ * with specializations for the scalar case (n_components == 1) and for the
+ * vector-valued case (n_components == dim).
+ */
+ void
+ submit_hessian(const hessian_type hessian_in, const unsigned int q_point);
+
/**
* Return the Hessian of a finite element function at quadrature point
* number @p q_point after a call to
*/
bool gradients_quad_submitted;
+ /**
+ * Debug information to track whether hessians on quadrature points have
+ * been submitted for integration before the integration is actually stared.
+ * Used to control exceptions when uninitialized data is used.
+ */
+ bool hessians_quad_submitted;
+
/**
* For a FiniteElement with more than one base element, select at which
* component this data structure should start.
using number_type = Number;
using value_type = VectorizedArrayType;
using gradient_type = Tensor<1, dim, VectorizedArrayType>;
+ using hessian_type = Tensor<2, dim, VectorizedArrayType>;
static constexpr unsigned int dimension = dim;
using BaseClass =
FEEvaluationBase<dim, 1, Number, is_face, VectorizedArrayType>;
/**
* @copydoc FEEvaluationBase<dim,1,Number,is_face>::get_hessian()
*/
- Tensor<2, dim, VectorizedArrayType>
+ hessian_type
get_hessian(unsigned int q_point) const;
/**
gradient_type
get_hessian_diagonal(const unsigned int q_point) const;
+ /**
+ * @copydoc FEEvaluationBase<dim,1,Number,is_face>::submit_hessian()
+ */
+ void
+ submit_hessian(const hessian_type hessian_in, const unsigned int q_point);
+
/**
* @copydoc FEEvaluationBase<dim,1,Number,is_face>::get_laplacian()
*/
using number_type = Number;
using value_type = VectorizedArrayType;
using gradient_type = Tensor<1, 1, VectorizedArrayType>;
+ using hessian_type = Tensor<2, 1, VectorizedArrayType>;
static constexpr unsigned int dimension = 1;
using BaseClass =
FEEvaluationBase<1, 1, Number, is_face, VectorizedArrayType>;
/**
* @copydoc FEEvaluationBase<1,1,Number,is_face>::get_hessian()
*/
- Tensor<2, 1, VectorizedArrayType>
+ hessian_type
get_hessian(unsigned int q_point) const;
/**
gradient_type
get_hessian_diagonal(const unsigned int q_point) const;
+ /**
+ * @copydoc FEEvaluationBase<1,1,Number,is_face>::submit_hessian()
+ */
+ void
+ submit_hessian(const hessian_type hessian_in, const unsigned int q_point);
+
/**
* @copydoc FEEvaluationBase<1,1,Number,is_face>::get_laplacian()
*/
FEEvaluationBase<dim, n_components_, Number, is_face, VectorizedArrayType>::
get_hessian(const unsigned int q_point) const
{
- Assert(!is_face, ExcNotImplemented());
# ifdef DEBUG
Assert(this->hessians_quad_initialized == true,
internal::ExcAccessToUninitializedField());
constexpr unsigned int hdim = (dim * (dim + 1)) / 2;
// Cartesian cell
- if (this->cell_type == internal::MatrixFreeFunctions::cartesian)
+ if (!is_face && this->cell_type == internal::MatrixFreeFunctions::cartesian)
{
for (unsigned int comp = 0; comp < n_components; comp++)
{
}
}
// cell with general Jacobian, but constant within the cell
- else if (this->cell_type == internal::MatrixFreeFunctions::affine)
+ else if (this->cell_type <= internal::MatrixFreeFunctions::affine)
{
for (unsigned int comp = 0; comp < n_components; comp++)
{
+template <int dim,
+ int n_components_,
+ typename Number,
+ bool is_face,
+ typename VectorizedArrayType>
+inline DEAL_II_ALWAYS_INLINE void
+FEEvaluationBase<dim, n_components_, Number, is_face, VectorizedArrayType>::
+ submit_hessian(
+ const Tensor<1, n_components_, Tensor<2, dim, VectorizedArrayType>>
+ hessian_in,
+ const unsigned int q_point)
+{
+ Assert(this->cell != numbers::invalid_unsigned_int, ExcNotInitialized());
+ AssertIndexRange(q_point, this->n_quadrature_points);
+ Assert(this->J_value != nullptr,
+ internal::ExcMatrixFreeAccessToUninitializedMappingField(
+ "update_hessians"));
+ Assert(this->jacobian != nullptr,
+ internal::ExcMatrixFreeAccessToUninitializedMappingField(
+ "update_hessians"));
+# ifdef DEBUG
+ this->hessians_quad_submitted = true;
+# endif
+
+ // compute hessian_unit = J^T * hessian_in(u) * J
+ const std::size_t nqp = this->n_quadrature_points;
+ constexpr unsigned int hdim = (dim * (dim + 1)) / 2;
+ if (!is_face && this->cell_type == internal::MatrixFreeFunctions::cartesian)
+ {
+ const VectorizedArrayType JxW =
+ this->J_value[0] * this->quadrature_weights[q_point];
+
+ // diagonal part
+ for (unsigned int d = 0; d < dim; ++d)
+ {
+ const auto jac_d = this->jacobian[0][d][d];
+ const VectorizedArrayType factor = jac_d * jac_d * JxW;
+ for (unsigned int comp = 0; comp < n_components; comp++)
+ hessians_quad[(comp * hdim + d) * nqp + q_point] =
+ hessian_in[comp][d][d] * factor;
+ }
+
+ // off diagonal part (use symmetry)
+ for (unsigned int d = 1, off_dia = dim; d < dim; ++d)
+ for (unsigned int e = 0; e < d; ++e, ++off_dia)
+ {
+ const auto jac_d = this->jacobian[0][d][d];
+ const auto jac_e = this->jacobian[0][e][e];
+ const VectorizedArrayType factor = jac_d * jac_e * JxW;
+ for (unsigned int comp = 0; comp < n_components; comp++)
+ hessians_quad[(comp * hdim + off_dia) * nqp + q_point] =
+ 2.0 * hessian_in[comp][d][e] * factor;
+ }
+ }
+ // cell with general Jacobian, but constant within the cell
+ else if (this->cell_type <= internal::MatrixFreeFunctions::affine)
+ {
+ const Tensor<2, dim, VectorizedArrayType> jac = this->jacobian[0];
+ const VectorizedArrayType JxW =
+ this->J_value[0] * this->quadrature_weights[q_point];
+ for (unsigned int comp = 0; comp < n_components; ++comp)
+ {
+ // 1. tmp = hessian_in(u) * J
+ VectorizedArrayType tmp[dim][dim];
+ for (unsigned int i = 0; i < dim; ++i)
+ for (unsigned int j = 0; j < dim; ++j)
+ tmp[i][j] = 0.;
+
+ for (unsigned int i = 0; i < dim; ++i)
+ for (unsigned int j = 0; j < dim; ++j)
+ for (unsigned int k = 0; k < dim; ++k)
+ tmp[i][j] += hessian_in[comp][i][k] * jac[k][j] * JxW;
+
+ // 2. hessian_unit = J^T * tmp
+
+ // diagonal part
+ for (unsigned int d = 0; d < dim; ++d)
+ {
+ hessians_quad[(comp * hdim + d) * nqp + q_point] = 0;
+ for (unsigned int i = 0; i < dim; ++i)
+ hessians_quad[(comp * hdim + d) * nqp + q_point] +=
+ jac[i][d] * tmp[i][d];
+ }
+
+ // off diagonal part (use symmetry)
+ for (unsigned int d = 1, off_dia = dim; d < dim; ++d)
+ for (unsigned int e = 0; e < d; ++e, ++off_dia)
+ {
+ hessians_quad[(comp * hdim + off_dia) * nqp + q_point] = 0;
+ for (unsigned int i = 0; i < dim; ++i)
+ hessians_quad[(comp * hdim + off_dia) * nqp + q_point] +=
+ 2.0 * jac[i][e] * tmp[i][d];
+ }
+ }
+ }
+ else
+ {
+ // TODO: implement general path submit_hessian
+ AssertThrow(false, ExcNotImplemented());
+ }
+}
+
+
+
template <int dim,
int n_components_,
typename Number,
+template <int dim, typename Number, bool is_face, typename VectorizedArrayType>
+inline DEAL_II_ALWAYS_INLINE void
+FEEvaluationAccess<dim, 1, Number, is_face, VectorizedArrayType>::
+ submit_hessian(const Tensor<2, dim, VectorizedArrayType> hessian_in,
+ const unsigned int q_point)
+{
+ Tensor<1, 1, Tensor<2, dim, VectorizedArrayType>> hessian;
+ hessian[0] = hessian_in;
+ BaseClass::submit_hessian(hessian, q_point);
+}
+
+
+
template <int dim, typename Number, bool is_face, typename VectorizedArrayType>
inline VectorizedArrayType
FEEvaluationAccess<dim, 1, Number, is_face, VectorizedArrayType>::
}
+template <typename Number, bool is_face, typename VectorizedArrayType>
+inline DEAL_II_ALWAYS_INLINE void
+FEEvaluationAccess<1, 1, Number, is_face, VectorizedArrayType>::submit_hessian(
+ const Tensor<2, 1, VectorizedArrayType> hessian_in,
+ const unsigned int q_point)
+{
+ Tensor<1, 1, Tensor<2, 1, VectorizedArrayType>> hessian;
+ hessian[0] = hessian_in;
+ BaseClass::submit_hessian(hessian, q_point);
+}
+
template <typename Number, bool is_face, typename VectorizedArrayType>
inline VectorizedArrayType
const internal::MatrixFreeFunctions::DoFInfo *dof_info,
VectorizedArrayType * values_quad,
VectorizedArrayType * gradients_quad,
+ VectorizedArrayType * hessians_quad,
VectorizedArrayType * scratch_data,
const internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType> *data,
const EvaluationFlags::EvaluationFlags integration_flag)
vec_values,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
true);
else
vec_values,
values_quad,
gradients_quad,
+ hessians_quad,
scratch_data,
true);
const VectorizedArrayType *,
const VectorizedArrayType *,
const VectorizedArrayType *,
+ const VectorizedArrayType *,
const internal::MatrixFreeFunctions::ShapeInfo<VectorizedArrayType> *,
const EvaluationFlags::EvaluationFlags)
{
ExcNotInitialized());
Assert(
- (integration_flag &
- ~(EvaluationFlags::values | EvaluationFlags::gradients)) == 0,
- ExcMessage(
- "Only EvaluationFlags::values and EvaluationFlags::gradients are supported."));
+ (integration_flag & ~(EvaluationFlags::values | EvaluationFlags::gradients |
+ EvaluationFlags::hessians)) == 0,
+ ExcMessage("Only EvaluationFlags::values, EvaluationFlags::gradients, and "
+ "EvaluationFlags::hessians are supported."));
if (fe_degree > -1)
SelectEvaluator<dim, fe_degree, n_q_points_1d, VectorizedArrayType>::
values_array,
this->values_quad,
this->gradients_quad,
+ this->hessians_quad,
this->scratch_data,
false);
else
values_array,
this->values_quad,
this->gradients_quad,
+ this->hessians_quad,
this->scratch_data,
false);
this->dof_info,
this->values_quad,
this->gradients_quad,
+ this->hessians_quad,
this->scratch_data,
this->data,
integration_flag) == false)
const EvaluationFlags::EvaluationFlags evaluation_flag)
{
Assert(
- (evaluation_flag &
- ~(EvaluationFlags::values | EvaluationFlags::gradients)) == 0,
+ (evaluation_flag & ~(EvaluationFlags::values | EvaluationFlags::gradients |
+ EvaluationFlags::hessians)) == 0,
ExcMessage(
- "Only EvaluationFlags::values and EvaluationFlags::gradients are supported."));
+ "Only EvaluationFlags::values, EvaluationFlags::gradients, and EvaluationFlags::hessians are supported."));
if (!(evaluation_flag & EvaluationFlags::values) &&
- !(evaluation_flag & EvaluationFlags::gradients))
+ !(evaluation_flag & EvaluationFlags::gradients) &&
+ !(evaluation_flag & EvaluationFlags::hessians))
return;
if (this->dof_access_index ==
values_array,
this->begin_values(),
this->begin_gradients(),
+ this->begin_hessians(),
this->scratch_data,
evaluation_flag & EvaluationFlags::values,
evaluation_flag & EvaluationFlags::gradients,
+ evaluation_flag & EvaluationFlags::hessians,
face_nos[0],
this->subface_index,
face_orientations[0],
values_array,
this->begin_values(),
this->begin_gradients(),
+ this->begin_hessians(),
this->scratch_data,
evaluation_flag & EvaluationFlags::values,
evaluation_flag & EvaluationFlags::gradients,
+ evaluation_flag & EvaluationFlags::hessians,
this->face_no,
this->subface_index,
this->face_orientation,
values_array,
this->begin_values(),
this->begin_gradients(),
+ this->begin_hessians(),
this->scratch_data,
evaluation_flag & EvaluationFlags::values,
evaluation_flag & EvaluationFlags::gradients,
+ evaluation_flag & EvaluationFlags::hessians,
this->face_no,
this->subface_index,
this->face_orientation,
this->values_quad_initialized = true;
if (evaluation_flag & EvaluationFlags::gradients)
this->gradients_quad_initialized = true;
+ if (evaluation_flag & EvaluationFlags::hessians)
+ this->hessians_quad_initialized = true;
# endif
}
VectorizedArrayType * values_array)
{
Assert(
- (evaluation_flag &
- ~(EvaluationFlags::values | EvaluationFlags::gradients)) == 0,
+ (evaluation_flag & ~(EvaluationFlags::values | EvaluationFlags::gradients |
+ EvaluationFlags::hessians)) == 0,
ExcMessage(
- "Only EvaluationFlags::values and EvaluationFlags::gradients are supported."));
+ "Only EvaluationFlags::values, EvaluationFlags::gradients, and EvaluationFlags::hessians are supported."));
if (!(evaluation_flag & EvaluationFlags::values) &&
- !(evaluation_flag & EvaluationFlags::gradients))
+ !(evaluation_flag & EvaluationFlags::gradients) &&
+ !(evaluation_flag & EvaluationFlags::hessians))
return;
if (fe_degree > -1)
values_array,
this->begin_values(),
this->begin_gradients(),
+ this->begin_hessians(),
this->scratch_data,
evaluation_flag & EvaluationFlags::values,
evaluation_flag & EvaluationFlags::gradients,
+ evaluation_flag & EvaluationFlags::hessians,
this->face_no,
this->subface_index,
this->face_orientation,
values_array,
this->begin_values(),
this->begin_gradients(),
+ this->begin_hessians(),
this->scratch_data,
evaluation_flag & EvaluationFlags::values,
evaluation_flag & EvaluationFlags::gradients,
+ evaluation_flag & EvaluationFlags::hessians,
this->face_no,
this->subface_index,
this->face_orientation,
const EvaluationFlags::EvaluationFlags evaluation_flag)
{
Assert(
- (evaluation_flag &
- ~(EvaluationFlags::values | EvaluationFlags::gradients)) == 0,
+ (evaluation_flag & ~(EvaluationFlags::values | EvaluationFlags::gradients |
+ EvaluationFlags::hessians)) == 0,
ExcMessage(
- "Only EvaluationFlags::values and EvaluationFlags::gradients are supported."));
+ "Only EvaluationFlags::values, EvaluationFlags::gradients, and EvaluationFlags::hessians are supported."));
const auto fu = [&]() {
const auto shared_vector_data = internal::get_shared_vector_data(
*this->dof_info,
this->begin_values(),
this->begin_gradients(),
+ this->begin_hessians(),
this->scratch_data,
evaluation_flag & EvaluationFlags::values,
evaluation_flag & EvaluationFlags::gradients,
+ evaluation_flag & EvaluationFlags::hessians,
this->active_fe_index,
this->first_selected_component,
cells,
*this->dof_info,
this->begin_values(),
this->begin_gradients(),
+ this->begin_hessians(),
this->scratch_data,
evaluation_flag & EvaluationFlags::values,
evaluation_flag & EvaluationFlags::gradients,
+ evaluation_flag & EvaluationFlags::hessians,
this->active_fe_index,
this->first_selected_component,
cells_,
*this->dof_info,
this->begin_values(),
this->begin_gradients(),
+ this->begin_hessians(),
this->scratch_data,
evaluation_flag & EvaluationFlags::values,
evaluation_flag & EvaluationFlags::gradients,
+ evaluation_flag & EvaluationFlags::hessians,
this->active_fe_index,
this->first_selected_component,
cells_,
this->values_quad_initialized = true;
if (evaluation_flag & EvaluationFlags::gradients)
this->gradients_quad_initialized = true;
+ if (evaluation_flag & EvaluationFlags::hessians)
+ this->hessians_quad_initialized = true;
# endif
}
this->begin_dof_values(),
this->begin_values(),
this->begin_gradients(),
+ this->begin_hessians(),
this->scratch_data,
evaluation_flag & EvaluationFlags::values,
evaluation_flag & EvaluationFlags::gradients,
+ evaluation_flag & EvaluationFlags::hessians,
this->active_fe_index,
this->first_selected_component,
cells_,
this->begin_dof_values(),
this->begin_values(),
this->begin_gradients(),
+ this->begin_hessians(),
this->scratch_data,
evaluation_flag & EvaluationFlags::values,
evaluation_flag & EvaluationFlags::gradients,
+ evaluation_flag & EvaluationFlags::hessians,
this->active_fe_index,
this->first_selected_component,
cells_,
const unsigned int n_q_points = my_data.descriptor[0].n_q_points;
const unsigned int n_mapping_points =
shape_info.dofs_per_component_on_cell;
+ constexpr unsigned int hess_dim = dim * (dim + 1) / 2;
AlignedVector<VectorizedDouble> cell_points(dim * n_mapping_points);
AlignedVector<VectorizedDouble> face_quads(dim * n_q_points);
AlignedVector<VectorizedDouble> face_grads(dim * dim * n_q_points);
+ AlignedVector<VectorizedDouble> face_grad_grads(dim * hess_dim *
+ n_q_points);
+
AlignedVector<VectorizedDouble> scratch_data(
dim * (2 * n_q_points + 3 * n_mapping_points));
cell_points.data(),
face_quads.data(),
face_grads.data(),
+ face_grad_grads.data(),
scratch_data.data(),
true,
true,
+ update_flags_faces & update_jacobian_grads,
face_no,
GeometryInfo<dim>::max_children_per_cell,
faces[face].face_orientation > 8 ?
my_data.jacobians[0][offset + q][d][e]);
}
+ if (update_flags_faces & update_jacobian_grads)
+ {
+ // TODO: implement jacobian grads for general path
+ AssertThrow(false, ExcNotImplemented());
+ }
+
std::array<Tensor<1, dim, VectorizedDouble>, dim - 1>
tangential_vectors;
for (unsigned int d = 0; d != dim - 1; ++d)
cell_points.data(),
face_quads.data(),
face_grads.data(),
+ face_grad_grads.data(),
scratch_data.data(),
false,
true,
+ update_flags_faces & update_jacobian_grads,
faces[face].exterior_face_no,
faces[face].subface_index,
faces[face].face_orientation < 8 ?
vv,
my_data.jacobians[1][offset + q][d][e]);
}
+
+ if (update_flags_faces & update_jacobian_grads)
+ {
+ // TODO: implement jacobian grads for general path
+ AssertThrow(false, ExcNotImplemented());
+ }
+
my_data.normals_times_jacobians[1][offset + q] =
my_data.normal_vectors[offset + q] *
my_data.jacobians[1][offset + q];
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:1d::Test without templates on FEEvaluation
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:1d::Testing FE_Q<1>(2)
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:1d::Testing FE_DGQ<1>(2)
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:1d::Test without templates on FEEvaluation
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:2d::Testing FE_DGQ<2>(0)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_Q<2>(1)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(1)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_Q<2>(2)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(2)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_Q<2>(3)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(3)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_Q<2>(4)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(4)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:3d::Testing FE_Q<3>(1)
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Testing FE_DGQ<3>(1)
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Test without templates on FEEvaluation
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Testing FE_Q<3>(2)
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Testing FE_DGQ<3>(2)
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Test without templates on FEEvaluation
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(3)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_Q<2>(4)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(4)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
{
fe_evalm.reinit(face);
fe_evalm.read_dof_values(src);
- fe_evalm.evaluate(EvaluationFlags::values | EvaluationFlags::gradients);
+ fe_evalm.evaluate(EvaluationFlags::values | EvaluationFlags::gradients |
+ EvaluationFlags::hessians);
fe_evalp.reinit(face);
fe_evalp.read_dof_values(src);
- fe_evalp.evaluate(EvaluationFlags::values | EvaluationFlags::gradients);
+ fe_evalp.evaluate(EvaluationFlags::values | EvaluationFlags::gradients |
+ EvaluationFlags::hessians);
for (unsigned int j = 0; j < VectorizedArray<Number>::size(); ++j)
{
++facem_times;
++facep_times;
Point<dim> p;
+
+ // interior face
for (unsigned int d = 0; d < dim; ++d)
p[d] = fe_evalm.quadrature_point(q)[d][j];
facem_errors[0] +=
std::abs(fe_evalm.get_value(q)[j] - function.value(p, 0));
for (unsigned int d = 0; d < dim; ++d)
- {
- // std::cout << fe_evalm.get_gradient(q)[d][j] << " ";
- facem_errors[1] += std::abs(fe_evalm.get_gradient(q)[d][j] -
- function.gradient(p, 0)[d]);
- }
+ facem_errors[1] += std::abs(fe_evalm.get_gradient(q)[d][j] -
+ function.gradient(p, 0)[d]);
+ for (unsigned int d = 0; d < dim; ++d)
+ for (unsigned int e = 0; e < dim; ++e)
+ facem_errors[2] +=
+ std::abs(fe_evalm.get_hessian(q)[d][e][j] -
+ function.hessian(p, 0)[d][e]);
double normal_derivative = 0;
for (unsigned int d = 0; d < dim; ++d)
normal_derivative += function.gradient(p, 0)[d] *
fe_evalm.get_normal_vector(q)[d][j];
facem_errors[3] += std::abs(
fe_evalm.get_normal_derivative(q)[j] - normal_derivative);
+
+ // exterior face
facep_errors[0] +=
std::abs(fe_evalp.get_value(q)[j] - function.value(p, 0));
for (unsigned int d = 0; d < dim; ++d)
facep_errors[1] += std::abs(fe_evalp.get_gradient(q)[d][j] -
function.gradient(p, 0)[d]);
+ for (unsigned int d = 0; d < dim; ++d)
+ for (unsigned int e = 0; e < dim; ++e)
+ facep_errors[2] +=
+ std::abs(fe_evalp.get_hessian(q)[d][e][j] -
+ function.hessian(p, 0)[d][e]);
facep_errors[3] += std::abs(
fe_evalp.get_normal_derivative(q)[j] - normal_derivative);
- // hessians not currently implemented in FEFaceEvaluation
- // for (unsigned int d=0; d<dim; ++d)
- // for (unsigned int e=0; e<dim; ++e)
- // facem_errors[2] +=
- // std::abs(fe_evalm.get_hessian(q)[d][e][j]-
- // function.hessian(q)[d][e]);
}
- // std::cout << std::endl;
}
}
}
{
fe_evalm.reinit(face);
fe_evalm.read_dof_values(src);
- fe_evalm.evaluate(EvaluationFlags::values | EvaluationFlags::gradients);
+ fe_evalm.evaluate(EvaluationFlags::values | EvaluationFlags::gradients |
+ EvaluationFlags::hessians);
for (unsigned int j = 0; j < VectorizedArray<Number>::size(); ++j)
{
boundary_errors[1] +=
std::abs(fe_evalm.get_gradient(q)[d][j] -
function.gradient(p, 0)[d]);
+ for (unsigned int d = 0; d < dim; ++d)
+ for (unsigned int e = 0; e < dim; ++e)
+ boundary_errors[2] +=
+ std::abs(fe_evalm.get_hessian(q)[d][e][j] -
+ function.hessian(p, 0)[d][e]);
double normal_derivative = 0;
for (unsigned int d = 0; d < dim; ++d)
normal_derivative += function.gradient(p, 0)[d] *
fe_evalm.get_normal_vector(q)[d][j];
boundary_errors[3] += std::abs(
fe_evalm.get_normal_derivative(q)[j] - normal_derivative);
- // hessians not currently implemented in FEFaceEvaluation
- // for (unsigned int d=0; d<dim; ++d)
- // for (unsigned int e=0; e<dim; ++e)
- // boundary_errors[2] +=
- // std::abs(fe_evalm.get_hessian(q)[d][e][j]-
- // function.hessian(q)[d][e]);
}
}
}
<< std::endl;
deallog << "Error face- gradients: " << facem_errors[1] / facem_times
<< std::endl;
- // deallog << "Error face- Hessians: " << facem_errors[2]/facem_times <<
- // std::endl;
+ deallog << "Error face- Hessians: " << facem_errors[2] / facem_times
+ << std::endl;
deallog << "Error face- grad*normal:" << facem_errors[3] / facem_times
<< std::endl;
deallog << "Error face+ values: " << facep_errors[0] / facep_times
<< std::endl;
deallog << "Error face+ gradients: " << facep_errors[1] / facep_times
<< std::endl;
- // deallog << "Error face+ Hessians: " << facep_errors[2]/facep_times <<
- // std::endl;
- deallog << "Error face+ grad*normal:" << facem_errors[3] / facep_times
+ deallog << "Error face+ Hessians: " << facep_errors[2] / facep_times
+ << std::endl;
+ deallog << "Error face+ grad*normal:" << facep_errors[3] / facep_times
<< std::endl;
deallog << "Error face b values: " << boundary_errors[0] / boundary_times
<< std::endl;
deallog << "Error face b gradients: " << boundary_errors[1] / boundary_times
<< std::endl;
- // deallog << "Error face b Hessians: " <<
- // boundary_errors[2]/boundary_times << std::endl;
+ deallog << "Error face b Hessians: " << boundary_errors[2] / boundary_times
+ << std::endl;
deallog << "Error face b grad*norm: " << boundary_errors[3] / boundary_times
<< std::endl;
deallog << std::endl;
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:1d::Testing FE_DGQ<1>(2)
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:1d::Test without templates on FEEvaluation
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:2d::Testing FE_Q<2>(1)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(1)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_Q<2>(2)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(2)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_Q<2>(3)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(3)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_Q<2>(4)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(4)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:3d::Testing FE_Q<3>(1)
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Testing FE_DGQ<3>(1)
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Test without templates on FEEvaluation
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Testing FE_Q<3>(2)
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Testing FE_DGQ<3>(2)
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Test without templates on FEEvaluation
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:1d::Testing FE_DGQ<1>(2)
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:1d::Test without templates on FEEvaluation
DEAL:1d::Error cell Hessians: 0
DEAL:1d::Error face- values: 0
DEAL:1d::Error face- gradients: 0
+DEAL:1d::Error face- Hessians: 0
DEAL:1d::Error face- grad*normal:0
DEAL:1d::Error face+ values: 0
DEAL:1d::Error face+ gradients: 0
+DEAL:1d::Error face+ Hessians: 0
DEAL:1d::Error face+ grad*normal:0
DEAL:1d::Error face b values: 0
DEAL:1d::Error face b gradients: 0
+DEAL:1d::Error face b Hessians: 0
DEAL:1d::Error face b grad*norm: 0
DEAL:1d::
DEAL:2d::Testing FE_Q<2>(2)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(2)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_Q<2>(3)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(3)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_Q<2>(4)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Testing FE_DGQ<2>(4)
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:2d::Test without templates on FEEvaluation
DEAL:2d::Error cell Hessians: 0
DEAL:2d::Error face- values: 0
DEAL:2d::Error face- gradients: 0
+DEAL:2d::Error face- Hessians: 0
DEAL:2d::Error face- grad*normal:0
DEAL:2d::Error face+ values: 0
DEAL:2d::Error face+ gradients: 0
+DEAL:2d::Error face+ Hessians: 0
DEAL:2d::Error face+ grad*normal:0
DEAL:2d::Error face b values: 0
DEAL:2d::Error face b gradients: 0
+DEAL:2d::Error face b Hessians: 0
DEAL:2d::Error face b grad*norm: 0
DEAL:2d::
DEAL:3d::Testing FE_Q<3>(2)
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Testing FE_DGQ<3>(2)
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
DEAL:3d::Test without templates on FEEvaluation
DEAL:3d::Error cell Hessians: 0
DEAL:3d::Error face- values: 0
DEAL:3d::Error face- gradients: 0
+DEAL:3d::Error face- Hessians: 0
DEAL:3d::Error face- grad*normal:0
DEAL:3d::Error face+ values: 0
DEAL:3d::Error face+ gradients: 0
+DEAL:3d::Error face+ Hessians: 0
DEAL:3d::Error face+ grad*normal:0
DEAL:3d::Error face b values: 0
DEAL:3d::Error face b gradients: 0
+DEAL:3d::Error face b Hessians: 0
DEAL:3d::Error face b grad*norm: 0
DEAL:3d::
phi.begin_values(),
temp.data(),
false,
+ false,
face);
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2013 - 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+#include <deal.II/base/logstream.h>
+
+#include <deal.II/distributed/tria.h>
+
+#include <deal.II/dofs/dof_handler.h>
+
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_values.h>
+#include <deal.II/fe/mapping_q_generic.h>
+
+#include <deal.II/grid/grid_generator.h>
+
+#include <deal.II/lac/affine_constraints.h>
+#include <deal.II/lac/la_parallel_vector.h>
+
+#include <deal.II/matrix_free/fe_evaluation.h>
+
+#include <deal.II/numerics/vector_tools.h>
+
+#include "../tests.h"
+
+// Tests if applying a matrix vector product with a matrix containing hessians
+// on faces produces the same result with FEFaceEvaluation and FEFaceValues.
+// This is checked for different combinations of EvaluationFlags, FE types and
+// polynomial degrees.
+
+template <int dim>
+void
+test_hessians(const unsigned int degree,
+ const dealii::FE_Poly<dim> & fe,
+ const dealii::EvaluationFlags::EvaluationFlags evaluation_flags)
+{
+ using namespace dealii;
+ using VectorizedArrayType = VectorizedArray<double>;
+
+ Triangulation<dim> tria;
+ GridGenerator::hyper_cube(tria);
+ tria.refine_global(3);
+
+ // TODO: make Hessians work for adaptive refinement
+
+ // unsigned int counter = 0;
+ // for (auto cell = tria.begin_active(); cell != tria.end(); ++cell,
+ // ++counter)
+ // if (cell->is_locally_owned() && counter % 3 == 0)
+ // cell->set_refine_flag();
+ // tria.execute_coarsening_and_refinement();
+
+ DoFHandler<dim> dof_handler(tria);
+ dof_handler.distribute_dofs(fe);
+ MappingQGeneric<dim> mapping(1);
+
+ AffineConstraints<double> constraints;
+ VectorTools::interpolate_boundary_values(
+ mapping, dof_handler, 0, Functions::ZeroFunction<dim>(), constraints);
+ constraints.close();
+
+ // FEFaceEvaluation
+ typename MatrixFree<dim, double, VectorizedArrayType>::AdditionalData
+ additional_data;
+ additional_data.mapping_update_flags_inner_faces =
+ update_values | update_gradients | update_hessians;
+
+ MatrixFree<dim, double, VectorizedArrayType> matrix_free;
+ matrix_free.reinit(
+ mapping, dof_handler, constraints, QGauss<1>(degree + 1), additional_data);
+
+ if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
+ deallog << "Working with " << fe.get_name() << " and "
+ << dof_handler.n_dofs() << " dofs" << std::endl;
+
+ LinearAlgebra::distributed::Vector<double> src, dst, dst2;
+ matrix_free.initialize_dof_vector(src);
+ for (auto &v : src)
+ v = random_value<double>();
+
+ matrix_free.initialize_dof_vector(dst);
+ matrix_free.initialize_dof_vector(dst2);
+
+ // Setup FEFaceValues
+ QGauss<dim - 1> quad(degree + 1);
+ FEFaceValues<dim> fe_face_values_m(mapping,
+ fe,
+ quad,
+ update_values | update_gradients |
+ update_hessians | update_JxW_values);
+ FEFaceValues<dim> fe_face_values_p(mapping,
+ fe,
+ quad,
+ update_values | update_gradients |
+ update_hessians | update_JxW_values);
+ Vector<double> solution_values_m(fe.dofs_per_cell);
+ Vector<double> solution_values_p(fe.dofs_per_cell);
+ std::vector<Tensor<2, dim>> solution_hessians(quad.size());
+ std::vector<Tensor<1, dim>> solution_gradients(quad.size());
+ std::vector<double> solution_values(quad.size());
+ std::vector<types::global_dof_index> dof_indices_m(fe.dofs_per_cell);
+ std::vector<types::global_dof_index> dof_indices_p(fe.dofs_per_cell);
+ src.update_ghost_values();
+ dst2 = 0;
+
+
+
+ matrix_free.template loop<LinearAlgebra::distributed::Vector<double>,
+ LinearAlgebra::distributed::Vector<double>>(
+ [&](
+ const auto &matrix_free, auto &dst, const auto &src, const auto &range) {
+ (void)matrix_free;
+ (void)dst;
+ (void)src;
+ (void)range;
+ },
+ [&](
+ const auto &matrix_free, auto &dst, const auto &src, const auto &range) {
+ FEFaceEvaluation<dim, -1, 0, 1, double, VectorizedArrayType> fe_eval_m(
+ matrix_free, true);
+ FEFaceEvaluation<dim, -1, 0, 1, double, VectorizedArrayType> fe_eval_p(
+ matrix_free, false);
+ for (unsigned int face = range.first; face < range.second; ++face)
+ {
+ // FEFaceEvaluation
+ fe_eval_m.reinit(face);
+ fe_eval_m.gather_evaluate(src, evaluation_flags);
+ for (unsigned int q = 0; q < fe_eval_m.n_q_points; ++q)
+ {
+ if (evaluation_flags & EvaluationFlags::hessians)
+ fe_eval_m.submit_hessian(fe_eval_m.get_hessian(q), q);
+ if (evaluation_flags & EvaluationFlags::gradients)
+ fe_eval_m.submit_gradient(fe_eval_m.get_gradient(q), q);
+ if (evaluation_flags & EvaluationFlags::values)
+ fe_eval_m.submit_value(fe_eval_m.get_value(q), q);
+ }
+ fe_eval_m.integrate(evaluation_flags);
+ fe_eval_m.distribute_local_to_global(dst);
+
+ fe_eval_p.reinit(face);
+ fe_eval_p.gather_evaluate(src, evaluation_flags);
+ for (unsigned int q = 0; q < fe_eval_p.n_q_points; ++q)
+ {
+ if (evaluation_flags & EvaluationFlags::hessians)
+ fe_eval_p.submit_hessian(fe_eval_p.get_hessian(q), q);
+ if (evaluation_flags & EvaluationFlags::gradients)
+ fe_eval_p.submit_gradient(fe_eval_p.get_gradient(q), q);
+ if (evaluation_flags & EvaluationFlags::values)
+ fe_eval_p.submit_value(fe_eval_p.get_value(q), q);
+ }
+ fe_eval_p.integrate(evaluation_flags);
+ fe_eval_p.distribute_local_to_global(dst);
+
+ // FEFaceValues
+ for (unsigned int v = 0;
+ v < matrix_free.n_active_entries_per_face_batch(face);
+ ++v)
+ {
+ const auto face_accessor_inside =
+ matrix_free.get_face_iterator(face, v, true);
+ const auto face_accessor_outside =
+ matrix_free.get_face_iterator(face, v, false);
+
+ fe_face_values_m.reinit(face_accessor_inside.first,
+ face_accessor_inside.second);
+
+ fe_face_values_p.reinit(face_accessor_outside.first,
+ face_accessor_outside.second);
+
+ face_accessor_inside.first->get_dof_indices(dof_indices_m);
+ face_accessor_outside.first->get_dof_indices(dof_indices_p);
+ constraints.get_dof_values(src,
+ dof_indices_m.begin(),
+ solution_values_m.begin(),
+ solution_values_m.end());
+
+ constraints.get_dof_values(src,
+ dof_indices_p.begin(),
+ solution_values_p.begin(),
+ solution_values_p.end());
+
+ for (unsigned int q = 0; q < quad.size(); ++q)
+ {
+ double values = 0.;
+ Tensor<1, dim> gradients;
+ Tensor<2, dim> hessians;
+ for (unsigned int i = 0; i < fe.dofs_per_cell; ++i)
+ {
+ if (evaluation_flags & EvaluationFlags::hessians)
+ hessians += solution_values_m(i) *
+ fe_face_values_m.shape_hessian(i, q);
+ if (evaluation_flags & EvaluationFlags::gradients)
+ gradients += solution_values_m(i) *
+ fe_face_values_m.shape_grad(i, q);
+ if (evaluation_flags & EvaluationFlags::values)
+ values += solution_values_m(i) *
+ fe_face_values_m.shape_value(i, q);
+ }
+ solution_hessians[q] = hessians * fe_face_values_m.JxW(q);
+ solution_gradients[q] = gradients * fe_face_values_m.JxW(q);
+ solution_values[q] = values * fe_face_values_m.JxW(q);
+ }
+ for (unsigned int i = 0; i < fe.dofs_per_cell; ++i)
+ {
+ double sum_hessians = 0.;
+ double sum_gradients = 0.;
+ double sum_values = 0.;
+ for (unsigned int q = 0; q < quad.size(); ++q)
+ {
+ if (evaluation_flags & EvaluationFlags::hessians)
+ sum_hessians += double_contract<0, 0, 1, 1>(
+ solution_hessians[q],
+ fe_face_values_m.shape_hessian(i, q));
+ if (evaluation_flags & EvaluationFlags::gradients)
+ sum_gradients += solution_gradients[q] *
+ fe_face_values_m.shape_grad(i, q);
+ if (evaluation_flags & EvaluationFlags::values)
+ sum_values += solution_values[q] *
+ fe_face_values_m.shape_value(i, q);
+ }
+ solution_values_m(i) =
+ sum_hessians + sum_gradients + sum_values;
+ }
+ constraints.distribute_local_to_global(solution_values_m,
+ dof_indices_m,
+ dst2);
+
+ for (unsigned int q = 0; q < quad.size(); ++q)
+ {
+ double values = 0.;
+ Tensor<1, dim> gradients;
+ Tensor<2, dim> hessians;
+ for (unsigned int i = 0; i < fe.dofs_per_cell; ++i)
+ {
+ if (evaluation_flags & EvaluationFlags::hessians)
+ hessians += solution_values_p(i) *
+ fe_face_values_p.shape_hessian(i, q);
+ if (evaluation_flags & EvaluationFlags::gradients)
+ gradients += solution_values_p(i) *
+ fe_face_values_p.shape_grad(i, q);
+ if (evaluation_flags & EvaluationFlags::values)
+ values += solution_values_p(i) *
+ fe_face_values_p.shape_value(i, q);
+ }
+ solution_hessians[q] = hessians * fe_face_values_p.JxW(q);
+ solution_gradients[q] = gradients * fe_face_values_p.JxW(q);
+ solution_values[q] = values * fe_face_values_p.JxW(q);
+ }
+ for (unsigned int i = 0; i < fe.dofs_per_cell; ++i)
+ {
+ double sum_hessians = 0.;
+ double sum_gradients = 0.;
+ double sum_values = 0.;
+ for (unsigned int q = 0; q < quad.size(); ++q)
+ {
+ if (evaluation_flags & EvaluationFlags::hessians)
+ sum_hessians += double_contract<0, 0, 1, 1>(
+ solution_hessians[q],
+ fe_face_values_p.shape_hessian(i, q));
+ if (evaluation_flags & EvaluationFlags::gradients)
+ sum_gradients += solution_gradients[q] *
+ fe_face_values_p.shape_grad(i, q);
+ if (evaluation_flags & EvaluationFlags::values)
+ sum_values += solution_values[q] *
+ fe_face_values_p.shape_value(i, q);
+ }
+ solution_values_p(i) =
+ sum_hessians + sum_gradients + sum_values;
+ }
+ constraints.distribute_local_to_global(solution_values_p,
+ dof_indices_p,
+ dst2);
+ }
+ }
+ },
+ [&](
+ const auto &matrix_free, auto &dst, const auto &src, const auto &range) {
+ (void)matrix_free;
+ (void)dst;
+ (void)src;
+ (void)range;
+ },
+ dst,
+ src,
+ true);
+
+ const unsigned int end_of_print_dst =
+ dof_handler.n_dofs() > 9 ? 9 : dof_handler.n_dofs();
+
+ deallog << "dst FEE: ";
+ for (unsigned int i = 0; i < end_of_print_dst; ++i)
+ deallog << dst[i] << " ";
+ deallog << std::endl;
+
+ deallog << "dst FEV: ";
+ for (unsigned int i = 0; i < end_of_print_dst; ++i)
+ deallog << dst2[i] << " ";
+ deallog << std::endl;
+
+ // compare solutions of matrix vector product
+ {
+ dst2 -= dst;
+
+ double error = 0.;
+ if (dst.l2_norm() > 0)
+ error = dst2.l2_norm() / dst.l2_norm();
+ else
+ error = dst2.l2_norm();
+
+ if (Utilities::MPI::this_mpi_process(MPI_COMM_WORLD) == 0)
+ deallog << "FEValues verification: " << error << std::endl << std::endl;
+ }
+}
+
+int
+main(int argc, char **argv)
+{
+ dealii::Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1);
+
+ initlog();
+ deallog << std::setprecision(10);
+
+ {
+ dealii::EvaluationFlags::EvaluationFlags evaluation_flags =
+ EvaluationFlags::hessians;
+ deallog << "test_hessians_only" << std::endl;
+ for (unsigned int i = 1; i < 4; ++i)
+ {
+ test_hessians<1>(i, dealii::FE_Q<1>(i), evaluation_flags);
+ test_hessians<2>(i, dealii::FE_Q<2>(i), evaluation_flags);
+ test_hessians<3>(i, dealii::FE_Q<3>(i), evaluation_flags);
+ test_hessians<1>(i, dealii::FE_DGQ<1>(i), evaluation_flags);
+ test_hessians<2>(i, dealii::FE_DGQ<2>(i), evaluation_flags);
+ test_hessians<3>(i, dealii::FE_DGQ<3>(i), evaluation_flags);
+ }
+ }
+
+ {
+ dealii::EvaluationFlags::EvaluationFlags evaluation_flags =
+ EvaluationFlags::values | EvaluationFlags::hessians;
+ deallog << "test_hessians_with_values" << std::endl;
+ for (unsigned int i = 1; i < 4; ++i)
+ {
+ test_hessians<1>(i, dealii::FE_Q<1>(i), evaluation_flags);
+ test_hessians<2>(i, dealii::FE_Q<2>(i), evaluation_flags);
+ test_hessians<3>(i, dealii::FE_Q<3>(i), evaluation_flags);
+ test_hessians<1>(i, dealii::FE_DGQ<1>(i), evaluation_flags);
+ test_hessians<2>(i, dealii::FE_DGQ<2>(i), evaluation_flags);
+ test_hessians<3>(i, dealii::FE_DGQ<3>(i), evaluation_flags);
+ }
+ }
+
+ {
+ dealii::EvaluationFlags::EvaluationFlags evaluation_flags =
+ EvaluationFlags::gradients | EvaluationFlags::hessians;
+ deallog << "test_hessians_with_gradients" << std::endl;
+ for (unsigned int i = 1; i < 4; ++i)
+ {
+ test_hessians<1>(i, dealii::FE_Q<1>(i), evaluation_flags);
+ test_hessians<2>(i, dealii::FE_Q<2>(i), evaluation_flags);
+ test_hessians<3>(i, dealii::FE_Q<3>(i), evaluation_flags);
+ test_hessians<1>(i, dealii::FE_DGQ<1>(i), evaluation_flags);
+ test_hessians<2>(i, dealii::FE_DGQ<2>(i), evaluation_flags);
+ test_hessians<3>(i, dealii::FE_DGQ<3>(i), evaluation_flags);
+ }
+ }
+
+ {
+ dealii::EvaluationFlags::EvaluationFlags evaluation_flags =
+ EvaluationFlags::values | EvaluationFlags::gradients |
+ EvaluationFlags::hessians;
+ deallog << "test_hessians_with_gradients_and_values" << std::endl;
+ for (unsigned int i = 1; i < 4; ++i)
+ {
+ test_hessians<1>(i, dealii::FE_Q<1>(i), evaluation_flags);
+ test_hessians<2>(i, dealii::FE_Q<2>(i), evaluation_flags);
+ test_hessians<3>(i, dealii::FE_Q<3>(i), evaluation_flags);
+ test_hessians<1>(i, dealii::FE_DGQ<1>(i), evaluation_flags);
+ test_hessians<2>(i, dealii::FE_DGQ<2>(i), evaluation_flags);
+ test_hessians<3>(i, dealii::FE_DGQ<3>(i), evaluation_flags);
+ }
+ }
+}
--- /dev/null
+
+DEAL::test_hessians_only
+DEAL::Working with FE_Q<1>(1) and 9 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(1) and 81 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 -3961.844511 0.000000000 3975.516664 0.000000000 2180.676992 4279.873135
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 -3961.844511 0.000000000 3975.516664 0.000000000 2180.676992 4279.873135
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(1) and 729 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 2671.158976 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 2671.158976 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(1) and 16 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(1) and 256 dofs
+DEAL::dst FEE: 186.7674628 -186.7674628 -186.7674628 186.7674628 1095.845132 -1095.845132 -1095.845132 1095.845132 -1770.726878
+DEAL::dst FEV: 186.7674628 -186.7674628 -186.7674628 186.7674628 1095.845132 -1095.845132 -1095.845132 1095.845132 -1770.726878
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(1) and 4096 dofs
+DEAL::dst FEE: -41.82729979 116.7876889 -72.84484080 -2.115548314 24.66611142 -99.62650053 90.00602917 -15.04564006 -135.6029421
+DEAL::dst FEV: -41.82729979 116.7876889 -72.84484080 -2.115548314 24.66611142 -99.62650053 90.00602917 -15.04564006 -135.6029421
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<1>(2) and 17 dofs
+DEAL::dst FEE: 0.000000000 67338.98941 16221.91627 174612.5300 -150899.8951 93741.61943 -198325.1649 100444.2774 10841.92603
+DEAL::dst FEV: 0.000000000 67338.98941 16221.91627 174612.5300 -150899.8951 93741.61943 -198325.1649 100444.2774 10841.92603
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(2) and 289 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 -104920.8906 0.000000000 -49521.70596 0.000000000 150903.6170 6452.694464
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 -104920.8906 0.000000000 -49521.70596 0.000000000 150903.6170 6452.694464
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(2) and 4913 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 31698.05993 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 31698.05993 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(2) and 24 dofs
+DEAL::dst FEE: -57370.44455 114740.8891 -57370.44455 -130179.7286 260359.4573 -130179.7286 -132118.6579 264237.3158 -132118.6579
+DEAL::dst FEV: -57370.44455 114740.8891 -57370.44455 -130179.7286 260359.4573 -130179.7286 -132118.6579 264237.3158 -132118.6579
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(2) and 576 dofs
+DEAL::dst FEE: 2621.950491 -1575.677282 -3332.535102 -31645.71115 63719.74128 -27501.50635 24104.29342 -52305.12954 25914.57422
+DEAL::dst FEV: 2621.950491 -1575.677282 -3332.535102 -31645.71115 63719.74128 -27501.50635 24104.29342 -52305.12954 25914.57422
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(2) and 13824 dofs
+DEAL::dst FEE: 50.88230764 -830.0219522 -818.2722402 -1019.544361 3222.190645 -744.7435660 811.4995692 1670.530726 -893.6889450
+DEAL::dst FEV: 50.88230764 -830.0219522 -818.2722402 -1019.544361 3222.190645 -744.7435660 811.4995692 1670.530726 -893.6889450
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<1>(3) and 25 dofs
+DEAL::dst FEE: 0.000000000 -737280.1891 -22080.92721 29817.48585 -1201734.837 1302772.885 -799181.1739 -1342610.859 3103193.420
+DEAL::dst FEV: 0.000000000 -737280.1891 -22080.92721 29817.48585 -1201734.837 1302772.885 -799181.1739 -1342610.859 3103193.420
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(3) and 625 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 -1547335.199 0.000000000 0.000000000 368454.5924 67296.69792 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 -1547335.199 0.000000000 0.000000000 368454.5924 67296.69792 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(3) and 15625 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 70210.21641 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 70210.21641 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(3) and 32 dofs
+DEAL::dst FEE: -254241.5072 725631.1855 -979872.6927 508483.0143 369638.5526 -633813.2719 321979.8414 -57805.12216 -296675.7776
+DEAL::dst FEV: -254241.5072 725631.1855 -979872.6927 508483.0143 369638.5526 -633813.2719 321979.8414 -57805.12216 -296675.7776
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(3) and 1024 dofs
+DEAL::dst FEE: 19404.30284 11401.17609 53539.65109 -156046.9045 -56073.49671 -37352.39201 -103984.8509 359118.5689 51864.00926
+DEAL::dst FEV: 19404.30284 11401.17609 53539.65109 -156046.9045 -56073.49671 -37352.39201 -103984.8509 359118.5689 51864.00926
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(3) and 32768 dofs
+DEAL::dst FEE: 118.4885551 -4107.177606 1812.382519 -363.6809581 -2068.812140 9700.899062 -6925.304917 4256.899099 -2307.998977
+DEAL::dst FEV: 118.4885551 -4107.177606 1812.382519 -363.6809581 -2068.812140 9700.899062 -6925.304917 4256.899099 -2307.998977
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::test_hessians_with_values
+DEAL::Working with FE_Q<1>(1) and 9 dofs
+DEAL::dst FEE: 0.000000000 1.485241358 0.4717695818 1.336048384 0.09186878712 0.5541558529 0.6274125085 0.1848661360 0.000000000
+DEAL::dst FEV: 0.000000000 1.485241358 0.4717695818 1.336048384 0.09186878712 0.5541558529 0.6274125085 0.1848661360 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(1) and 81 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 1070.969383 0.000000000 2878.096119 0.000000000 -3268.391675 -412.5247158
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 1070.969383 0.000000000 2878.096119 0.000000000 -3268.391675 -412.5247158
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(1) and 729 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 3500.106392 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 3500.106392 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(1) and 16 dofs
+DEAL::dst FEE: 0.000000000 0.9178876015 0.3258156191 0.6197056345 0.1537631360 0.2945842553 0.3555839413 0.1329097907 0.7070026150
+DEAL::dst FEV: 0.000000000 0.9178876015 0.3258156191 0.6197056345 0.1537631360 0.2945842553 0.3555839413 0.1329097907 0.7070026150
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(1) and 256 dofs
+DEAL::dst FEE: -569.5871632 569.6347435 569.6092785 -569.5287350 -1511.130148 1511.164480 1511.202041 -1511.108327 -2579.882528
+DEAL::dst FEV: -569.5871632 569.6347435 569.6092785 -569.5287350 -1511.130148 1511.164480 1511.202041 -1511.108327 -2579.882528
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(1) and 4096 dofs
+DEAL::dst FEE: 7.424074176 241.6432657 -274.2813513 25.22374027 79.00171305 -328.0612417 187.8648553 61.21168368 78.30753564
+DEAL::dst FEV: 7.424074176 241.6432657 -274.2813513 25.22374027 79.00171305 -328.0612417 187.8648553 61.21168368 78.30753564
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<1>(2) and 17 dofs
+DEAL::dst FEE: 0.000000000 -167135.6577 229538.9869 -97797.09182 104732.8852 -39369.24028 90862.39763 -109611.4362 -12122.77601
+DEAL::dst FEV: 0.000000000 -167135.6577 229538.9869 -97797.09182 104732.8852 -39369.24028 90862.39763 -109611.4362 -12122.77601
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(2) and 289 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 60209.63375 0.000000000 69766.67093 0.000000000 -58791.04965 -39753.16339
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 60209.63375 0.000000000 69766.67093 0.000000000 -58791.04965 -39753.16339
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(2) and 4913 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 12306.73607 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 12306.73607 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(2) and 24 dofs
+DEAL::dst FEE: 47813.68795 -95627.37590 47814.42813 -96189.12898 192379.2740 -96189.51817 55225.52064 -110450.7480 55226.32803
+DEAL::dst FEV: 47813.68795 -95627.37590 47814.42813 -96189.12898 192379.2740 -96189.51817 55225.52064 -110450.7480 55226.32803
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(2) and 576 dofs
+DEAL::dst FEE: -5287.684869 6572.726856 -20184.16775 24229.94350 -48365.22252 61933.61384 -25037.39476 53982.83575 -47844.48685
+DEAL::dst FEV: -5287.684869 6572.726856 -20184.16775 24229.94350 -48365.22252 61933.61384 -25037.39476 53982.83575 -47844.48685
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(2) and 13824 dofs
+DEAL::dst FEE: 370.0566403 -292.3081868 -720.4052507 -622.7558219 238.6508182 1862.634403 281.8888045 2091.765602 -1871.063632
+DEAL::dst FEV: 370.0566403 -292.3081868 -720.4052507 -622.7558219 238.6508182 1862.634403 281.8888045 2091.765602 -1871.063632
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<1>(3) and 25 dofs
+DEAL::dst FEE: 0.000000000 299377.0807 289587.9028 -391051.6578 -3248169.708 -134872.5658 4149.691300 2108415.055 7239160.344
+DEAL::dst FEV: 0.000000000 299377.0807 289587.9028 -391051.6578 -3248169.708 -134872.5658 4149.691300 2108415.055 7239160.344
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(3) and 625 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 -293073.3740 0.000000000 0.000000000 675561.0838 -427952.8498 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 -293073.3740 0.000000000 0.000000000 675561.0838 -427952.8498 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(3) and 15625 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 -85712.86767 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 -85712.86767 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(3) and 32 dofs
+DEAL::dst FEE: -252775.6199 721447.3938 -974223.0137 505551.7303 1611958.374 -3761443.564 4015469.154 -1865982.968 -51596.77459
+DEAL::dst FEV: -252775.6199 721447.3938 -974223.0137 505551.7303 1611958.374 -3761443.564 4015469.154 -1865982.968 -51596.77459
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(3) and 1024 dofs
+DEAL::dst FEE: -10040.20263 53844.53594 -182876.7764 97992.24016 38492.66295 -176621.9667 403072.9022 -164030.7460 150191.7481
+DEAL::dst FEV: -10040.20263 53844.53594 -182876.7764 97992.24016 38492.66295 -176621.9667 403072.9022 -164030.7460 150191.7481
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(3) and 32768 dofs
+DEAL::dst FEE: -674.3169640 814.0460283 -4216.899024 5057.280398 2850.106346 -2615.447288 14600.65907 -8650.675187 -3053.718113
+DEAL::dst FEV: -674.3169640 814.0460283 -4216.899024 5057.280398 2850.106346 -2615.447288 14600.65907 -8650.675187 -3053.718113
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::test_hessians_with_gradients
+DEAL::Working with FE_Q<1>(1) and 9 dofs
+DEAL::dst FEE: 0.000000000 64.15179569 36.90356960 -40.28644325 -68.37967620 39.09859355 46.88549312 -62.19046491 31.72474919
+DEAL::dst FEV: 0.000000000 64.15179569 36.90356960 -40.28644325 -68.37967620 39.09859355 46.88549312 -62.19046491 31.72474919
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(1) and 81 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 442.0107683 0.000000000 3542.943855 0.000000000 -9575.588132 7499.031588
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 442.0107683 0.000000000 3542.943855 0.000000000 -9575.588132 7499.031588
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(1) and 729 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 -2821.551619 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 -2821.551619 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(1) and 16 dofs
+DEAL::dst FEE: 10.06237412 -10.06237412 -27.76560022 27.76560022 34.04355665 -34.04355665 -21.65914739 21.65914739 -7.379041318
+DEAL::dst FEV: 10.06237412 -10.06237412 -27.76560022 27.76560022 34.04355665 -34.04355665 -21.65914739 21.65914739 -7.379041318
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(1) and 256 dofs
+DEAL::dst FEE: -1170.185020 1177.791704 1171.983635 -1179.590319 -370.5322257 366.0001104 361.6813287 -357.1492134 373.2952228
+DEAL::dst FEV: -1170.185020 1177.791704 1171.983635 -1179.590319 -370.5322257 366.0001104 361.6813287 -357.1492134 373.2952228
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(1) and 4096 dofs
+DEAL::dst FEE: -21.58586802 226.0914318 181.3343315 -385.5024982 28.95060815 -234.0415895 -189.5067703 394.2603544 -219.6450550
+DEAL::dst FEV: -21.58586802 226.0914318 181.3343315 -385.5024982 28.95060815 -234.0415895 -189.5067703 394.2603544 -219.6450550
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<1>(2) and 17 dofs
+DEAL::dst FEE: 0.000000000 22439.52331 189626.7649 11978.51334 -234916.2278 -265869.0386 211155.8797 -262006.1754 320433.1750
+DEAL::dst FEV: 0.000000000 22439.52331 189626.7649 11978.51334 -234916.2278 -265869.0386 211155.8797 -262006.1754 320433.1750
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(2) and 289 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 -10856.65148 0.000000000 -25321.75472 0.000000000 95280.37568 4511.212957
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 -10856.65148 0.000000000 -25321.75472 0.000000000 95280.37568 4511.212957
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(2) and 4913 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 9080.423198 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 9080.423198 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(2) and 24 dofs
+DEAL::dst FEE: 9696.728295 -19428.05746 9731.329162 -57259.37641 114443.5555 -57184.17906 159527.5899 -319143.7449 159616.1550
+DEAL::dst FEV: 9696.728295 -19428.05746 9731.329162 -57259.37641 114443.5555 -57184.17906 159527.5899 -319143.7449 159616.1550
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(2) and 576 dofs
+DEAL::dst FEE: 10194.28422 -35453.44115 24293.82892 -45055.92628 130503.0063 -83489.42850 30711.39427 -86727.69612 55023.97831
+DEAL::dst FEV: 10194.28422 -35453.44115 24293.82892 -45055.92628 130503.0063 -83489.42850 30711.39427 -86727.69612 55023.97831
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(2) and 13824 dofs
+DEAL::dst FEE: 448.1446413 -510.9189963 1634.494590 755.7615764 -6223.772600 800.6586105 529.4850620 1959.380310 711.9024984
+DEAL::dst FEV: 448.1446413 -510.9189963 1634.494590 755.7615764 -6223.772600 800.6586105 529.4850620 1959.380310 711.9024984
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<1>(3) and 25 dofs
+DEAL::dst FEE: 0.000000000 2252292.265 1573930.091 -2126197.154 -174516.4842 -2501369.733 2393718.162 -640706.0941 -2054230.489
+DEAL::dst FEV: 0.000000000 2252292.265 1573930.091 -2126197.154 -174516.4842 -2501369.733 2393718.162 -640706.0941 -2054230.489
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(3) and 625 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 -1317837.470 0.000000000 0.000000000 -264992.6854 1082118.987 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 -1317837.470 0.000000000 0.000000000 -264992.6854 1082118.987 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(3) and 15625 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 -65742.12001 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 -65742.12001 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(3) and 32 dofs
+DEAL::dst FEE: 607362.1400 -1733480.886 2340968.548 -1214849.802 -1047359.729 2045725.963 -1566917.253 568551.0199 -3864497.016
+DEAL::dst FEV: 607362.1400 -1733480.886 2340968.548 -1214849.802 -1047359.729 2045725.963 -1566917.253 568551.0199 -3864497.016
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(3) and 1024 dofs
+DEAL::dst FEE: -27040.61278 3573.267393 57805.44291 58432.19141 166809.4545 -257530.0470 225343.1275 -358312.5442 -334186.2219
+DEAL::dst FEV: -27040.61278 3573.267393 57805.44291 58432.19141 166809.4545 -257530.0470 225343.1275 -358312.5442 -334186.2219
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(3) and 32768 dofs
+DEAL::dst FEE: -464.3635240 2028.142271 763.8205053 -549.6966376 327.0786350 -596.5268953 -20821.36159 10310.01212 -817.8506617
+DEAL::dst FEV: -464.3635240 2028.142271 763.8205053 -549.6966376 327.0786350 -596.5268953 -20821.36159 10310.01212 -817.8506617
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::test_hessians_with_gradients_and_values
+DEAL::Working with FE_Q<1>(1) and 9 dofs
+DEAL::dst FEE: 0.000000000 -37.00478557 48.59625945 0.8581540048 -7.323365868 59.63530424 -26.02881961 -15.72262317 -0.3217377366
+DEAL::dst FEV: 0.000000000 -37.00478557 48.59625945 0.8581540048 -7.323365868 59.63530424 -26.02881961 -15.72262317 -0.3217377366
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(1) and 81 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 9507.556791 0.000000000 -7321.070916 0.000000000 -12431.59166 12502.87908
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 9507.556791 0.000000000 -7321.070916 0.000000000 -12431.59166 12502.87908
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(1) and 729 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 280.4491881 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 280.4491881 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(1) and 16 dofs
+DEAL::dst FEE: 39.29021196 -39.01623275 -9.963182797 11.39805324 -15.84908361 16.43351670 80.59747247 -79.58483093 53.75457114
+DEAL::dst FEV: 39.29021196 -39.01623275 -9.963182797 11.39805324 -15.84908361 16.43351670 80.59747247 -79.58483093 53.75457114
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(1) and 256 dofs
+DEAL::dst FEE: 1516.712105 -1527.764366 -1523.902913 1535.086538 -1666.971584 1682.840844 1665.583768 -1681.331729 -76.91218429
+DEAL::dst FEV: 1516.712105 -1527.764366 -1523.902913 1535.086538 -1666.971584 1682.840844 1665.583768 -1681.331729 -76.91218429
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(1) and 4096 dofs
+DEAL::dst FEE: -1.189912404 234.8477478 180.7860225 -413.6249483 -121.5440256 -112.1298409 -57.51853159 290.3907440 289.0256656
+DEAL::dst FEV: -1.189912404 234.8477478 180.7860225 -413.6249483 -121.5440256 -112.1298409 -57.51853159 290.3907440 289.0256656
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<1>(2) and 17 dofs
+DEAL::dst FEE: 0.000000000 29141.70579 65351.23242 -37541.71163 -123874.9385 15975.67678 199194.1660 179480.0205 -231343.3646
+DEAL::dst FEV: 0.000000000 29141.70579 65351.23242 -37541.71163 -123874.9385 15975.67678 199194.1660 179480.0205 -231343.3646
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(2) and 289 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 5622.410746 0.000000000 -124785.4497 0.000000000 -28083.28907 140221.5951
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 5622.410746 0.000000000 -124785.4497 0.000000000 -28083.28907 140221.5951
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(2) and 4913 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 4030.831352 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 4030.831352 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(2) and 24 dofs
+DEAL::dst FEE: 1260.253162 -2590.873853 1331.260025 -26516.06494 53162.96332 -26645.70214 192196.9013 -384333.9650 192138.6030
+DEAL::dst FEV: 1260.253162 -2590.873853 1331.260025 -26516.06494 53162.96332 -26645.70214 192196.9013 -384333.9650 192138.6030
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(2) and 576 dofs
+DEAL::dst FEE: -5169.333357 29063.34663 -23105.45307 24704.96903 -97773.79532 71498.06891 -4767.981572 39154.65283 -33604.35418
+DEAL::dst FEV: -5169.333357 29063.34663 -23105.45307 24704.96903 -97773.79532 71498.06891 -4767.981572 39154.65283 -33604.35418
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(2) and 13824 dofs
+DEAL::dst FEE: -95.58676656 1276.489009 -1334.575175 274.9405999 -4042.776356 2530.457071 1137.743263 1405.963407 -1676.022846
+DEAL::dst FEV: -95.58676656 1276.489009 -1334.575175 274.9405999 -4042.776356 2530.457071 1137.743263 1405.963407 -1676.022846
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<1>(3) and 25 dofs
+DEAL::dst FEE: 0.000000000 -3677129.436 -674168.6624 910870.7211 167797.2252 7072904.503 -6925910.014 1212429.567 6407996.850
+DEAL::dst FEV: 0.000000000 -3677129.436 -674168.6624 910870.7211 167797.2252 7072904.503 -6925910.014 1212429.567 6407996.850
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<2>(3) and 625 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 901522.2209 0.000000000 0.000000000 497122.6409 -783742.6748 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 901522.2209 0.000000000 0.000000000 497122.6409 -783742.6748 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_Q<3>(3) and 15625 dofs
+DEAL::dst FEE: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 61931.57932 0.000000000
+DEAL::dst FEV: 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 0.000000000 61931.57932 0.000000000
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<1>(3) and 32 dofs
+DEAL::dst FEE: -561890.2156 1603751.370 -2166648.089 1124787.675 -928689.6365 1915272.283 -1654708.155 668125.8408 -2711210.037
+DEAL::dst FEV: -561890.2156 1603751.370 -2166648.089 1124787.675 -928689.6365 1915272.283 -1654708.155 668125.8408 -2711210.037
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<2>(3) and 1024 dofs
+DEAL::dst FEE: -2674.684464 105342.7462 -38524.91232 -48556.52335 48969.79796 -415911.0458 256625.4211 38926.74483 -131829.3972
+DEAL::dst FEV: -2674.684464 105342.7462 -38524.91232 -48556.52335 48969.79796 -415911.0458 256625.4211 38926.74483 -131829.3972
+DEAL::FEValues verification: 0.000000000
+DEAL::
+DEAL::Working with FE_DGQ<3>(3) and 32768 dofs
+DEAL::dst FEE: 692.2851052 972.7257414 820.4660691 1018.320908 -949.9325546 -10318.08673 14538.53642 -19371.16378 5741.320571
+DEAL::dst FEV: 692.2851052 972.7257414 820.4660691 1018.320908 -949.9325546 -10318.08673 14538.53642 -19371.16378 5741.320571
+DEAL::FEValues verification: 0.000000000
+DEAL::