LandIce: fixes for issue #1092

sandialabs · Nov 21, 2024 · 2f779ae · 2f779ae
1 parent ef72472
commit 2f779ae
Show file tree

Hide file tree

Showing 3 changed files with 10 additions and 8 deletions.
diff --git a/src/landIce/evaluators/LandIce_BasalFrictionCoefficient_Def.hpp b/src/landIce/evaluators/LandIce_BasalFrictionCoefficient_Def.hpp
@@ -186,6 +186,7 @@ BasalFrictionCoefficient (const Teuchos::ParameterList& p,
     } else if (effectivePressureType == "HYDROSTATIC") {
       effectivePressure_type = EFFECTIVE_PRESSURE_TYPE::HYDROSTATIC;
       use_pressurized_bed = beta_list.get<bool>("Use Pressurized Bed Above Sea Level", false);
+      save_pressure_field = p.isParameter("Effective Pressure Output Variable Name");
       if(save_pressure_field && nodal) {
         outN = PHX::MDField<EffPressureST>(p.get<std::string> ("Effective Pressure Output Variable Name"), nodal_layout);
         this->addEvaluatedField (outN);
@@ -397,6 +398,9 @@ operator() (const BasalFrictionCoefficient_Tag&, const int& cell) const {
   ParamScalarT muValue = 1.0;
   typename Albany::StrongestScalarType<EffPressureST,MeshScalarT>::type NVal = N_val;
 
+  TEUCHOS_TEST_FOR_EXCEPTION ((save_pressure_field && !std::is_constructible<EffPressureST,MeshScalarT>::value), std::logic_error,
+                                "Error! BasalFrictionCoefficient: Trying to convert a FAD type (NVal) into a double (outN).\n");
+
   if(beta_type != BETA_TYPE::CONSTANT) {
     for (int ipt=0; ipt<dim; ++ipt) {
 
@@ -435,7 +439,7 @@ operator() (const BasalFrictionCoefficient_Tag&, const int& cell) const {
                     thickness_field(cell,ipt)*f_p) + (1.0 - f_p)*
                     KU::max(-1.0 * rho_w*bed_topo_field(cell,ipt),0.0) ),0.0);
           if(save_pressure_field) {
-            outN(cell,ipt) = static_cast<EffPressureST>(NVal);
+            outN(cell,ipt) = Albany::convertScalar<EffPressureST>(NVal);
           }
 	} else {
           MeshScalarT thickness(0), bed_topo(0);
@@ -456,7 +460,7 @@ operator() (const BasalFrictionCoefficient_Tag&, const int& cell) const {
                     thickness_field(cell,ipt)*f_p) + (1.0 - f_p)*
                     KU::max(-1.0 * rho_w*bed_topo_field(cell,ipt),0.0) ),0.0);
           if(save_pressure_field) {
-            outN(cell,ipt) = static_cast<EffPressureST>(NVal);
+            outN(cell,ipt) = Albany::convertScalar<EffPressureST>(NVal);
           }
 	} else {
           NVal = 0;

diff --git a/src/landIce/evaluators/LandIce_Gather2DField_Def.hpp b/src/landIce/evaluators/LandIce_Gather2DField_Def.hpp
@@ -299,8 +299,6 @@ evaluateFields (typename PHALTraits::EvalData workset)
   const auto& elem_dof_lids = dof_mgr->elem_dof_lids().host();
   const auto& node_dof_mgr  = workset.disc->getNodeDOFManager();
 
-  const int   neq = dof_mgr->getNumFields();
-
   if (extruded) {
 #ifdef ALBANY_DEBUG
     check_topology(dof_mgr->get_topology());

diff --git a/src/landIce/problems/LandIce_StokesFOBase.cpp b/src/landIce/problems/LandIce_StokesFOBase.cpp
@@ -563,10 +563,10 @@ void StokesFOBase::setFieldsProperties ()
   setSingleFieldProperties(ice_thickness_name,                FRT::Scalar, FST::MeshScalar);
   setSingleFieldProperties(surface_height_name,               FRT::Scalar, FST::MeshScalar);
   setSingleFieldProperties(vertically_averaged_velocity_name, FRT::Vector, FST::Scalar);
-  setSingleFieldProperties(corrected_temperature_name,        FRT::Scalar);
+  setSingleFieldProperties(corrected_temperature_name,        FRT::Scalar, viscosity_use_corrected_temperature ? FST::MeshScalar : FST::Real);
   setSingleFieldProperties(bed_topography_name,               FRT::Scalar, FST::MeshScalar);
   setSingleFieldProperties(body_force_name,                   FRT::Vector);
-  setSingleFieldProperties(flow_factor_name,                  FRT::Scalar);
+  setSingleFieldProperties(flow_factor_name,                  FRT::Scalar, viscosity_use_corrected_temperature ? FST::MeshScalar : FST::Real);
   setSingleFieldProperties(flux_divergence_name,              FRT::Scalar, FST::Scalar);
 
   setSingleFieldProperties(Albany::coord_vec_name,            FRT::Vector, FST::MeshScalar);
@@ -621,7 +621,7 @@ void StokesFOBase::setupEvaluatorRequests ()
       ss_build_interp_ev[ssName][effective_pressure_name][IReq::GRAD_QP_VAL ] = true;
     }
     ss_build_interp_ev[ssName][flow_factor_name][IReq::CELL_TO_SIDE] = true;
-    setSingleFieldProperties(effective_pressure_name, FRT::Scalar);
+    setSingleFieldProperties(effective_pressure_name, FRT::Scalar, FST::MeshScalar);
     setSingleFieldProperties(sliding_velocity_name,   FRT::Scalar, FST::Scalar);
 
     auto& bfc = it->sublist("Basal Friction Coefficient");
@@ -656,7 +656,7 @@ void StokesFOBase::setupEvaluatorRequests ()
             ss_build_interp_ev[ssName][fname][IReq::CELL_TO_SIDE] = true;
             ss_build_interp_ev[ssName][fname][IReq::QP_VAL] = true;
           }
-          setSingleFieldProperties(fname, FRT::Scalar);
+          setSingleFieldProperties(fname, FRT::Scalar, FST::ParamScalar);
         }
       }
     }