pkg/seaice/seaice_init_varia.F

809c36b928 Patr*0001 #include "SEAICE_OPTIONS.h"
fc004572cd Jean*0002 #ifdef ALLOW_OBCS
                0003 # include "OBCS_OPTIONS.h"
                0004 #endif
809c36b928 Patr*0005
                0006 CStartOfInterface
7e4f7741f6 Patr*0007       SUBROUTINE SEAICE_INIT_VARIA( myThid )
e5b8ebbabe Jean*0008 C     *==========================================================*
7e4f7741f6 Patr*0009 C     | SUBROUTINE SEAICE_INIT_VARIA                             |
809c36b928 Patr*0010 C     | o Initialization of sea ice model.                       |
e5b8ebbabe Jean*0011 C     *==========================================================*
                0012 C     *==========================================================*
809c36b928 Patr*0013       IMPLICIT NONE
f4cb8ca7db Jean*0014
809c36b928 Patr*0015 C     === Global variables ===
                0016 #include "SIZE.h"
                0017 #include "EEPARAMS.h"
                0018 #include "PARAMS.h"
                0019 #include "GRID.h"
36c84b38c7 Dimi*0020 #include "DYNVARS.h"
f3ce416a61 Jean*0021 #include "FFIELDS.h"
ccaa3c61f4 Patr*0022 #include "SEAICE_SIZE.h"
8c50aa8796 Jean*0023 #include "SEAICE_PARAMS.h"
                0024 #include "SEAICE.h"
ccaa3c61f4 Patr*0025 #include "SEAICE_TRACER.h"
fc004572cd Jean*0026 #ifdef OBCS_UVICE_OLD
a9eb030de2 Jean*0027 # include "OBCS_GRID.h"
ae1fb66b64 Dimi*0028 #endif
809c36b928 Patr*0029
                0030 C     === Routine arguments ===
425e8efc36 Jean*0031 C     myThid :: Thread no. that called this routine.
809c36b928 Patr*0032       INTEGER myThid
                0033 CEndOfInterface
f4cb8ca7db Jean*0034
809c36b928 Patr*0035 C     === Local variables ===
425e8efc36 Jean*0036 C     i,j,k,bi,bj :: Loop counters
809c36b928 Patr*0037
edfdf5fa1d Jean*0038       INTEGER i, j, bi, bj
0320e25227 Mart*0039       INTEGER kSrf
8c50aa8796 Jean*0040       INTEGER k
f50f58ec54 Gael*0041 #ifdef ALLOW_SITRACER
                0042       INTEGER iTr, jTh
                0043 #endif
fc004572cd Jean*0044 #ifdef OBCS_UVICE_OLD
dc8f395c29 Dimi*0045       INTEGER I_obc, J_obc
                0046 #endif /* ALLOW_OBCS */
0dcda34642 Jean*0047 #ifdef SEAICE_CGRID
45315406aa Mart*0048       _RS mask_uice
2f5e8addfd Mart*0049       _RL recip_tensilDepth
0dcda34642 Jean*0050 #endif
486cba1d02 Mart*0051
0320e25227 Mart*0052       IF ( usingPCoords ) THEN
                0053        kSrf = Nr
f4cb8ca7db Jean*0054       ELSE
0320e25227 Mart*0055        kSrf = 1
f4cb8ca7db Jean*0056       ENDIF
a6d2313c04 Dimi*0057
f4cb8ca7db Jean*0058 C--   Initialise all variables in common blocks:
7109a141b2 Patr*0059       DO bj=myByLo(myThid),myByHi(myThid)
                0060        DO bi=myBxLo(myThid),myBxHi(myThid)
04016f2c47 Mart*0061         DO j=1-OLy,sNy+OLy
                0062          DO i=1-OLx,sNx+OLx
772590b63c Mart*0063           HEFF(i,j,bi,bj)=0. _d 0
                0064           AREA(i,j,bi,bj)=0. _d 0
45315406aa Mart*0065           HSNOW(i,j,bi,bj)   = 0. _d 0
86b84a92fc Patr*0066 #ifdef SEAICE_ITD
                0067           DO k=1,nITD
                0068            AREAITD(i,j,k,bi,bj) =0. _d 0
                0069            HEFFITD(i,j,k,bi,bj) =0. _d 0
