From: Maximilian Bergbauer Date: Tue, 24 Aug 2021 07:40:34 +0000 (+0200) Subject: Hessains on general cells and faces X-Git-Tag: v9.4.0-rc1~991^2 X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=refs%2Fpull%2F12744%2Fhead;p=dealii.git Hessains on general cells and faces --- diff --git a/doc/news/changes/minor/20210728MaximilianBergbauer b/doc/news/changes/major/20210906MaximilianBergbauer similarity index 64% rename from doc/news/changes/minor/20210728MaximilianBergbauer rename to doc/news/changes/major/20210906MaximilianBergbauer index 1ccf8c26f8..22dbc6a9aa 100644 --- a/doc/news/changes/minor/20210728MaximilianBergbauer +++ b/doc/news/changes/major/20210906MaximilianBergbauer @@ -1,4 +1,4 @@ Added: Evaluation and integration of hessians both on cells and on faces is now available in the matrix-free framework.
-(Maximilian Bergbauer, 2021/07/28) +(Maximilian Bergbauer, Martin Kronbichler, Peter Munch, 2021/09/06) diff --git a/include/deal.II/matrix_free/fe_evaluation.h b/include/deal.II/matrix_free/fe_evaluation.h index 50e61feb1f..972e122f09 100644 --- a/include/deal.II/matrix_free/fe_evaluation.h +++ b/include/deal.II/matrix_free/fe_evaluation.h @@ -1315,6 +1315,21 @@ protected: */ VectorizedArrayType *gradients_quad; + /** + * This field stores the gradients of the finite element function on + * quadrature points after applying unit cell transformations or before + * integrating. The methods get_hessian() and submit_hessian() (as well as + * some specializations like get_hessian_diagonal() or get_laplacian()) + * access this field for general cell/face types. + * + * The values of this array are stored in the start section of + * @p scratch_data_array. Due to its access as a thread local memory, the + * memory can get reused between different calls. As opposed to requesting + * memory on the stack, this approach allows for very large polynomial + * degrees. + */ + VectorizedArrayType *gradients_from_hessians_quad; + /** * This field stores the Hessians of the finite element function on * quadrature points after applying unit cell transformations. The methods @@ -4330,7 +4345,8 @@ FEEvaluationBase:: 2 * n_quadrature_points; const unsigned int allocated_size = shift + n_components_ * dofs_per_component + - (n_components_ * ((dim * (dim + 1)) / 2 + dim + 1) * n_quadrature_points); + (n_components_ * ((dim * (dim + 1)) / 2 + 2 * dim + 1) * + n_quadrature_points); this->scratch_data_array->resize_fast(allocated_size); // set the pointers to the correct position in the data array @@ -4343,12 +4359,16 @@ FEEvaluationBase:: this->scratch_data_array->begin() + n_components * dofs_per_component; gradients_quad = this->scratch_data_array->begin() + n_components * (dofs_per_component + n_quadrature_points); - hessians_quad = + gradients_from_hessians_quad = this->scratch_data_array->begin() + n_components * (dofs_per_component + (dim + 1) * n_quadrature_points); - this->scratch_data = - this->scratch_data_array->begin() + n_components_ * dofs_per_component + - (n_components_ * ((dim * (dim + 1)) / 2 + dim + 1) * n_quadrature_points); + hessians_quad = + this->scratch_data_array->begin() + + n_components * (dofs_per_component + (2 * dim + 1) * n_quadrature_points); + this->scratch_data = this->scratch_data_array->begin() + + n_components_ * dofs_per_component + + (n_components_ * ((dim * (dim + 1)) / 2 + 2 * dim + 1) * + n_quadrature_points); } @@ -6019,8 +6039,6 @@ FEEvaluationBase:: [this->mapping_data->data_index_offsets[this->cell] + q_point]; for (unsigned int comp = 0; comp < n_components; ++comp) { - // compute laplacian before the gradient because it needs to access - // unscaled gradient data VectorizedArrayType tmp[dim][dim]; internal::hessian_unit_times_jac( jac, hessians_quad + comp * hdim * nqp + q_point, nqp, tmp); @@ -6039,7 +6057,7 @@ FEEvaluationBase:: for (unsigned int e = 0; e < dim; ++e) hessian_out[comp][d][d] += jac_grad[d][e] * - gradients_quad[(comp * dim + e) * nqp + q_point]; + gradients_from_hessians_quad[(comp * dim + e) * nqp + q_point]; // add off-diagonal part of J' * grad(u) for (unsigned int d = 0, count = dim; d < dim; ++d) @@ -6047,7 +6065,8 @@ FEEvaluationBase:: for (unsigned int f = 0; f < dim; ++f) hessian_out[comp][d][e] += jac_grad[count][f] * - gradients_quad[(comp * dim + f) * nqp + q_point]; + gradients_from_hessians_quad[(comp * dim + f) * nqp + + q_point]; // take symmetric part for (unsigned int d = 0; d < dim; ++d) @@ -6144,7 +6163,7 @@ FEEvaluationBase:: for (unsigned int e = 0; e < dim; ++e) hessian_out[comp][d] += jac_grad[d][e] * - gradients_quad[(comp * dim + e) * nqp + q_point]; + gradients_from_hessians_quad[(comp * dim + e) * nqp + q_point]; } } return hessian_out; @@ -6389,7 +6408,7 @@ FEEvaluationBase:: hessian_in[comp][d][d] * factor; } - // off diagonal part (use symmetry) + // off diagonal part for (unsigned int d = 1, off_dia = dim; d < dim; ++d) for (unsigned int e = 0; e < d; ++e, ++off_dia) { @@ -6398,7 +6417,7 @@ FEEvaluationBase:: 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; + (hessian_in[comp][d][e] + hessian_in[comp][e][d]) * factor; } } // cell with general Jacobian, but constant within the cell @@ -6413,39 +6432,87 @@ FEEvaluationBase:: VectorizedArrayType tmp[dim][dim]; for (unsigned int i = 0; i < dim; ++i) for (unsigned int j = 0; j < dim; ++j) - tmp[i][j] = 0.; + { + tmp[i][j] = hessian_in[comp][i][0] * jac[0][j]; + for (unsigned int k = 1; k < dim; ++k) + tmp[i][j] += hessian_in[comp][i][k] * jac[k][j]; + } + // 2. hessian_unit = J^T * tmp + VectorizedArrayType tmp2[dim][dim]; 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 + { + tmp2[i][j] = jac[0][i] * tmp[0][j]; + for (unsigned int k = 1; k < dim; ++k) + tmp2[i][j] += jac[k][i] * tmp[k][j]; + } // 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]; - } + hessians_quad[(comp * hdim + d) * nqp + q_point] = tmp2[d][d] * JxW; - // 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]; - } + // off diagonal part + for (unsigned int d = 0, off_diag = dim; d < dim; ++d) + for (unsigned int e = d + 1; e < dim; ++e, ++off_diag) + hessians_quad[(comp * hdim + off_diag) * nqp + q_point] = + (tmp2[d][e] + tmp2[e][d]) * JxW; } } else { - // TODO: implement general path submit_hessian - AssertThrow(false, ExcNotImplemented()); + const Tensor<2, dim, VectorizedArrayType> jac = this->jacobian[q_point]; + const VectorizedArrayType JxW = this->J_value[q_point]; + const auto & jac_grad = + this->mapping_data->jacobian_gradients + [1 - this->is_interior_face] + [this->mapping_data->data_index_offsets[this->cell] + 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] = hessian_in[comp][i][0] * jac[0][j]; + for (unsigned int k = 1; k < dim; ++k) + tmp[i][j] += hessian_in[comp][i][k] * jac[k][j]; + } + + // 2. hessian_unit = J^T * tmp + VectorizedArrayType tmp2[dim][dim]; + for (unsigned int i = 0; i < dim; ++i) + for (unsigned int j = 0; j < dim; ++j) + { + tmp2[i][j] = jac[0][i] * tmp[0][j]; + for (unsigned int k = 1; k < dim; ++k) + tmp2[i][j] += jac[k][i] * tmp[k][j]; + } + + // diagonal part + for (unsigned int d = 0; d < dim; ++d) + hessians_quad[(comp * hdim + d) * nqp + q_point] = tmp2[d][d] * JxW; + + // off diagonal part + for (unsigned int d = 0, off_diag = dim; d < dim; ++d) + for (unsigned int e = d + 1; e < dim; ++e, ++off_diag) + hessians_quad[(comp * hdim + off_diag) * nqp + q_point] = + (tmp2[d][e] + tmp2[e][d]) * JxW; + + // 3. gradient_unit = J' ** hessian_in + for (unsigned int d = 0; d < dim; ++d) + { + VectorizedArrayType sum = 0; + for (unsigned int e = 0; e < dim; ++e) + sum += hessian_in[comp][e][e] * jac_grad[e][d]; + for (unsigned int e = 0, count = dim; e < dim; ++e) + for (unsigned int f = e + 1; f < dim; ++f, ++count) + sum += (hessian_in[comp][e][f] + hessian_in[comp][f][e]) * + jac_grad[count][d]; + gradients_from_hessians_quad[(comp * dim + d) * nqp + q_point] = + sum * JxW; + } + } } } @@ -8281,12 +8348,19 @@ FEEvaluation:: evaluate(const VectorizedArrayType * values_array, - const EvaluationFlags::EvaluationFlags evaluation_flags) + const EvaluationFlags::EvaluationFlags evaluation_flag) { + const bool hessians_on_general_cells = + evaluation_flag & EvaluationFlags::hessians && + (this->cell_type > internal::MatrixFreeFunctions::affine); + EvaluationFlags::EvaluationFlags evaluation_flag_actual = evaluation_flag; + if (hessians_on_general_cells) + evaluation_flag_actual |= EvaluationFlags::gradients; + if (fe_degree > -1) SelectEvaluator:: evaluate(n_components, - evaluation_flags, + evaluation_flag_actual, *this->data, const_cast(values_array), this->values_quad, @@ -8296,7 +8370,7 @@ FEEvaluation::evaluate( n_components, - evaluation_flags, + evaluation_flag_actual, *this->data, const_cast(values_array), this->values_quad, @@ -8304,12 +8378,16 @@ FEEvaluationhessians_quad, this->scratch_data); + if (hessians_on_general_cells) + for (unsigned int i = 0; i < n_components * dim * n_q_points; ++i) + this->gradients_from_hessians_quad[i] = this->gradients_quad[i]; + # ifdef DEBUG - if (evaluation_flags & EvaluationFlags::values) + if (evaluation_flag_actual & EvaluationFlags::values) this->values_quad_initialized = true; - if (evaluation_flags & EvaluationFlags::gradients) + if (evaluation_flag_actual & EvaluationFlags::gradients) this->gradients_quad_initialized = true; - if (evaluation_flags & EvaluationFlags::hessians) + if (evaluation_flag_actual & EvaluationFlags::hessians) this->hessians_quad_initialized = true; # endif } @@ -8678,6 +8756,9 @@ FEEvaluationgradients_quad_submitted == true, internal::ExcAccessToUninitializedField()); + if (integration_flag & EvaluationFlags::hessians) + Assert(this->hessians_quad_submitted == true, + internal::ExcAccessToUninitializedField()); # endif Assert(this->matrix_info != nullptr || this->mapped_geometry->is_initialized(), @@ -8689,10 +8770,28 @@ FEEvaluationcell_type > internal::MatrixFreeFunctions::affine)) + { + unsigned int size = n_components * dim * n_q_points; + if (integration_flag & EvaluationFlags::gradients) + { + for (unsigned int i = 0; i < size; ++i) + this->gradients_quad[i] += this->gradients_from_hessians_quad[i]; + } + else + { + for (unsigned int i = 0; i < size; ++i) + this->gradients_quad[i] = this->gradients_from_hessians_quad[i]; + integration_flag_actual |= EvaluationFlags::gradients; + } + } + if (fe_degree > -1) SelectEvaluator:: integrate(n_components, - integration_flag, + integration_flag_actual, *this->data, values_array, this->values_quad, @@ -8703,7 +8802,7 @@ FEEvaluation::integrate( n_components, - integration_flag, + integration_flag_actual, *this->data, values_array, this->values_quad, @@ -9112,6 +9211,13 @@ FEFaceEvaluationcell_type > internal::MatrixFreeFunctions::affine); + EvaluationFlags::EvaluationFlags evaluation_flag_actual = evaluation_flag; + if (hessians_on_general_cells) + evaluation_flag_actual |= EvaluationFlags::gradients; + if (this->dof_access_index == internal::MatrixFreeFunctions::DoFInfo::dof_access_cell && this->is_interior_face == false) @@ -9140,9 +9246,9 @@ FEFaceEvaluationbegin_gradients(), this->begin_hessians(), this->scratch_data, - evaluation_flag & EvaluationFlags::values, - evaluation_flag & EvaluationFlags::gradients, - evaluation_flag & EvaluationFlags::hessians, + evaluation_flag_actual & EvaluationFlags::values, + evaluation_flag_actual & EvaluationFlags::gradients, + evaluation_flag_actual & EvaluationFlags::hessians, face_nos[0], this->subface_index, face_orientations[0], @@ -9161,9 +9267,9 @@ FEFaceEvaluationbegin_gradients(), this->begin_hessians(), this->scratch_data, - evaluation_flag & EvaluationFlags::values, - evaluation_flag & EvaluationFlags::gradients, - evaluation_flag & EvaluationFlags::hessians, + evaluation_flag_actual & EvaluationFlags::values, + evaluation_flag_actual & EvaluationFlags::gradients, + evaluation_flag_actual & EvaluationFlags::hessians, this->face_no, this->subface_index, this->face_orientation, @@ -9177,21 +9283,28 @@ FEFaceEvaluationbegin_gradients(), this->begin_hessians(), this->scratch_data, - evaluation_flag & EvaluationFlags::values, - evaluation_flag & EvaluationFlags::gradients, - evaluation_flag & EvaluationFlags::hessians, + evaluation_flag_actual & EvaluationFlags::values, + evaluation_flag_actual & EvaluationFlags::gradients, + evaluation_flag_actual & EvaluationFlags::hessians, this->face_no, this->subface_index, this->face_orientation, this->descriptor->face_orientations); } + if (hessians_on_general_cells) + { + unsigned int size = n_components * dim * n_q_points; + for (unsigned int i = 0; i < size; ++i) + this->gradients_from_hessians_quad[i] = this->gradients_quad[i]; + } + # ifdef DEBUG - if (evaluation_flag & EvaluationFlags::values) + if (evaluation_flag_actual & EvaluationFlags::values) this->values_quad_initialized = true; - if (evaluation_flag & EvaluationFlags::gradients) + if (evaluation_flag_actual & EvaluationFlags::gradients) this->gradients_quad_initialized = true; - if (evaluation_flag & EvaluationFlags::hessians) + if (evaluation_flag_actual & EvaluationFlags::hessians) this->hessians_quad_initialized = true; # endif } @@ -9300,6 +9413,24 @@ FEFaceEvaluationcell_type > internal::MatrixFreeFunctions::affine)) + { + unsigned int size = n_components * dim * n_q_points; + if (evaluation_flag & EvaluationFlags::gradients) + { + for (unsigned int i = 0; i < size; ++i) + this->gradients_quad[i] += this->gradients_from_hessians_quad[i]; + } + else + { + for (unsigned int i = 0; i < size; ++i) + this->gradients_quad[i] = this->gradients_from_hessians_quad[i]; + evaluation_flag_actual |= EvaluationFlags::gradients; + } + } + if (fe_degree > -1) internal::FEFaceEvaluationImplIntegrateSelector:: template run( @@ -9310,9 +9441,9 @@ FEFaceEvaluationbegin_gradients(), this->begin_hessians(), this->scratch_data, - evaluation_flag & EvaluationFlags::values, - evaluation_flag & EvaluationFlags::gradients, - evaluation_flag & EvaluationFlags::hessians, + evaluation_flag_actual & EvaluationFlags::values, + evaluation_flag_actual & EvaluationFlags::gradients, + evaluation_flag_actual & EvaluationFlags::hessians, this->face_no, this->subface_index, this->face_orientation, @@ -9326,9 +9457,9 @@ FEFaceEvaluationbegin_gradients(), this->begin_hessians(), this->scratch_data, - evaluation_flag & EvaluationFlags::values, - evaluation_flag & EvaluationFlags::gradients, - evaluation_flag & EvaluationFlags::hessians, + evaluation_flag_actual & EvaluationFlags::values, + evaluation_flag_actual & EvaluationFlags::gradients, + evaluation_flag_actual & EvaluationFlags::hessians, this->face_no, this->subface_index, this->face_orientation, diff --git a/include/deal.II/matrix_free/mapping_info.templates.h b/include/deal.II/matrix_free/mapping_info.templates.h index 2e3e1ecf17..2c3072008f 100644 --- a/include/deal.II/matrix_free/mapping_info.templates.h +++ b/include/deal.II/matrix_free/mapping_info.templates.h @@ -358,6 +358,10 @@ namespace internal update_quadrature_points ? update_quadrature_points : update_default) | + ((update_flags_inner_faces | update_flags_boundary_faces) & + (update_jacobian_grads | update_hessians) ? + update_jacobian_grads : + update_default) | update_normal_vectors | update_JxW_values | update_jacobians; this->update_flags_inner_faces = this->update_flags_boundary_faces; this->update_flags_faces_by_cells = update_flags_faces_by_cells; @@ -428,7 +432,7 @@ namespace internal const auto quad_face = quad[my_q][hpq].get_tensor_basis()[0]; face_data[my_q].descriptor[hpq * scale].initialize( - quad_face, update_flags_boundary_faces); + quad_face, this->update_flags_boundary_faces); face_data[my_q].q_collection[hpq] = dealii::hp::QCollection(quad_face); face_data_by_cells[my_q].descriptor[hpq * scale].initialize( @@ -496,8 +500,6 @@ namespace internal - /* ------------------------- initialization of cells ------------------- */ - // Copy a vectorized array of one type to another type template inline DEAL_II_ALWAYS_INLINE void @@ -516,6 +518,76 @@ namespace internal + // For second derivatives on the real cell, we need the gradient of the + // inverse Jacobian J. This involves some calculus and is done + // vectorized. If L is the gradient of the jacobian on the unit cell, + // the gradient of the inverse is given by (multidimensional calculus) - + // J * (J * L) * J (the third J is because we need to transform the + // gradient L from the unit to the real cell, and then apply the inverse + // Jacobian). Compare this with 1D with j(x) = 1/k(phi(x)), where j = + // phi' is the inverse of the jacobian and k is the derivative of the + // jacobian on the unit cell. Then j' = phi' k'/k^2 = j k' j^2. + template + Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>> + process_jacobian_gradient(const Tensor<2, dim, Number> &inv_jac_permut, + const Tensor<2, dim, Number> &inv_jac, + const Tensor<3, dim, Number> &jac_grad) + { + Number inv_jac_grad[dim][dim][dim]; + + // compute: inv_jac_grad = inv_jac_permut * grad_unit(jac) + for (unsigned int d = 0; d < dim; ++d) + for (unsigned int e = 0; e < dim; ++e) + for (unsigned int f = 0; f < dim; ++f) + { + inv_jac_grad[f][e][d] = + (inv_jac_permut[f][0] * jac_grad[d][e][0]); + for (unsigned int g = 1; g < dim; ++g) + inv_jac_grad[f][e][d] += + (inv_jac_permut[f][g] * jac_grad[d][e][g]); + } + + // compute: transpose (-inv_jac_permut * inv_jac_grad[d] * inv_jac) + Number tmp[dim]; + Number grad_jac_inv[dim][dim][dim]; + for (unsigned int d = 0; d < dim; ++d) + for (unsigned int e = 0; e < dim; ++e) + { + for (unsigned int f = 0; f < dim; ++f) + { + tmp[f] = -inv_jac_grad[d][f][0] * inv_jac[0][e]; + for (unsigned int g = 1; g < dim; ++g) + tmp[f] -= inv_jac_grad[d][f][g] * inv_jac[g][e]; + } + + // needed for non-diagonal part of Jacobian grad + for (unsigned int f = 0; f < dim; ++f) + { + grad_jac_inv[f][d][e] = inv_jac_permut[f][0] * tmp[0]; + for (unsigned int g = 1; g < dim; ++g) + grad_jac_inv[f][d][e] += inv_jac_permut[f][g] * tmp[g]; + } + } + + Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>> result; + + // the diagonal part of Jacobian gradient comes first + for (unsigned int d = 0; d < dim; ++d) + for (unsigned int f = 0; f < dim; ++f) + result[d][f] = grad_jac_inv[d][d][f]; + + // then the upper-diagonal part + for (unsigned int d = 0, count = 0; d < dim; ++d) + for (unsigned int e = d + 1; e < dim; ++e, ++count) + for (unsigned int f = 0; f < dim; ++f) + result[dim + count][f] = grad_jac_inv[d][e][f]; + return result; + } + + + + /* ------------------------- initialization of cells ------------------- */ + // Namespace with implementation of extraction of values on cell // range namespace ExtractCellHelper @@ -589,68 +661,7 @@ namespace internal } } - // For second derivatives on the real cell, we need the gradient of the - // inverse Jacobian J. This involves some calculus and is done - // vectorized. If L is the gradient of the jacobian on the unit cell, - // the gradient of the inverse is given by (multidimensional calculus) - - // J * (J * L) * J (the third J is because we need to transform the - // gradient L from the unit to the real cell, and then apply the inverse - // Jacobian). Compare this with 1D with j(x) = 1/k(phi(x)), where j = - // phi' is the inverse of the jacobian and k is the derivative of the - // jacobian on the unit cell. Then j' = phi' k'/k^2 = j k' j^2. - template - Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>> - process_jacobian_gradient(const Tensor<2, dim, Number> &inv_jac, - const Tensor<3, dim, Number> &jac_grad) - { - Number inv_jac_grad[dim][dim][dim]; - - // compute: inv_jac_grad = J*grad_unit(J^-1) - for (unsigned int d = 0; d < dim; ++d) - for (unsigned int e = 0; e < dim; ++e) - for (unsigned int f = 0; f < dim; ++f) - { - inv_jac_grad[f][e][d] = (inv_jac[f][0] * jac_grad[d][e][0]); - for (unsigned int g = 1; g < dim; ++g) - inv_jac_grad[f][e][d] += (inv_jac[f][g] * jac_grad[d][e][g]); - } - - // compute: transpose (-jac * jac_grad[d] * jac) - Number tmp[dim]; - Number grad_jac_inv[dim][dim][dim]; - for (unsigned int d = 0; d < dim; ++d) - for (unsigned int e = 0; e < dim; ++e) - { - for (unsigned int f = 0; f < dim; ++f) - { - tmp[f] = Number(); - for (unsigned int g = 0; g < dim; ++g) - tmp[f] -= inv_jac_grad[d][f][g] * inv_jac[g][e]; - } - // needed for non-diagonal part of Jacobian grad - for (unsigned int f = 0; f < dim; ++f) - { - grad_jac_inv[f][d][e] = inv_jac[f][0] * tmp[0]; - for (unsigned int g = 1; g < dim; ++g) - grad_jac_inv[f][d][e] += inv_jac[f][g] * tmp[g]; - } - } - - Tensor<1, dim *(dim + 1) / 2, Tensor<1, dim, Number>> result; - - // the diagonal part of Jacobian gradient comes first - for (unsigned int d = 0; d < dim; ++d) - for (unsigned int e = 0; e < dim; ++e) - result[d][e] = grad_jac_inv[d][d][e]; - - // then the upper-diagonal part - for (unsigned int d = 0, count = 0; d < dim; ++d) - for (unsigned int e = d + 1; e < dim; ++e, ++count) - for (unsigned int f = 0; f < dim; ++f) - result[dim + count][f] = grad_jac_inv[d][e][f]; - return result; - } /** * Helper function called internally during the initialize function. @@ -1049,7 +1060,9 @@ namespace internal if (update_flags & update_jacobian_grads) data.first[my_q].jacobian_gradients[0].push_back( - process_jacobian_gradient(inv_jac, jacobian_grad)); + process_jacobian_gradient(inv_jac, + inv_jac, + jacobian_grad)); } } @@ -1504,7 +1517,9 @@ namespace internal cell_grad_grads[q + (d * hess_dim + c) * n_q_points]; const auto inv_jac_grad = - process_jacobian_gradient(inv_jac, jac_grad); + process_jacobian_gradient(inv_jac, + inv_jac, + jac_grad); for (unsigned int d = 0; d < hess_dim; ++d) for (unsigned int e = 0; e < dim; ++e) store_vectorized_array( @@ -2355,8 +2370,56 @@ namespace internal if (process_face[face] == false) continue; - // go through the faces and fill the result const unsigned int offset = my_data.data_index_offsets[face]; + + const auto compute_jacobian_grad = + [&](unsigned int face_no, + int is_exterior, + unsigned int q, + Tensor<2, dim, VectorizedDouble> inv_jac) { + Tensor<2, dim, VectorizedDouble> inv_transp_jac_permut; + for (unsigned int d = 0; d < dim; ++d) + for (unsigned int e = 0; e < dim; ++e) + { + const unsigned int ee = + ExtractFaceHelper::reorder_face_derivative_indices< + dim>(face_no, e); + inv_transp_jac_permut[d][e] = inv_jac[ee][d]; + } + Tensor<2, dim, VectorizedDouble> jacobi; + for (unsigned int e = 0; e < dim; ++e) + for (unsigned int d = 0; d < dim; ++d) + jacobi[d][e] = face_grads[(d * dim + e) * n_q_points + q]; + Tensor<2, dim, VectorizedDouble> inv_transp_jac = + transpose(invert(jacobi)); + Tensor<3, dim, VectorizedDouble> jac_grad; + for (unsigned int d = 0; d < dim; ++d) + { + for (unsigned int e = 0; e < dim; ++e) + jac_grad[d][e][e] = + face_grad_grads[q + (d * hess_dim + e) * n_q_points]; + for (unsigned int c = dim, e = 0; e < dim; ++e) + for (unsigned int f = e + 1; f < dim; ++f, ++c) + jac_grad[d][e][f] = jac_grad[d][f][e] = + face_grad_grads[q + + (d * hess_dim + c) * n_q_points]; + const auto inv_jac_grad = + process_jacobian_gradient(inv_transp_jac_permut, + inv_transp_jac, + jac_grad); + for (unsigned int e = 0; e < dim; ++e) + { + for (unsigned int d = 0; d < hess_dim; ++d) + store_vectorized_array( + inv_jac_grad[d][e], + vv, + my_data.jacobian_gradients[is_exterior] + [offset + q][d][e]); + } + } + }; + + // go through the faces and fill the result const unsigned int n_points_compute = face_type[face] <= affine ? 1 : n_q_points; for (unsigned int q = 0; q < n_points_compute; ++q) @@ -2385,8 +2448,7 @@ namespace internal if (update_flags_faces & update_jacobian_grads) { - // TODO: implement jacobian grads for general path - AssertThrow(false, ExcNotImplemented()); + compute_jacobian_grad(face_no, 0, q, inv_jac); } std::array, dim - 1> @@ -2494,8 +2556,10 @@ namespace internal if (update_flags_faces & update_jacobian_grads) { - // TODO: implement jacobian grads for general path - AssertThrow(false, ExcNotImplemented()); + compute_jacobian_grad(faces[face].exterior_face_no, + 1, + q, + inv_jac); } my_data.normals_times_jacobians[1][offset + q] = diff --git a/tests/matrix_free/integrate_functions_hessians.cc b/tests/matrix_free/integrate_functions_hessians.cc new file mode 100644 index 0000000000..11363a5720 --- /dev/null +++ b/tests/matrix_free/integrate_functions_hessians.cc @@ -0,0 +1,263 @@ +// --------------------------------------------------------------------- +// +// 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. +// +// --------------------------------------------------------------------- + + + +// this function tests the correctness of the implementation of matrix free +// operations in integrating functions, gradients and hessians on a hyperball +// mesh with adaptive refinement. + +#include + +#include +#include + +#include +#include + +#include +#include +#include + +#include +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + + + +template +class MatrixFreeTest +{ +public: + using VectorType = std::vector *>; + + MatrixFreeTest(const MatrixFree &data_in) + : data(data_in) + , fe_val(data.get_dof_handler().get_fe(), + Quadrature(data.get_quadrature(0)), + update_values | update_gradients | update_hessians | + update_JxW_values){}; + + void + operator()(const MatrixFree & data, + VectorType & dst, + const VectorType & src, + const std::pair &cell_range) const; + + void + test_functions(Vector &dst, Vector &dst_deal) const + { + dst = 0; + dst_deal = 0; + VectorType dst_data(2); + dst_data[0] = &dst; + dst_data[1] = &dst_deal; + VectorType src_dummy; + data.cell_loop(&MatrixFreeTest::operator(), + this, + dst_data, + src_dummy); + }; + +private: + const MatrixFree &data; + mutable FEValues fe_val; +}; + + + +template +void +MatrixFreeTest::operator()( + const MatrixFree &data, + std::vector *> &dst, + const std::vector *> &, + const std::pair &cell_range) const +{ + FEEvaluation fe_eval(data); + const unsigned int n_q_points = fe_eval.n_q_points; + const unsigned int dofs_per_cell = fe_eval.dofs_per_cell; + AlignedVector> values(n_q_points); + AlignedVector> gradients(dim * n_q_points); + AlignedVector> hessians(dim * dim * n_q_points); + std::vector dof_indices(dofs_per_cell); + for (unsigned int cell = cell_range.first; cell < cell_range.second; ++cell) + { + fe_eval.reinit(cell); + // compare values with the ones the FEValues + // gives us. Those are seen as reference + for (unsigned int j = 0; j < data.n_components_filled(cell); ++j) + { + // generate random numbers at quadrature + // points and test them with basis functions + // and their gradients + for (unsigned int q = 0; q < n_q_points; ++q) + { + values[q][j] = random_value(); + for (unsigned int d = 0; d < dim; ++d) + gradients[q * dim + d][j] = -1. + 2. * (random_value()); + for (unsigned int d = 0; d < dim; ++d) + for (unsigned int e = 0; e < dim; ++e) + hessians[q * dim * dim + d * dim + e][j] = + -1. + 2. * (random_value()); + } + fe_val.reinit(data.get_cell_iterator(cell, j)); + data.get_cell_iterator(cell, j)->get_dof_indices(dof_indices); + + for (unsigned int i = 0; i < dofs_per_cell; ++i) + { + double sum = 0.; + for (unsigned int q = 0; q < n_q_points; ++q) + { + sum += + values[q][j] * fe_val.shape_value(i, q) * fe_val.JxW(q); + for (unsigned int d = 0; d < dim; ++d) + sum += (gradients[q * dim + d][j] * + fe_val.shape_grad(i, q)[d] * fe_val.JxW(q)); + for (unsigned int d = 0; d < dim; ++d) + for (unsigned int e = 0; e < dim; ++e) + sum += (hessians[q * dim * dim + d * dim + e][j] * + fe_val.shape_hessian(i, q)[d][e] * fe_val.JxW(q)); + } + (*dst[1])(dof_indices[i]) += sum; + } + } + for (unsigned int q = 0; q < n_q_points; ++q) + { + fe_eval.submit_value(values[q], q); + Tensor<1, dim, VectorizedArray> submit; + for (unsigned int d = 0; d < dim; ++d) + submit[d] = gradients[q * dim + d]; + fe_eval.submit_gradient(submit, q); + Tensor<2, dim, VectorizedArray> submit2; + for (unsigned int d = 0; d < dim; ++d) + for (unsigned int e = 0; e < dim; ++e) + submit2[d][e] = hessians[q * dim * dim + d * dim + e]; + fe_eval.submit_hessian(submit2, q); + } + fe_eval.integrate(EvaluationFlags::values | EvaluationFlags::gradients | + EvaluationFlags::hessians); + fe_eval.distribute_local_to_global(*dst[0]); + } +} + + + +template +void +test() +{ + using number = double; + const SphericalManifold manifold; + Triangulation tria; + GridGenerator::hyper_ball(tria); + typename Triangulation::active_cell_iterator cell = tria.begin_active(), + endc = tria.end(); + for (; cell != endc; ++cell) + for (const unsigned int f : GeometryInfo::face_indices()) + if (cell->at_boundary(f)) + cell->face(f)->set_all_manifold_ids(0); + tria.set_manifold(0, manifold); + + cell = tria.begin_active(); + for (; cell != endc; ++cell) + if (cell->center().norm() < 1e-8) + cell->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + cell = tria.begin_active(); + for (; cell != endc; ++cell) + if (cell->center().norm() < 0.2) + cell->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + if (dim < 3 || fe_degree < 2) + tria.refine_global(1); + tria.begin(tria.n_levels() - 1)->set_refine_flag(); + tria.last()->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + cell = tria.begin_active(); + for (unsigned int i = 0; i < 7 - 2 * dim; ++i) + { + cell = tria.begin_active(); + unsigned int counter = 0; + for (; cell != endc; ++cell, ++counter) + if (counter % (7 - i) == 0) + cell->set_refine_flag(); + tria.execute_coarsening_and_refinement(); + } + + FE_Q fe(fe_degree); + DoFHandler dof(tria); + dof.distribute_dofs(fe); + deallog << "Testing " << fe.get_name() << std::endl; + // std::cout << "Number of cells: " << tria.n_active_cells() << std::endl; + // std::cout << "Number of degrees of freedom: " << dof.n_dofs() << std::endl; + + AffineConstraints constraints; + DoFTools::make_hanging_node_constraints(dof, constraints); + constraints.close(); + + MatrixFree mf_data; + { + typename MatrixFree::AdditionalData data; + data.tasks_parallel_scheme = MatrixFree::AdditionalData::none; + data.mapping_update_flags |= update_hessians; + const QGauss<1> quad(fe_degree + 1); + mf_data.reinit(dof, constraints, quad, data); + } + + MatrixFreeTest mf(mf_data); + Vector solution(dof.n_dofs()); + Vector solution_dist(dof.n_dofs()); + + mf.test_functions(solution_dist, solution); + + constraints.condense(solution); + + Vector compare(solution_dist); + compare -= solution; + const double diff_norm = compare.linfty_norm(); + + deallog << "Norm of difference: " << diff_norm << std::endl << std::endl; +} + + +int +main() +{ + initlog(); + deallog << std::setprecision(3); + + { + deallog.push("2d"); + test<2, 1>(); + test<2, 2>(); + test<2, 3>(); + test<2, 4>(); + deallog.pop(); + deallog.push("3d"); + test<3, 1>(); + test<3, 2>(); + test<3, 3>(); + deallog.pop(); + } +} diff --git a/tests/matrix_free/integrate_functions_hessians.output b/tests/matrix_free/integrate_functions_hessians.output new file mode 100644 index 0000000000..3829c8ffdb --- /dev/null +++ b/tests/matrix_free/integrate_functions_hessians.output @@ -0,0 +1,22 @@ + +DEAL:2d::Testing FE_Q<2>(1) +DEAL:2d::Norm of difference: 0 +DEAL:2d:: +DEAL:2d::Testing FE_Q<2>(2) +DEAL:2d::Norm of difference: 0 +DEAL:2d:: +DEAL:2d::Testing FE_Q<2>(3) +DEAL:2d::Norm of difference: 0 +DEAL:2d:: +DEAL:2d::Testing FE_Q<2>(4) +DEAL:2d::Norm of difference: 0 +DEAL:2d:: +DEAL:3d::Testing FE_Q<3>(1) +DEAL:3d::Norm of difference: 0 +DEAL:3d:: +DEAL:3d::Testing FE_Q<3>(2) +DEAL:3d::Norm of difference: 0 +DEAL:3d:: +DEAL:3d::Testing FE_Q<3>(3) +DEAL:3d::Norm of difference: 0 +DEAL:3d:: diff --git a/tests/matrix_free/matrix_vector_hessians_general_cells.cc b/tests/matrix_free/matrix_vector_hessians_general_cells.cc new file mode 100644 index 0000000000..d9171bca55 --- /dev/null +++ b/tests/matrix_free/matrix_vector_hessians_general_cells.cc @@ -0,0 +1,310 @@ +// --------------------------------------------------------------------- +// +// 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 + +#include + +#include + +#include +#include +#include +#include + +#include + +#include +#include + +#include + +#include + +#include "../tests.h" + +// Tests if applying a matrix vector product with a matrix containing hessians +// on faces produces the same result with FEEvaluation and FEValues. +// This is checked for different combinations of EvaluationFlags, FE types and +// polynomial degrees. + +template +void +test_hessians(const dealii::FE_Poly & fe, + const dealii::Quadrature & quad, + const dealii::EvaluationFlags::EvaluationFlags evaluation_flags) +{ + using namespace dealii; + using VectorizedArrayType = VectorizedArray; + + Triangulation tria; + GridGenerator::hyper_ball(tria); + tria.refine_global(1); + + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(fe); + MappingQGeneric mapping(1); + + AffineConstraints constraints; + VectorTools::interpolate_boundary_values( + mapping, dof_handler, 0, Functions::ZeroFunction(), constraints); + constraints.close(); + + // FEEvaluation + typename MatrixFree::AdditionalData + additional_data; + additional_data.mapping_update_flags = + update_values | update_gradients | update_hessians; + + MatrixFree matrix_free; + matrix_free.reinit(mapping, dof_handler, constraints, quad, 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 src, dst, dst2; + matrix_free.initialize_dof_vector(src); + for (auto &v : src) + v = random_value(); + + matrix_free.initialize_dof_vector(dst); + matrix_free.initialize_dof_vector(dst2); + + // Setup FEFaceValues + FEValues fe_values(mapping, + fe, + quad, + update_values | update_gradients | update_hessians | + update_JxW_values); + Vector solution_values_local(fe.dofs_per_cell); + std::vector> solution_hessians(quad.size()); + std::vector> solution_gradients(quad.size()); + std::vector solution_values(quad.size()); + std::vector dof_indices(fe.dofs_per_cell); + src.update_ghost_values(); + dst2 = 0; + + + + matrix_free.template loop, + LinearAlgebra::distributed::Vector>( + [&]( + const auto &matrix_free, auto &dst, const auto &src, const auto &range) { + FEEvaluation fe_eval( + matrix_free); + for (unsigned int cell = range.first; cell < range.second; ++cell) + { + // FEEvaluation + fe_eval.reinit(cell); + fe_eval.gather_evaluate(src, evaluation_flags); + for (unsigned int q = 0; q < fe_eval.n_q_points; ++q) + { + if (evaluation_flags & EvaluationFlags::hessians) + fe_eval.submit_hessian(fe_eval.get_hessian(q), q); + if (evaluation_flags & EvaluationFlags::gradients) + fe_eval.submit_gradient(fe_eval.get_gradient(q), q); + if (evaluation_flags & EvaluationFlags::values) + fe_eval.submit_value(fe_eval.get_value(q), q); + } + fe_eval.integrate(evaluation_flags); + fe_eval.distribute_local_to_global(dst); + + // FEValues + for (unsigned int v = 0; + v < matrix_free.n_active_entries_per_cell_batch(cell); + ++v) + { + const auto cell_iterator = matrix_free.get_cell_iterator(cell, v); + + fe_values.reinit(cell_iterator); + + cell_iterator->get_dof_indices(dof_indices); + constraints.get_dof_values(src, + dof_indices.begin(), + solution_values_local.begin(), + solution_values_local.end()); + + std::vector> hessians_function(quad.size()); + fe_values.get_function_hessians(src, hessians_function); + + 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_local(i) * + fe_values.shape_hessian(i, q); + if (evaluation_flags & EvaluationFlags::gradients) + gradients += + solution_values_local(i) * fe_values.shape_grad(i, q); + if (evaluation_flags & EvaluationFlags::values) + values += solution_values_local(i) * + fe_values.shape_value(i, q); + } + solution_hessians[q] = hessians * fe_values.JxW(q); + solution_gradients[q] = gradients * fe_values.JxW(q); + solution_values[q] = values * fe_values.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_values.shape_hessian(i, q)); + if (evaluation_flags & EvaluationFlags::gradients) + sum_gradients += + solution_gradients[q] * fe_values.shape_grad(i, q); + if (evaluation_flags & EvaluationFlags::values) + sum_values += + solution_values[q] * fe_values.shape_value(i, q); + } + solution_values_local(i) = + sum_hessians + sum_gradients + sum_values; + } + constraints.distribute_local_to_global(solution_values_local, + dof_indices, + dst2); + } + } + }, + [&]( + 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) { + (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; + } +} + +void +test_qiterated(dealii::EvaluationFlags::EvaluationFlags evaluation_flags) +{ + using namespace dealii; + deallog << "test_qiterated" << std::endl; + for (unsigned int i = 1; i < 4; ++i) + { + QIterated<2> quad_2(QGauss<1>(i + 1), 2); + QIterated<3> quad_3(QGauss<1>(i + 1), 2); + test_hessians<2>(FE_Q<2>(i), quad_2, evaluation_flags); + test_hessians<3>(FE_Q<3>(i), quad_3, evaluation_flags); + test_hessians<2>(FE_DGQ<2>(i), quad_2, evaluation_flags); + test_hessians<3>(FE_DGQ<3>(i), quad_3, evaluation_flags); + } +} + +void +test_qgauss(dealii::EvaluationFlags::EvaluationFlags evaluation_flags) +{ + using namespace dealii; + deallog << "test_qgauss" << std::endl; + for (unsigned int i = 1; i < 4; ++i) + { + QGauss<2> quad_2(i + 1); + QGauss<3> quad_3(i + 1); + test_hessians<2>(FE_Q<2>(i), quad_2, evaluation_flags); + test_hessians<3>(FE_Q<3>(i), quad_3, evaluation_flags); + test_hessians<2>(FE_DGQ<2>(i), quad_2, evaluation_flags); + test_hessians<3>(FE_DGQ<3>(i), quad_3, evaluation_flags); + } +} + +int +main(int argc, char **argv) +{ + using namespace dealii; + Utilities::MPI::MPI_InitFinalize mpi(argc, argv, 1); + + initlog(); + deallog << std::setprecision(10); + + { + EvaluationFlags::EvaluationFlags evaluation_flags = + EvaluationFlags::hessians; + deallog << "test_hessians_only" << std::endl; + test_qgauss(evaluation_flags); + test_qiterated(evaluation_flags); + } + + { + EvaluationFlags::EvaluationFlags evaluation_flags = + EvaluationFlags::values | EvaluationFlags::hessians; + deallog << "test_hessians_with_values" << std::endl; + test_qgauss(evaluation_flags); + test_qiterated(evaluation_flags); + } + + { + EvaluationFlags::EvaluationFlags evaluation_flags = + EvaluationFlags::gradients | EvaluationFlags::hessians; + deallog << "test_hessians_with_gradients" << std::endl; + test_qgauss(evaluation_flags); + test_qiterated(evaluation_flags); + } + + { + EvaluationFlags::EvaluationFlags evaluation_flags = + EvaluationFlags::values | EvaluationFlags::gradients | + EvaluationFlags::hessians; + deallog << "test_hessians_with_gradients_and_values" << std::endl; + test_qgauss(evaluation_flags); + test_qiterated(evaluation_flags); + } +} diff --git a/tests/matrix_free/matrix_vector_hessians_general_cells.output b/tests/matrix_free/matrix_vector_hessians_general_cells.output new file mode 100644 index 0000000000..412bae749c --- /dev/null +++ b/tests/matrix_free/matrix_vector_hessians_general_cells.output @@ -0,0 +1,493 @@ + +DEAL::test_hessians_only +DEAL::test_qgauss +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 45.22219818 -7.037353151 0.000000000 11.88687292 -61.49576085 37.57202054 -35.68314307 +DEAL::dst FEV: 0.000000000 0.000000000 45.22219818 -7.037353151 0.000000000 11.88687292 -61.49576085 37.57202054 -35.68314307 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: 33.10696022 -22.71237546 -68.25761235 64.86438841 -52.49601275 16.86513406 30.83643936 -7.704940388 -34.30613074 +DEAL::dst FEV: 33.10696022 -22.71237546 -68.25761235 64.86438841 -52.49601275 16.86513406 30.83643936 -7.704940388 -34.30613074 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: -10.10812224 7.655480511 14.29504356 -11.84240183 -2.477974473 3.271861099 3.833223715 -4.627110340 8.081683648 +DEAL::dst FEV: -10.10812224 7.655480511 14.29504356 -11.84240183 -2.477974473 3.271861099 3.833223715 -4.627110340 8.081683648 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: 6.235165197 -7.858320343 0.1523344538 1.470820692 -4.214851267 5.838006413 -2.172648384 0.5494932382 3.923280641 +DEAL::dst FEV: 6.235165197 -7.858320343 0.1523344538 1.470820692 -4.214851267 5.838006413 -2.172648384 0.5494932382 3.923280641 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -207.3090085 -399.8368570 365.2209748 349.7535808 0.000000000 438.9124523 -467.2428521 +DEAL::dst FEV: 0.000000000 0.000000000 -207.3090085 -399.8368570 365.2209748 349.7535808 0.000000000 438.9124523 -467.2428521 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: -89.44696014 33.15396135 43.60705670 122.9701105 202.3600767 30.62778421 36.41906605 26.08882904 -14.27490495 +DEAL::dst FEV: -89.44696014 33.15396135 43.60705670 122.9701105 202.3600767 30.62778421 36.41906605 26.08882904 -14.27490495 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: 66.71712651 251.4306103 -12.27594710 -170.3467774 -483.4850601 35.48250021 117.8761478 214.8381899 -20.23679013 +DEAL::dst FEV: 66.71712651 251.4306103 -12.27594710 -170.3467774 -483.4850601 35.48250021 117.8761478 214.8381899 -20.23679013 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -9.588563680 19.11134992 -2.739954593 46.74620058 -92.01359251 28.89086228 -5.744434228 -16.14481508 19.53643851 +DEAL::dst FEV: -9.588563680 19.11134992 -2.739954593 46.74620058 -92.01359251 28.89086228 -5.744434228 -16.14481508 19.53643851 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 115.2792512 1514.207777 129.6997082 -385.2037949 562.9695734 -1541.447662 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 115.2792512 1514.207777 129.6997082 -385.2037949 562.9695734 -1541.447662 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: 182.6670278 -173.2496744 42.37257073 233.4142098 294.1859475 286.3887194 -175.4075400 -146.4736019 -864.9517464 +DEAL::dst FEV: 182.6670278 -173.2496744 42.37257073 233.4142098 294.1859475 286.3887194 -175.4075400 -146.4736019 -864.9517464 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: -19.97069668 -2263.033643 1394.341009 -180.1207746 205.5919186 4843.256322 -3031.471695 363.4727737 -584.4861111 +DEAL::dst FEV: -19.97069668 -2263.033643 1394.341009 -180.1207746 205.5919186 4843.256322 -3031.471695 363.4727737 -584.4861111 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 5.120375759 -102.7942743 35.13091124 26.35855211 -3.524129155 270.9555244 -269.4815327 0.9612496484 122.2229695 +DEAL::dst FEV: 5.120375759 -102.7942743 35.13091124 26.35855211 -3.524129155 270.9555244 -269.4815327 0.9612496484 122.2229695 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::test_qiterated +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -78.31750939 8.028509910 0.000000000 31.94914913 112.2418132 -8.404529286 -57.17545763 +DEAL::dst FEV: 0.000000000 0.000000000 -78.31750939 8.028509910 0.000000000 31.94914913 112.2418132 -8.404529286 -57.17545763 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: -28.62894665 12.83395912 -3.022244288 11.75298450 54.81383740 -40.08721361 -25.26188941 -8.528413431 28.07903454 +DEAL::dst FEV: -28.62894665 12.83395912 -3.022244288 11.75298450 54.81383740 -40.08721361 -25.26188941 -8.528413431 28.07903454 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: 8.563495216 -6.485642852 -12.11061108 10.03275871 -0.4091739206 0.5402639324 0.6329585687 -0.7640485805 -0.7596303223 +DEAL::dst FEV: 8.563495216 -6.485642852 -12.11061108 10.03275871 -0.4091739206 0.5402639324 0.6329585687 -0.7640485805 -0.7596303223 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: -0.04925434952 4.129829593 5.129818715 -9.210393958 -5.559736708 1.479161465 0.4791723432 3.601402900 -0.9781090982 +DEAL::dst FEV: -0.04925434952 4.129829593 5.129818715 -9.210393958 -5.559736708 1.479161465 0.4791723432 3.601402900 -0.9781090982 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -360.3224836 -234.0666285 579.7961020 132.2256446 0.000000000 313.0605472 -3.884425931 +DEAL::dst FEV: 0.000000000 0.000000000 -360.3224836 -234.0666285 579.7961020 132.2256446 0.000000000 313.0605472 -3.884425931 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: 139.4934997 117.6177543 38.83313679 49.94720211 -56.48270646 -7.179565374 38.65311026 -21.47267519 -149.4955229 +DEAL::dst FEV: 139.4934997 117.6177543 38.83313679 49.94720211 -56.48270646 -7.179565374 38.65311026 -21.47267519 -149.4955229 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: 26.94269068 24.79455169 -50.15544477 -9.626629392 -173.3905474 175.5729219 -18.80078412 158.8729847 -134.2097432 +DEAL::dst FEV: 26.94269068 24.79455169 -50.15544477 -9.626629392 -173.3905474 175.5729219 -18.80078412 158.8729847 -134.2097432 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -22.35910006 48.84059752 -5.609505384 36.26733584 -116.6366681 34.78761915 -11.86986656 34.27145167 -13.50944896 +DEAL::dst FEV: -22.35910006 48.84059752 -5.609505384 36.26733584 -116.6366681 34.78761915 -11.86986656 34.27145167 -13.50944896 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 2281.203835 53.14339166 799.5213508 -1660.550741 -61.70284234 218.9584292 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 2281.203835 53.14339166 799.5213508 -1660.550741 -61.70284234 218.9584292 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: 411.0414196 19.69629401 -45.09261305 -160.0937108 157.9595318 -63.70495743 -285.2088913 -88.92577089 854.5629205 +DEAL::dst FEV: 411.0414196 19.69629401 -45.09261305 -160.0937108 157.9595318 -63.70495743 -285.2088913 -88.92577089 854.5629205 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: -255.5427784 812.0180402 956.7833695 -196.0315102 455.8205352 -1387.185608 -2097.876988 411.3732530 -14.36300381 +DEAL::dst FEV: -255.5427784 812.0180402 956.7833695 -196.0315102 455.8205352 -1387.185608 -2097.876988 411.3732530 -14.36300381 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 17.32821876 -95.95926190 145.4413654 -15.35003562 -43.17407779 208.8850676 -131.1611595 -40.83909124 3.768804606 +DEAL::dst FEV: 17.32821876 -95.95926190 145.4413654 -15.35003562 -43.17407779 208.8850676 -131.1611595 -40.83909124 3.768804606 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::test_hessians_with_values +DEAL::test_qgauss +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -16.58940782 -4.782648286 0.000000000 -19.53954235 59.54922800 -66.05480814 64.73504286 +DEAL::dst FEV: 0.000000000 0.000000000 -16.58940782 -4.782648286 0.000000000 -19.53954235 59.54922800 -66.05480814 64.73504286 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: 32.15282885 21.33825418 -46.88409574 -5.642916990 4.268878289 21.02923008 -6.835454147 1.347064398 -35.37416089 +DEAL::dst FEV: 32.15282885 21.33825418 -46.88409574 -5.642916990 4.268878289 21.02923008 -6.835454147 1.347064398 -35.37416089 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: 10.79722351 -8.139267187 -15.21835329 12.64734630 -1.429251592 1.946800889 2.271259133 -2.693814615 18.10526726 +DEAL::dst FEV: 10.79722351 -8.139267187 -15.21835329 12.64734630 -1.429251592 1.946800889 2.271259133 -2.693814615 18.10526726 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: -3.738241635 8.395409955 0.6254841474 -5.280834462 -0.9752822717 -3.679925410 4.089943235 0.5673075883 0.3087058789 +DEAL::dst FEV: -3.738241635 8.395409955 0.6254841474 -5.280834462 -0.9752822717 -3.679925410 4.089943235 0.5673075883 0.3087058789 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 380.5559204 -105.0922063 -145.8738327 65.76715372 0.000000000 -153.4360487 450.3778700 +DEAL::dst FEV: 0.000000000 0.000000000 380.5559204 -105.0922063 -145.8738327 65.76715372 0.000000000 -153.4360487 450.3778700 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: 2.937489230 132.6971713 50.52850669 -15.20168086 -139.5233581 6.756647264 39.08008100 26.94906549 304.6782079 +DEAL::dst FEV: 2.937489230 132.6971713 50.52850669 -15.20168086 -139.5233581 6.756647264 39.08008100 26.94906549 304.6782079 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: -81.08771912 334.1346990 13.63891270 221.2008497 -725.0794917 -21.98122432 -146.3265523 390.8131029 14.77884105 +DEAL::dst FEV: -81.08771912 334.1346990 13.63891270 221.2008497 -725.0794917 -21.98122432 -146.3265523 390.8131029 14.77884105 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -4.494825168 -37.24295230 3.658353095 46.00761389 17.77100136 -4.702653614 -19.70483155 37.19542383 -11.33928896 +DEAL::dst FEV: -4.494825168 -37.24295230 3.658353095 46.00761389 17.77100136 -4.702653614 -19.70483155 37.19542383 -11.33928896 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 1015.004669 418.3204606 890.1073227 -1057.084816 752.3299712 -759.4644896 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 1015.004669 418.3204606 890.1073227 -1057.084816 752.3299712 -759.4644896 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: -161.7750270 -180.4296561 -16.05309565 -18.79545541 34.84814822 -168.2536254 -18.17675410 -29.88875304 618.8759516 +DEAL::dst FEV: -161.7750270 -180.4296561 -16.05309565 -18.79545541 34.84814822 -168.2536254 -18.17675410 -29.88875304 618.8759516 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: -129.9833408 232.5435605 -1068.363295 159.0051914 261.6113130 -257.6833985 2175.493847 -287.8736806 -200.4477312 +DEAL::dst FEV: -129.9833408 232.5435605 -1068.363295 159.0051914 261.6113130 -257.6833985 2175.493847 -287.8736806 -200.4477312 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 8.573820572 -70.82155289 -103.2850267 17.61450370 -51.13396496 124.0614029 312.1554282 -84.46207607 67.21754006 +DEAL::dst FEV: 8.573820572 -70.82155289 -103.2850267 17.61450370 -51.13396496 124.0614029 312.1554282 -84.46207607 67.21754006 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::test_qiterated +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 18.53584598 -36.66169783 0.000000000 84.22806902 -40.72609295 109.2107485 -146.1168399 +DEAL::dst FEV: 0.000000000 0.000000000 18.53584598 -36.66169783 0.000000000 84.22806902 -40.72609295 109.2107485 -146.1168399 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: 75.52280024 -48.98798905 -17.99142928 17.81212936 -39.28533289 -43.32180002 2.565804398 42.49911425 43.87285076 +DEAL::dst FEV: 75.52280024 -48.98798905 -17.99142928 17.81212936 -39.28533289 -43.32180002 2.565804398 42.49911425 43.87285076 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: 3.331492342 -2.476604953 -4.673138931 3.912736086 -1.153920507 1.551152426 1.809397922 -2.168905382 22.77362202 +DEAL::dst FEV: 3.331492342 -2.476604953 -4.673138931 3.912736086 -1.153920507 1.551152426 1.809397922 -2.168905382 22.77362202 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: -1.084466294 -1.172066883 6.039579603 -3.780049703 -2.346661305 4.606206981 -2.605353929 0.3489988567 -0.4208673164 +DEAL::dst FEV: -1.084466294 -1.172066883 6.039579603 -3.780049703 -2.346661305 4.606206981 -2.605353929 0.3489988567 -0.4208673164 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -309.5249641 -65.11577579 378.5689837 81.96287877 0.000000000 -18.49549227 -95.28776221 +DEAL::dst FEV: 0.000000000 0.000000000 -309.5249641 -65.11577579 378.5689837 81.96287877 0.000000000 -18.49549227 -95.28776221 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: 92.07855370 37.90950393 32.56690227 -101.3425164 -77.46971424 73.09994955 -38.69427253 1.202702857 -86.04224867 +DEAL::dst FEV: 92.07855370 37.90950393 32.56690227 -101.3425164 -77.46971424 73.09994955 -38.69427253 1.202702857 -86.04224867 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: -47.82882231 -102.6420815 76.20983134 145.0526142 166.6021274 -155.7349457 -123.9589010 -22.78484976 65.20942197 +DEAL::dst FEV: -47.82882231 -102.6420815 76.20983134 145.0526142 166.6021274 -155.7349457 -123.9589010 -22.78484976 65.20942197 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -15.97522003 7.821566253 17.51495995 45.64355144 -65.09372864 -6.523935492 3.012415787 -1.810943197 9.600318192 +DEAL::dst FEV: -15.97522003 7.821566253 17.51495995 45.64355144 -65.09372864 -6.523935492 3.012415787 -1.810943197 9.600318192 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -209.0416960 -346.9818449 564.0182222 -1099.889775 575.2820845 -212.5628742 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 -209.0416960 -346.9818449 564.0182222 -1099.889775 575.2820845 -212.5628742 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: 456.2392295 -118.0476821 229.4245288 -73.52012482 282.7625847 70.52133886 -142.1717202 163.4847811 271.8552957 +DEAL::dst FEV: 456.2392295 -118.0476821 229.4245288 -73.52012482 282.7625847 70.52133886 -142.1717202 163.4847811 271.8552957 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: 243.7688278 -2179.037732 70.00926913 93.64438209 -357.9827311 4397.708618 149.2573164 -206.5012884 -35.52626013 +DEAL::dst FEV: 243.7688278 -2179.037732 70.00926913 93.64438209 -357.9827311 4397.708618 149.2573164 -206.5012884 -35.52626013 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 42.60953463 1.905386931 4.105747246 -32.68370149 -66.55720248 -119.7933249 135.5493554 -15.93352778 22.40456232 +DEAL::dst FEV: 42.60953463 1.905386931 4.105747246 -32.68370149 -66.55720248 -119.7933249 135.5493554 -15.93352778 22.40456232 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::test_hessians_with_gradients +DEAL::test_qgauss +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 7.311269713 -5.156295237 0.000000000 -19.59465764 65.07757006 -96.82039038 98.14454960 +DEAL::dst FEV: 0.000000000 0.000000000 7.311269713 -5.156295237 0.000000000 -19.59465764 65.07757006 -96.82039038 98.14454960 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: -1.505294846 -24.57206121 52.11921332 39.11055377 -9.219780077 13.43074516 -61.01667614 9.190914949 -12.82165271 +DEAL::dst FEV: -1.505294846 -24.57206121 52.11921332 39.11055377 -9.219780077 13.43074516 -61.01667614 9.190914949 -12.82165271 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: -6.317046236 4.302693679 8.925812636 -6.911460079 -0.3150298140 0.1815725597 0.4627341489 -0.3292768946 -1.150968081 +DEAL::dst FEV: -6.317046236 4.302693679 8.925812636 -6.911460079 -0.3150298140 0.1815725597 0.4627341489 -0.3292768946 -1.150968081 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: 0.1232404125 0.01022319883 1.416411963 -1.568321407 2.617554109 -2.833395117 -4.137893678 4.372180519 -1.115248686 +DEAL::dst FEV: 0.1232404125 0.01022319883 1.416411963 -1.568321407 2.617554109 -2.833395117 -4.137893678 4.372180519 -1.115248686 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -205.1606811 -109.2816798 169.2814957 25.66427265 0.000000000 -115.4425259 402.3574467 +DEAL::dst FEV: 0.000000000 0.000000000 -205.1606811 -109.2816798 169.2814957 25.66427265 0.000000000 -115.4425259 402.3574467 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: -229.0756197 -169.2001752 0.4466035661 -74.58939161 -142.9520136 -29.63421128 68.81675020 57.23240132 -54.66379827 +DEAL::dst FEV: -229.0756197 -169.2001752 0.4466035661 -74.58939161 -142.9520136 -29.63421128 68.81675020 57.23240132 -54.66379827 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: 56.25534771 319.1648241 -91.58165669 -170.0871302 -603.7836113 196.7016090 151.6182738 225.7414787 -84.02913509 +DEAL::dst FEV: 56.25534771 319.1648241 -91.58165669 -170.0871302 -603.7836113 196.7016090 151.6182738 225.7414787 -84.02913509 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -15.84166619 0.2353443988 -2.115338897 23.85811429 49.84048831 -7.585819571 -14.10981493 -20.17993682 15.77644128 +DEAL::dst FEV: -15.84166619 0.2353443988 -2.115338897 23.85811429 49.84048831 -7.585819571 -14.10981493 -20.17993682 15.77644128 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 257.9275170 938.3914844 1444.117736 -1047.901673 1340.863608 -1542.857351 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 257.9275170 938.3914844 1444.117736 -1047.901673 1340.863608 -1542.857351 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: -272.6279247 29.31795285 -149.4439137 17.98294068 82.28977435 -333.7417511 170.4438834 50.90367540 571.9652347 +DEAL::dst FEV: -272.6279247 29.31795285 -149.4439137 17.98294068 82.28977435 -333.7417511 170.4438834 50.90367540 571.9652347 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: 167.8239843 937.4290293 373.7473509 -19.98309967 -453.9251491 -1886.074462 -704.8832500 -71.02701366 558.1818119 +DEAL::dst FEV: 167.8239843 937.4290293 373.7473509 -19.98309967 -453.9251491 -1886.074462 -704.8832500 -71.02701366 558.1818119 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: -37.20119468 142.7420109 -3.007994596 1.347129280 58.28259413 -229.2925423 20.73287757 31.08558645 -49.78341272 +DEAL::dst FEV: -37.20119468 142.7420109 -3.007994596 1.347129280 58.28259413 -229.2925423 20.73287757 31.08558645 -49.78341272 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::test_qiterated +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -39.34040403 53.15927296 0.000000000 -83.88672176 36.18288001 -85.82630612 133.8837425 +DEAL::dst FEV: 0.000000000 0.000000000 -39.34040403 53.15927296 0.000000000 -83.88672176 36.18288001 -85.82630612 133.8837425 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: -6.763525875 -29.84983569 8.787229829 29.57510339 35.20802355 -6.515494560 3.199548103 -17.09671493 61.28472734 +DEAL::dst FEV: -6.763525875 -29.84983569 8.787229829 29.57510339 35.20802355 -6.515494560 3.199548103 -17.09671493 61.28472734 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: -9.799798542 6.710971666 13.84374919 -10.75492232 -1.428274410 1.045904167 1.947944821 -1.565574578 5.993067307 +DEAL::dst FEV: -9.799798542 6.710971666 13.84374919 -10.75492232 -1.428274410 1.045904167 1.947944821 -1.565574578 5.993067307 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: -0.8390714818 3.507305214 2.477527164 -5.202526796 3.371341804 -5.961446157 -5.043747453 7.690617706 -1.205573925 +DEAL::dst FEV: -0.8390714818 3.507305214 2.477527164 -5.202526796 3.371341804 -5.961446157 -5.043747453 7.690617706 -1.205573925 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 76.24995757 -151.4652529 139.3454674 112.8200482 0.000000000 -86.44907039 268.3337342 +DEAL::dst FEV: 0.000000000 0.000000000 76.24995757 -151.4652529 139.3454674 112.8200482 0.000000000 -86.44907039 268.3337342 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: -109.3389581 -55.45400574 -15.87490632 -43.83124134 -209.9211798 40.69143970 -52.74249772 -18.77328792 -148.3336084 +DEAL::dst FEV: -109.3389581 -55.45400574 -15.87490632 -43.83124134 -209.9211798 40.69143970 -52.74249772 -18.77328792 -148.3336084 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: 59.12248949 145.3263730 20.39413228 -145.1793073 -240.0765308 -62.93699888 99.98252044 67.73515150 55.63217032 +DEAL::dst FEV: 59.12248949 145.3263730 20.39413228 -145.1793073 -240.0765308 -62.93699888 99.98252044 67.73515150 55.63217032 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: 0.1855879691 22.48428433 4.315303570 -54.04857952 108.0103722 -48.99698879 23.31684388 -30.68686048 11.60957481 +DEAL::dst FEV: 0.1855879691 22.48428433 4.315303570 -54.04857952 108.0103722 -48.99698879 23.31684388 -30.68686048 11.60957481 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 158.5095907 -1232.437308 572.2153320 -244.6069909 -498.6954353 994.6613382 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 158.5095907 -1232.437308 572.2153320 -244.6069909 -498.6954353 994.6613382 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: -43.97463300 -73.92611758 -241.0763625 16.24844181 -292.3228855 -131.0871170 71.84204246 -200.0407478 -502.1313285 +DEAL::dst FEV: -43.97463300 -73.92611758 -241.0763625 16.24844181 -292.3228855 -131.0871170 71.84204246 -200.0407478 -502.1313285 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: -79.52218037 -253.6335126 -828.4329925 -48.67475436 315.0517120 165.4890443 1618.475137 295.2068982 -521.4704040 +DEAL::dst FEV: -79.52218037 -253.6335126 -828.4329925 -48.67475436 315.0517120 165.4890443 1618.475137 295.2068982 -521.4704040 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 5.976291794 41.32606645 28.53127386 15.47828625 29.83484294 -321.4991292 -201.2675788 50.72178853 -45.30966870 +DEAL::dst FEV: 5.976291794 41.32606645 28.53127386 15.47828625 29.83484294 -321.4991292 -201.2675788 50.72178853 -45.30966870 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::test_hessians_with_gradients_and_values +DEAL::test_qgauss +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 10.65897096 -39.52678335 0.000000000 103.0744784 -17.64340333 94.38631066 -164.5032332 +DEAL::dst FEV: 0.000000000 0.000000000 10.65897096 -39.52678335 0.000000000 103.0744784 -17.64340333 94.38631066 -164.5032332 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: 1.102627435 -23.74365223 12.79134451 15.14153610 20.35616429 11.08290824 -9.664669182 14.19970889 -61.95383325 +DEAL::dst FEV: 1.102627435 -23.74365223 12.79134451 15.14153610 20.35616429 11.08290824 -9.664669182 14.19970889 -61.95383325 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: -5.658488958 4.454226542 8.341036788 -7.037923219 3.277233078 -3.216068598 -4.875453948 4.916510798 -38.29668774 +DEAL::dst FEV: -5.658488958 4.454226542 8.341036788 -7.037923219 3.277233078 -3.216068598 -4.875453948 4.916510798 -38.29668774 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: 3.537497829 -7.088058162 3.891919744 -0.3425928421 -1.788199416 5.329676319 -5.720480567 2.184775324 -1.566814782 +DEAL::dst FEV: 3.537497829 -7.088058162 3.891919744 -0.3425928421 -1.788199416 5.329676319 -5.720480567 2.184775324 -1.566814782 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 227.5797473 350.2646716 14.46622536 -299.5647869 0.000000000 -246.4604494 81.46044522 +DEAL::dst FEV: 0.000000000 0.000000000 227.5797473 350.2646716 14.46622536 -299.5647869 0.000000000 -246.4604494 81.46044522 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: -179.0434278 23.06902656 82.29275118 15.02013198 71.55602698 -124.1388389 48.89168750 20.79615558 269.8927219 +DEAL::dst FEV: -179.0434278 23.06902656 82.29275118 15.02013198 71.55602698 -124.1388389 48.89168750 20.79615558 269.8927219 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: 39.15418328 -29.11997591 67.98689134 -159.4974460 182.2488961 -183.0088824 126.0094564 -155.9890929 112.3426688 +DEAL::dst FEV: 39.15418328 -29.11997591 67.98689134 -159.4974460 182.2488961 -183.0088824 126.0094564 -155.9890929 112.3426688 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -1.274756990 -43.23619053 2.299787190 -40.17458947 170.9257619 8.426704364 11.03438425 -43.08065331 0.6067802486 +DEAL::dst FEV: -1.274756990 -43.23619053 2.299787190 -40.17458947 170.9257619 8.426704364 11.03438425 -43.08065331 0.6067802486 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 390.6713500 -209.1815116 453.4456649 -417.3813654 133.4270439 -135.8576152 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 390.6713500 -209.1815116 453.4456649 -417.3813654 133.4270439 -135.8576152 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: 154.9190172 -360.4883318 -356.1558831 157.2464180 117.8532397 252.9031712 -66.01918249 96.95024383 287.9441329 +DEAL::dst FEV: 154.9190172 -360.4883318 -356.1558831 157.2464180 117.8532397 252.9031712 -66.01918249 96.95024383 287.9441329 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: -89.45060520 259.6951529 848.1451912 -41.23094175 178.0046639 -291.4637914 -1901.434532 0.9182705078 -205.8252742 +DEAL::dst FEV: -89.45060520 259.6951529 848.1451912 -41.23094175 178.0046639 -291.4637914 -1901.434532 0.9182705078 -205.8252742 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 16.88952936 -26.19076464 166.9424606 -5.798811776 44.16364356 -140.5368456 -273.2726993 49.83102002 -17.72949070 +DEAL::dst FEV: 16.88952936 -26.19076464 166.9424606 -5.798811776 44.16364356 -140.5368456 -273.2726993 49.83102002 -17.72949070 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::test_qiterated +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 11.52339170 33.12164009 0.000000000 -31.44757979 31.42806080 -91.11960114 74.83230770 +DEAL::dst FEV: 0.000000000 0.000000000 11.52339170 33.12164009 0.000000000 -31.44757979 31.42806080 -91.11960114 74.83230770 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: -26.72770100 16.57344973 38.86689629 -30.27824765 33.97005374 -24.12586571 -42.87121363 43.97939426 -35.86348976 +DEAL::dst FEV: -26.72770100 16.57344973 38.86689629 -30.27824765 33.97005374 -24.12586571 -42.87121363 43.97939426 -35.86348976 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: 3.371302359 -2.385095954 -5.105225396 4.213858047 -1.101697459 1.730125280 1.611434656 -2.168999319 23.58480366 +DEAL::dst FEV: 3.371302359 -2.385095954 -5.105225396 4.213858047 -1.101697459 1.730125280 1.611434656 -2.168999319 23.58480366 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: 1.207718911 -3.743044786 4.895166874 -2.351876459 -3.005052078 5.530366961 -3.058534796 0.5306904929 0.1734889179 +DEAL::dst FEV: 1.207718911 -3.743044786 4.895166874 -2.351876459 -3.005052078 5.530366961 -3.058534796 0.5306904929 0.1734889179 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 304.2620687 -12.11556339 17.19536074 -102.9230813 0.000000000 -102.4623557 -65.22429151 +DEAL::dst FEV: 0.000000000 0.000000000 304.2620687 -12.11556339 17.19536074 -102.9230813 0.000000000 -102.4623557 -65.22429151 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: 176.7401543 -83.32894758 71.39944018 88.63852540 -190.9847564 -135.9080833 -40.88743431 53.54002808 -70.69405833 +DEAL::dst FEV: 176.7401543 -83.32894758 71.39944018 88.63852540 -190.9847564 -135.9080833 -40.88743431 53.54002808 -70.69405833 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: -88.91974786 39.17050152 34.41914498 239.8242491 -151.0085069 -51.10455761 -177.7689307 155.8642570 -0.3653264004 +DEAL::dst FEV: -88.91974786 39.17050152 34.41914498 239.8242491 -151.0085069 -51.10455761 -177.7689307 155.8642570 -0.3653264004 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -16.00930361 29.32089522 -4.691924057 13.64430384 -67.03254929 35.16094568 -2.368318354 0.8062879409 -2.621673503 +DEAL::dst FEV: -16.00930361 29.32089522 -4.691924057 13.64430384 -67.03254929 35.16094568 -2.368318354 0.8062879409 -2.621673503 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 602.3311513 349.8347114 865.7547210 -1258.503940 -599.7204024 190.5161394 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 602.3311513 349.8347114 865.7547210 -1258.503940 -599.7204024 190.5161394 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: 891.3558785 -64.82585347 144.7351978 -144.2971717 348.2630789 -131.1675690 -350.8844751 -25.03595897 145.5246669 +DEAL::dst FEV: 891.3558785 -64.82585347 144.7351978 -144.2971717 348.2630789 -131.1675690 -350.8844751 -25.03595897 145.5246669 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: 75.49607652 271.9395872 206.8530002 51.99120011 -332.7199127 -226.1663066 -472.3855211 -160.5810108 800.9167395 +DEAL::dst FEV: 75.49607652 271.9395872 206.8530002 51.99120011 -332.7199127 -226.1663066 -472.3855211 -160.5810108 800.9167395 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 10.13419666 26.50225800 -89.47868398 15.45030973 44.40065804 -53.19843171 279.5017995 -63.22517727 -2.347590852 +DEAL::dst FEV: 10.13419666 26.50225800 -89.47868398 15.45030973 44.40065804 -53.19843171 279.5017995 -63.22517727 -2.347590852 +DEAL::FEValues verification: 0.0 +DEAL:: diff --git a/tests/matrix_free/matrix_vector_hessians_general_faces.cc b/tests/matrix_free/matrix_vector_hessians_general_faces.cc new file mode 100644 index 0000000000..d327a0f594 --- /dev/null +++ b/tests/matrix_free/matrix_vector_hessians_general_faces.cc @@ -0,0 +1,383 @@ +// --------------------------------------------------------------------- +// +// 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 + +#include + +#include + +#include +#include +#include +#include + +#include +#include + +#include +#include + +#include + +#include +#include + +#include + +#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 +void +test_hessians(const unsigned int degree, + const dealii::FE_Poly & fe, + const dealii::EvaluationFlags::EvaluationFlags evaluation_flags) +{ + using namespace dealii; + using VectorizedArrayType = VectorizedArray; + + Triangulation tria; + GridGenerator::hyper_ball(tria, Point(), 1., true); + tria.refine_global(1); + + DoFHandler dof_handler(tria); + dof_handler.distribute_dofs(fe); + MappingQGeneric mapping(2); + + AffineConstraints constraints; + VectorTools::interpolate_boundary_values( + mapping, dof_handler, 0, Functions::ZeroFunction(), constraints); + constraints.close(); + + // FEFaceEvaluation + typename MatrixFree::AdditionalData + additional_data; + additional_data.mapping_update_flags_inner_faces = + update_values | update_gradients | update_hessians; + + MatrixFree 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 src, dst, dst2; + matrix_free.initialize_dof_vector(src); + for (auto &v : src) + v = random_value(); + + matrix_free.initialize_dof_vector(dst); + matrix_free.initialize_dof_vector(dst2); + + // Setup FEFaceValues + QGauss quad(degree + 1); + FEFaceValues fe_face_values_m(mapping, + fe, + quad, + update_values | update_gradients | + update_hessians | update_JxW_values); + FEFaceValues fe_face_values_p(mapping, + fe, + quad, + update_values | update_gradients | + update_hessians | update_JxW_values); + Vector solution_values_m(fe.dofs_per_cell); + Vector solution_values_p(fe.dofs_per_cell); + std::vector> solution_hessians(quad.size()); + std::vector> solution_gradients(quad.size()); + std::vector solution_values(quad.size()); + std::vector dof_indices_m(fe.dofs_per_cell); + std::vector dof_indices_p(fe.dofs_per_cell); + src.update_ghost_values(); + dst2 = 0; + + matrix_free.template loop, + LinearAlgebra::distributed::Vector>( + [&]( + 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 fe_eval_m( + matrix_free, true); + FEFaceEvaluation 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.read_dof_values(src); + fe_eval_m.evaluate(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) +{ + using namespace dealii; + 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<2>(i, dealii::FE_Q<2>(i), evaluation_flags); + test_hessians<3>(i, dealii::FE_Q<3>(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<2>(i, dealii::FE_Q<2>(i), evaluation_flags); + test_hessians<3>(i, dealii::FE_Q<3>(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<2>(i, dealii::FE_Q<2>(i), evaluation_flags); + test_hessians<3>(i, dealii::FE_Q<3>(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<2>(i, dealii::FE_Q<2>(i), evaluation_flags); + test_hessians<3>(i, dealii::FE_Q<3>(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); + } + } +} diff --git a/tests/matrix_free/matrix_vector_hessians_general_faces.output b/tests/matrix_free/matrix_vector_hessians_general_faces.output new file mode 100644 index 0000000000..90ff70b992 --- /dev/null +++ b/tests/matrix_free/matrix_vector_hessians_general_faces.output @@ -0,0 +1,245 @@ + +DEAL::test_hessians_only +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 1023.282938 -241.4191116 0.000000000 404.0338925 -1610.938705 830.1574688 -815.6478478 +DEAL::dst FEV: 0.000000000 0.000000000 1023.282938 -241.4191116 0.000000000 404.0338925 -1610.938705 830.1574688 -815.6478478 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: 2085.339062 -984.5969128 -3156.855020 2414.250954 -2457.242215 914.7438039 1549.567722 -598.1729479 -1625.597093 +DEAL::dst FEV: 2085.339062 -984.5969128 -3156.855020 2414.250954 -2457.242215 914.7438039 1549.567722 -598.1729479 -1625.597093 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: -70.25497193 56.55633625 99.11541602 -85.41678034 -18.84008428 22.15098884 28.29206961 -31.60297417 224.7350853 +DEAL::dst FEV: -70.25497193 56.55633625 99.11541602 -85.41678034 -18.84008428 22.15098884 28.29206961 -31.60297417 224.7350853 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: 230.8404882 -276.9256039 -49.48465092 95.56976663 -173.4791192 219.5642349 -7.876718073 -38.20839764 154.6373788 +DEAL::dst FEV: 230.8404882 -276.9256039 -49.48465092 95.56976663 -173.4791192 219.5642349 -7.876718073 -38.20839764 154.6373788 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 1995.727331 -5872.174659 4254.122091 4893.906481 0.000000000 5978.307265 -5884.314781 +DEAL::dst FEV: 0.000000000 0.000000000 1995.727331 -5872.174659 4254.122091 4893.906481 0.000000000 5978.307265 -5884.314781 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: -18232.67381 -63.84155745 8315.900756 10959.22675 18566.85721 1642.297186 5567.506956 3577.708281 3777.934353 +DEAL::dst FEV: -18232.67381 -63.84155745 8315.900756 10959.22675 18566.85721 1642.297186 5567.506956 3577.708281 3777.934353 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: 781.6486510 699.7352769 -40.96207937 -1563.480011 -1942.133396 521.4080939 724.0336595 1483.639308 -663.8895026 +DEAL::dst FEV: 781.6486510 699.7352769 -40.96207937 -1563.480011 -1942.133396 521.4080939 724.0336595 1483.639308 -663.8895026 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -1145.622366 1317.428972 -661.2759029 3000.099316 -4245.714206 1889.713020 -394.8445158 -1418.351983 1342.575424 +DEAL::dst FEV: -1145.622366 1317.428972 -661.2759029 3000.099316 -4245.714206 1889.713020 -394.8445158 -1418.351983 1342.575424 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -40598.10535 46880.86712 3405.824484 -6942.530759 19670.93356 -35093.42245 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 -40598.10535 46880.86712 3405.824484 -6942.530759 19670.93356 -35093.42245 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: 76024.91722 -21521.48567 7429.238626 34892.37890 61041.37781 39385.92760 -23643.19008 -28773.40415 -133970.1244 +DEAL::dst FEV: 76024.91722 -21521.48567 7429.238626 34892.37890 61041.37781 39385.92760 -23643.19008 -28773.40415 -133970.1244 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: 1588.303554 -14029.09500 10487.91417 -4673.201800 -2725.353375 37066.63526 -27593.29159 11153.52368 5951.323359 +DEAL::dst FEV: 1588.303554 -14029.09500 10487.91417 -4673.201800 -2725.353375 37066.63526 -27593.29159 11153.52368 5951.323359 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 374.2238240 -10778.38656 5315.827310 3061.915154 -3526.094877 27036.97679 -25348.35949 3956.883602 14680.05201 +DEAL::dst FEV: 374.2238240 -10778.38656 5315.827310 3061.915154 -3526.094877 27036.97679 -25348.35949 3956.883602 14680.05201 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::test_hessians_with_values +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -1584.165210 65.43910208 0.000000000 872.7146628 2609.017889 -88.44154865 -1485.472063 +DEAL::dst FEV: 0.000000000 0.000000000 -1584.165210 65.43910208 0.000000000 872.7146628 2609.017889 -88.44154865 -1485.472063 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: -1621.730760 643.4195686 293.0631459 279.0288608 2304.285927 -1409.468004 -1137.774174 -131.8455981 984.9798012 +DEAL::dst FEV: -1621.730760 643.4195686 293.0631459 279.0288608 2304.285927 -1409.468004 -1137.774174 -131.8455981 984.9798012 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: 57.37233989 -52.08768998 -81.52364724 76.90652431 -1.864777112 3.864526751 3.668178727 -5.283659601 -41.91915588 +DEAL::dst FEV: 57.37233989 -52.08768998 -81.52364724 76.90652431 -1.864777112 3.864526751 3.668178727 -5.283659601 -41.91915588 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: -35.93990419 151.8225756 180.2088970 -296.0414848 -123.2946863 7.466972238 -20.93161579 136.8238756 -43.26334548 +DEAL::dst FEV: -35.93990419 151.8225756 180.2088970 -296.0414848 -123.2946863 7.466972238 -20.93161579 136.8238756 -43.26334548 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 880.6007238 -2516.066873 7387.748162 1913.628818 0.000000000 2864.036480 -2737.954009 +DEAL::dst FEV: 0.000000000 0.000000000 880.6007238 -2516.066873 7387.748162 1913.628818 0.000000000 2864.036480 -2737.954009 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: 18596.12613 15528.57481 6884.355281 5777.026592 -2674.374463 -1007.329718 4432.749714 -2911.678883 -13831.23738 +DEAL::dst FEV: 18596.12613 15528.57481 6884.355281 5777.026592 -2674.374463 -1007.329718 4432.749714 -2911.678883 -13831.23738 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: 103.2620949 532.2380623 -633.7928408 419.9454288 -2675.339147 2218.294506 -622.4386240 2469.622503 -1811.077935 +DEAL::dst FEV: 103.2620949 532.2380623 -633.7928408 419.9454288 -2675.339147 2218.294506 -622.4386240 2469.622503 -1811.077935 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -2094.840483 3896.320685 -974.7402212 3024.687552 -6977.400539 2671.522120 -1333.566959 3011.785402 -1481.004655 +DEAL::dst FEV: -2094.840483 3896.320685 -974.7402212 3024.687552 -6977.400539 2671.522120 -1333.566959 3011.785402 -1481.004655 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 126464.1026 -1087.849656 22905.39195 -37835.72538 -736.0791996 9064.904645 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 126464.1026 -1087.849656 22905.39195 -37835.72538 -736.0791996 9064.904645 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: -11637.00801 -886.8694441 472.9210085 -28012.52627 36260.35838 -16255.65596 -35494.55425 -2416.417224 123004.4958 +DEAL::dst FEV: -11637.00801 -886.8694441 472.9210085 -28012.52627 36260.35838 -16255.65596 -35494.55425 -2416.417224 123004.4958 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: -3827.887005 -560.1149871 6336.511188 -2424.587270 7291.112492 3998.764339 -19339.34113 6836.645423 613.0280050 +DEAL::dst FEV: -3827.887005 -560.1149871 6336.511188 -2424.587270 7291.112492 3998.764339 -19339.34113 6836.645423 613.0280050 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 3267.074706 -11654.92306 15702.51697 -2120.107249 -3708.542275 13910.79441 -12561.65555 -4917.218117 -643.8540446 +DEAL::dst FEV: 3267.074706 -11654.92306 15702.51697 -2120.107249 -3708.542275 13910.79441 -12561.65555 -4917.218117 -643.8540446 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::test_hessians_with_gradients +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -621.4111919 114.1857906 0.000000000 -564.6157944 1395.687107 -1160.831541 1354.790422 +DEAL::dst FEV: 0.000000000 0.000000000 -621.4111919 114.1857906 0.000000000 -564.6157944 1395.687107 -1160.831541 1354.790422 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: 1278.032634 776.1433800 -1937.270825 -166.3396816 164.8182756 646.8512138 73.41373899 27.59426419 -1571.717176 +DEAL::dst FEV: 1278.032634 776.1433800 -1937.270825 -166.3396816 164.8182756 646.8512138 73.41373899 27.59426419 -1571.717176 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: 75.36615516 -64.43827144 -107.2029444 96.27506071 -10.15977218 14.65805710 15.33219542 -19.83048034 485.1973589 +DEAL::dst FEV: 75.36615516 -64.43827144 -107.2029444 96.27506071 -10.15977218 14.65805710 15.33219542 -19.83048034 485.1973589 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: -147.6510186 280.3568712 58.16675185 -191.7314131 12.97777322 -145.4542511 76.29975751 57.03552894 21.57976763 +DEAL::dst FEV: -147.6510186 280.3568712 58.16675185 -191.7314131 12.97777322 -145.4542511 76.29975751 57.03552894 21.57976763 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 22159.26719 2260.566599 -2765.305051 714.0319473 0.000000000 -6768.530153 2345.620782 +DEAL::dst FEV: 0.000000000 0.000000000 22159.26719 2260.566599 -2765.305051 714.0319473 0.000000000 -6768.530153 2345.620782 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: -3899.489149 15520.20968 5914.037728 -1363.044631 -15353.90033 2244.139771 4155.919950 3503.216608 22003.29974 +DEAL::dst FEV: -3899.489149 15520.20968 5914.037728 -1363.044631 -15353.90033 2244.139771 4155.919950 3503.216608 22003.29974 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: -1292.549827 1556.133239 6.674065417 3834.116898 -4965.432205 581.6592615 -2903.715710 4050.659028 -867.5447504 +DEAL::dst FEV: -1292.549827 1556.133239 6.674065417 3834.116898 -4965.432205 581.6592615 -2903.715710 4050.659028 -867.5447504 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -169.5401972 -2238.676603 552.5699787 2580.643343 483.3065112 -585.8533834 -1693.719654 1992.806308 -614.3341662 +DEAL::dst FEV: -169.5401972 -2238.676603 552.5699787 2580.643343 483.3065112 -585.8533834 -1693.719654 1992.806308 -614.3341662 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 56775.09128 12613.03416 12887.52510 -18193.82181 17692.30415 -22298.47264 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 56775.09128 12613.03416 12887.52510 -18193.82181 17692.30415 -22298.47264 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: -42336.21123 -26120.26166 -2297.564253 795.3646393 195.9688350 -26011.43162 -6710.972876 -1437.849217 111350.4393 +DEAL::dst FEV: -42336.21123 -26120.26166 -2297.564253 795.3646393 195.9688350 -26011.43162 -6710.972876 -1437.849217 111350.4393 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: -5791.865257 1722.236271 -6101.199466 3515.520543 12454.33999 -1254.417975 14202.82726 -7529.176527 -12023.46778 +DEAL::dst FEV: -5791.865257 1722.236271 -6101.199466 3515.520543 12454.33999 -1254.417975 14202.82726 -7529.176527 -12023.46778 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 873.2174032 -6121.609196 -13551.64405 4309.630048 -4242.986432 9017.066410 28639.88440 -10457.78950 5235.775871 +DEAL::dst FEV: 873.2174032 -6121.609196 -13551.64405 4309.630048 -4242.986432 9017.066410 28639.88440 -10457.78950 5235.775871 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::test_hessians_with_gradients_and_values +DEAL::Working with FE_Q<2>(1) and 25 dofs +DEAL::dst FEE: 0.000000000 0.000000000 520.1325322 -690.2874624 0.000000000 1816.271572 -973.5840476 2143.792127 -3254.583688 +DEAL::dst FEV: 0.000000000 0.000000000 520.1325322 -690.2874624 0.000000000 1816.271572 -973.5840476 2143.792127 -3254.583688 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(1) and 79 dofs +DEAL::dst FEE: 3244.549979 -1911.946038 -996.0944826 666.3414459 -1630.451414 -1193.463312 174.5157649 1414.847115 1650.753381 +DEAL::dst FEV: 3244.549979 -1911.946038 -996.0944826 666.3414459 -1630.451414 -1193.463312 174.5157649 1414.847115 1650.753381 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(1) and 80 dofs +DEAL::dst FEE: 22.28659957 -15.47011943 -34.59406422 28.35461984 -4.581461453 13.39644349 8.875622598 -17.53301884 575.4344949 +DEAL::dst FEV: 22.28659957 -15.47011943 -34.59406422 28.35461984 -4.581461453 13.39644349 8.875622598 -17.53301884 575.4344949 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(1) and 448 dofs +DEAL::dst FEE: -44.84164597 -20.30585088 186.0485414 -121.5573779 -53.47487885 118.0871478 -87.61030991 23.83004676 -39.39444906 +DEAL::dst FEV: -44.84164597 -20.30585088 186.0485414 -121.5573779 -53.47487885 118.0871478 -87.61030991 23.83004676 -39.39444906 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(2) and 89 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -7209.671714 -3603.056480 4376.848314 994.6112263 0.000000000 2021.933328 -981.7698190 +DEAL::dst FEV: 0.000000000 0.000000000 -7209.671714 -3603.056480 4376.848314 994.6112263 0.000000000 2021.933328 -981.7698190 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(2) and 517 dofs +DEAL::dst FEE: 1367.885593 4424.107523 5513.989531 -7434.495343 -3908.076539 6480.955358 -330.9300163 -1136.937005 -6414.069171 +DEAL::dst FEV: 1367.885593 4424.107523 5513.989531 -7434.495343 -3908.076539 6480.955358 -330.9300163 -1136.937005 -6414.069171 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(2) and 180 dofs +DEAL::dst FEE: -437.8541486 -496.2779202 809.0422323 1332.294523 719.5400451 -1784.743121 -1085.021345 29.27757678 914.5814941 +DEAL::dst FEV: -437.8541486 -496.2779202 809.0422323 1332.294523 719.5400451 -1784.743121 -1085.021345 29.27757678 914.5814941 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(2) and 1512 dofs +DEAL::dst FEE: -1566.061648 375.3459780 1252.028307 3338.351514 -3056.749505 -440.7156811 -382.1943741 -181.2283054 682.1878856 +DEAL::dst FEV: -1566.061648 375.3459780 1252.028307 3338.351514 -3056.749505 -440.7156811 -382.1943741 -181.2283054 682.1878856 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<2>(3) and 193 dofs +DEAL::dst FEE: 0.000000000 0.000000000 -72513.09849 2440.045444 5936.767657 366.5967395 13588.30105 -15481.35114 0.000000000 +DEAL::dst FEV: 0.000000000 0.000000000 -72513.09849 2440.045444 5936.767657 366.5967395 13588.30105 -15481.35114 0.000000000 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_Q<3>(3) and 1651 dofs +DEAL::dst FEE: 33810.56033 -17394.35348 34025.74741 -12344.52484 48066.48418 4634.340963 -24021.08665 26673.43063 30975.17620 +DEAL::dst FEV: 33810.56033 -17394.35348 34025.74741 -12344.52484 48066.48418 4634.340963 -24021.08665 26673.43063 30975.17620 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<2>(3) and 320 dofs +DEAL::dst FEE: 5527.094242 -10178.50933 242.7247396 -696.3185404 -11066.69999 24621.98310 2200.859720 1533.485743 11062.88761 +DEAL::dst FEV: 5527.094242 -10178.50933 242.7247396 -696.3185404 -11066.69999 24621.98310 2200.859720 1533.485743 11062.88761 +DEAL::FEValues verification: 0.0 +DEAL:: +DEAL::Working with FE_DGQ<3>(3) and 3584 dofs +DEAL::dst FEE: 7134.467624 -3170.261114 -745.5916238 -4621.556435 -6399.737030 -7674.914630 14942.24049 -3133.337521 5780.416539 +DEAL::dst FEV: 7134.467624 -3170.261114 -745.5916238 -4621.556435 -6399.737030 -7674.914630 14942.24049 -3133.337521 5780.416539 +DEAL::FEValues verification: 0.0 +DEAL::