File indexing completed on 2023-09-21 05:10:53 UTC
view on githubraw file Latest commit 96b00645 on 2023-09-20 15:15:14 UTC
5ca83cd8f7 Dani*0001 #include "STREAMICE_OPTIONS.h"
9d8b0af494 Jean*0002 #ifdef ALLOW_AUTODIFF
0003 # include "AUTODIFF_OPTIONS.h"
0004 #endif
5ca83cd8f7 Dani*0005
0006
0007
0008
0009 SUBROUTINE STREAMICE_VISC_BETA_HYBRID ( myThid )
0010
0011
9d8b0af494 Jean*0012
5ca83cd8f7 Dani*0013
0014
0015
0016
0017 IMPLICIT NONE
0018
0019
0020 #include "SIZE.h"
0021 #include "GRID.h"
0022 #include "EEPARAMS.h"
0023 #include "PARAMS.h"
0024 #include "STREAMICE.h"
0025 #include "STREAMICE_CG.h"
0026 #ifdef ALLOW_AUTODIFF_TAMC
0027 # include "tamc.h"
0028 #endif
0029
0030
0031 INTEGER myThid
07e785229e dngo*0032
0033
5ca83cd8f7 Dani*0034
0035 #ifdef ALLOW_STREAMICE
0036 #ifdef STREAMICE_HYBRID_STRESS
0037
0038
0039 INTEGER i, j, bi, bj, k, l, m
0040 _RL ux, uy, vx, vy, exx, eyy, exy, unorm, second_inv
0041 _RL ub, vb, fb, mean_u_shear, mean_v_shear, umid, vmid
0042 _RL omega_temp (Nr+1), u_shear(Nr+1), v_shear(Nr+1)
0fbff46b46 dngo*0043 _RL C_fric_val, Bglen_val
1cb54b8236 Dani*0044 #ifdef STREAMICE_FLOWLINE_BUTTRESS
0045 _RL buttr_param, pwr
0046 #endif
07e785229e dngo*0047 #ifdef STREAMICE_COULOMB_SLIDING
0048 _RL effective_stress, hf
0049 _RL ETA_GL_STREAMICE
0fbff46b46 dngo*0050 _RL i_nbasalfric
07e785229e dngo*0051 EXTERNAL ETA_GL_STREAMICE
0052 #endif
0053 #ifdef ALLOW_AUTODIFF_TAMC
0054 INTEGER ikey_1
0055 #endif
7493c1d2cc Patr*0056 _RL STREAMICE_BSTRESS_EXP
07e785229e dngo*0057
0058 EXTERNAL STREAMICE_BSTRESS_EXP
5ca83cd8f7 Dani*0059
1cb54b8236 Dani*0060 #ifdef STREAMICE_FLOWLINE_BUTTRESS
0061 buttr_param = 5**(1./3.) / (streamice_buttr_width/2.)**(4./3.)
0062 pwr = 1./n_glen
0063 #endif
0064
0fbff46b46 dngo*0065 #ifdef STREAMICE_COULOMB_SLIDING
0066 i_nbasalfric = 1. / n_basal_friction
0067 #endif
0068
5ca83cd8f7 Dani*0069 DO bj=myByLo(myThid),myByHi(myThid)
0070 DO bi=myBxLo(myThid),myBxHi(myThid)
0071 DO j=1,sNy
0072 DO i=1,sNx
0073 IF (STREAMICE_hmask(i,j,bi,bj).eq.1) THEN
9d8b0af494 Jean*0074
5ca83cd8f7 Dani*0075 umid = 0
0076 vmid = 0
0077
0078 DO k=0,1
0079 DO l=0,1
9d8b0af494 Jean*0080 umid = umid + 0.25 *
0081 & dxG(i,j+l,bi,bj)*dyG(i+k,j,bi,bj) *
5ca83cd8f7 Dani*0082 & recip_rA(i,j,bi,bj) *
9d8b0af494 Jean*0083 & U_streamice(i+k,j+l,bi,bj)
0084 vmid = vmid + 0.25 *
0085 & dxG(i,j+l,bi,bj)*dyG(i+k,j,bi,bj) *
5ca83cd8f7 Dani*0086 & recip_rA(i,j,bi,bj) *
0087 & V_streamice(i+k,j+l,bi,bj)
0088 ENDDO
0089 ENDDO
9d8b0af494 Jean*0090
0091 ux = (U_streamice(i+1,j+1,bi,bj) +
5ca83cd8f7 Dani*0092 & U_streamice(i+1,j,bi,bj) -
0093 & U_streamice(i,j+1,bi,bj) -
0094 & U_streamice(i,j,bi,bj)) /
0095 & (2. * dxF(i,j,bi,bj))
9d8b0af494 Jean*0096 vx = (V_streamice(i+1,j+1,bi,bj) +
5ca83cd8f7 Dani*0097 & V_streamice(i+1,j,bi,bj) -
0098 & V_streamice(i,j+1,bi,bj) -
0099 & V_streamice(i,j,bi,bj)) /
0100 & (2. * dxF(i,j,bi,bj))
9d8b0af494 Jean*0101 uy = (U_streamice(i+1,j+1,bi,bj) -
5ca83cd8f7 Dani*0102 & U_streamice(i+1,j,bi,bj) +
0103 & U_streamice(i,j+1,bi,bj) -
0104 & U_streamice(i,j,bi,bj)) /
0105 & (2. * dyF(i,j,bi,bj))
9d8b0af494 Jean*0106 vy = (V_streamice(i+1,j+1,bi,bj) -
5ca83cd8f7 Dani*0107 & V_streamice(i+1,j,bi,bj) +
0108 & V_streamice(i,j+1,bi,bj) -
0109 & V_streamice(i,j,bi,bj)) /
0110 & (2. * dyF(i,j,bi,bj))
0111
0112 exx = ux + k1AtC_str(i,j,bi,bj)*vmid
0113 eyy = vy + k2AtC_str(i,j,bi,bj)*umid
9d8b0af494 Jean*0114 exy = .5*(uy+vx) +
5ca83cd8f7 Dani*0115 & k1AtC_str(i,j,bi,bj)*umid + k2AtC_str(i,j,bi,bj)*vmid
0116
0117 visc_streamice (i,j,bi,bj) = 0.0
0118 streamice_omega(i,j,bi,bj) = 0.0
0119 omega_temp (Nr+1) = 0.0
0120 u_shear(Nr+1) = 0.0
0121 v_shear(Nr+1) = 0.0
0122
0123 DO m=Nr,1,-1
0124
0125 #ifdef ALLOW_AUTODIFF_TAMC
20dee61641 Mart*0126 ikey_1 = bi + (bj-1)*nSx + (ikey_dynamics-1)*nSx*nSy
0127 ikey_1 = m + ((i-1) + (j-1)*sNx + (ikey_1-1)*sNx*sNy)*Nr
5ca83cd8f7 Dani*0128
0129
0130 #endif
0131
9d8b0af494 Jean*0132 streamice_vert_shear_uz (m) = streamice_taubx(i,j,bi,bj) /
5ca83cd8f7 Dani*0133 & visc_streamice_full(i,j,m,bi,bj)
0134 & * streamice_sigma_coord(m)
0135
9d8b0af494 Jean*0136 streamice_vert_shear_vz (m) = streamice_tauby(i,j,bi,bj) /
5ca83cd8f7 Dani*0137 & visc_streamice_full(i,j,m,bi,bj)
0138 & * streamice_sigma_coord(m)
0139
9d8b0af494 Jean*0140 second_inv =
5ca83cd8f7 Dani*0141 & sqrt(exx**2+eyy**2+exx*eyy+exy**2+eps_glen_min**2+
0142 & 0.25 * streamice_vert_shear_uz(m)**2 +
0143 & 0.25 * streamice_vert_shear_vz(m)**2)
0144
bdd8102d3e Dani*0145 #ifdef STREAMICE_3D_GLEN_CONST
0fbff46b46 dngo*0146 IF (.not.STREAMICE_use_log_ctrl) THEN
96b006450c dngo*0147 Bglen_val = (B_glen(i,j,m,bi,bj))**2
0fbff46b46 dngo*0148 ELSE
0149 Bglen_val = exp(B_glen(i,j,m,bi,bj))
0150 ENDIF
95afe7199b Dani*0151 #else
0fbff46b46 dngo*0152 IF (.not.STREAMICE_use_log_ctrl) THEN
96b006450c dngo*0153 Bglen_val = (B_glen(i,j,bi,bj))**2
0fbff46b46 dngo*0154 ELSE
0155 Bglen_val = exp(B_glen(i,j,bi,bj))
0156 ENDIF
95afe7199b Dani*0157 #endif
0fbff46b46 dngo*0158
351fd6b6a4 Dani*0159 #if (defined (ALLOW_STREAMICE_OAD_FP))
95afe7199b Dani*0160 visc_full_new_si(i,j,m,bi,bj) =
0161 #else
9d8b0af494 Jean*0162 visc_streamice_full(i,j,m,bi,bj) =
95afe7199b Dani*0163 #endif
0fbff46b46 dngo*0164 & .5 * Bglen_val *
9d8b0af494 Jean*0165 & second_inv**((1-n_glen)/n_glen)
5ca83cd8f7 Dani*0166
0167 visc_streamice (i,j,bi,bj) = visc_streamice (i,j,bi,bj) +
9d8b0af494 Jean*0168 & H_streamice(i,j,bi,bj) * streamice_delsigma (m) *
351fd6b6a4 Dani*0169 #if (defined (ALLOW_STREAMICE_OAD_FP))
95afe7199b Dani*0170 & visc_full_new_si(i,j,m,bi,bj)
0171 #else
5ca83cd8f7 Dani*0172 & visc_streamice_full(i,j,m,bi,bj)
95afe7199b Dani*0173 #endif
5ca83cd8f7 Dani*0174
9d8b0af494 Jean*0175 omega_temp (m) = omega_temp(m+1) +
5ca83cd8f7 Dani*0176 & streamice_sigma_coord(m) * streamice_delsigma(m) /
351fd6b6a4 Dani*0177 #if (defined (ALLOW_STREAMICE_OAD_FP))
95afe7199b Dani*0178 & visc_full_new_si(i,j,m,bi,bj)
0179 #else
5ca83cd8f7 Dani*0180 & visc_streamice_full(i,j,m,bi,bj)
95afe7199b Dani*0181 #endif
5ca83cd8f7 Dani*0182
9d8b0af494 Jean*0183 u_shear (m) = u_shear (m+1) +
5ca83cd8f7 Dani*0184 & streamice_vert_shear_uz (m) * streamice_delsigma (m) *
0185 & H_streamice(i,j,bi,bj)
0186
9d8b0af494 Jean*0187 v_shear (m) = v_shear (m+1) +
5ca83cd8f7 Dani*0188 & streamice_vert_shear_vz (m) * streamice_delsigma (m) *
0189 & H_streamice(i,j,bi,bj)
0190
0191 ENDDO
0192
0193 mean_u_shear = 0.0
0194 mean_v_shear = 0.0
0195
0196 DO m=Nr,1,-1
0197
9d8b0af494 Jean*0198 streamice_omega(i,j,bi,bj) = streamice_omega(i,j,bi,bj) +
5ca83cd8f7 Dani*0199 & streamice_delsigma(m)*(omega_temp(m)+omega_temp(m+1))*.5
0200 & * H_streamice(i,j,bi,bj)**2
0201
9d8b0af494 Jean*0202 mean_u_shear = mean_u_shear +
5ca83cd8f7 Dani*0203 & streamice_delsigma(m)*(u_shear(m)+u_shear(m+1))*.5
0204
9d8b0af494 Jean*0205 mean_v_shear = mean_v_shear +
5ca83cd8f7 Dani*0206 & streamice_delsigma(m)*(v_shear(m)+v_shear(m+1))*.5
0207
0208 ENDDO
0209
9d8b0af494 Jean*0210 streamice_u_surf(i,j,bi,bj) =
5ca83cd8f7 Dani*0211 & u_shear(1) + umid - mean_u_shear
0212
9d8b0af494 Jean*0213 streamice_v_surf(i,j,bi,bj) =
5ca83cd8f7 Dani*0214 & v_shear(1) + vmid - mean_v_shear
0215
0216 ub = umid - streamice_taubx(i,j,bi,bj) *
0217 & streamice_omega(i,j,bi,bj) / H_streamice(i,j,bi,bj)
9d8b0af494 Jean*0218
cceca0b33c Dani*0219 streamice_u_bed (i,j,bi,bj) = ub
5ca83cd8f7 Dani*0220
0221 vb = vmid - streamice_tauby(i,j,bi,bj) *
0222 & streamice_omega(i,j,bi,bj) / H_streamice(i,j,bi,bj)
0223
cceca0b33c Dani*0224 streamice_v_bed (i,j,bi,bj) = vb
0225
5ca83cd8f7 Dani*0226 unorm = sqrt(ub**2+vb**2+eps_u_min**2)
0227
0fbff46b46 dngo*0228 if (.not.STREAMICE_use_log_ctrl) THEN
96b006450c dngo*0229 C_fric_val =
0230 & (C_basal_friction(i,j,bi,bj))**2
0fbff46b46 dngo*0231 else
96b006450c dngo*0232 C_fric_val =
0233 & exp(C_basal_friction(i,j,bi,bj))
0fbff46b46 dngo*0234 endif
0235
07e785229e dngo*0236 #ifdef STREAMICE_COULOMB_SLIDING
0fbff46b46 dngo*0237 IF (.not.streamice_allow_reg_coulomb) THEN
07e785229e dngo*0238 #endif
0fbff46b46 dngo*0239 fb = C_fric_val *
7493c1d2cc Patr*0240 & STREAMICE_BSTRESS_EXP (unorm,n_basal_friction) *
5ca83cd8f7 Dani*0241 & streamice_basal_geom(i,j,bi,bj) *
0242 & float_frac_streamice(i,j,bi,bj)
0243
07e785229e dngo*0244 #ifdef STREAMICE_COULOMB_SLIDING
0fbff46b46 dngo*0245 ELSE
07e785229e dngo*0246
0247 hf = max (0.0, -1.0 *streamice_density_ocean_avg
0248 & / streamice_density * R_low_si(i,j,bi,bj))
0249
0250 effective_stress = streamice_density * gravity *
0fbff46b46 dngo*0251 & (H_streamice(i,j,bi,bj)-hf)
07e785229e dngo*0252 effective_stress= max(0.0, effective_stress)
0fbff46b46 dngo*0253 fb = C_fric_val *
07e785229e dngo*0254 & unorm ** n_basal_friction *
0255 & .5 * effective_stress /
0fbff46b46 dngo*0256 & (C_fric_val**(i_nbasalfric)*unorm+
0257 & (0.5 * effective_stress) ** (i_nbasalfric)) **
07e785229e dngo*0258 & n_basal_friction *
0259 & streamice_basal_geom(i,j,bi,bj) *
0260 & float_frac_streamice(i,j,bi,bj) / unorm
0fbff46b46 dngo*0261 ENDIF
07e785229e dngo*0262
0263 #endif
0264
9d8b0af494 Jean*0265 tau_beta_eff_streamice(i,j,bi,bj) =
5ca83cd8f7 Dani*0266 & fb /
0267 & (1+fb*streamice_omega(i,j,bi,bj)/H_streamice(i,j,bi,bj))
1cb54b8236 Dani*0268 #ifdef STREAMICE_FLOWLINE_BUTTRESS
0269 if (usestreamiceflowlineButtr) then
07e785229e dngo*0270 tau_beta_eff_streamice(i,j,bi,bj) =
0271 & tau_beta_eff_streamice(i,j,bi,bj) +
0272 & buttr_param*B_glen(i,j,bi,bj)**2 *
1cb54b8236 Dani*0273 & STREAMICE_BSTRESS_EXP (unorm,pwr)*
0274 & H_streamice(i,j,bi,bj)
0275 endif
0276 #endif
0277
5ca83cd8f7 Dani*0278 ENDIF
0279 ENDDO
0280 ENDDO
0281 ENDDO
0282 ENDDO
891d3745f3 Dani*0283 #ifdef ALLOW_OPENAD
0284 #ifdef STREAMICE_HYBRID_STRESS
07e785229e dngo*0285 tau_beta_eff_streamice(1,1,1,1) =
891d3745f3 Dani*0286 & tau_beta_eff_streamice(1,1,1,1) +
07e785229e dngo*0287 & 0. * streamice_u_surf(1,1,1,1) +
0288 & 0. * streamice_v_surf(1,1,1,1)
891d3745f3 Dani*0289 #endif
0290 #endif
5ca83cd8f7 Dani*0291
0292 #endif
0293 #endif
0294 RETURN
0295 END