45315406aa Mart*0070            HSNOWITD(i,j,k,bi,bj)=0. _d 0
86b84a92fc Patr*0071           ENDDO
                0072 #endif
772590b63c Mart*0073           UICE(i,j,bi,bj)=0. _d 0
                0074           VICE(i,j,bi,bj)=0. _d 0
                0075 C
5fe78992ba Mart*0076           uIceNm1(i,j,bi,bj) = 0. _d 0
                0077           vIceNm1(i,j,bi,bj) = 0. _d 0
45315406aa Mart*0078           DWATN  (i,j,bi,bj) = 0. _d 0
                0079 #ifdef SEAICE_CGRID
                0080           stressDivergenceX(i,j,bi,bj) = 0. _d 0
                0081           stressDivergenceY(i,j,bi,bj) = 0. _d 0
                0082 # ifdef SEAICE_ALLOW_EVP
                0083           seaice_sigma1 (i,j,bi,bj) = 0. _d 0
                0084           seaice_sigma2 (i,j,bi,bj) = 0. _d 0
                0085           seaice_sigma12(i,j,bi,bj) = 0. _d 0
                0086 # endif /* SEAICE_ALLOW_EVP */
                0087 #endif
                0088 #if ( defined SEAICE_CGRID || defined SEAICE_BGRID_DYNAMICS )
5fe78992ba Mart*0089           e11    (i,j,bi,bj) = 0. _d 0
                0090           e22    (i,j,bi,bj) = 0. _d 0
                0091           e12    (i,j,bi,bj) = 0. _d 0
                0092           deltaC (i,j,bi,bj) = 0. _d 0
45315406aa Mart*0093           PRESS  (i,j,bi,bj) = 0. _d 0
5fe78992ba Mart*0094           ETA    (i,j,bi,bj) = 0. _d 0
                0095           etaZ   (i,j,bi,bj) = 0. _d 0
                0096           ZETA   (i,j,bi,bj) = 0. _d 0
                0097           FORCEX (i,j,bi,bj) = 0. _d 0
                0098           FORCEY (i,j,bi,bj) = 0. _d 0
45315406aa Mart*0099           tensileStrFac(i,j,bi,bj) = 0. _d 0
                0100           PRESS0(i,j,bi,bj)  = 0. _d 0
                0101           FORCEX0(i,j,bi,bj) = 0. _d 0
                0102           FORCEY0(i,j,bi,bj) = 0. _d 0
                0103           SEAICE_zMax(i,j,bi,bj) = 0. _d 0
                0104           SEAICE_zMin(i,j,bi,bj) = 0. _d 0
                0105 #endif
f4cb8ca7db Jean*0106 #ifdef SEAICE_CGRID
                0107           seaiceMassC(i,j,bi,bj)=0. _d 0
                0108           seaiceMassU(i,j,bi,bj)=0. _d 0
                0109           seaiceMassV(i,j,bi,bj)=0. _d 0
45315406aa Mart*0110 # ifdef SEAICE_ALLOW_FREEDRIFT
                0111           uice_fd(i,j,bi,bj) = 0. _d 0
                0112           vice_fd(i,j,bi,bj) = 0. _d 0
                0113 # endif
                0114 # ifdef SEAICE_ALLOW_BOTTOMDRAG
                0115           CbotC(i,j,bi,bj)   = 0. _d 0
                0116 # endif /* SEAICE_ALLOW_BOTTOMDRAG */
5bb179ddc2 Mart*0117 # ifdef SEAICE_ALLOW_SIDEDRAG
                0118           sideDragU(i,j,bi,bj) = 0. _d 0
                0119           sideDragV(i,j,bi,bj) = 0. _d 0
                0120 # endif /* SEAICE_ALLOW_SIDEDRAG */
45315406aa Mart*0121 #endif /* SEAICE_CGRID */
                0122 #ifdef SEAICE_BGRID_DYNAMICS
                0123           uIceB(i,j,bi,bj)   = 0. _d 0
                0124           vIceB(i,j,bi,bj)   = 0. _d 0
5fe78992ba Mart*0125           AMASS(i,j,bi,bj)   = 0. _d 0
                0126           DAIRN(i,j,bi,bj)   = 0. _d 0
                0127           WINDX(i,j,bi,bj)   = 0. _d 0
                0128           WINDY(i,j,bi,bj)   = 0. _d 0
                0129           GWATX(i,j,bi,bj)   = 0. _d 0
                0130           GWATY(i,j,bi,bj)   = 0. _d 0
