From 9ee8c56a5f0b3e8f076dd13210a9be8172fd3754 Mon Sep 17 00:00:00 2001
From: =?UTF-8?q?Einar=20=C3=96rn=20=C3=93lason?= <einar.olason@nersc.no>
Date: Fri, 6 Dec 2024 06:18:16 +0100
Subject: [PATCH] Fixes to getIceOceanStress functions

They didn't work properly because of changes in the structure of params.
---
 dynamics/src/include/BrittleCGDynamicsKernel.hpp |  6 ++++--
 dynamics/src/include/FreeDriftDynamicsKernel.hpp | 10 ++++++----
 dynamics/src/include/VPCGDynamicsKernel.hpp      |  8 +++++---
 3 files changed, 15 insertions(+), 9 deletions(-)

diff --git a/dynamics/src/include/BrittleCGDynamicsKernel.hpp b/dynamics/src/include/BrittleCGDynamicsKernel.hpp
index a09e2965e..d366dc1c7 100644
--- a/dynamics/src/include/BrittleCGDynamicsKernel.hpp
+++ b/dynamics/src/include/BrittleCGDynamicsKernel.hpp
@@ -1,7 +1,7 @@
 /*!
  * @file BrittleCGDynamicsKernel.hpp
  *
- * @date 19 Nov 2024
+ * @date 06 Dec 2024
  * @author Tim Spain <timothy.spain@nersc.no>
  * @author Einar Ólason <einar.olason@nersc.no>
  */
@@ -169,11 +169,13 @@ template <int DGadvection> class BrittleCGDynamicsKernel : public CGDynamicsKern
         taux.resizeLike(avgU);
         tauy.resizeLike(avgV);
 
+        const double FOcean = params.COcean * params.rhoOcean;
+
 #pragma omp parallel for
         for (int i = 0; i < taux.rows(); ++i) {
             const double uOceanRel = uOcean(i) - avgU(i);
             const double vOceanRel = vOcean(i) - avgV(i);
-            const double cPrime = params.F_ocean * std::hypot(uOceanRel, vOceanRel);
+            const double cPrime = FOcean * std::hypot(uOceanRel, vOceanRel);
 
             taux(i) = cPrime * (uOceanRel * cosOceanAngle - vOceanRel * sinOceanAngle);
             tauy(i) = cPrime * (vOceanRel * cosOceanAngle + uOceanRel * sinOceanAngle);
diff --git a/dynamics/src/include/FreeDriftDynamicsKernel.hpp b/dynamics/src/include/FreeDriftDynamicsKernel.hpp
index 2ebd7d160..b7972fa09 100644
--- a/dynamics/src/include/FreeDriftDynamicsKernel.hpp
+++ b/dynamics/src/include/FreeDriftDynamicsKernel.hpp
@@ -4,7 +4,7 @@
  * Implementation of "classic free drift", where we ignore all \rho h terms in the momentum
  * equation. This is equivalent to assuming that the ice is very thin.
  *
- * @date 19 Nov 2024
+ * @date 06 Dec 2024
  * @author Tim Spain <timothy.spain@nersc.no>
  * @author Einar Ólason <einar.olason@nersc.no>
  */
@@ -76,11 +76,13 @@ template <int DGadvection> class FreeDriftDynamicsKernel : public CGDynamicsKern
         taux.resizeLike(u);
         tauy.resizeLike(v);
 
+        const double FOcean = params.COcean * params.rhoOcean;
+
 #pragma omp parallel for
         for (int i = 0; i < taux.rows(); ++i) {
             const double uOceanRel = uOcean(i) - u(i);
             const double vOceanRel = vOcean(i) - v(i);
-            const double cPrime = params.F_ocean * std::hypot(uOceanRel, vOceanRel);
+            const double cPrime = FOcean * std::hypot(uOceanRel, vOceanRel);
 
             taux(i) = cPrime * (uOceanRel * cosOceanAngle - vOceanRel * sinOceanAngle);
             tauy(i) = cPrime * (vOceanRel * cosOceanAngle + uOceanRel * sinOceanAngle);
@@ -91,8 +93,8 @@ template <int DGadvection> class FreeDriftDynamicsKernel : public CGDynamicsKern
         } else if (name == vIOStressName) {
             return tauy;
         } else {
-            throw std::logic_error(std::string(__func__) + " called with an unknown argument " + name
-                                   + ". Only " + uIOStressName + " and " + vIOStressName + " are supported\n");
+            throw std::logic_error(std::string(__func__) + " called with an unknown argument "
+                + name + ". Only " + uIOStressName + " and " + vIOStressName + " are supported\n");
         }
     }
 };
diff --git a/dynamics/src/include/VPCGDynamicsKernel.hpp b/dynamics/src/include/VPCGDynamicsKernel.hpp
index e67465d44..2db4ac630 100644
--- a/dynamics/src/include/VPCGDynamicsKernel.hpp
+++ b/dynamics/src/include/VPCGDynamicsKernel.hpp
@@ -1,7 +1,7 @@
 /*!
  * @file VPCGDynamicsKernel.hpp
  *
- * @date 19 Nov 2024
+ * @date 06 Dec 2024
  * @author Tim Spain <timothy.spain@nersc.no>
  */
 
@@ -95,14 +95,16 @@ template <int DGadvection> class VPCGDynamicsKernel : public CGDynamicsKernel<DG
         taux.resizeLike(u);
         tauy.resizeLike(v);
 
+        const double FOcean = params.COcean * params.rhoOcean;
+
 #pragma omp parallel for
         for (int i = 0; i < taux.rows(); ++i) {
             double uOcnRel = u(i) - uOcean(i);
             double vOcnRel = v(i) - vOcean(i);
             double absocn = sqrt(SQR(uOcnRel) + SQR(vOcnRel));
 
-            taux(i) = params.F_ocean * absocn * uOcnRel;
-            tauy(i) = params.F_ocean * absocn * vOcnRel;
+            taux(i) = FOcean * absocn * uOcnRel;
+            tauy(i) = FOcean * absocn * vOcnRel;
         }
 
         if (name == uIOStressName) {