From 7223a399063f2e212a1ee0328e66eb583cb7c98c Mon Sep 17 00:00:00 2001 From: Jiaqi Zhang Date: Mon, 11 Apr 2022 11:27:01 -0400 Subject: [PATCH] feiv::get_*_function_values --- include/deal.II/fe/fe_interface_values.h | 298 +++++++++++++++++++++++ 1 file changed, 298 insertions(+) diff --git a/include/deal.II/fe/fe_interface_values.h b/include/deal.II/fe/fe_interface_values.h index 3cd4763398..543a8c8a44 100644 --- a/include/deal.II/fe/fe_interface_values.h +++ b/include/deal.II/fe/fe_interface_values.h @@ -1954,6 +1954,126 @@ public: * @} */ + + + /** + * @name Access to jumps in the function values and derivatives + * @{ + */ + + /** + * Return the jump in the values of the + * finite element function characterized by fe_function at the + * quadrature points of the cell interface selected the last time + * the reinit function of the FEInterfaceValues object was called. + * + * @dealiiRequiresUpdateFlags{update_values} + */ + template + void + get_jump_in_function_values( + const InputVector & fe_function, + std::vector &values) const; + + /** + * Return the jump in the gradients of the + * finite element function characterized by fe_function at the + * quadrature points of the cell interface selected the last time + * the reinit function of the FEInterfaceValues object was called. + * + * @dealiiRequiresUpdateFlags{update_gradients} + */ + template + void + get_jump_in_function_gradients( + const InputVector &fe_function, + std::vector> + &gradients) const; + + /** + * Return the jump in the Hessians of the + * finite element function characterized by fe_function at the + * quadrature points of the cell interface selected the last time + * the reinit function of the FEInterfaceValues object was called. + * @dealiiRequiresUpdateFlags{update_hessians} + */ + template + void + get_jump_in_function_hessians( + const InputVector &fe_function, + std::vector> + &hessians) const; + + /** + * Return the jump in the third derivatives of the + * the finite element function characterized by fe_function at + * the quadrature points of the cell interface selected the last time + * the reinit function of the FEInterfaceValues object was called. + * + * @dealiiRequiresUpdateFlags{update_third_derivatives} + */ + template + void + get_jump_in_function_third_derivatives( + const InputVector &fe_function, + std::vector> + &third_derivatives) const; + + //@} + + /** + * @name Access to the average of the function values and derivatives + */ + //@{ + + /** + * Return the average of the values of the + * finite element function characterized by fe_function at the + * quadrature points of the cell interface selected the last time + * the reinit function of the FEInterfaceValues object was called. + * + * @dealiiRequiresUpdateFlags{update_values} + */ + template + void + get_average_of_function_values( + const InputVector & fe_function, + std::vector &values) const; + + /** + * Return the average of the gradients of the + * the finite element function characterized by fe_function at the + * quadrature points of the cell interface selected the last time + * the reinit function of the FEInterfaceValues object was called. + * @dealiiRequiresUpdateFlags{update_gradients} + */ + template + void + get_average_of_function_gradients( + const InputVector &fe_function, + std::vector> + &gradients) const; + + /** + * Return the average of the Hessians of the + * the finite element function characterized by fe_function at the + * quadrature points of the cell interface selected the last time + * the reinit function of the FEInterfaceValues object was called. + * @dealiiRequiresUpdateFlags{update_hessians} + */ + template + void + get_average_of_function_hessians( + const InputVector &fe_function, + std::vector> + &hessians) const; + + /** + * @} + */ + + + /** * @name Extractors Methods to extract individual components * @{ @@ -2715,6 +2835,184 @@ FEInterfaceValues::jump_3rd_derivative( +template +template +void +FEInterfaceValues::get_jump_in_function_values( + const InputVector & fe_function, + std::vector &values) const +{ + AssertDimension(values.size(), n_quadrature_points); + + std::vector values_0(n_quadrature_points); + std::vector values_1(n_quadrature_points); + + get_fe_face_values(0).get_function_values(fe_function, values_0); + get_fe_face_values(1).get_function_values(fe_function, values_1); + + for (unsigned int q = 0; q < n_quadrature_points; ++q) + { + values[q] = values_0[q] - values_1[q]; + } +} + + + +template +template +void +FEInterfaceValues::get_jump_in_function_gradients( + const InputVector &fe_function, + std::vector> &gradients) + const +{ + AssertDimension(gradients.size(), n_quadrature_points); + + std::vector> + gradients_0(n_quadrature_points); + std::vector> + gradients_1(n_quadrature_points); + + get_fe_face_values(0).get_function_gradients(fe_function, gradients_0); + get_fe_face_values(1).get_function_gradients(fe_function, gradients_1); + + for (unsigned int q = 0; q < n_quadrature_points; ++q) + { + gradients[q] = gradients_0[q] - gradients_1[q]; + } +} + + + +template +template +void +FEInterfaceValues::get_jump_in_function_hessians( + const InputVector &fe_function, + std::vector> &hessians) + const +{ + AssertDimension(hessians.size(), n_quadrature_points); + + std::vector> hessians_0( + n_quadrature_points); + std::vector> hessians_1( + n_quadrature_points); + + get_fe_face_values(0).get_function_hessians(fe_function, hessians_0); + get_fe_face_values(1).get_function_hessians(fe_function, hessians_1); + + for (unsigned int q = 0; q < n_quadrature_points; ++q) + { + hessians[q] = (hessians_0[q] - hessians_1[q]); + } +} + + + +template +template +void +FEInterfaceValues::get_jump_in_function_third_derivatives( + const InputVector &fe_function, + std::vector> + &third_derivatives) const +{ + AssertDimension(third_derivatives.size(), n_quadrature_points); + + std::vector> + third_derivatives_0(n_quadrature_points); + std::vector> + third_derivatives_1(n_quadrature_points); + + get_fe_face_values(0).get_function_third_derivatives(fe_function, + third_derivatives_0); + get_fe_face_values(1).get_function_third_derivatives(fe_function, + third_derivatives_1); + + for (unsigned int q = 0; q < n_quadrature_points; ++q) + { + third_derivatives[q] = (third_derivatives_0[q] - third_derivatives_1[q]); + } +} + + + +template +template +void +FEInterfaceValues::get_average_of_function_values( + const InputVector & fe_function, + std::vector &values) const +{ + AssertDimension(values.size(), n_quadrature_points); + + std::vector values_0(n_quadrature_points); + std::vector values_1(n_quadrature_points); + + get_fe_face_values(0).get_function_values(fe_function, values_0); + get_fe_face_values(1).get_function_values(fe_function, values_1); + + for (unsigned int q = 0; q < n_quadrature_points; ++q) + { + values[q] = (values_0[q] + values_1[q]) * 0.5; + } +} + + + +template +template +void +FEInterfaceValues::get_average_of_function_gradients( + const InputVector &fe_function, + std::vector> &gradients) + const +{ + AssertDimension(gradients.size(), n_quadrature_points); + + std::vector> + gradients_0(n_quadrature_points); + std::vector> + gradients_1(n_quadrature_points); + + get_fe_face_values(0).get_function_gradients(fe_function, gradients_0); + get_fe_face_values(1).get_function_gradients(fe_function, gradients_1); + + for (unsigned int q = 0; q < n_quadrature_points; ++q) + { + gradients[q] = (gradients_0[q] + gradients_1[q]) * 0.5; + } +} + + + +template +template +void +FEInterfaceValues::get_average_of_function_hessians( + const InputVector &fe_function, + std::vector> &hessians) + const +{ + AssertDimension(hessians.size(), n_quadrature_points); + + std::vector> hessians_0( + n_quadrature_points); + std::vector> hessians_1( + n_quadrature_points); + + get_fe_face_values(0).get_function_hessians(fe_function, hessians_0); + get_fe_face_values(1).get_function_hessians(fe_function, hessians_1); + + for (unsigned int q = 0; q < n_quadrature_points; ++q) + { + hessians[q] = (hessians_0[q] + hessians_1[q]) * 0.5; + } +} + + + /*------------ Inline functions: FEInterfaceValues------------*/ template inline const FEInterfaceViews::Scalar -- 2.39.5