45315406aa Mart*0131 #endif /* SEAICE_BGRID_DYNAMICS */
a98c4b8072 Ian *0132 #ifdef SEAICE_VARIABLE_SALINITY
f4cb8ca7db Jean*0133           HSALT(i,j,bi,bj)  = 0. _d 0
f681b7f5d4 Dimi*0134 #endif
f50f58ec54 Gael*0135 #ifdef ALLOW_SITRACER
                0136           DO iTr = 1, SItrMaxNum
                0137            SItracer(i,j,bi,bj,iTr) = 0. _d 0
                0138            SItrBucket(i,j,bi,bj,iTr) = 0. _d 0
                0139 c "ice concentration" tracer that should remain .EQ.1.
                0140            if (SItrName(iTr).EQ.'one') SItracer(i,j,bi,bj,iTr)=1. _d 0
                0141           ENDDO
                0142           DO jTh = 1, 5
                0143            SItrHEFF (i,j,bi,bj,jTh) = 0. _d 0
                0144           ENDDO
bb24b8a3e6 Gael*0145           DO jTh = 1, 3
                0146            SItrAREA (i,j,bi,bj,jTh) = 0. _d 0
                0147           ENDDO
f4cb8ca7db Jean*0148 #endif
f913c5a485 Mart*0149           DO k=1,nITD
f5282c5b03 Gael*0150             TICES(i,j,k,bi,bj)=0. _d 0
                0151           ENDDO
1cfef3b6ad Patr*0152           saltWtrIce(i,j,bi,bj) = 0. _d 0
                0153           frWtrIce(i,j,bi,bj)   = 0. _d 0
7109a141b2 Patr*0154          ENDDO
                0155         ENDDO
                0156        ENDDO
                0157       ENDDO
                0158
c69ccc91fb antn*0159 #ifdef ALLOW_AUTODIFF
                0160 C- Note: To simplify dependency & avoid recomputations, when compiling
                0161 C        pkg/autodiff, we always call SEAICE_INIT_VARIA to initialise control
                0162 C        variables (as done above) without condition on useSEAICE.
                0163 C        Therefore, in this case, the "If useSEAICE" is added back here:
                0164       IF ( useSEAICE ) THEN
                0165 #endif
                0166
84cb676f82 Mart*0167 C--   Initialize (variable) grid info. As long as we allow masking of
                0168 C--   velocities outside of ice covered areas (in seaice_dynsolver)
b7b61e618a Mart*0169 C--   we need to re-initialize seaiceMaskU/V here for TAF
809c36b928 Patr*0170       DO bj=myByLo(myThid),myByHi(myThid)
                0171        DO bi=myBxLo(myThid),myBxHi(myThid)
45315406aa Mart*0172 #ifdef SEAICE_CGRID
f4cb8ca7db Jean*0173         DO j=1-OLy+1,sNy+OLy
                0174          DO i=1-OLx+1,sNx+OLx
                0175           seaiceMaskU(i,j,bi,bj)=   0.0 _d 0
                0176           seaiceMaskV(i,j,bi,bj)=   0.0 _d 0
                0177           mask_uice=HEFFM(i,j,bi,bj)+HEFFM(i-1,j  ,bi,bj)
a808b72b4b Patr*0178           IF(mask_uice.GT.1.5 _d 0) seaiceMaskU(i,j,bi,bj)=1.0 _d 0
f4cb8ca7db Jean*0179           mask_uice=HEFFM(i,j,bi,bj)+HEFFM(i  ,j-1,bi,bj)
a808b72b4b Patr*0180           IF(mask_uice.GT.1.5 _d 0) seaiceMaskV(i,j,bi,bj)=1.0 _d 0
ab37184352 Mart*0181          ENDDO
                0182         ENDDO
45315406aa Mart*0183 # ifdef OBCS_UVICE_OLD
eb0337bb06 Dimi*0184         IF (useOBCS) THEN
b48f786726 Dimi*0185 C--   If OBCS is turned on, close southern and western boundaries
48474b8e08 Jean*0186          DO i=1-OLx,sNx+OLx
a6d2313c04 Dimi*0187 C Southern boundary
dc8f395c29 Dimi*0188           J_obc = OB_Js(i,bi,bj)
d6145066cc Jean*0189           IF ( J_obc.NE.OB_indexNone ) THEN
dc8f395c29 Dimi*0190            seaiceMaskU(i,J_obc,bi,bj)=   0.0 _d 0
                0191            seaiceMaskV(i,J_obc,bi,bj)=   0.0 _d 0
