File indexing completed on 2022-04-25 05:09:44 UTC
view on githubraw file Latest commit 2e7aec99 on 2022-04-25 03:54:25 UTC
0e09621e3e Patr*0001 #include "PACKAGES_CONFIG.h"
0002 #include "CPP_OPTIONS.h"
52d51b46b8 Jean*0003 #ifdef ALLOW_AUTODIFF
0004 # include "AUTODIFF_OPTIONS.h"
0005 #endif
987a76ae74 Jean*0006 #ifdef ALLOW_SALT_PLUME
0007 # include "SALT_PLUME_OPTIONS.h"
0008 #endif
0e09621e3e Patr*0009 #undef CHECK_OVERLAP_FORCING
0010
0011
0012
0013
0014 SUBROUTINE EXTERNAL_FORCING_SURF(
0015 I iMin, iMax, jMin, jMax,
0016 I myTime, myIter, myThid )
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026 IMPLICIT NONE
0027
0028 #include "SIZE.h"
0029 #include "EEPARAMS.h"
0030 #include "PARAMS.h"
0031 #include "FFIELDS.h"
0032 #include "DYNVARS.h"
0033 #include "GRID.h"
0034 #include "SURFACE.h"
7c50f07931 Mart*0035 #ifdef ALLOW_AUTODIFF_TAMC
0e09621e3e Patr*0036 # include "tamc.h"
0037 #endif
0038
0039
0040
0041
0042
0043
0044
0045 INTEGER iMin, iMax
0046 INTEGER jMin, jMax
0047 _RL myTime
0048 INTEGER myIter
0049 INTEGER myThid
0050
0051
0052
0053
0054
0055
0056 INTEGER bi,bj
0057 INTEGER i,j
0058 INTEGER ks
2e7aec9951 dngo*0059 #ifdef ALLOW_AUTODIFF_TAMC
0060 INTEGER act1, act2, act3, act4
0061 INTEGER max1, max2, max3
0062 INTEGER ikey
0063 #endif
52d51b46b8 Jean*0064 _RL recip_Cp
0e09621e3e Patr*0065 #ifdef ALLOW_BALANCE_FLUXES
0066 _RS tmpVar(1)
0067 #endif
0068 #ifdef CHECK_OVERLAP_FORCING
0069 _RS fixVal
0070 #endif
0071
0072
0073 IF ( usingPCoords ) THEN
0074 ks = Nr
0075 ELSE
0076 ks = 1
0077 ENDIF
52d51b46b8 Jean*0078 recip_Cp = 1. _d 0 / HeatCapacity_Cp
0e09621e3e Patr*0079
0080
0081
0082
0083
0084
0085 #ifdef ALLOW_BALANCE_FLUXES
0086
2e7aec9951 dngo*0087 # ifdef ALLOW_AUTODIFF
0e09621e3e Patr*0088 tmpVar(1) = oneRS
2e7aec9951 dngo*0089 # endif
0090 IF ( selectBalanceEmPmR.GE.1 .AND.
0091 & (.NOT.useSeaice .OR. useThSIce) ) THEN
0092 IF ( selectBalanceEmPmR.EQ.1 ) THEN
0093 tmpVar(1) = oneRS
0094 CALL REMOVE_MEAN_RS( 1, EmPmR, maskInC, maskInC, rA,
0095 & tmpVar, 'EmPmR', myTime, myIter, myThid )
0096 ELSEIF ( selectBalanceEmPmR.EQ.2 ) THEN
0097 tmpVar(1) = -oneRS
0098 CALL REMOVE_MEAN_RS( 1, EmPmR, weight2BalanceFlx, maskInC, rA,
0099 & tmpVar, 'EmPmR', myTime, myIter, myThid )
0100 ENDIF
0e09621e3e Patr*0101 ENDIF
0102 IF ( balanceQnet .AND. (.NOT.useSeaice .OR. useThSIce) ) THEN
2e7aec9951 dngo*0103 tmpVar(1) = oneRS
0104 CALL REMOVE_MEAN_RS( 1, Qnet, maskInC, maskInC, rA,
0105 & tmpVar, 'Qnet ', myTime, myIter, myThid )
0e09621e3e Patr*0106 ENDIF
0107 #endif /* ALLOW_BALANCE_FLUXES */
0108
0109
0110
0111 #ifdef CHECK_OVERLAP_FORCING
0112
0113
0114 fixVal = 1.
0115 CALL RESET_HALO_RS ( EmPmR, fixVal, 1, myThid )
0116 fixVal = 400.
0117 CALL RESET_HALO_RS ( Qnet, fixVal, 1, myThid )
0118 fixVal = -200.
0119 CALL RESET_HALO_RS ( Qsw, fixVal, 1, myThid )
0120 fixVal = 40.
0121 CALL RESET_HALO_RS ( saltFlux, fixVal, 1, myThid )
0122
0123 #endif /* CHECK_OVERLAP_FORCING */
0124
0125
0126
0127 #ifdef EXACT_CONSERV
0128
0129
0130 # ifdef ALLOW_AUTODIFF_TAMC
0131
52d51b46b8 Jean*0132
0e09621e3e Patr*0133 # endif
0134 DO bj=myByLo(myThid),myByHi(myThid)
0135 DO bi=myBxLo(myThid),myBxHi(myThid)
0136 IF ( staggerTimeStep ) THEN
0137 DO j=1-OLy,sNy+OLy
0138 DO i=1-OLx,sNx+OLx
0139 PmEpR(i,j,bi,bj) = -EmPmR(i,j,bi,bj)
0140 ENDDO
0141 ENDDO
0142 ENDIF
0143 ENDDO
0144 ENDDO
0145 #endif /* EXACT_CONSERV */
0146
0147
0148 DO bj=myByLo(myThid),myByHi(myThid)
0149 DO bi=myBxLo(myThid),myBxHi(myThid)
0150 DO j=1-OLy,sNy+OLy
0151 DO i=1-OLx,sNx+OLx
0152 surfaceForcingT(i,j,bi,bj) = 0. _d 0
0153 surfaceForcingS(i,j,bi,bj) = 0. _d 0
0154 ENDDO
0155 ENDDO
0156 ENDDO
0157 ENDDO
0158
0159
0160
0161
0162 IF ( doThetaClimRelax .OR. doSaltClimRelax ) THEN
0163 CALL FORCING_SURF_RELAX(
0164 I iMin, iMax, jMin, jMax,
0165 I myTime, myIter, myThid )
0166 ENDIF
0167
0168 #ifdef ALLOW_PTRACERS
0169
0170 #ifdef ALLOW_AUTODIFF_TAMC
0171
0172
2e7aec9951 dngo*0173 #ifdef ALLOW_BLING
0174
0175 #endif
0e09621e3e Patr*0176 #endif
0177 IF ( usePTRACERS ) THEN
0178 DO bj=myByLo(myThid),myByHi(myThid)
0179 DO bi=myBxLo(myThid),myBxHi(myThid)
0180 CALL PTRACERS_FORCING_SURF(
0181 I surfaceForcingS(1-OLx,1-OLy,bi,bj),
0182 I bi, bj, iMin, iMax, jMin, jMax,
52d51b46b8 Jean*0183 I myTime, myIter, myThid )
0e09621e3e Patr*0184 ENDDO
0185 ENDDO
0186 ENDIF
0187 #endif /* ALLOW_PTRACERS */
0188
0189
0190
0191
0192
0193
0194 DO bj=myByLo(myThid),myByHi(myThid)
0195 DO bi=myBxLo(myThid),myBxHi(myThid)
0196
52d51b46b8 Jean*0197 #ifdef ALLOW_AUTODIFF_TAMC
7c50f07931 Mart*0198 act1 = bi - myBxLo(myThid)
0199 max1 = myBxHi(myThid) - myBxLo(myThid) + 1
0200 act2 = bj - myByLo(myThid)
0201 max2 = myByHi(myThid) - myByLo(myThid) + 1
0202 act3 = myThid - 1
0203 max3 = nTx*nTy
0204 act4 = ikey_dynamics - 1
0205 ikey = (act1 + 1) + act2*max1
0206 & + act3*max1*max2
0207 & + act4*max1*max2*max3
52d51b46b8 Jean*0208 #endif /* ALLOW_AUTODIFF_TAMC */
0209
0210 #ifdef ALLOW_AUTODIFF_TAMC
0211
0212 #ifdef EXACT_CONSERV
0213
0214 #endif
0215 #endif /* ALLOW_AUTODIFF_TAMC */
0216
0e09621e3e Patr*0217
0218 DO j = jMin, jMax
0219 DO i = iMin, iMax
0220
0221
0222 surfaceForcingU(i,j,bi,bj) = fu(i,j,bi,bj)*mass2rUnit
0223
0224 surfaceForcingV(i,j,bi,bj) = fv(i,j,bi,bj)*mass2rUnit
0225
0226 surfaceForcingT(i,j,bi,bj) = surfaceForcingT(i,j,bi,bj)
0227 & - ( Qnet(i,j,bi,bj)
0228 #ifdef SHORTWAVE_HEATING
0229 & -Qsw(i,j,bi,bj)
0230 #endif
2e7aec9951 dngo*0231
0e09621e3e Patr*0232 & +Qnetm(i,j,bi,bj)
0233 & ) *recip_Cp*mass2rUnit
0234
0235 surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)
0236 & -saltFlux(i,j,bi,bj)*mass2rUnit
0237
0238 ENDDO
0239 ENDDO
0240
0241 #ifdef ALLOW_SALT_PLUME
0242
0243
0244
987a76ae74 Jean*0245
0246
0247
0248 #ifndef SALT_PLUME_VOLUME
0e09621e3e Patr*0249 IF ( useSALT_PLUME ) THEN
0250 CALL SALT_PLUME_FORCING_SURF(
0251 I bi, bj, iMin, iMax, jMin, jMax,
52d51b46b8 Jean*0252 I myTime, myIter, myThid )
0e09621e3e Patr*0253 ENDIF
987a76ae74 Jean*0254 #endif /* SALT_PLUME_VOLUME */
0e09621e3e Patr*0255 #endif /* ALLOW_SALT_PLUME */
0256
0257
0258
0259
0260
0261
0262
0263 #ifdef EXACT_CONSERV
0264 IF ( (nonlinFreeSurf.GT.0 .OR. usingPCoords)
0265 & .AND. useRealFreshWaterFlux ) THEN
0266
0267
0268
0269
0270
0271 IF (temp_EvPrRn.NE.UNSET_RL) THEN
0272 DO j = jMin, jMax
0273 DO i = iMin, iMax
0274 surfaceForcingT(i,j,bi,bj) = surfaceForcingT(i,j,bi,bj)
0275 & + PmEpR(i,j,bi,bj)
0276 & *( temp_EvPrRn - theta(i,j,ks,bi,bj) )
0277 & *mass2rUnit
0278 ENDDO
0279 ENDDO
0280 ENDIF
0281
0282 IF (salt_EvPrRn.NE.UNSET_RL) THEN
0283 DO j = jMin, jMax
0284 DO i = iMin, iMax
0285 surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)
0286 & + PmEpR(i,j,bi,bj)
0287 & *( salt_EvPrRn - salt(i,j,ks,bi,bj) )
0288 & *mass2rUnit
0289 ENDDO
0290 ENDDO
0291 ENDIF
0292
0293
0294 ELSE
0295 #else /* EXACT_CONSERV */
0296 IF (.TRUE.) THEN
0297 #endif /* EXACT_CONSERV */
0298
0299
0300
0301
0302 IF (convertFW2Salt .EQ. -1.) THEN
0303
0304
0305 IF (temp_EvPrRn.NE.UNSET_RL) THEN
0306
0307 DO j = jMin, jMax
0308 DO i = iMin, iMax
0309 surfaceForcingT(i,j,bi,bj) = surfaceForcingT(i,j,bi,bj)
0310 & + EmPmR(i,j,bi,bj)
0311 & *( theta(i,j,ks,bi,bj) - temp_EvPrRn )
0312 & *mass2rUnit
0313 ENDDO
0314 ENDDO
0315 ENDIF
0316 IF (salt_EvPrRn.NE.UNSET_RL) THEN
0317
0318 DO j = jMin, jMax
0319 DO i = iMin, iMax
0320 surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)
0321 & + EmPmR(i,j,bi,bj)
0322 & *( salt(i,j,ks,bi,bj) - salt_EvPrRn )
0323 & *mass2rUnit
0324 ENDDO
0325 ENDDO
0326 ENDIF
0327
0328 ELSE
0329
0330
0331 IF (temp_EvPrRn.NE.UNSET_RL) THEN
0332
0333 DO j = jMin, jMax
0334 DO i = iMin, iMax
0335 surfaceForcingT(i,j,bi,bj) = surfaceForcingT(i,j,bi,bj)
0336 & + EmPmR(i,j,bi,bj)
0337 & *( tRef(ks) - temp_EvPrRn )
0338 & *mass2rUnit
0339 ENDDO
0340 ENDDO
0341 ENDIF
0342 IF (salt_EvPrRn.NE.UNSET_RL) THEN
0343
0344 DO j = jMin, jMax
0345 DO i = iMin, iMax
0346 surfaceForcingS(i,j,bi,bj) = surfaceForcingS(i,j,bi,bj)
0347 & + EmPmR(i,j,bi,bj)
0348 & *( convertFW2Salt - salt_EvPrRn )
0349 & *mass2rUnit
0350 ENDDO
0351 ENDDO
0352 ENDIF
0353
0354
0355 ENDIF
0356
0357 ENDIF
0358
0359
0360
0361 #ifdef ATMOSPHERIC_LOADING
0362
0363
0364
0365
0366
0367
0368
0369
0370 IF ( usingZCoords ) THEN
0371
0372
0373 IF ( useRealFreshWaterFlux ) THEN
0374 DO j = jMin, jMax
0375 DO i = iMin, iMax
0376 phi0surf(i,j,bi,bj) = ( pLoad(i,j,bi,bj)
2e7aec9951 dngo*0377 & +sIceLoad(i,j,bi,bj)*gravity*sIceLoadFac
0e09621e3e Patr*0378 & )*recip_rhoConst
0379 ENDDO
0380 ENDDO
0381 ELSE
0382 DO j = jMin, jMax
0383 DO i = iMin, iMax
0384 phi0surf(i,j,bi,bj) = pLoad(i,j,bi,bj)*recip_rhoConst
0385 ENDDO
0386 ENDDO
0387 ENDIF
0388
0389
0390
0391
0392
0393
0394
0395
0396
0397 ENDIF
0398 #endif /* ATMOSPHERIC_LOADING */
0399
0400 #ifdef ALLOW_SHELFICE
52d51b46b8 Jean*0401 IF ( useSHELFICE) THEN
0402 CALL SHELFICE_FORCING_SURF(
0403 I bi, bj, iMin, iMax, jMin, jMax,
0404 I myTime, myIter, myThid )
0e09621e3e Patr*0405 ENDIF
0406 #endif /* ALLOW_SHELFICE */
0407
0408
0409 ENDDO
0410 ENDDO
0411
0412 RETURN
0413 END