From 6787636517dd0191e2607a07d1eaf7b7baf57ebe Mon Sep 17 00:00:00 2001 From: Alex Lindsay Date: Tue, 23 May 2023 22:14:55 -0700 Subject: [PATCH] Fix adCurvatures method --- framework/include/base/Assembly.h | 6 +----- framework/src/base/Assembly.C | 13 +++++++++++++ .../test/include/bcs/CrazyKCPlantFitsBoundary.h | 2 ++ .../test/src/bcs/CrazyKCPlantFitsBoundary.C | 4 ++++ 4 files changed, 20 insertions(+), 5 deletions(-) diff --git a/framework/include/base/Assembly.h b/framework/include/base/Assembly.h index ab236d0d22a7..a72022487f04 100644 --- a/framework/include/base/Assembly.h +++ b/framework/include/base/Assembly.h @@ -243,11 +243,7 @@ class Assembly const MooseArray & adJxWFace() const { return _ad_JxW_face; } - const MooseArray & adCurvatures() const - { - _calculate_curvatures = true; - return _ad_curvatures; - } + const MooseArray & adCurvatures() const; /** * Returns the reference to the coordinate transformation coefficients diff --git a/framework/src/base/Assembly.C b/framework/src/base/Assembly.C index c5a8e00ba190..e539034b32a9 100644 --- a/framework/src/base/Assembly.C +++ b/framework/src/base/Assembly.C @@ -4736,6 +4736,19 @@ Assembly::processResidualsAndJacobian(const std::vector & residuals, cacheJacobian(row_indices[i], column_indices[j], element_matrix(i, j), matrix_tags); } +const MooseArray & +Assembly::adCurvatures() const +{ + _calculate_curvatures = true; + const Order helper_order = _mesh.hasSecondOrderElements() ? SECOND : FIRST; + const FEType helper_type(helper_order, LAGRANGE); + // Must prerequest the second derivatives. Sadly because there is only one _need_second_derivative + // map for both volumetric and face FE objects we must request both here + feSecondPhi(helper_type); + feSecondPhiFace(helper_type); + return _ad_curvatures; +} + template void coordTransformFactor(const SubProblem & s, SubdomainID sub_id, const Point & point, diff --git a/modules/navier_stokes/test/include/bcs/CrazyKCPlantFitsBoundary.h b/modules/navier_stokes/test/include/bcs/CrazyKCPlantFitsBoundary.h index 538b3d250751..3454776f1d6e 100644 --- a/modules/navier_stokes/test/include/bcs/CrazyKCPlantFitsBoundary.h +++ b/modules/navier_stokes/test/include/bcs/CrazyKCPlantFitsBoundary.h @@ -47,6 +47,8 @@ class CrazyKCPlantFitsBoundary : public ADMaterial ADMaterialProperty & _surface_tension; ADMaterialProperty & _grad_surface_tension; const MooseArray & _ad_normals; + const MooseArray & _ad_curvatures; + ADMaterialProperty & _surface_term_curvature; ADMaterialProperty & _surface_term_gradient1; ADMaterialProperty & _surface_term_gradient2; diff --git a/modules/navier_stokes/test/src/bcs/CrazyKCPlantFitsBoundary.C b/modules/navier_stokes/test/src/bcs/CrazyKCPlantFitsBoundary.C index 547a1676ae81..70b332efedf1 100644 --- a/modules/navier_stokes/test/src/bcs/CrazyKCPlantFitsBoundary.C +++ b/modules/navier_stokes/test/src/bcs/CrazyKCPlantFitsBoundary.C @@ -82,6 +82,8 @@ CrazyKCPlantFitsBoundary::CrazyKCPlantFitsBoundary(const InputParameters & param _grad_surface_tension(declareADProperty( getParam("grad_surface_tension_name"))), _ad_normals(_assembly.adNormals()), + _ad_curvatures(_assembly.adCurvatures()), + _surface_term_curvature(declareADProperty("surface_term_curvature")), _surface_term_gradient1(declareADProperty("surface_term_gradient1")), _surface_term_gradient2(declareADProperty("surface_term_gradient2")), _length_units_per_meter(1. / std::pow(10, getParam("length_unit_exponent"))), @@ -115,6 +117,8 @@ CrazyKCPlantFitsBoundary::computeQpProperties() _grad_surface_tension[_qp] = _alpha * _mass_units_per_kilogram / (_time_units_per_second * _time_units_per_second) / _temperature_units_per_kelvin * _grad_temperature[_qp]; + _surface_term_curvature[_qp] = + -2. * _ad_curvatures[_qp] * _surface_tension[_qp] * _ad_normals[_qp]; _surface_term_gradient1[_qp] = -_grad_surface_tension[_qp]; _surface_term_gradient2[_qp] = _ad_normals[_qp] * (_ad_normals[_qp] * _grad_surface_tension[_qp]); }