eb0337bb06 Dimi*0192           ENDIF
                0193          ENDDO
48474b8e08 Jean*0194          DO j=1-OLy,sNy+OLy
a6d2313c04 Dimi*0195 C Western boundary
d6145066cc Jean*0196           I_obc = OB_Iw(j,bi,bj)
                0197           IF ( I_obc.NE.OB_indexNone ) THEN
dc8f395c29 Dimi*0198            seaiceMaskU(I_obc,j,bi,bj)=   0.0 _d 0
                0199            seaiceMaskV(I_obc,j,bi,bj)=   0.0 _d 0
b48f786726 Dimi*0200           ENDIF
809c36b928 Patr*0201          ENDDO
b48f786726 Dimi*0202         ENDIF
45315406aa Mart*0203 # endif /* OBCS_UVICE_OLD */
a6d2313c04 Dimi*0204 #endif /* SEAICE_CGRID */
69c991a1ec Dimi*0205
809c36b928 Patr*0206         DO j=1-OLy,sNy+OLy
                0207          DO i=1-OLx,sNx+OLx
f913c5a485 Mart*0208           DO k=1,nITD
f4cb8ca7db Jean*0209            TICES(i,j,k,bi,bj)=273.0 _d 0
09510da3bb Dimi*0210           ENDDO
45315406aa Mart*0211 #ifdef SEAICE_CGRID
f4cb8ca7db Jean*0212           seaiceMassC(i,j,bi,bj)=1000.0 _d 0
                0213           seaiceMassU(i,j,bi,bj)=1000.0 _d 0
                0214           seaiceMassV(i,j,bi,bj)=1000.0 _d 0
45315406aa Mart*0215 #endif
                0216 #ifdef SEAICE_BGRID_DYNAMICS
                0217           AMASS      (i,j,bi,bj)=1000.0 _d 0
53092bcb42 Mart*0218 #endif
809c36b928 Patr*0219          ENDDO
                0220         ENDDO
d5bc48d0a9 Dimi*0221
809c36b928 Patr*0222        ENDDO
                0223       ENDDO
                0224
                0225 C--   Update overlap regions
3ed8349f04 Mart*0226 #ifdef SEAICE_CGRID
b4f21203f7 Jean*0227       CALL EXCH_UV_XY_RL(seaiceMaskU,seaiceMaskV,.FALSE.,myThid)
45315406aa Mart*0228 #endif
                0229 #ifdef SEAICE_BGRID_DYNAMICS
48474b8e08 Jean*0230       _EXCH_XY_RS(UVM, myThid)
3ed8349f04 Mart*0231 #endif
809c36b928 Patr*0232
a9a8db28c6 Alis*0233 C--   Now lets look at all these beasts
a24915ab1a Jean*0234       IF ( plotLevel.GE.debLevC ) THEN
3e0af136db Dimi*0235          CALL PLOT_FIELD_XYRL( HEFFM   , 'Current HEFFM   ' ,
23142459d0 Jean*0236      &        nIter0, myThid )
3ed8349f04 Mart*0237 #ifdef SEAICE_CGRID
                0238          CALL PLOT_FIELD_XYRL( seaiceMaskU, 'Current seaiceMaskU',
23142459d0 Jean*0239      &        nIter0, myThid )
3ed8349f04 Mart*0240          CALL PLOT_FIELD_XYRL( seaiceMaskV, 'Current seaiceMaskV',
23142459d0 Jean*0241      &        nIter0, myThid )
45315406aa Mart*0242 #endif
                0243 #ifdef SEAICE_BGRID_DYNAMICS
48474b8e08 Jean*0244          CALL PLOT_FIELD_XYRS( UVM     , 'Current UVM     ' ,
23142459d0 Jean*0245      &        nIter0, myThid )
3ed8349f04 Mart*0246 #endif
3e0af136db Dimi*0247       ENDIF
09510da3bb Dimi*0248
809c36b928 Patr*0249 C--   Set model variables to initial/restart conditions
2a682dc462 Mart*0250       IF ( .NOT. ( startTime .EQ. baseTime .AND.  nIter0 .EQ. 0
                0251      &     .AND. pickupSuff .EQ. ' ') ) THEN
