From ed8cee95f32a4be38cbf81d74effc31b951ca256 Mon Sep 17 00:00:00 2001
From: Philippe Blain <philippe.blain@canada.ca>
Date: Mon, 15 Nov 2021 16:09:16 -0500
Subject: [PATCH] ice_dyn_shared: add 'viscous_coeffs_and_rep_pressure_T'
 subroutine

Add a subroutine mimicking what 'viscous_coeffs_and_rep_pressure' does,
but at a single location.

Name it '*_T' since it's going to be used to compute the viscous
coefficients and replacement pressure at the T point.
---
 cicecore/cicedynB/dynamics/ice_dyn_shared.F90 | 55 +++++++++++++++++++
 1 file changed, 55 insertions(+)

diff --git a/cicecore/cicedynB/dynamics/ice_dyn_shared.F90 b/cicecore/cicedynB/dynamics/ice_dyn_shared.F90
index ab6466ee2..9f1790df0 100755
--- a/cicecore/cicedynB/dynamics/ice_dyn_shared.F90
+++ b/cicecore/cicedynB/dynamics/ice_dyn_shared.F90
@@ -1448,6 +1448,61 @@ subroutine viscous_coeffs_and_rep_pressure (strength,  tinyarea, &
       
        end subroutine viscous_coeffs_and_rep_pressure
       
+
+ !=======================================================================
+ ! Computes viscous coefficients and replacement pressure for stress 
+ ! calculations. Note that tensile strength is included here.
+ !
+ ! Hibler, W. D. (1979). A dynamic thermodynamic sea ice model. J. Phys.
+ ! Oceanogr., 9, 817-846.
+ !
+ ! Konig Beatty, C. and Holland, D. M.  (2010). Modeling landfast ice by
+ ! adding tensile strength. J. Phys. Oceanogr. 40, 185-198.
+ !
+ ! Lemieux, J. F. et al. (2016). Improving the simulation of landfast ice
+ ! by combining tensile strength and a parameterization for grounded ridges.
+ ! J. Geophys. Res. Oceans, 121, 7354-7368.
+
+      subroutine viscous_coeffs_and_rep_pressure_T (strength, tinyarea, &
+                                                    Delta   , zetax2  , &
+                                                    etax2   , rep_prs , &
+                                                    capping)
+
+      real (kind=dbl_kind), intent(in)::  &
+        strength, tinyarea
+
+      real (kind=dbl_kind), intent(in)::  &
+        Delta
+
+      logical, intent(in):: capping
+
+      real (kind=dbl_kind), intent(out):: &
+        zetax2, etax2, rep_prs ! 2 x visous coeffs, replacement pressure
+
+      ! local variables
+      real (kind=dbl_kind) :: &
+        tmpcalc
+
+      ! NOTE: for comp. efficiency 2 x zeta and 2 x eta are used in the code
+
+!      if (trim(yield_curve) == 'ellipse') then
+
+      if (capping) then
+         tmpcalc = strength/max(Delta,tinyarea)
+      else
+         tmpcalc = strength/(Delta + tinyarea)
+      endif
+
+      zetax2 = (c1+Ktens)*tmpcalc
+      rep_prs = (c1-Ktens)*tmpcalc*Delta
+      etax2 = epp2i*zetax2
+
+!      else
+
+!      endif
+
+       end subroutine viscous_coeffs_and_rep_pressure_T
+
 !=======================================================================
 
 ! Load velocity components into array for boundary updates