pkg/streamice/streamice_driving_stress.F

5ca83cd8f7 Dani*0001 #include "STREAMICE_OPTIONS.h"
                0002
                0003 C---+----1----+----2----+----3----+----4----+----5----+----6----+----7-|--+----|
                0004
                0005 CBOP
96b006450c dngo*0006       SUBROUTINE STREAMICE_DRIVING_STRESS( myThid )
                0007 c      O taudx,
                0008 c      O taudy )
5ca83cd8f7 Dani*0009
                0010 C     /============================================================\
96b006450c dngo*0011 C     | SUBROUTINE                                                 |
5ca83cd8f7 Dani*0012 C     | o                                                          |
                0013 C     |============================================================|
                0014 C     |                                                            |
                0015 C     \============================================================/
                0016       IMPLICIT NONE
                0017
                0018 C     === Global variables ===
                0019 #include "SIZE.h"
                0020 #include "EEPARAMS.h"
                0021 #include "PARAMS.h"
                0022 #include "GRID.h"
                0023 #include "STREAMICE.h"
                0024 #include "STREAMICE_CG.h"
                0025
                0026 C     !INPUT/OUTPUT ARGUMENTS
                0027       INTEGER myThid
96b006450c dngo*0028 c       _RL taudx (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0029 c       _RL taudx (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
5ca83cd8f7 Dani*0030
                0031 #ifdef ALLOW_STREAMICE
                0032 C     LOCAL VARIABLES
96b006450c dngo*0033       INTEGER i, j, bi, bj, Gi, Gj
                0034       LOGICAL at_west_bdry, at_east_bdry,
5ca83cd8f7 Dani*0035      &        at_north_bdry, at_south_bdry
050473d034 Dani*0036       _RL grd_below_sl
96b006450c dngo*0037       _RL sx, sy, diffx, diffy
                0038 c     _RL geom_fac
e4cfce0a6c Dani*0039       _RL i_rhow
                0040       _RL avg_density
                0041      & (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
                0042 #ifdef STREAMICE_FIRN_CORRECTION
                0043       _RL firn_depth, h
                0044 #endif
96b006450c dngo*0045
5ca83cd8f7 Dani*0046       IF (myXGlobalLo.eq.1) at_west_bdry = .true.
                0047       IF (myYGlobalLo.eq.1) at_south_bdry = .true.
96b006450c dngo*0048       IF (myXGlobalLo-1+sNx*nSx.eq.Nx)
5ca83cd8f7 Dani*0049      & at_east_bdry = .false.
96b006450c dngo*0050       IF (myYGlobalLo-1+sNy*nSy.eq.Ny)
5ca83cd8f7 Dani*0051      & at_north_bdry = .false.
e4cfce0a6c Dani*0052 #ifdef STREAMICE_FIRN_CORRECTION
96b006450c dngo*0053       firn_depth = streamice_density *
e4cfce0a6c Dani*0054      &    streamice_firn_correction
                0055      & / (streamice_density-streamice_density_firn)
                0056 #endif
                0057       i_rhow = 1./streamice_density_ocean_avg
5ca83cd8f7 Dani*0058
                0059       DO bj = myByLo(myThid), myByHi(myThid)
                0060        DO bi = myBxLo(myThid), myBxHi(myThid)
bdd8102d3e Dani*0061         DO j=1-OLy+1,sNy+OLy-1
                0062          DO i=1-OLx+1,sNx+OLx-1
5ca83cd8f7 Dani*0063           taudx_SI(i,j,bi,bj) = 0. _d 0
                0064           taudy_SI(i,j,bi,bj) = 0. _d 0
e4cfce0a6c Dani*0065 #ifdef STREAMICE_FIRN_CORRECTION
                0066           if (STREAMICE_apply_firn_correction) then
                0067            if (streamice_hmask(i,j,bi,bj).eq.1) then
                0068             h = h_streamice(i,j,bi,bj)
                0069             if (h.lt.firn_depth) then
                0070              avg_density(i,j,bi,bj) = streamice_density_firn
                0071             else
96b006450c dngo*0072              avg_density(i,j,bi,bj) = streamice_density *
e4cfce0a6c Dani*0073      &        (h - streamice_firn_correction) / h
                0074             endif
                0075            endif
                0076           else
                0077 #endif
                0078            avg_density(i,j,bi,bj) = streamice_density
                0079 #ifdef STREAMICE_FIRN_CORRECTION
                0080           endif
                0081 #endif
5ca83cd8f7 Dani*0082          ENDDO
                0083         ENDDO
                0084        ENDDO
                0085       ENDDO
                0086
                0087       DO bj = myByLo(myThid), myByHi(myThid)
                0088        DO bi = myBxLo(myThid), myBxHi(myThid)
96b006450c dngo*0089
6a651fd056 Dani*0090         DO i=1,sNx
                0091          DO j=1,sNy
5ca83cd8f7 Dani*0092
                0093           diffx = 0. _d 0
                0094           diffy = 0. _d 0
                0095           sx = 0. _d 0
                0096           sy = 0. _d 0
                0097
                0098           Gi = (myXGlobalLo-1)+(bi-1)*sNx+i
                0099           Gj = (myYGlobalLo-1)+(bj-1)*sNy+j
                0100
6a651fd056 Dani*0101           IF (streamice_umask(i,j,bi,bj).eq.1.0) THEN
5ca83cd8f7 Dani*0102
6a651fd056 Dani*0103            IF (streamice_hmask(i-1,j,bi,bj).eq.1.0.AND.
                0104      &      streamice_hmask(i,j,bi,bj).eq.1.0) THEN
5ca83cd8f7 Dani*0105
96b006450c dngo*0106 c             geom_fac = sqrt(rA(i-1,j,bi,bj)*recip_rA(i,j,bi,bj)*
                0107 c      &        dxF(i,j,bi,bj)*recip_dxF(i-1,j,bi,bj))
e4cfce0a6c Dani*0108
96b006450c dngo*0109             taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
26481ac222 Dani*0110      &       0.25 * dyG(i,j,bi,bj) *
e4cfce0a6c Dani*0111      &       gravity *
                0112      &       (H_streamice(i,j,bi,bj)*avg_density(i,j,bi,bj)+
                0113      &        H_streamice(i-1,j,bi,bj)*avg_density(i-1,j,bi,bj)) *
6a651fd056 Dani*0114      &       (surf_el_streamice(i,j,bi,bj)-
26481ac222 Dani*0115      &              surf_el_streamice(i-1,j,bi,bj))
5ca83cd8f7 Dani*0116
e4cfce0a6c Dani*0117 CCC
4a973b81ed Dani*0118             taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
96b006450c dngo*0119      &       streamice_density * gravity *
4a973b81ed Dani*0120      &       streamice_bg_surf_slope_x * .25 * rA(i,j,bi,bj) *
                0121      &       (H_streamice(i-1,j,bi,bj) + H_streamice(i,j,bi,bj))
e4cfce0a6c Dani*0122 CCC
4a973b81ed Dani*0123
6a651fd056 Dani*0124            ELSE IF (streamice_hmask(i-1,j,bi,bj).eq.1.0) THEN
5ca83cd8f7 Dani*0125
6a651fd056 Dani*0126             IF (float_frac_streamice(i-1,j,bi,bj) .eq. 1.0) THEN
5ca83cd8f7 Dani*0127
050473d034 Dani*0128 #ifdef USE_ALT_RLOW
                0129              IF (R_low_si(i-1,j,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0130 #else
050473d034 Dani*0131              IF (R_low(i-1,j,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0132 #endif
050473d034 Dani*0133               grd_below_sl = 1. _d 0
                0134              else
                0135               grd_below_sl = 0. _d 0
                0136              endif
96b006450c dngo*0137
6a651fd056 Dani*0138              taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) +
26481ac222 Dani*0139      &        0.25 * dyG(i,j,bi,bj) *
96b006450c dngo*0140      &        gravity *
e4cfce0a6c Dani*0141      &        (avg_density(i-1,j,bi,bj) * H_streamice(i-1,j,bi,bj)**2-
6a651fd056 Dani*0142 #ifdef USE_ALT_RLOW
96b006450c dngo*0143      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0144      &         R_low_si(i-1,j,bi,bj)**2)
96b006450c dngo*0145 #else
                0146      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0147      &         R_low(i-1,j,bi,bj)**2)
6a651fd056 Dani*0148 #endif
5ca83cd8f7 Dani*0149
                0150             ELSE
                0151
e4cfce0a6c Dani*0152 #ifdef STREAMICE_FIRN_CORRECTION
96b006450c dngo*0153
e4cfce0a6c Dani*0154              if (STREAMICE_apply_firn_correction) then
                0155              if (H_streamice(i-1,j,bi,bj).lt.firn_depth) then
                0156               taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) +
                0157      &         0.25 * dyG(i,j,bi,bj) *
                0158      &         streamice_density_firn * gravity *
                0159      &         (1-streamice_density_firn*i_rhow) *
96b006450c dngo*0160      &         H_streamice(i-1,j,bi,bj)**2
e4cfce0a6c Dani*0161              else
                0162               taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) +
                0163      &         0.25 * dyG(i,j,bi,bj) * gravity * (
96b006450c dngo*0164      &         streamice_density_firn * firn_depth**2 +
e4cfce0a6c Dani*0165      &         (h_streamice(i-1,j,bi,bj)-firn_depth) *
                0166      &          (streamice_density_firn*firn_depth+streamice_density*
                0167      &          (h_streamice(i-1,j,bi,bj)-streamice_firn_correction)) -
                0168      &         streamice_density**2*i_rhow*
                0169      &          (h_streamice(i-1,j,bi,bj)-streamice_firn_correction)**2
96b006450c dngo*0170      &         )
e4cfce0a6c Dani*0171              endif
                0172              else
                0173
                0174 #endif
6a651fd056 Dani*0175              taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) +
26481ac222 Dani*0176      &        0.25 * dyG(i,j,bi,bj) *
6a651fd056 Dani*0177      &        streamice_density * gravity *
e4cfce0a6c Dani*0178      &        (1-streamice_density*i_rhow) *
96b006450c dngo*0179      &         H_streamice(i-1,j,bi,bj)**2
5ca83cd8f7 Dani*0180
e4cfce0a6c Dani*0181 #ifdef STREAMICE_FIRN_CORRECTION
                0182              endif
                0183 #endif
                0184
5ca83cd8f7 Dani*0185             ENDIF
                0186
6a651fd056 Dani*0187            ELSE IF (streamice_hmask(i,j,bi,bj).eq.1.0) THEN
5ca83cd8f7 Dani*0188
6a651fd056 Dani*0189             IF (float_frac_streamice(i,j,bi,bj) .eq. 1.0) THEN
5ca83cd8f7 Dani*0190
050473d034 Dani*0191 #ifdef USE_ALT_RLOW
                0192              IF (R_low_si(i,j,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0193 #else
050473d034 Dani*0194              IF (R_low(i,j,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0195 #endif
                0196               grd_below_sl = 1. _d 0
050473d034 Dani*0197              else
                0198               grd_below_sl = 0. _d 0
                0199              endif
                0200
6a651fd056 Dani*0201              taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
26481ac222 Dani*0202      &        0.25 * dyG(i,j,bi,bj) *
6a651fd056 Dani*0203      &        gravity *
e4cfce0a6c Dani*0204      &        (avg_density(i,j,bi,bj) * H_streamice(i,j,bi,bj)**2 -
6a651fd056 Dani*0205 #ifdef USE_ALT_RLOW
96b006450c dngo*0206      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0207      &         R_low_si(i,j,bi,bj)**2)
96b006450c dngo*0208 #else
                0209      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0210      &         R_low(i,j,bi,bj)**2)
6a651fd056 Dani*0211 #endif
5ca83cd8f7 Dani*0212
6a651fd056 Dani*0213             ELSE
                0214
e4cfce0a6c Dani*0215 #ifdef STREAMICE_FIRN_CORRECTION
                0216
                0217              if (STREAMICE_apply_firn_correction) then
                0218              if (H_streamice(i,j,bi,bj).lt.firn_depth) then
e0987c9b93 Dani*0219               taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
e4cfce0a6c Dani*0220      &         0.25 * dyG(i,j,bi,bj) *
                0221      &         streamice_density_firn * gravity *
                0222      &         (1-streamice_density_firn*i_rhow) *
96b006450c dngo*0223      &         H_streamice(i,j,bi,bj)**2
e4cfce0a6c Dani*0224              else
e0987c9b93 Dani*0225               taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
e4cfce0a6c Dani*0226      &         0.25 * dyG(i,j,bi,bj) * gravity * (
96b006450c dngo*0227      &         streamice_density_firn *  firn_depth**2 +
e4cfce0a6c Dani*0228      &         (h_streamice(i,j,bi,bj)-firn_depth) *
                0229      &          (streamice_density_firn*firn_depth+streamice_density*
                0230      &          (h_streamice(i,j,bi,bj)-streamice_firn_correction)) -
                0231      &         streamice_density**2*i_rhow*
                0232      &          (h_streamice(i,j,bi,bj)-streamice_firn_correction)**2
96b006450c dngo*0233      &         )
e4cfce0a6c Dani*0234              endif
                0235              else
                0236
                0237 #endif
e0987c9b93 Dani*0238              taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
26481ac222 Dani*0239      &        0.25 * dyG(i,j,bi,bj) *
6a651fd056 Dani*0240      &        streamice_density * gravity *
e4cfce0a6c Dani*0241      &        (1-streamice_density*i_rhow) *
96b006450c dngo*0242      &         H_streamice(i,j,bi,bj)**2
e4cfce0a6c Dani*0243 #ifdef STREAMICE_FIRN_CORRECTION
                0244              endif
                0245 #endif
                0246
6a651fd056 Dani*0247             END IF
                0248            END IF
                0249
96b006450c dngo*0250 c cells below
                0251
6a651fd056 Dani*0252            IF (streamice_hmask(i-1,j-1,bi,bj).eq.1.0.AND.
                0253      &      streamice_hmask(i,j-1,bi,bj).eq.1.0) THEN
                0254
96b006450c dngo*0255 c             geom_fac = sqrt(rA(i-1,j-1,bi,bj)*recip_rA(i,j-1,bi,bj)*
                0256 c      &             dxF(i,j-1,bi,bj)*recip_dxF(i-1,j-1,bi,bj))
                0257
                0258             taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
26481ac222 Dani*0259      &       0.25 * dyg(i,j-1,bi,bj) *
e4cfce0a6c Dani*0260      &       gravity *
                0261      &       (H_streamice(i,j-1,bi,bj)*avg_density(i,j-1,bi,bj)+
                0262      &        H_streamice(i-1,j-1,bi,bj)*avg_density(i-1,j-1,bi,bj)) *
6a651fd056 Dani*0263      &       (surf_el_streamice(i,j-1,bi,bj)-
26481ac222 Dani*0264      &              surf_el_streamice(i-1,j-1,bi,bj))
6a651fd056 Dani*0265
4a973b81ed Dani*0266             taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
                0267      &       streamice_density * gravity *
96b006450c dngo*0268      &       streamice_bg_surf_slope_x * .25 * rA(i,j-1,bi,bj) *
4a973b81ed Dani*0269      &       (H_streamice(i-1,j-1,bi,bj) + H_streamice(i,j-1,bi,bj))
                0270
6a651fd056 Dani*0271            ELSE IF (streamice_hmask(i-1,j-1,bi,bj).eq.1.0) THEN
                0272
                0273             IF (float_frac_streamice(i-1,j-1,bi,bj) .eq. 1.0) THEN
                0274
050473d034 Dani*0275 #ifdef USE_ALT_RLOW
                0276              IF (R_low_si(i-1,j-1,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0277 #else
050473d034 Dani*0278              IF (R_low(i-1,j-1,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0279 #endif
050473d034 Dani*0280               grd_below_sl = 1. _d 0
                0281              else
                0282               grd_below_sl = 0. _d 0
                0283              endif
                0284
6a651fd056 Dani*0285              taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) +
26481ac222 Dani*0286      &        0.25 * dyg(i,j-1,bi,bj) *
96b006450c dngo*0287      &        gravity *
e4cfce0a6c Dani*0288      &        (avg_density(i-1,j-1,bi,bj)*H_streamice(i-1,j-1,bi,bj)**2-
6a651fd056 Dani*0289 #ifdef USE_ALT_RLOW
96b006450c dngo*0290      &         streamice_density_ocean_avg*grd_below_sl *
050473d034 Dani*0291      &         R_low_si(i-1,j-1,bi,bj)**2)
96b006450c dngo*0292 #else
                0293      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0294      &         R_low(i-1,j-1,bi,bj)**2)
6a651fd056 Dani*0295 #endif
5ca83cd8f7 Dani*0296
                0297             ELSE
                0298
e4cfce0a6c Dani*0299 #ifdef STREAMICE_FIRN_CORRECTION
                0300
                0301              if (STREAMICE_apply_firn_correction) then
                0302              if (H_streamice(i-1,j-1,bi,bj).lt.firn_depth) then
                0303               taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) +
                0304      &         0.25 * dyG(i,j-1,bi,bj) *
                0305      &         streamice_density_firn * gravity *
                0306      &         (1-streamice_density_firn*i_rhow) *
96b006450c dngo*0307      &         H_streamice(i-1,j-1,bi,bj)**2
e4cfce0a6c Dani*0308              else
                0309               taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) +
                0310      &         0.25 * dyG(i,j-1,bi,bj) * gravity * (
96b006450c dngo*0311      &         streamice_density_firn * firn_depth**2 +
e4cfce0a6c Dani*0312      &         (h_streamice(i-1,j-1,bi,bj)-firn_depth) *
                0313      &          (streamice_density_firn*firn_depth+streamice_density*
                0314      &          (h_streamice(i-1,j-1,bi,bj)-streamice_firn_correction))-
                0315      &         streamice_density**2*i_rhow*
                0316      &         (h_streamice(i-1,j-1,bi,bj)-streamice_firn_correction)**2
96b006450c dngo*0317      &         )
e4cfce0a6c Dani*0318              endif
                0319              else
                0320
                0321 #endif
6a651fd056 Dani*0322              taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) +
e4cfce0a6c Dani*0323      &        0.25 * dyG(i,j-1,bi,bj) *
6a651fd056 Dani*0324      &        streamice_density * gravity *
e4cfce0a6c Dani*0325      &        (1-streamice_density*i_rhow) *
96b006450c dngo*0326      &         H_streamice(i-1,j-1,bi,bj)**2
e4cfce0a6c Dani*0327 #ifdef STREAMICE_FIRN_CORRECTION
                0328              endif
                0329 #endif
                0330
5ca83cd8f7 Dani*0331             ENDIF
                0332
6a651fd056 Dani*0333            ELSE IF (streamice_hmask(i,j-1,bi,bj).eq.1.0) THEN
5ca83cd8f7 Dani*0334
6a651fd056 Dani*0335             IF (float_frac_streamice(i,j-1,bi,bj) .eq. 1.0) THEN
5ca83cd8f7 Dani*0336
050473d034 Dani*0337 #ifdef USE_ALT_RLOW
                0338              IF (R_low_si(i,j-1,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0339 #else
050473d034 Dani*0340              IF (R_low(i,j-1,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0341 #endif
050473d034 Dani*0342               grd_below_sl = 1. _d 0
                0343              else
                0344               grd_below_sl = 0. _d 0
                0345              endif
                0346
6a651fd056 Dani*0347              taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
26481ac222 Dani*0348      &        0.25 * dyg(i,j-1,bi,bj) *
6a651fd056 Dani*0349      &        gravity *
e4cfce0a6c Dani*0350      &        (avg_density(i,j-1,bi,bj) * H_streamice(i,j-1,bi,bj)**2 -
6a651fd056 Dani*0351 #ifdef USE_ALT_RLOW
96b006450c dngo*0352      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0353      &         R_low_si(i,j-1,bi,bj)**2)
96b006450c dngo*0354 #else
                0355      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0356      &         R_low(i,j-1,bi,bj)**2)
6a651fd056 Dani*0357 #endif
5ca83cd8f7 Dani*0358
                0359             ELSE
                0360
e4cfce0a6c Dani*0361 #ifdef STREAMICE_FIRN_CORRECTION
                0362              if (STREAMICE_apply_firn_correction) then
                0363              if (H_streamice(i,j-1,bi,bj).lt.firn_depth) then
e0987c9b93 Dani*0364               taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
e4cfce0a6c Dani*0365      &         0.25 * dyG(i,j-1,bi,bj) *
                0366      &         streamice_density_firn * gravity *
                0367      &         (1-streamice_density_firn*i_rhow) *
96b006450c dngo*0368      &         H_streamice(i,j-1,bi,bj)**2
e4cfce0a6c Dani*0369              else
e0987c9b93 Dani*0370               taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
e4cfce0a6c Dani*0371      &         0.25 * dyG(i,j-1,bi,bj) * gravity * (
96b006450c dngo*0372      &         streamice_density_firn * firn_depth**2 +
e4cfce0a6c Dani*0373      &         (h_streamice(i,j-1,bi,bj)-firn_depth) *
                0374      &          (streamice_density_firn*firn_depth+streamice_density*
                0375      &          (h_streamice(i,j-1,bi,bj)-streamice_firn_correction))-
                0376      &         streamice_density**2*i_rhow*
                0377      &         (h_streamice(i,j-1,bi,bj)-streamice_firn_correction)**2
96b006450c dngo*0378      &         )
e4cfce0a6c Dani*0379              endif
                0380              else
                0381
                0382 #endif
e0987c9b93 Dani*0383              taudx_si(i,j,bi,bj) = taudx_si(i,j,bi,bj) -
e4cfce0a6c Dani*0384      &        0.25 * dyG(i,j-1,bi,bj) *
6a651fd056 Dani*0385      &        streamice_density * gravity *
e4cfce0a6c Dani*0386      &        (1-streamice_density*i_rhow) *
96b006450c dngo*0387      &         H_streamice(i,j-1,bi,bj)**2
e4cfce0a6c Dani*0388 #ifdef STREAMICE_FIRN_CORRECTION
                0389              endif
                0390 #endif
                0391
96b006450c dngo*0392             END IF
6a651fd056 Dani*0393            END IF
96b006450c dngo*0394           END IF       ! if umask==1
5ca83cd8f7 Dani*0395
96b006450c dngo*0396 c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
5ca83cd8f7 Dani*0397
6a651fd056 Dani*0398           IF (streamice_vmask(i,j,bi,bj).eq.1.0) THEN
                0399
                0400            IF (streamice_hmask(i,j-1,bi,bj).eq.1.0.AND.
                0401      &      streamice_hmask(i,j,bi,bj).eq.1.0) THEN
                0402
96b006450c dngo*0403 c             geom_fac = sqrt(rA(i,j-1,bi,bj)*recip_rA(i,j,bi,bj)*
                0404 c      &        dxF(i,j,bi,bj)*recip_dyF(i,j-1,bi,bj))
                0405
                0406             taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
26481ac222 Dani*0407      &       0.25 * dxG(i,j,bi,bj) *
e4cfce0a6c Dani*0408      &       gravity *
                0409      &       (H_streamice(i,j,bi,bj)*avg_density(i,j,bi,bj)+
                0410      &        H_streamice(i,j-1,bi,bj)*avg_density(i,j-1,bi,bj))*
6a651fd056 Dani*0411      &       (surf_el_streamice(i,j,bi,bj)-
26481ac222 Dani*0412      &              surf_el_streamice(i,j-1,bi,bj))
6a651fd056 Dani*0413
4a973b81ed Dani*0414             taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
                0415      &       streamice_density * gravity *
96b006450c dngo*0416      &       streamice_bg_surf_slope_y * .25 * rA(i,j,bi,bj) *
4a973b81ed Dani*0417      &       (H_streamice(i,j-1,bi,bj) + H_streamice(i,j,bi,bj))
                0418
6a651fd056 Dani*0419            ELSE IF (streamice_hmask(i,j-1,bi,bj).eq.1.0) THEN
                0420
                0421             IF (float_frac_streamice(i,j-1,bi,bj) .eq. 1.0) THEN
                0422
050473d034 Dani*0423 #ifdef USE_ALT_RLOW
                0424              IF (R_low_si(i,j-1,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0425 #else
050473d034 Dani*0426              IF (R_low(i,j-1,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0427 #endif
050473d034 Dani*0428               grd_below_sl = 1. _d 0
                0429              else
                0430               grd_below_sl = 0. _d 0
                0431              endif
                0432
6a651fd056 Dani*0433              taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) +
26481ac222 Dani*0434      &        0.25 * dxG(i,j,bi,bj) *
96b006450c dngo*0435      &        gravity *
e4cfce0a6c Dani*0436      &        (avg_density(i,j-1,bi,bj) * H_streamice(i,j-1,bi,bj)**2 -
6a651fd056 Dani*0437 #ifdef USE_ALT_RLOW
96b006450c dngo*0438      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0439      &         R_low_si(i,j-1,bi,bj)**2)
96b006450c dngo*0440 #else
                0441      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0442      &         R_low(i,j-1,bi,bj)**2)
6a651fd056 Dani*0443 #endif
5ca83cd8f7 Dani*0444
                0445             ELSE
6a651fd056 Dani*0446
e4cfce0a6c Dani*0447 #ifdef STREAMICE_FIRN_CORRECTION
                0448              if (STREAMICE_apply_firn_correction) then
                0449              if (H_streamice(i,j-1,bi,bj).lt.firn_depth) then
751a06e847 Dani*0450               taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) +
c04e36832e Dani*0451      &         0.25 * dxG(i,j,bi,bj) *
e4cfce0a6c Dani*0452      &         streamice_density_firn * gravity *
                0453      &         (1-streamice_density_firn*i_rhow) *
96b006450c dngo*0454      &         H_streamice(i,j-1,bi,bj)**2
e4cfce0a6c Dani*0455              else
751a06e847 Dani*0456               taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) +
c04e36832e Dani*0457      &         0.25 * dxG(i,j,bi,bj) * gravity * (
96b006450c dngo*0458      &         streamice_density_firn * firn_depth**2 +
e4cfce0a6c Dani*0459      &         (h_streamice(i,j-1,bi,bj)-firn_depth) *
                0460      &          (streamice_density_firn*firn_depth+streamice_density*
                0461      &          (h_streamice(i,j-1,bi,bj)-streamice_firn_correction))-
                0462      &         streamice_density**2*i_rhow*
                0463      &         (h_streamice(i,j-1,bi,bj)-streamice_firn_correction)**2
96b006450c dngo*0464      &         )
e4cfce0a6c Dani*0465              endif
                0466              else
                0467
                0468 #endif
751a06e847 Dani*0469              taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) +
c04e36832e Dani*0470      &        0.25 * dxG(i,j,bi,bj) *
6a651fd056 Dani*0471      &        streamice_density * gravity *
e4cfce0a6c Dani*0472      &        (1-streamice_density*i_rhow) *
96b006450c dngo*0473      &         H_streamice(i,j-1,bi,bj)**2
e4cfce0a6c Dani*0474 #ifdef STREAMICE_FIRN_CORRECTION
                0475              endif
                0476 #endif
                0477
5ca83cd8f7 Dani*0478             ENDIF
6a651fd056 Dani*0479
                0480            ELSE IF (streamice_hmask(i,j,bi,bj).eq.1.0) THEN
                0481
                0482             IF (float_frac_streamice(i,j,bi,bj) .eq. 1.0) THEN
                0483
050473d034 Dani*0484 #ifdef USE_ALT_RLOW
                0485              IF (R_low_si(i,j,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0486 #else
050473d034 Dani*0487              IF (R_low(i,j,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0488 #endif
050473d034 Dani*0489               grd_below_sl = 1. _d 0
                0490              else
                0491               grd_below_sl = 0. _d 0
                0492              endif
                0493
6a651fd056 Dani*0494              taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
26481ac222 Dani*0495      &        0.25 * dxG(i,j,bi,bj) *
6a651fd056 Dani*0496      &        gravity *
e4cfce0a6c Dani*0497      &        (avg_density(i,j,bi,bj) * H_streamice(i,j,bi,bj)**2 -
6a651fd056 Dani*0498 #ifdef USE_ALT_RLOW
96b006450c dngo*0499      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0500      &         R_low_si(i,j,bi,bj)**2)
96b006450c dngo*0501 #else
                0502      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0503      &         R_low(i,j,bi,bj)**2)
6a651fd056 Dani*0504 #endif
                0505
5ca83cd8f7 Dani*0506             ELSE
6a651fd056 Dani*0507
e4cfce0a6c Dani*0508 #ifdef STREAMICE_FIRN_CORRECTION
                0509              if (STREAMICE_apply_firn_correction) then
                0510              if (H_streamice(i,j,bi,bj).lt.firn_depth) then
751a06e847 Dani*0511               taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
c04e36832e Dani*0512      &         0.25 * dxG(i,j,bi,bj) *
e4cfce0a6c Dani*0513      &         streamice_density_firn * gravity *
                0514      &         (1-streamice_density_firn*i_rhow) *
96b006450c dngo*0515      &         H_streamice(i,j,bi,bj)**2
e4cfce0a6c Dani*0516              else
751a06e847 Dani*0517               taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
c04e36832e Dani*0518      &         0.25 * dxG(i,j,bi,bj) * gravity * (
96b006450c dngo*0519      &         streamice_density_firn * firn_depth**2 +
e4cfce0a6c Dani*0520      &         (h_streamice(i,j,bi,bj)-firn_depth) *
                0521      &          (streamice_density_firn*firn_depth+streamice_density*
                0522      &          (h_streamice(i,j,bi,bj)-streamice_firn_correction))-
                0523      &         streamice_density**2*i_rhow*
                0524      &         (h_streamice(i,j,bi,bj)-streamice_firn_correction)**2
96b006450c dngo*0525      &         )
e4cfce0a6c Dani*0526              endif
                0527              else
                0528
                0529 #endif
751a06e847 Dani*0530              taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
c04e36832e Dani*0531      &        0.25 * dxG(i,j,bi,bj) *
6a651fd056 Dani*0532      &        streamice_density * gravity *
e4cfce0a6c Dani*0533      &        (1-streamice_density*i_rhow) *
96b006450c dngo*0534      &         H_streamice(i,j,bi,bj)**2
e4cfce0a6c Dani*0535 #ifdef STREAMICE_FIRN_CORRECTION
                0536              endif
                0537 #endif
                0538
6a651fd056 Dani*0539             END IF
                0540            END IF
                0541
96b006450c dngo*0542 c cells to left
                0543
6a651fd056 Dani*0544            IF (streamice_hmask(i-1,j-1,bi,bj).eq.1.0.AND.
                0545      &      streamice_hmask(i-1,j,bi,bj).eq.1.0) THEN
96b006450c dngo*0546
                0547 c             geom_fac = sqrt(rA(i-1,j-1,bi,bj)*recip_rA(i-1,j,bi,bj)*
                0548 c      &        dxF(i-1,j,bi,bj)*recip_dxF(i-1,j-1,bi,bj))
6a651fd056 Dani*0549
                0550             taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
26481ac222 Dani*0551      &       0.25 * dxG(i-1,j,bi,bj) *
e4cfce0a6c Dani*0552      &       gravity *
                0553      &       (H_streamice(i-1,j,bi,bj)*avg_density(i-1,j,bi,bj)+
                0554      &        H_streamice(i-1,j-1,bi,bj)*avg_density(i-1,j-1,bi,bj))*
6a651fd056 Dani*0555      &       (surf_el_streamice(i-1,j,bi,bj)-
26481ac222 Dani*0556      &              surf_el_streamice(i-1,j-1,bi,bj))
6a651fd056 Dani*0557
4a973b81ed Dani*0558             taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
                0559      &       streamice_density * gravity *
26481ac222 Dani*0560      &       streamice_bg_surf_slope_y * .25 * rA(i-1,j,bi,bj) *
4a973b81ed Dani*0561      &       (H_streamice(i-1,j-1,bi,bj) + H_streamice(i-1,j,bi,bj))
                0562
6a651fd056 Dani*0563            ELSE IF (streamice_hmask(i-1,j-1,bi,bj).eq.1.0) THEN
                0564
                0565             IF (float_frac_streamice(i-1,j-1,bi,bj) .eq. 1.0) THEN
                0566
050473d034 Dani*0567 #ifdef USE_ALT_RLOW
                0568              IF (R_low_si(i-1,j-1,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0569 #else
050473d034 Dani*0570              IF (R_low(i-1,j-1,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0571 #endif
050473d034 Dani*0572               grd_below_sl = 1. _d 0
                0573              else
                0574               grd_below_sl = 0. _d 0
                0575              endif
                0576
6a651fd056 Dani*0577              taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) +
26481ac222 Dani*0578      &        0.25 * dxG(i-1,j,bi,bj) *
96b006450c dngo*0579      &        gravity *
e4cfce0a6c Dani*0580      &        (avg_density(i-1,j-1,bi,bj)*H_streamice(i-1,j-1,bi,bj)**2-
6a651fd056 Dani*0581 #ifdef USE_ALT_RLOW
96b006450c dngo*0582      &         streamice_density_ocean_avg*grd_below_sl *
050473d034 Dani*0583      &         R_low_si(i-1,j-1,bi,bj)**2)
96b006450c dngo*0584 #else
                0585      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0586      &         R_low(i-1,j-1,bi,bj)**2)
6a651fd056 Dani*0587 #endif
                0588
5ca83cd8f7 Dani*0589             ELSE
6a651fd056 Dani*0590
e4cfce0a6c Dani*0591 #ifdef STREAMICE_FIRN_CORRECTION
                0592              if (STREAMICE_apply_firn_correction) then
                0593              if (H_streamice(i-1,j-1,bi,bj).lt.firn_depth) then
751a06e847 Dani*0594               taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) +
c04e36832e Dani*0595      &         0.25 * dxG(i-1,j,bi,bj) *
e4cfce0a6c Dani*0596      &         streamice_density_firn * gravity *
                0597      &         (1-streamice_density_firn*i_rhow) *
96b006450c dngo*0598      &         H_streamice(i-1,j-1,bi,bj)**2
e4cfce0a6c Dani*0599              else
751a06e847 Dani*0600               taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) +
c04e36832e Dani*0601      &         0.25 * dxG(i-1,j,bi,bj) * gravity * (
96b006450c dngo*0602      &         streamice_density_firn * firn_depth**2 +
e4cfce0a6c Dani*0603      &         (h_streamice(i-1,j-1,bi,bj)-firn_depth) *
                0604      &          (streamice_density_firn*firn_depth+streamice_density*
                0605      &          (h_streamice(i-1,j-1,bi,bj)-streamice_firn_correction))-
                0606      &         streamice_density**2*i_rhow*
                0607      &        (h_streamice(i-1,j-1,bi,bj)-streamice_firn_correction)**2
96b006450c dngo*0608      &         )
e4cfce0a6c Dani*0609              endif
                0610              else
                0611 #endif
751a06e847 Dani*0612              taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) +
c04e36832e Dani*0613      &        0.25 * dxG(i-1,j,bi,bj) *
6a651fd056 Dani*0614      &        streamice_density * gravity *
e4cfce0a6c Dani*0615      &        (1-streamice_density*i_rhow) *
96b006450c dngo*0616      &         H_streamice(i-1,j-1,bi,bj)**2
e4cfce0a6c Dani*0617 #ifdef STREAMICE_FIRN_CORRECTION
                0618              endif
                0619 #endif
                0620
5ca83cd8f7 Dani*0621             ENDIF
6a651fd056 Dani*0622
                0623            ELSE IF (streamice_hmask(i-1,j,bi,bj).eq.1.0) THEN
                0624
                0625             IF (float_frac_streamice(i-1,j,bi,bj) .eq. 1.0) THEN
                0626
050473d034 Dani*0627 #ifdef USE_ALT_RLOW
                0628              IF (R_low_si(i-1,j,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0629 #else
050473d034 Dani*0630              IF (R_low(i-1,j,bi,bj) .lt. 0. _d 0) then
96b006450c dngo*0631 #endif
050473d034 Dani*0632               grd_below_sl = 1. _d 0
                0633              else
                0634               grd_below_sl = 0. _d 0
                0635              endif
                0636
6a651fd056 Dani*0637              taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
26481ac222 Dani*0638      &        0.25 * dxG(i-1,j,bi,bj) *
6a651fd056 Dani*0639      &        gravity *
e4cfce0a6c Dani*0640      &        (avg_density(i-1,j,bi,bj) * H_streamice(i-1,j,bi,bj)**2 -
5ca83cd8f7 Dani*0641 #ifdef USE_ALT_RLOW
96b006450c dngo*0642      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0643      &         R_low_si(i-1,j,bi,bj)**2)
96b006450c dngo*0644 #else
                0645      &         streamice_density_ocean_avg * grd_below_sl *
050473d034 Dani*0646      &         R_low(i-1,j,bi,bj)**2)
5ca83cd8f7 Dani*0647 #endif
6a651fd056 Dani*0648
                0649             ELSE
                0650
e4cfce0a6c Dani*0651 #ifdef STREAMICE_FIRN_CORRECTION
                0652              if (STREAMICE_apply_firn_correction) then
                0653              if (H_streamice(i-1,j,bi,bj).lt.firn_depth) then
751a06e847 Dani*0654               taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
c04e36832e Dani*0655      &         0.25 * dxG(i-1,j,bi,bj) *
e4cfce0a6c Dani*0656      &         streamice_density_firn * gravity *
                0657      &         (1-streamice_density_firn*i_rhow) *
96b006450c dngo*0658      &         H_streamice(i-1,j,bi,bj)**2
e4cfce0a6c Dani*0659              else
751a06e847 Dani*0660               taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
c04e36832e Dani*0661      &         0.25 * dxG(i-1,j,bi,bj) * gravity * (
96b006450c dngo*0662      &         streamice_density_firn * firn_depth**2 +
e4cfce0a6c Dani*0663      &         (h_streamice(i-1,j,bi,bj)-firn_depth) *
                0664      &          (streamice_density_firn*firn_depth+streamice_density*
                0665      &          (h_streamice(i-1,j,bi,bj)-streamice_firn_correction))-
                0666      &         streamice_density**2*i_rhow*
                0667      &        (h_streamice(i-1,j,bi,bj)-streamice_firn_correction)**2
96b006450c dngo*0668      &         )
e4cfce0a6c Dani*0669              endif
                0670              else
                0671 #endif
751a06e847 Dani*0672              taudy_si(i,j,bi,bj) = taudy_si(i,j,bi,bj) -
c04e36832e Dani*0673      &        0.25 * dxG(i-1,j,bi,bj) *
6a651fd056 Dani*0674      &        streamice_density * gravity *
e4cfce0a6c Dani*0675      &        (1-streamice_density*i_rhow) *
96b006450c dngo*0676      &         H_streamice(i-1,j,bi,bj)**2
e4cfce0a6c Dani*0677 #ifdef STREAMICE_FIRN_CORRECTION
                0678              endif
                0679 #endif
                0680
96b006450c dngo*0681             END IF
6a651fd056 Dani*0682            END IF
                0683           END IF      ! if vmask ==1
                0684
5ca83cd8f7 Dani*0685          ENDDO
                0686         ENDDO
                0687        ENDDO
                0688       ENDDO
                0689
96b006450c dngo*0690 c      taudx_SI (1,1,1,1) = taudx_SI (1,1,1,1) +
                0691 c     & streamice_v_normal_pert (1,1,1,1)
                0692 c       call write_fld_xy_rl("TAUDX_SI","",taudx_si,0,myThid)
                0693 c       call write_fld_xy_rl("TAUDY_SI","",taudy_si,0,myThid)
5ca83cd8f7 Dani*0694
                0695 #endif
                0696       RETURN
                0697       END