36f6faaf99 Dimi*0252
6060ec2938 Dimi*0253          CALL SEAICE_READ_PICKUP ( myThid )
36f6faaf99 Dimi*0254
6060ec2938 Dimi*0255       ELSE
36f6faaf99 Dimi*0256
e7c33124ce Dimi*0257        DO bj=myByLo(myThid),myByHi(myThid)
                0258         DO bi=myBxLo(myThid),myBxHi(myThid)
                0259          DO j=1-OLy,sNy+OLy
                0260           DO i=1-OLx,sNx+OLx
772590b63c Mart*0261            HEFF(i,j,bi,bj)=SEAICE_initialHEFF*HEFFM(i,j,bi,bj)
                0262            UICE(i,j,bi,bj)=ZERO
                0263            VICE(i,j,bi,bj)=ZERO
e7c33124ce Dimi*0264           ENDDO
809c36b928 Patr*0265          ENDDO
e7c33124ce Dimi*0266         ENDDO
                0267        ENDDO
6060ec2938 Dimi*0268
425e8efc36 Jean*0269 C--   Read initial sea-ice velocity from file (if available)
                0270        IF ( uIceFile .NE. ' ' )
                0271      &  CALL READ_FLD_XY_RL( uIceFile, ' ', uIce, 0, myThid )
                0272        IF ( vIceFile .NE. ' ' )
                0273      &  CALL READ_FLD_XY_RL( vIceFile, ' ', vIce, 0, myThid )
                0274        IF ( uIceFile .NE. ' ' .OR. vIceFile .NE. ' ' ) THEN
                0275 #ifdef SEAICE_CGRID
                0276         DO bj=myByLo(myThid),myByHi(myThid)
                0277          DO bi=myBxLo(myThid),myBxHi(myThid)
                0278           DO j=1-OLy,sNy+OLy
                0279            DO i=1-OLx,sNx+OLx
                0280             uIce(i,j,bi,bj) = uIce(i,j,bi,bj)*seaiceMaskU(i,j,bi,bj)
                0281             vIce(i,j,bi,bj) = vIce(i,j,bi,bj)*seaiceMaskV(i,j,bi,bj)
                0282            ENDDO
                0283           ENDDO
                0284          ENDDO
                0285         ENDDO
                0286 #endif /* SEAICE_CGRID */
                0287         CALL EXCH_UV_XY_RL( uIce, vIce, .TRUE., myThid )
                0288        ENDIF
                0289
6060ec2938 Dimi*0290 C--   Read initial sea-ice thickness from file if available.
e7c33124ce Dimi*0291        IF ( HeffFile .NE. ' ' ) THEN
772590b63c Mart*0292         CALL READ_FLD_XY_RL( HeffFile, ' ', HEFF, 0, myThid )
                0293         _EXCH_XY_RL(HEFF,myThid)
e7c33124ce Dimi*0294         DO bj=myByLo(myThid),myByHi(myThid)
                0295          DO bi=myBxLo(myThid),myBxHi(myThid)
                0296           DO j=1-OLy,sNy+OLy
                0297            DO i=1-OLx,sNx+OLx
772590b63c Mart*0298             HEFF(i,j,bi,bj) = MAX(HEFF(i,j,bi,bj),ZERO)
e7c33124ce Dimi*0299            ENDDO
                0300           ENDDO
                0301          ENDDO
                0302         ENDDO
                0303        ENDIF
                0304
                0305        DO bj=myByLo(myThid),myByHi(myThid)
                0306         DO bi=myBxLo(myThid),myBxHi(myThid)
                0307          DO j=1-OLy,sNy+OLy
                0308           DO i=1-OLx,sNx+OLx
772590b63c Mart*0309            IF(HEFF(i,j,bi,bj).GT.ZERO) AREA(i,j,bi,bj)=ONE
e7c33124ce Dimi*0310           ENDDO
                0311          ENDDO
                0312         ENDDO
                0313        ENDDO
f4cb8ca7db Jean*0314
ff73878394 Dimi*0315 C--   Read initial sea-ice area from file if available.
09066b09cb Mart*0316        IF ( AreaFile .NE. ' ' ) THEN
772590b63c Mart*0317         CALL READ_FLD_XY_RL( AreaFile, ' ', AREA, 0, myThid )
                0318         _EXCH_XY_RL(AREA,myThid)
09066b09cb Mart*0319         DO bj=myByLo(myThid),myByHi(myThid)
                0320          DO bi=myBxLo(myThid),myBxHi(myThid)
                0321           DO j=1-OLy,sNy+OLy
                0322            DO i=1-OLx,sNx+OLx
772590b63c Mart*0323             AREA(i,j,bi,bj) = MAX(AREA(i,j,bi,bj),ZERO)
                0324             AREA(i,j,bi,bj) = MIN(AREA(i,j,bi,bj),ONE)
                0325             IF ( AREA(i,j,bi,bj) .LE. ZERO ) HEFF(i,j,bi,bj) = ZERO
                0326             IF ( HEFF(i,j,bi,bj) .LE. ZERO ) AREA(i,j,bi,bj) = ZERO
09066b09cb Mart*0327            ENDDO
                0328           ENDDO
                0329          ENDDO
                0330         ENDDO
                0331        ENDIF
                0332
                0333        DO bj=myByLo(myThid),myByHi(myThid)
                0334         DO bi=myBxLo(myThid),myBxHi(myThid)
                0335          DO j=1-OLy,sNy+OLy
                0336           DO i=1-OLx,sNx+OLx
e5b8ebbabe Jean*0337            HSNOW(i,j,bi,bj) = 0.2 _d 0 * AREA(i,j,bi,bj)
09066b09cb Mart*0338           ENDDO
                0339          ENDDO
                0340         ENDDO
                0341        ENDDO
de31ea8481 Dimi*0342
                0343 C--   Read initial snow thickness from file if available.
                0344        IF ( HsnowFile .NE. ' ' ) THEN
772590b63c Mart*0345         CALL READ_FLD_XY_RL( HsnowFile, ' ', HSNOW, 0, myThid )
                0346         _EXCH_XY_RL(HSNOW,myThid)
de31ea8481 Dimi*0347         DO bj=myByLo(myThid),myByHi(myThid)
                0348          DO bi=myBxLo(myThid),myBxHi(myThid)
                0349           DO j=1-OLy,sNy+OLy
                0350            DO i=1-OLx,sNx+OLx
772590b63c Mart*0351             HSNOW(i,j,bi,bj) = MAX(HSNOW(i,j,bi,bj),ZERO)
de31ea8481 Dimi*0352            ENDDO
                0353           ENDDO
                0354          ENDDO
                0355         ENDDO
                0356        ENDIF
fdfa8e151f Dimi*0357
86b84a92fc Patr*0358 #ifdef SEAICE_ITD
                0359        DO bj=myByLo(myThid),myByHi(myThid)
                0360         DO bi=myBxLo(myThid),myBxHi(myThid)
                0361          DO j=1-OLy,sNy+OLy
                0362           DO i=1-OLx,sNx+OLx
5fe78992ba Mart*0363            AREAITD(i,j,1,bi,bj)   = AREA(i,j,bi,bj)
                0364            HEFFITD(i,j,1,bi,bj)   = HEFF(i,j,bi,bj)
                0365            HSNOWITD(i,j,1,bi,bj)  = HSNOW(i,j,bi,bj)
                0366            opnWtrFrac(i,j,bi,bj)  = 1. _d 0 - AREA(i,j,bi,bj)
                0367            fw2ObyRidge(i,j,bi,bj) = 0. _d 0
86b84a92fc Patr*0368           ENDDO
                0369          ENDDO
fd2abd71b8 Mart*0370          CALL SEAICE_ITD_REDIST(bi, bj, baseTime, nIter0, myThid)
86b84a92fc Patr*0371         ENDDO
                0372        ENDDO
                0373 #endif
                0374
a98c4b8072 Ian *0375 #ifdef SEAICE_VARIABLE_SALINITY
3183247f03 Dimi*0376        DO bj=myByLo(myThid),myByHi(myThid)
                0377         DO bi=myBxLo(myThid),myBxHi(myThid)
f4cb8ca7db Jean*0378          DO j=1-OLy,sNy+OLy
                0379           DO i=1-OLx,sNx+OLx
0320e25227 Mart*0380            HSALT(i,j,bi,bj)=HEFF(i,j,bi,bj)*salt(i,j,kSrf,bi,bj)*
4eb4a54cba Jean*0381      &            SEAICE_rhoIce*SEAICE_saltFrac
                0382 cif     &            ICE2WATR*rhoConstFresh*SEAICE_saltFrac
a98c4b8072 Ian *0383
3183247f03 Dimi*0384           ENDDO
                0385          ENDDO
                0386         ENDDO
                0387        ENDDO
                0388
fdfa8e151f Dimi*0389 C--   Read initial sea ice salinity from file if available.
                0390        IF ( HsaltFile .NE. ' ' ) THEN
772590b63c Mart*0391         CALL READ_FLD_XY_RL( HsaltFile, ' ', HSALT, 0, myThid )
                0392         _EXCH_XY_RL(HSALT,myThid)
fdfa8e151f Dimi*0393        ENDIF
a98c4b8072 Ian *0394 #endif /* SEAICE_VARIABLE_SALINITY */
f4cb8ca7db Jean*0395
e54fe3e1f9 Gael*0396 #ifdef ALLOW_SITRACER
f681b7f5d4 Dimi*0397 C--   Read initial sea ice age from file if available.
e54fe3e1f9 Gael*0398        DO iTr = 1, SItrMaxNum
                0399         IF ( SItrFile(iTr) .NE. ' ' ) THEN
                0400         CALL READ_FLD_XY_RL( siTrFile(iTr), ' ',
                0401      &   SItracer(1-OLx,1-OLy,1,1,iTr), 0, myThid )
                0402         _EXCH_XY_RL(SItracer(1-OLx,1-OLy,1,1,iTr),myThid)
ccaa3c61f4 Patr*0403         ENDIF
                0404        ENDDO
e54fe3e1f9 Gael*0405 #endif /* ALLOW_SITRACER */
f681b7f5d4 Dimi*0406
809c36b928 Patr*0407       ENDIF
                0408
4f6e464f38 Jean*0409 #ifdef ALLOW_OBCS
3352db77f4 Jean*0410 C--   In case we use scheme with a large stencil that extends into overlap:
                0411 C     no longer needed with the right masking in advection & diffusion S/R.
                0412 c     IF ( useOBCS ) THEN
                0413 c       DO bj=myByLo(myThid),myByHi(myThid)
                0414 c        DO bi=myBxLo(myThid),myBxHi(myThid)
48474b8e08 Jean*0415 c          CALL OBCS_COPY_TRACER( HEFF(1-OLx,1-OLy,bi,bj),
3352db77f4 Jean*0416 c    I                            1, bi, bj, myThid )
48474b8e08 Jean*0417 c          CALL OBCS_COPY_TRACER( AREA(1-OLx,1-OLy,bi,bj),
3352db77f4 Jean*0418 c    I                            1, bi, bj, myThid )
48474b8e08 Jean*0419 c          CALL OBCS_COPY_TRACER( HSNOW(1-OLx,1-OLy,bi,bj),
3352db77f4 Jean*0420 c    I                            1, bi, bj, myThid )
a98c4b8072 Ian *0421 #ifdef SEAICE_VARIABLE_SALINITY
48474b8e08 Jean*0422 c          CALL OBCS_COPY_TRACER( HSALT(1-OLx,1-OLy,bi,bj),
3352db77f4 Jean*0423 c    I                            1, bi, bj, myThid )
4f6e464f38 Jean*0424 #endif
3352db77f4 Jean*0425 c        ENDDO
                0426 c       ENDDO
                0427 c     ENDIF
4f6e464f38 Jean*0428 #endif /* ALLOW_OBCS */
                0429
809c36b928 Patr*0430 C---  Complete initialization
                0431       DO bj=myByLo(myThid),myByHi(myThid)
                0432        DO bi=myBxLo(myThid),myBxHi(myThid)
45315406aa Mart*0433 #if ( defined SEAICE_CGRID || defined SEAICE_BGRID_DYNAMICS )
809c36b928 Patr*0434         DO j=1-OLy,sNy+OLy
                0435          DO i=1-OLx,sNx+OLx
8e32c48b8f Mart*0436           ZETA(i,j,bi,bj)        = HEFF(i,j,bi,bj)*(1.0 _d 11)
                0437           ETA(i,j,bi,bj)         = ZETA(i,j,bi,bj)/SEAICE_eccen**2
                0438           PRESS0(i,j,bi,bj)      = SEAICE_strength*HEFF(i,j,bi,bj)
ba6cfc5714 Mart*0439      &         *EXP(-SEAICE_cStar*(ONE-AREA(i,j,bi,bj)))
8e32c48b8f Mart*0440           SEAICE_zMax(i,j,bi,bj) = SEAICE_zetaMaxFac*PRESS0(i,j,bi,bj)
                0441           SEAICE_zMin(i,j,bi,bj) = SEAICE_zetaMin
                0442           PRESS0(i,j,bi,bj)      = PRESS0(i,j,bi,bj)*HEFFM(i,j,bi,bj)
809c36b928 Patr*0443          ENDDO
                0444         ENDDO
45315406aa Mart*0445 #endif
f3ce416a61 Jean*0446         IF ( useRealFreshWaterFlux .AND. .NOT.useThSIce ) THEN
                0447          DO j=1-OLy,sNy+OLy
                0448           DO i=1-OLx,sNx+OLx
772590b63c Mart*0449            sIceLoad(i,j,bi,bj) = HEFF(i,j,bi,bj)*SEAICE_rhoIce
13533b2b42 Mart*0450      &                         + HSNOW(i,j,bi,bj)*SEAICE_rhoSnow
f4cb8ca7db Jean*0451
f3ce416a61 Jean*0452           ENDDO
                0453          ENDDO
                0454         ENDIF
809c36b928 Patr*0455        ENDDO
                0456       ENDDO
2f5e8addfd Mart*0457 #ifdef SEAICE_CGRID
                0458 C     compute tensile strength factor k: tensileStrength = k*PRESS
0dcda34642 Jean*0459 C     can be done in initialisation phase as long as it depends only
2f5e8addfd Mart*0460 C     on depth
                0461       IF ( SEAICE_tensilFac .NE. 0. _d 0 ) THEN
                0462        recip_tensilDepth = 0. _d 0
0dcda34642 Jean*0463        IF ( SEAICE_tensilDepth .GT. 0. _d 0 )
2f5e8addfd Mart*0464      &      recip_tensilDepth = 1. _d 0 / SEAICE_tensilDepth
                0465        DO bj=myByLo(myThid),myByHi(myThid)
                0466         DO bi=myBxLo(myThid),myBxHi(myThid)
                0467          DO j=1-OLy,sNy+OLy
                0468           DO i=1-OLx,sNx+OLx
                0469            tensileStrFac(i,j,bi,bj) = SEAICE_tensilFac*HEFFM(i,j,bi,bj)
                0470      &          *exp(-ABS(R_low(i,j,bi,bj))*recip_tensilDepth)
                0471           ENDDO
                0472          ENDDO
                0473         ENDDO
                0474        ENDDO
                0475       ENDIF
45315406aa Mart*0476 # if ( defined SEAICE_ALLOW_JFNK || defined SEAICE_ALLOW_KRYLOV )
                0477 C     Computing this metric cannot be done in S/R SEAICE_INIT_FIXED
                0478 C     where it belongs, because globalArea is only defined later after
                0479 C     S/R PACKAGES_INIT_FIXED, so we move this computation here.
                0480       CALL SEAICE_MAP_RS2VEC( nVec, rAw, rAs,
                0481      &     scalarProductMetric, .TRUE., myThid )
                0482       DO bj=myByLo(myThid),myByHi(myThid)
                0483        DO bi=myBxLo(myThid),myBxHi(myThid)
                0484         DO i=1,nVec
                0485          scalarProductMetric(i,1,bi,bj) =
                0486      &        scalarProductMetric(i,1,bi,bj)/globalArea
                0487         ENDDO
                0488        ENDDO
                0489       ENDDO
                0490 # endif /* SEAICE_ALLOW_JFNK or KRYLOV */
                0491
2f5e8addfd Mart*0492 #endif /* SEAICE_CGRID */
809c36b928 Patr*0493
c69ccc91fb antn*0494 #ifdef ALLOW_AUTODIFF
                0495 C-    end if useSEAICE block
                0496       ENDIF
                0497 #endif
                0498
809c36b928 Patr*0499       RETURN
                0500       END