- Timestamp:
- 2019-12-09T13:55:34+01:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_AGRIF-01-05_merged/src/NST/agrif_oce_interp.F90
r10068 r12123 33 33 USE agrif_oce_sponge 34 34 USE lib_mpp 35 USE vremap 35 36 36 37 IMPLICIT NONE 37 38 PRIVATE 38 39 39 PUBLIC Agrif_dyn, Agrif_ssh, Agrif_dyn_ts, Agrif_ ssh_ts, Agrif_dta_ts40 PUBLIC Agrif_dyn, Agrif_ssh, Agrif_dyn_ts, Agrif_dyn_ts_flux, Agrif_ssh_ts, Agrif_dta_ts 40 41 PUBLIC Agrif_tra, Agrif_avm 41 42 PUBLIC interpun , interpvn 42 43 PUBLIC interptsn, interpsshn, interpavm 43 44 PUBLIC interpunb, interpvnb , interpub2b, interpvb2b 44 PUBLIC interpe3t, interpumsk, interpvmsk 45 45 PUBLIC interpe3t 46 #if defined key_vertical 47 PUBLIC interpht0, interpmbkt 48 # endif 46 49 INTEGER :: bdy_tinterp = 0 47 50 … … 78 81 ! 79 82 INTEGER :: ji, jj, jk ! dummy loop indices 80 INTEGER :: j1, j2, i1, i281 83 INTEGER :: ibdy1, jbdy1, ibdy2, jbdy2 82 84 REAL(wp), DIMENSION(jpi,jpj) :: zub, zvb … … 93 95 Agrif_UseSpecialValue = .FALSE. 94 96 ! 95 ! prevent smoothing in ghost cells96 i1 = 1 ; i2 = nlci97 j1 = 1 ; j2 = nlcj98 IF( nbondj == -1 .OR. nbondj == 2 ) j1 = 2 + nbghostcells99 IF( nbondj == +1 .OR. nbondj == 2 ) j2 = nlcj - nbghostcells - 1100 IF( nbondi == -1 .OR. nbondi == 2 ) i1 = 2 + nbghostcells101 IF( nbondi == +1 .OR. nbondi == 2 ) i2 = nlci - nbghostcells - 1102 103 97 ! --- West --- ! 104 IF( nbondi == -1 .OR. nbondi == 2 ) THEN 105 ibdy1 = 2 106 ibdy2 = 1+nbghostcells 107 ! 108 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport 109 ua_b(ibdy1:ibdy2,:) = 0._wp 98 ibdy1 = 2 99 ibdy2 = 1+nbghostcells 100 ! 101 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport 102 DO ji = mi0(ibdy1), mi1(ibdy2) 103 ua_b(ji,:) = 0._wp 104 110 105 DO jk = 1, jpkm1 111 106 DO jj = 1, jpj 112 ua_b( ibdy1:ibdy2,jj) = ua_b(ibdy1:ibdy2,jj) &113 & + e3u_a(ibdy1:ibdy2,jj,jk) * ua(ibdy1:ibdy2,jj,jk) * umask(ibdy1:ibdy2,jj,jk)114 115 END DO 107 ua_b(ji,jj) = ua_b(ji,jj) + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 108 END DO 109 END DO 110 116 111 DO jj = 1, jpj 117 ua_b(ibdy1:ibdy2,jj) = ua_b(ibdy1:ibdy2,jj) * r1_hu_a(ibdy1:ibdy2,jj) 118 END DO 119 ENDIF 120 ! 121 IF( .NOT.lk_agrif_clp ) THEN 122 DO jk=1,jpkm1 ! Smooth 123 DO jj=j1,j2 124 ua(ibdy2,jj,jk) = 0.25_wp*(ua(ibdy2-1,jj,jk)+2._wp*ua(ibdy2,jj,jk)+ua(ibdy2+1,jj,jk)) 125 END DO 126 END DO 127 ENDIF 128 ! 129 zub(ibdy1:ibdy2,:) = 0._wp ! Correct transport 112 ua_b(ji,jj) = ua_b(ji,jj) * r1_hu_a(ji,jj) 113 END DO 114 END DO 115 ENDIF 116 ! 117 DO ji = mi0(ibdy1), mi1(ibdy2) 118 zub(ji,:) = 0._wp ! Correct transport 130 119 DO jk = 1, jpkm1 131 120 DO jj = 1, jpj 132 zub( ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) &133 & + e3u_a( ibdy1:ibdy2,jj,jk) * ua(ibdy1:ibdy2,jj,jk)*umask(ibdy1:ibdy2,jj,jk)121 zub(ji,jj) = zub(ji,jj) & 122 & + e3u_a(ji,jj,jk) * ua(ji,jj,jk)*umask(ji,jj,jk) 134 123 END DO 135 124 END DO 136 125 DO jj=1,jpj 137 zub( ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu_a(ibdy1:ibdy2,jj)126 zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj) 138 127 END DO 139 128 140 129 DO jk = 1, jpkm1 141 130 DO jj = 1, jpj 142 ua( ibdy1:ibdy2,jj,jk) = ( ua(ibdy1:ibdy2,jj,jk) &143 & + ua_b(ibdy1:ibdy2,jj)-zub(ibdy1:ibdy2,jj)) * umask(ibdy1:ibdy2,jj,jk)144 145 131 ua(ji,jj,jk) = ( ua(ji,jj,jk) + ua_b(ji,jj)-zub(ji,jj)) * umask(ji,jj,jk) 132 END DO 133 END DO 134 END DO 146 135 147 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 148 zvb(ibdy1:ibdy2,:) = 0._wp 136 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 137 DO ji = mi0(ibdy1), mi1(ibdy2) 138 zvb(ji,:) = 0._wp 149 139 DO jk = 1, jpkm1 150 140 DO jj = 1, jpj 151 zvb(ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) & 152 & + e3v_a(ibdy1:ibdy2,jj,jk) * va(ibdy1:ibdy2,jj,jk) * vmask(ibdy1:ibdy2,jj,jk) 141 zvb(ji,jj) = zvb(ji,jj) + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 153 142 END DO 154 143 END DO 155 144 DO jj = 1, jpj 156 zvb( ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv_a(ibdy1:ibdy2,jj)145 zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj) 157 146 END DO 158 147 DO jk = 1, jpkm1 159 148 DO jj = 1, jpj 160 va(ibdy1:ibdy2,jj,jk) = ( va(ibdy1:ibdy2,jj,jk) & 161 & + va_b(ibdy1:ibdy2,jj)-zvb(ibdy1:ibdy2,jj))*vmask(ibdy1:ibdy2,jj,jk) 162 END DO 163 END DO 164 ENDIF 165 ! 166 DO jk = 1, jpkm1 ! Mask domain edges 167 DO jj = 1, jpj 168 ua(1,jj,jk) = 0._wp 169 va(1,jj,jk) = 0._wp 170 END DO 171 END DO 149 va(ji,jj,jk) = ( va(ji,jj,jk) + va_b(ji,jj)-zvb(ji,jj))*vmask(ji,jj,jk) 150 END DO 151 END DO 152 END DO 172 153 ENDIF 173 154 174 155 ! --- East --- ! 175 IF( nbondi == 1 .OR. nbondi == 2 ) THEN176 ibdy1 = nlci-1-nbghostcells177 ibdy2 = nlci-2178 !179 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport180 ua_b( ibdy1:ibdy2,:) = 0._wp156 ibdy1 = jpiglo-1-nbghostcells 157 ibdy2 = jpiglo-2 158 ! 159 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport 160 DO ji = mi0(ibdy1), mi1(ibdy2) 161 ua_b(ji,:) = 0._wp 181 162 DO jk = 1, jpkm1 182 163 DO jj = 1, jpj 183 ua_b( ibdy1:ibdy2,jj) = ua_b(ibdy1:ibdy2,jj) &184 & + e3u_a( ibdy1:ibdy2,jj,jk) * ua(ibdy1:ibdy2,jj,jk) * umask(ibdy1:ibdy2,jj,jk)164 ua_b(ji,jj) = ua_b(ji,jj) & 165 & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 185 166 END DO 186 167 END DO 187 168 DO jj = 1, jpj 188 ua_b(ibdy1:ibdy2,jj) = ua_b(ibdy1:ibdy2,jj) * r1_hu_a(ibdy1:ibdy2,jj) 189 END DO 190 ENDIF 191 ! 192 IF( .NOT.lk_agrif_clp ) THEN 193 DO jk=1,jpkm1 ! Smooth 194 DO jj=j1,j2 195 ua(ibdy1,jj,jk) = 0.25_wp*(ua(ibdy1-1,jj,jk)+2._wp*ua(ibdy1,jj,jk)+ua(ibdy1+1,jj,jk)) 196 END DO 197 END DO 198 ENDIF 199 ! 200 zub(ibdy1:ibdy2,:) = 0._wp ! Correct transport 169 ua_b(ji,jj) = ua_b(ji,jj) * r1_hu_a(ji,jj) 170 END DO 171 END DO 172 ENDIF 173 ! 174 DO ji = mi0(ibdy1), mi1(ibdy2) 175 zub(ji,:) = 0._wp ! Correct transport 201 176 DO jk = 1, jpkm1 202 177 DO jj = 1, jpj 203 zub( ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) &204 & + e3u_a( ibdy1:ibdy2,jj,jk) * ua(ibdy1:ibdy2,jj,jk) * umask(ibdy1:ibdy2,jj,jk)178 zub(ji,jj) = zub(ji,jj) & 179 & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 205 180 END DO 206 181 END DO 207 182 DO jj=1,jpj 208 zub( ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu_a(ibdy1:ibdy2,jj)183 zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj) 209 184 END DO 210 185 211 186 DO jk = 1, jpkm1 212 187 DO jj = 1, jpj 213 ua(ibdy1:ibdy2,jj,jk) = ( ua(ibdy1:ibdy2,jj,jk) & 214 & + ua_b(ibdy1:ibdy2,jj)-zub(ibdy1:ibdy2,jj))*umask(ibdy1:ibdy2,jj,jk) 215 END DO 216 END DO 188 ua(ji,jj,jk) = ( ua(ji,jj,jk) & 189 & + ua_b(ji,jj)-zub(ji,jj))*umask(ji,jj,jk) 190 END DO 191 END DO 192 END DO 217 193 218 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 219 ibdy1 = ibdy1 + 1 220 ibdy2 = ibdy2 + 1 221 zvb(ibdy1:ibdy2,:) = 0._wp 194 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 195 ibdy1 = jpiglo-nbghostcells 196 ibdy2 = jpiglo-1 197 DO ji = mi0(ibdy1), mi1(ibdy2) 198 zvb(ji,:) = 0._wp 222 199 DO jk = 1, jpkm1 223 200 DO jj = 1, jpj 224 zvb( ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) &225 & + e3v_a( ibdy1:ibdy2,jj,jk) * va(ibdy1:ibdy2,jj,jk) * vmask(ibdy1:ibdy2,jj,jk)201 zvb(ji,jj) = zvb(ji,jj) & 202 & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 226 203 END DO 227 204 END DO 228 205 DO jj = 1, jpj 229 zvb( ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv_a(ibdy1:ibdy2,jj)206 zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj) 230 207 END DO 231 208 DO jk = 1, jpkm1 232 209 DO jj = 1, jpj 233 va(ibdy1:ibdy2,jj,jk) = ( va(ibdy1:ibdy2,jj,jk) & 234 & + va_b(ibdy1:ibdy2,jj)-zvb(ibdy1:ibdy2,jj)) * vmask(ibdy1:ibdy2,jj,jk) 235 END DO 236 END DO 237 ENDIF 238 ! 239 DO jk = 1, jpkm1 ! Mask domain edges 240 DO jj = 1, jpj 241 ua(nlci-1,jj,jk) = 0._wp 242 va(nlci ,jj,jk) = 0._wp 243 END DO 244 END DO 210 va(ji,jj,jk) = ( va(ji,jj,jk) & 211 & + va_b(ji,jj)-zvb(ji,jj)) * vmask(ji,jj,jk) 212 END DO 213 END DO 214 END DO 245 215 ENDIF 246 216 247 217 ! --- South --- ! 248 IF( nbondj == -1 .OR. nbondj == 2 ) THEN249 jbdy1 = 2250 jbdy2 = 1+nbghostcells251 !252 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport253 va_b(:,j bdy1:jbdy2) = 0._wp218 jbdy1 = 2 219 jbdy2 = 1+nbghostcells 220 ! 221 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport 222 DO jj = mj0(jbdy1), mj1(jbdy2) 223 va_b(:,jj) = 0._wp 254 224 DO jk = 1, jpkm1 255 225 DO ji = 1, jpi 256 va_b(ji,j bdy1:jbdy2) = va_b(ji,jbdy1:jbdy2) &257 & + e3v_a(ji,j bdy1:jbdy2,jk) * va(ji,jbdy1:jbdy2,jk) * vmask(ji,jbdy1:jbdy2,jk)226 va_b(ji,jj) = va_b(ji,jj) & 227 & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 258 228 END DO 259 229 END DO 260 230 DO ji=1,jpi 261 va_b(ji,jbdy1:jbdy2) = va_b(ji,jbdy1:jbdy2) * r1_hv_a(ji,jbdy1:jbdy2) 262 END DO 263 ENDIF 264 ! 265 IF ( .NOT.lk_agrif_clp ) THEN 266 DO jk = 1, jpkm1 ! Smooth 267 DO ji = i1, i2 268 va(ji,jbdy2,jk) = 0.25_wp*(va(ji,jbdy2-1,jk)+2._wp*va(ji,jbdy2,jk)+va(ji,jbdy2+1,jk)) 269 END DO 270 END DO 271 ENDIF 272 ! 273 zvb(:,jbdy1:jbdy2) = 0._wp ! Correct transport 231 va_b(ji,jj) = va_b(ji,jj) * r1_hv_a(ji,jj) 232 END DO 233 END DO 234 ENDIF 235 ! 236 DO jj = mj0(jbdy1), mj1(jbdy2) 237 zvb(:,jj) = 0._wp ! Correct transport 274 238 DO jk=1,jpkm1 275 239 DO ji=1,jpi 276 zvb(ji,j bdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) &277 & + e3v_a(ji,j bdy1:jbdy2,jk) * va(ji,jbdy1:jbdy2,jk) * vmask(ji,jbdy1:jbdy2,jk)240 zvb(ji,jj) = zvb(ji,jj) & 241 & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 278 242 END DO 279 243 END DO 280 244 DO ji = 1, jpi 281 zvb(ji,j bdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv_a(ji,jbdy1:jbdy2)245 zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj) 282 246 END DO 283 247 284 248 DO jk = 1, jpkm1 285 249 DO ji = 1, jpi 286 va(ji,jbdy1:jbdy2,jk) = ( va(ji,jbdy1:jbdy2,jk) & 287 & + va_b(ji,jbdy1:jbdy2) - zvb(ji,jbdy1:jbdy2) ) * vmask(ji,jbdy1:jbdy2,jk) 288 END DO 289 END DO 250 va(ji,jj,jk) = ( va(ji,jj,jk) & 251 & + va_b(ji,jj) - zvb(ji,jj) ) * vmask(ji,jj,jk) 252 END DO 253 END DO 254 END DO 290 255 291 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 292 zub(:,jbdy1:jbdy2) = 0._wp 256 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 257 DO jj = mj0(jbdy1), mj1(jbdy2) 258 zub(:,jj) = 0._wp 293 259 DO jk = 1, jpkm1 294 260 DO ji = 1, jpi 295 zub(ji,j bdy1:jbdy2) = zub(ji,jbdy1:jbdy2) &296 & + e3u_a(ji,j bdy1:jbdy2,jk) * ua(ji,jbdy1:jbdy2,jk) * umask(ji,jbdy1:jbdy2,jk)261 zub(ji,jj) = zub(ji,jj) & 262 & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 297 263 END DO 298 264 END DO 299 265 DO ji = 1, jpi 300 zub(ji,j bdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu_a(ji,jbdy1:jbdy2)266 zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj) 301 267 END DO 302 268 303 269 DO jk = 1, jpkm1 304 270 DO ji = 1, jpi 305 ua(ji,jbdy1:jbdy2,jk) = ( ua(ji,jbdy1:jbdy2,jk) & 306 & + ua_b(ji,jbdy1:jbdy2) - zub(ji,jbdy1:jbdy2) ) * umask(ji,jbdy1:jbdy2,jk) 307 END DO 308 END DO 309 ENDIF 310 ! 311 DO jk = 1, jpkm1 ! Mask domain edges 312 DO ji = 1, jpi 313 ua(ji,1,jk) = 0._wp 314 va(ji,1,jk) = 0._wp 315 END DO 316 END DO 271 ua(ji,jj,jk) = ( ua(ji,jj,jk) & 272 & + ua_b(ji,jj) - zub(ji,jj) ) * umask(ji,jj,jk) 273 END DO 274 END DO 275 END DO 317 276 ENDIF 318 277 319 278 ! --- North --- ! 320 IF( nbondj == 1 .OR. nbondj == 2 ) THEN321 jbdy1 = nlcj-1-nbghostcells322 jbdy2 = nlcj-2323 !324 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport325 va_b(:,j bdy1:jbdy2) = 0._wp279 jbdy1 = jpjglo-1-nbghostcells 280 jbdy2 = jpjglo-2 281 ! 282 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport 283 DO jj = mj0(jbdy1), mj1(jbdy2) 284 va_b(:,jj) = 0._wp 326 285 DO jk = 1, jpkm1 327 286 DO ji = 1, jpi 328 va_b(ji,j bdy1:jbdy2) = va_b(ji,jbdy1:jbdy2) &329 & + e3v_a(ji,j bdy1:jbdy2,jk) * va(ji,jbdy1:jbdy2,jk) * vmask(ji,jbdy1:jbdy2,jk)287 va_b(ji,jj) = va_b(ji,jj) & 288 & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 330 289 END DO 331 290 END DO 332 291 DO ji=1,jpi 333 va_b(ji,jbdy1:jbdy2) = va_b(ji,jbdy1:jbdy2) * r1_hv_a(ji,jbdy1:jbdy2) 334 END DO 335 ENDIF 336 ! 337 IF ( .NOT.lk_agrif_clp ) THEN 338 DO jk = 1, jpkm1 ! Smooth 339 DO ji = i1, i2 340 va(ji,jbdy1,jk) = 0.25_wp*(va(ji,jbdy1-1,jk)+2._wp*va(ji,jbdy1,jk)+va(ji,jbdy1+1,jk)) 341 END DO 342 END DO 343 ENDIF 344 ! 345 zvb(:,jbdy1:jbdy2) = 0._wp ! Correct transport 292 va_b(ji,jj) = va_b(ji,jj) * r1_hv_a(ji,jj) 293 END DO 294 END DO 295 ENDIF 296 ! 297 DO jj = mj0(jbdy1), mj1(jbdy2) 298 zvb(:,jj) = 0._wp ! Correct transport 346 299 DO jk=1,jpkm1 347 300 DO ji=1,jpi 348 zvb(ji,j bdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) &349 & + e3v_a(ji,j bdy1:jbdy2,jk) * va(ji,jbdy1:jbdy2,jk) * vmask(ji,jbdy1:jbdy2,jk)301 zvb(ji,jj) = zvb(ji,jj) & 302 & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 350 303 END DO 351 304 END DO 352 305 DO ji = 1, jpi 353 zvb(ji,j bdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv_a(ji,jbdy1:jbdy2)306 zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj) 354 307 END DO 355 308 356 309 DO jk = 1, jpkm1 357 310 DO ji = 1, jpi 358 va(ji,jbdy1:jbdy2,jk) = ( va(ji,jbdy1:jbdy2,jk) & 359 & + va_b(ji,jbdy1:jbdy2) - zvb(ji,jbdy1:jbdy2) ) * vmask(ji,jbdy1:jbdy2,jk) 360 END DO 361 END DO 311 va(ji,jj,jk) = ( va(ji,jj,jk) & 312 & + va_b(ji,jj) - zvb(ji,jj) ) * vmask(ji,jj,jk) 313 END DO 314 END DO 315 END DO 362 316 363 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 364 jbdy1 = jbdy1 + 1 365 jbdy2 = jbdy2 + 1 366 zub(:,jbdy1:jbdy2) = 0._wp 317 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 318 jbdy1 = jpjglo-nbghostcells 319 jbdy2 = jpjglo-1 320 DO jj = mj0(jbdy1), mj1(jbdy2) 321 zub(:,jj) = 0._wp 367 322 DO jk = 1, jpkm1 368 323 DO ji = 1, jpi 369 zub(ji,j bdy1:jbdy2) = zub(ji,jbdy1:jbdy2) &370 & + e3u_a(ji,j bdy1:jbdy2,jk) * ua(ji,jbdy1:jbdy2,jk) * umask(ji,jbdy1:jbdy2,jk)324 zub(ji,jj) = zub(ji,jj) & 325 & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 371 326 END DO 372 327 END DO 373 328 DO ji = 1, jpi 374 zub(ji,j bdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu_a(ji,jbdy1:jbdy2)329 zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj) 375 330 END DO 376 331 377 332 DO jk = 1, jpkm1 378 333 DO ji = 1, jpi 379 ua(ji,jbdy1:jbdy2,jk) = ( ua(ji,jbdy1:jbdy2,jk) & 380 & + ua_b(ji,jbdy1:jbdy2) - zub(ji,jbdy1:jbdy2) ) * umask(ji,jbdy1:jbdy2,jk) 381 END DO 382 END DO 383 ENDIF 384 ! 385 DO jk = 1, jpkm1 ! Mask domain edges 386 DO ji = 1, jpi 387 ua(ji,nlcj ,jk) = 0._wp 388 va(ji,nlcj-1,jk) = 0._wp 389 END DO 390 END DO 334 ua(ji,jj,jk) = ( ua(ji,jj,jk) & 335 & + ua_b(ji,jj) - zub(ji,jj) ) * umask(ji,jj,jk) 336 END DO 337 END DO 338 END DO 391 339 ENDIF 392 340 ! … … 401 349 !! 402 350 INTEGER :: ji, jj 351 INTEGER :: istart, iend, jstart, jend 403 352 !!---------------------------------------------------------------------- 404 353 ! 405 354 IF( Agrif_Root() ) RETURN 406 355 ! 407 IF((nbondi == -1).OR.(nbondi == 2)) THEN 356 !--- West ---! 357 istart = 2 358 iend = nbghostcells+1 359 DO ji = mi0(istart), mi1(iend) 408 360 DO jj=1,jpj 409 va_e(2:nbghostcells+1,jj) = vbdy_w(1:nbghostcells,jj) * hvr_e(2:nbghostcells+1,jj) 410 ! Specified fluxes: 411 ua_e(2:nbghostcells+1,jj) = ubdy_w(1:nbghostcells,jj) * hur_e(2:nbghostcells+1,jj) 412 ! Characteristics method (only if ghostcells=1): 413 !alt ua_e(2,jj) = 0.5_wp * ( ubdy_w(jj) * hur_e(2,jj) + ua_e(3,jj) & 414 !alt & - sqrt(grav * hur_e(2,jj)) * (sshn_e(3,jj) - hbdy_w(jj)) ) 415 END DO 416 ENDIF 417 ! 418 IF((nbondi == 1).OR.(nbondi == 2)) THEN 361 va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 362 ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 363 END DO 364 END DO 365 ! 366 !--- East ---! 367 istart = jpiglo-nbghostcells 368 iend = jpiglo-1 369 DO ji = mi0(istart), mi1(iend) 419 370 DO jj=1,jpj 420 va_e(nlci-nbghostcells:nlci-1,jj) = vbdy_e(1:nbghostcells,jj) * hvr_e(nlci-nbghostcells:nlci-1,jj) 421 ! Specified fluxes: 422 ua_e(nlci-nbghostcells-1:nlci-2,jj) = ubdy_e(1:nbghostcells,jj) * hur_e(nlci-nbghostcells-1:nlci-2,jj) 423 ! Characteristics method (only if ghostcells=1): 424 !alt ua_e(nlci-2,jj) = 0.5_wp * ( ubdy_e(jj) * hur_e(nlci-2,jj) + ua_e(nlci-3,jj) & 425 !alt & + sqrt(grav * hur_e(nlci-2,jj)) * (sshn_e(nlci-2,jj) - hbdy_e(jj)) ) 426 END DO 427 ENDIF 428 ! 429 IF((nbondj == -1).OR.(nbondj == 2)) THEN 371 va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 372 END DO 373 END DO 374 istart = jpiglo-nbghostcells-1 375 iend = jpiglo-2 376 DO ji = mi0(istart), mi1(iend) 377 DO jj=1,jpj 378 ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 379 END DO 380 END DO 381 ! 382 !--- South ---! 383 jstart = 2 384 jend = nbghostcells+1 385 DO jj = mj0(jstart), mj1(jend) 430 386 DO ji=1,jpi 431 ua_e(ji,2:nbghostcells+1) = ubdy_s(ji,1:nbghostcells) * hur_e(ji,2:nbghostcells+1) 432 ! Specified fluxes: 433 va_e(ji,2:nbghostcells+1) = vbdy_s(ji,1:nbghostcells) * hvr_e(ji,2:nbghostcells+1) 434 ! Characteristics method (only if ghostcells=1): 435 !alt va_e(ji,2) = 0.5_wp * ( vbdy_s(ji) * hvr_e(ji,2) + va_e(ji,3) & 436 !alt & - sqrt(grav * hvr_e(ji,2)) * (sshn_e(ji,3) - hbdy_s(ji)) ) 437 END DO 438 ENDIF 439 ! 440 IF((nbondj == 1).OR.(nbondj == 2)) THEN 387 ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 388 va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 389 END DO 390 END DO 391 ! 392 !--- North ---! 393 jstart = jpjglo-nbghostcells 394 jend = jpjglo-1 395 DO jj = mj0(jstart), mj1(jend) 441 396 DO ji=1,jpi 442 ua_e(ji,nlcj-nbghostcells:nlcj-1) = ubdy_n(ji,1:nbghostcells) * hur_e(ji,nlcj-nbghostcells:nlcj-1) 443 ! Specified fluxes: 444 va_e(ji,nlcj-nbghostcells-1:nlcj-2) = vbdy_n(ji,1:nbghostcells) * hvr_e(ji,nlcj-nbghostcells-1:nlcj-2) 445 ! Characteristics method (only if ghostcells=1): 446 !alt va_e(ji,nlcj-2) = 0.5_wp * ( vbdy_n(ji) * hvr_e(ji,nlcj-2) + va_e(ji,nlcj-3) & 447 !alt & + sqrt(grav * hvr_e(ji,nlcj-2)) * (sshn_e(ji,nlcj-2) - hbdy_n(ji)) ) 448 END DO 449 ENDIF 397 ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 398 END DO 399 END DO 400 jstart = jpjglo-nbghostcells-1 401 jend = jpjglo-2 402 DO jj = mj0(jstart), mj1(jend) 403 DO ji=1,jpi 404 va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 405 END DO 406 END DO 450 407 ! 451 408 END SUBROUTINE Agrif_dyn_ts 452 409 410 SUBROUTINE Agrif_dyn_ts_flux( jn, zu, zv ) 411 !!---------------------------------------------------------------------- 412 !! *** ROUTINE Agrif_dyn_ts_flux *** 413 !!---------------------------------------------------------------------- 414 INTEGER, INTENT(in) :: jn 415 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: zu, zv 416 !! 417 INTEGER :: ji, jj 418 INTEGER :: istart, iend, jstart, jend 419 !!---------------------------------------------------------------------- 420 ! 421 IF( Agrif_Root() ) RETURN 422 ! 423 !--- West ---! 424 istart = 2 425 iend = nbghostcells+1 426 DO ji = mi0(istart), mi1(iend) 427 DO jj=1,jpj 428 zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 429 zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 430 END DO 431 END DO 432 ! 433 !--- East ---! 434 istart = jpiglo-nbghostcells 435 iend = jpiglo-1 436 DO ji = mi0(istart), mi1(iend) 437 DO jj=1,jpj 438 zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 439 END DO 440 END DO 441 istart = jpiglo-nbghostcells-1 442 iend = jpiglo-2 443 DO ji = mi0(istart), mi1(iend) 444 DO jj=1,jpj 445 zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 446 END DO 447 END DO 448 ! 449 !--- South ---! 450 jstart = 2 451 jend = nbghostcells+1 452 DO jj = mj0(jstart), mj1(jend) 453 DO ji=1,jpi 454 zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 455 zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 456 END DO 457 END DO 458 ! 459 !--- North ---! 460 jstart = jpjglo-nbghostcells 461 jend = jpjglo-1 462 DO jj = mj0(jstart), mj1(jend) 463 DO ji=1,jpi 464 zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 465 END DO 466 END DO 467 jstart = jpjglo-nbghostcells-1 468 jend = jpjglo-2 469 DO jj = mj0(jstart), mj1(jend) 470 DO ji=1,jpi 471 zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 472 END DO 473 END DO 474 ! 475 END SUBROUTINE Agrif_dyn_ts_flux 453 476 454 477 SUBROUTINE Agrif_dta_ts( kt ) … … 470 493 ! 471 494 ! Interpolate barotropic fluxes 472 Agrif_SpecialValue =0._wp495 Agrif_SpecialValue = 0._wp 473 496 Agrif_UseSpecialValue = ln_spc_dyn 497 ! 498 ! Set bdy time interpolation stage to 0 (latter incremented locally do deal with corners) 499 utint_stage(:,:) = 0 500 vtint_stage(:,:) = 0 474 501 ! 475 502 IF( ll_int_cons ) THEN ! Conservative interpolation 476 503 ! order matters here !!!!!! 477 504 CALL Agrif_Bc_variable( ub2b_interp_id, calledweight=1._wp, procname=interpub2b ) ! Time integrated 478 CALL Agrif_Bc_variable( vb2b_interp_id, calledweight=1._wp, procname=interpvb2b ) 505 CALL Agrif_Bc_variable( vb2b_interp_id, calledweight=1._wp, procname=interpvb2b ) 506 ! 479 507 bdy_tinterp = 1 480 508 CALL Agrif_Bc_variable( unb_id , calledweight=1._wp, procname=interpunb ) ! After 481 CALL Agrif_Bc_variable( vnb_id , calledweight=1._wp, procname=interpvnb ) 509 CALL Agrif_Bc_variable( vnb_id , calledweight=1._wp, procname=interpvnb ) 510 ! 482 511 bdy_tinterp = 2 483 512 CALL Agrif_Bc_variable( unb_id , calledweight=0._wp, procname=interpunb ) ! Before 484 CALL Agrif_Bc_variable( vnb_id , calledweight=0._wp, procname=interpvnb ) 513 CALL Agrif_Bc_variable( vnb_id , calledweight=0._wp, procname=interpvnb ) 485 514 ELSE ! Linear interpolation 486 bdy_tinterp = 0 487 ubdy_w(:,:) = 0._wp ; vbdy_w(:,:) = 0._wp 488 ubdy_e(:,:) = 0._wp ; vbdy_e(:,:) = 0._wp 489 ubdy_n(:,:) = 0._wp ; vbdy_n(:,:) = 0._wp 490 ubdy_s(:,:) = 0._wp ; vbdy_s(:,:) = 0._wp 515 ! 516 ubdy(:,:) = 0._wp ; vbdy(:,:) = 0._wp 491 517 CALL Agrif_Bc_variable( unb_id, procname=interpunb ) 492 518 CALL Agrif_Bc_variable( vnb_id, procname=interpvnb ) … … 503 529 INTEGER, INTENT(in) :: kt 504 530 ! 505 INTEGER :: ji, jj, indx, indy 531 INTEGER :: ji, jj 532 INTEGER :: istart, iend, jstart, jend 506 533 !!---------------------------------------------------------------------- 507 534 ! … … 516 543 ! 517 544 ! --- West --- ! 518 IF((nbondi == -1).OR.(nbondi == 2)) THEN 519 indx = 1+nbghostcells 545 istart = 2 546 iend = 1 + nbghostcells 547 DO ji = mi0(istart), mi1(iend) 520 548 DO jj = 1, jpj 521 DO ji = 2, indx 522 ssha(ji,jj) = hbdy_w(ji-1,jj) 523 ENDDO 549 ssha(ji,jj) = hbdy(ji,jj) 524 550 ENDDO 525 END IF551 ENDDO 526 552 ! 527 553 ! --- East --- ! 528 IF((nbondi == 1).OR.(nbondi == 2)) THEN 529 indx = nlci-nbghostcells 554 istart = jpiglo - nbghostcells 555 iend = jpiglo - 1 556 DO ji = mi0(istart), mi1(iend) 530 557 DO jj = 1, jpj 531 DO ji = indx, nlci-1 532 ssha(ji,jj) = hbdy_e(ji-indx+1,jj) 533 ENDDO 558 ssha(ji,jj) = hbdy(ji,jj) 534 559 ENDDO 535 END IF560 ENDDO 536 561 ! 537 562 ! --- South --- ! 538 IF((nbondj == -1).OR.(nbondj == 2)) THEN 539 indy = 1+nbghostcells 540 DO jj = 2, indy 541 DO ji = 1, jpi 542 ssha(ji,jj) = hbdy_s(ji,jj-1) 543 ENDDO 563 jstart = 2 564 jend = 1 + nbghostcells 565 DO jj = mj0(jstart), mj1(jend) 566 DO ji = 1, jpi 567 ssha(ji,jj) = hbdy(ji,jj) 544 568 ENDDO 545 END IF569 ENDDO 546 570 ! 547 571 ! --- North --- ! 548 IF((nbondj == 1).OR.(nbondj == 2)) THEN 549 indy = nlcj-nbghostcells 550 DO jj = indy, nlcj-1 551 DO ji = 1, jpi 552 ssha(ji,jj) = hbdy_n(ji,jj-indy+1) 553 ENDDO 572 jstart = jpjglo - nbghostcells 573 jend = jpjglo - 1 574 DO jj = mj0(jstart), mj1(jend) 575 DO ji = 1, jpi 576 ssha(ji,jj) = hbdy(ji,jj) 554 577 ENDDO 555 END IF578 ENDDO 556 579 ! 557 580 END SUBROUTINE Agrif_ssh … … 564 587 INTEGER, INTENT(in) :: jn 565 588 !! 566 INTEGER :: ji, jj , indx, indy567 !!----------------------------------------------------------------------568 !! clem ghost (starting at i,j=1 is important I think otherwise you introduce a grad(ssh)/=0 at point 2)589 INTEGER :: ji, jj 590 INTEGER :: istart, iend, jstart, jend 591 !!---------------------------------------------------------------------- 569 592 ! 570 593 IF( Agrif_Root() ) RETURN 571 594 ! 572 595 ! --- West --- ! 573 IF((nbondi == -1).OR.(nbondi == 2)) THEN 574 indx = 1+nbghostcells 596 istart = 2 597 iend = 1+nbghostcells 598 DO ji = mi0(istart), mi1(iend) 575 599 DO jj = 1, jpj 576 DO ji = 2, indx 577 ssha_e(ji,jj) = hbdy_w(ji-1,jj) 578 ENDDO 600 ssha_e(ji,jj) = hbdy(ji,jj) 579 601 ENDDO 580 END IF602 ENDDO 581 603 ! 582 604 ! --- East --- ! 583 IF((nbondi == 1).OR.(nbondi == 2)) THEN 584 indx = nlci-nbghostcells 605 istart = jpiglo - nbghostcells 606 iend = jpiglo - 1 607 DO ji = mi0(istart), mi1(iend) 585 608 DO jj = 1, jpj 586 DO ji = indx, nlci-1 587 ssha_e(ji,jj) = hbdy_e(ji-indx+1,jj) 588 ENDDO 609 ssha_e(ji,jj) = hbdy(ji,jj) 589 610 ENDDO 590 END IF611 ENDDO 591 612 ! 592 613 ! --- South --- ! 593 IF((nbondj == -1).OR.(nbondj == 2)) THEN 594 indy = 1+nbghostcells 595 DO jj = 2, indy 596 DO ji = 1, jpi 597 ssha_e(ji,jj) = hbdy_s(ji,jj-1) 598 ENDDO 614 jstart = 2 615 jend = 1+nbghostcells 616 DO jj = mj0(jstart), mj1(jend) 617 DO ji = 1, jpi 618 ssha_e(ji,jj) = hbdy(ji,jj) 599 619 ENDDO 600 END IF620 ENDDO 601 621 ! 602 622 ! --- North --- ! 603 IF((nbondj == 1).OR.(nbondj == 2)) THEN 604 indy = nlcj-nbghostcells 605 DO jj = indy, nlcj-1 606 DO ji = 1, jpi 607 ssha_e(ji,jj) = hbdy_n(ji,jj-indy+1) 608 ENDDO 623 jstart = jpjglo - nbghostcells 624 jend = jpjglo - 1 625 DO jj = mj0(jstart), mj1(jend) 626 DO ji = 1, jpi 627 ssha_e(ji,jj) = hbdy(ji,jj) 609 628 ENDDO 610 END IF629 ENDDO 611 630 ! 612 631 END SUBROUTINE Agrif_ssh_ts … … 634 653 635 654 636 SUBROUTINE interptsn( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before , nb, ndir)655 SUBROUTINE interptsn( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before ) 637 656 !!---------------------------------------------------------------------- 638 657 !! *** ROUTINE interptsn *** … … 641 660 INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2 642 661 LOGICAL , INTENT(in ) :: before 643 INTEGER , INTENT(in ) :: nb , ndir 644 ! 645 INTEGER :: ji, jj, jk, jn, iref, jref, ibdy, jbdy ! dummy loop indices 646 INTEGER :: imin, imax, jmin, jmax, N_in, N_out 647 REAL(wp) :: zrho, z1, z2, z3, z4, z5, z6, z7 648 LOGICAL :: western_side, eastern_side,northern_side,southern_side 662 ! 663 INTEGER :: ji, jj, jk, jn ! dummy loop indices 664 INTEGER :: N_in, N_out 649 665 ! vertical interpolation: 650 REAL(wp) , DIMENSION(i1:i2,j1:j2,1:jpk,n1:n2) :: ptab_child651 REAL(wp), DIMENSION(k1:k2, n1:n2-1) :: tabin666 REAL(wp) :: zhtot 667 REAL(wp), DIMENSION(k1:k2,1:jpts) :: tabin 652 668 REAL(wp), DIMENSION(k1:k2) :: h_in 653 669 REAL(wp), DIMENSION(1:jpk) :: h_out 654 REAL(wp) :: h_diff670 !!---------------------------------------------------------------------- 655 671 656 672 IF( before ) THEN … … 666 682 667 683 # if defined key_vertical 684 ! Interpolate thicknesses 685 ! Warning: these are masked, hence extrapolated prior interpolation. 668 686 DO jk=k1,k2 669 687 DO jj=j1,j2 670 688 DO ji=i1,i2 671 ptab(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t_n(ji,jj,jk)689 ptab(ji,jj,jk,jpts+1) = tmask(ji,jj,jk) * e3t_n(ji,jj,jk) 672 690 END DO 673 691 END DO 674 692 END DO 693 694 ! Extrapolate thicknesses in partial bottom cells: 695 ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 696 IF (ln_zps) THEN 697 DO jj=j1,j2 698 DO ji=i1,i2 699 jk = mbkt(ji,jj) 700 ptab(ji,jj,jk,jpts+1) = 0._wp 701 END DO 702 END DO 703 END IF 704 705 ! Save ssh at last level: 706 IF (.NOT.ln_linssh) THEN 707 ptab(i1:i2,j1:j2,k2,jpts+1) = sshn(i1:i2,j1:j2)*tmask(i1:i2,j1:j2,1) 708 ELSE 709 ptab(i1:i2,j1:j2,k2,jpts+1) = 0._wp 710 END IF 675 711 # endif 676 712 ELSE 677 713 678 western_side = (nb == 1).AND.(ndir == 1) ; eastern_side = (nb == 1).AND.(ndir == 2) 679 southern_side = (nb == 2).AND.(ndir == 1) ; northern_side = (nb == 2).AND.(ndir == 2) 680 681 # if defined key_vertical 714 # if defined key_vertical 715 IF (ln_linssh) ptab(i1:i2,j1:j2,k2,n2) = 0._wp 716 682 717 DO jj=j1,j2 683 718 DO ji=i1,i2 684 iref = ji685 jref = jj686 if(western_side) iref=MAX(2,ji)687 if(eastern_side) iref=MIN(nlci-1,ji)688 if(southern_side) jref=MAX(2,jj)689 if(northern_side) jref=MIN(nlcj-1,jj)690 N_in = 0691 DO jk=k1,k2 !k2 = jpk of parent grid692 IF (ptab(ji,jj,jk,n2) == 0) EXIT693 N_in = N_in + 1719 tsa(ji,jj,:,:) = 0._wp 720 N_in = mbkt_parent(ji,jj) 721 zhtot = 0._wp 722 DO jk=1,N_in !k2 = jpk of parent grid 723 IF (jk==N_in) THEN 724 h_in(jk) = ht0_parent(ji,jj) + ptab(ji,jj,k2,n2) - zhtot 725 ELSE 726 h_in(jk) = ptab(ji,jj,jk,n2) 727 ENDIF 728 zhtot = zhtot + h_in(jk) 694 729 tabin(jk,:) = ptab(ji,jj,jk,n1:n2-1) 695 h_in(N_in) = ptab(ji,jj,jk,n2)696 730 END DO 697 731 N_out = 0 698 732 DO jk=1,jpk ! jpk of child grid 699 IF (tmask( iref,jref,jk) == 0) EXIT733 IF (tmask(ji,jj,jk) == 0._wp) EXIT 700 734 N_out = N_out + 1 701 h_out(jk) = e3t_ n(iref,jref,jk)735 h_out(jk) = e3t_a(ji,jj,jk) 702 736 ENDDO 703 IF (N_in > 0) THEN 704 DO jn=1,jpts 705 call reconstructandremap(tabin(1:N_in,jn),h_in,ptab_child(ji,jj,1:N_out,jn),h_out,N_in,N_out) 706 ENDDO 737 IF (N_in*N_out > 0) THEN 738 CALL reconstructandremap(tabin(1:N_in,1:jpts),h_in(1:N_in),tsa(ji,jj,1:N_out,1:jpts),h_out(1:N_out),N_in,N_out,jpts) 707 739 ENDIF 708 740 ENDDO 709 741 ENDDO 710 742 # else 711 ptab_child(i1:i2,j1:j2,1:jpk,1:jpts) = ptab(i1:i2,j1:j2,1:jpk,1:jpts)712 # endif713 743 ! 714 744 DO jn=1, jpts 715 tsa(i1:i2,j1:j2,1:jpk,jn)=ptab_child(i1:i2,j1:j2,1:jpk,jn)*tmask(i1:i2,j1:j2,1:jpk) 716 END DO 717 718 IF ( .NOT.lk_agrif_clp ) THEN 719 ! 720 imin = i1 ; imax = i2 721 jmin = j1 ; jmax = j2 722 ! 723 ! Remove CORNERS 724 IF((nbondj == -1).OR.(nbondj == 2)) jmin = 2 + nbghostcells 725 IF((nbondj == +1).OR.(nbondj == 2)) jmax = nlcj - nbghostcells - 1 726 IF((nbondi == -1).OR.(nbondi == 2)) imin = 2 + nbghostcells 727 IF((nbondi == +1).OR.(nbondi == 2)) imax = nlci - nbghostcells - 1 728 ! 729 IF( eastern_side ) THEN 730 zrho = Agrif_Rhox() 731 z1 = ( zrho - 1._wp ) * 0.5_wp 732 z3 = ( zrho - 1._wp ) / ( zrho + 1._wp ) 733 z6 = 2._wp * ( zrho - 1._wp ) / ( zrho + 1._wp ) 734 z7 = - ( zrho - 1._wp ) / ( zrho + 3._wp ) 735 z2 = 1._wp - z1 ; z4 = 1._wp - z3 ; z5 = 1._wp - z6 - z7 736 ! 737 ibdy = nlci-nbghostcells 738 DO jn = 1, jpts 739 tsa(ibdy+1,jmin:jmax,1:jpkm1,jn) = z1 * ptab_child(ibdy+1,jmin:jmax,1:jpkm1,jn) + z2 * ptab_child(ibdy,jmin:jmax,1:jpkm1,jn) 740 DO jk = 1, jpkm1 741 DO jj = jmin,jmax 742 IF( umask(ibdy-1,jj,jk) == 0._wp ) THEN 743 tsa(ibdy,jj,jk,jn) = tsa(ibdy+1,jj,jk,jn) * tmask(ibdy,jj,jk) 744 ELSE 745 tsa(ibdy,jj,jk,jn)=(z4*tsa(ibdy+1,jj,jk,jn)+z3*tsa(ibdy-1,jj,jk,jn))*tmask(ibdy,jj,jk) 746 IF( un(ibdy-1,jj,jk) > 0._wp ) THEN 747 tsa(ibdy,jj,jk,jn)=( z6*tsa(ibdy-1,jj,jk,jn)+z5*tsa(ibdy+1,jj,jk,jn) & 748 + z7*tsa(ibdy-2,jj,jk,jn) ) * tmask(ibdy,jj,jk) 749 ENDIF 750 ENDIF 751 END DO 752 END DO 753 ! Restore ghost points: 754 tsa(ibdy+1,jmin:jmax,1:jpkm1,jn) = ptab_child(ibdy+1,jmin:jmax,1:jpkm1,jn) * tmask(ibdy+1,jmin:jmax,1:jpkm1) 755 END DO 756 ENDIF 757 ! 758 IF( northern_side ) THEN 759 zrho = Agrif_Rhoy() 760 z1 = ( zrho - 1._wp ) * 0.5_wp 761 z3 = ( zrho - 1._wp ) / ( zrho + 1._wp ) 762 z6 = 2._wp * ( zrho - 1._wp ) / ( zrho + 1._wp ) 763 z7 = - ( zrho - 1._wp ) / ( zrho + 3._wp ) 764 z2 = 1._wp - z1 ; z4 = 1._wp - z3 ; z5 = 1._wp - z6 - z7 765 ! 766 jbdy = nlcj-nbghostcells 767 DO jn = 1, jpts 768 tsa(imin:imax,jbdy+1,1:jpkm1,jn) = z1 * ptab_child(imin:imax,jbdy+1,1:jpkm1,jn) + z2 * ptab_child(imin:imax,jbdy,1:jpkm1,jn) 769 DO jk = 1, jpkm1 770 DO ji = imin,imax 771 IF( vmask(ji,jbdy-1,jk) == 0._wp ) THEN 772 tsa(ji,jbdy,jk,jn) = tsa(ji,jbdy+1,jk,jn) * tmask(ji,jbdy,jk) 773 ELSE 774 tsa(ji,jbdy,jk,jn)=(z4*tsa(ji,jbdy+1,jk,jn)+z3*tsa(ji,jbdy-1,jk,jn))*tmask(ji,jbdy,jk) 775 IF (vn(ji,jbdy-1,jk) > 0._wp ) THEN 776 tsa(ji,jbdy,jk,jn)=( z6*tsa(ji,jbdy-1,jk,jn)+z5*tsa(ji,jbdy+1,jk,jn) & 777 + z7*tsa(ji,jbdy-2,jk,jn) ) * tmask(ji,jbdy,jk) 778 ENDIF 779 ENDIF 780 END DO 781 END DO 782 ! Restore ghost points: 783 tsa(imin:imax,jbdy+1,1:jpkm1,jn) = ptab_child(imin:imax,jbdy+1,1:jpkm1,jn) * tmask(imin:imax,jbdy+1,1:jpkm1) 784 END DO 785 ENDIF 786 ! 787 IF( western_side ) THEN 788 zrho = Agrif_Rhox() 789 z1 = ( zrho - 1._wp ) * 0.5_wp 790 z3 = ( zrho - 1._wp ) / ( zrho + 1._wp ) 791 z6 = 2._wp * ( zrho - 1._wp ) / ( zrho + 1._wp ) 792 z7 = - ( zrho - 1._wp ) / ( zrho + 3._wp ) 793 z2 = 1._wp - z1 ; z4 = 1._wp - z3 ; z5 = 1._wp - z6 - z7 794 ! 795 ibdy = 1+nbghostcells 796 DO jn = 1, jpts 797 tsa(ibdy-1,jmin:jmax,1:jpkm1,jn) = z1 * ptab_child(ibdy-1,jmin:jmax,1:jpkm1,jn) + z2 * ptab_child(ibdy,jmin:jmax,1:jpkm1,jn) 798 DO jk = 1, jpkm1 799 DO jj = jmin,jmax 800 IF( umask(ibdy,jj,jk) == 0._wp ) THEN 801 tsa(ibdy,jj,jk,jn) = tsa(ibdy-1,jj,jk,jn) * tmask(ibdy,jj,jk) 802 ELSE 803 tsa(ibdy,jj,jk,jn)=(z4*tsa(ibdy-1,jj,jk,jn)+z3*tsa(ibdy+1,jj,jk,jn))*tmask(ibdy,jj,jk) 804 IF( un(ibdy,jj,jk) < 0._wp ) THEN 805 tsa(ibdy,jj,jk,jn)=( z6*tsa(ibdy+1,jj,jk,jn)+z5*tsa(ibdy-1,jj,jk,jn) & 806 + z7*tsa(ibdy+2,jj,jk,jn) ) * tmask(ibdy,jj,jk) 807 ENDIF 808 ENDIF 809 END DO 810 END DO 811 ! Restore ghost points: 812 tsa(ibdy-1,jmin:jmax,1:jpkm1,jn) = ptab_child(ibdy-1,jmin:jmax,1:jpkm1,jn) * tmask(ibdy-1,jmin:jmax,1:jpkm1) 813 END DO 814 ENDIF 815 ! 816 IF( southern_side ) THEN 817 zrho = Agrif_Rhoy() 818 z1 = ( zrho - 1._wp ) * 0.5_wp 819 z3 = ( zrho - 1._wp ) / ( zrho + 1._wp ) 820 z6 = 2._wp * ( zrho - 1._wp ) / ( zrho + 1._wp ) 821 z7 = - ( zrho - 1._wp ) / ( zrho + 3._wp ) 822 z2 = 1._wp - z1 ; z4 = 1._wp - z3 ; z5 = 1._wp - z6 - z7 823 ! 824 jbdy=1+nbghostcells 825 DO jn = 1, jpts 826 tsa(imin:imax,jbdy-1,1:jpkm1,jn) = z1 * ptab_child(imin:imax,jbdy-1,1:jpkm1,jn) + z2 * ptab_child(imin:imax,jbdy,1:jpkm1,jn) 827 DO jk = 1, jpkm1 828 DO ji = imin,imax 829 IF( vmask(ji,jbdy,jk) == 0._wp ) THEN 830 tsa(ji,jbdy,jk,jn)=tsa(ji,jbdy-1,jk,jn) * tmask(ji,jbdy,jk) 831 ELSE 832 tsa(ji,jbdy,jk,jn)=(z4*tsa(ji,jbdy-1,jk,jn)+z3*tsa(ji,jbdy+1,jk,jn))*tmask(ji,jbdy,jk) 833 IF( vn(ji,jbdy,jk) < 0._wp ) THEN 834 tsa(ji,jbdy,jk,jn)=( z6*tsa(ji,jbdy+1,jk,jn)+z5*tsa(ji,jbdy-1,jk,jn) & 835 + z7*tsa(ji,jbdy+2,jk,jn) ) * tmask(ji,jbdy,jk) 836 ENDIF 837 ENDIF 838 END DO 839 END DO 840 ! Restore ghost points: 841 tsa(imin:imax,jbdy-1,1:jpkm1,jn) = ptab_child(imin:imax,jbdy-1,1:jpkm1,jn) * tmask(imin:imax,jbdy-1,1:jpkm1) 842 END DO 843 ENDIF 844 ! 845 ENDIF 745 tsa(i1:i2,j1:j2,1:jpk,jn)=ptab(i1:i2,j1:j2,1:jpk,jn)*tmask(i1:i2,j1:j2,1:jpk) 746 END DO 747 # endif 748 846 749 ENDIF 847 750 ! 848 751 END SUBROUTINE interptsn 849 752 850 SUBROUTINE interpsshn( ptab, i1, i2, j1, j2, before , nb, ndir)753 SUBROUTINE interpsshn( ptab, i1, i2, j1, j2, before ) 851 754 !!---------------------------------------------------------------------- 852 755 !! *** ROUTINE interpsshn *** … … 855 758 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 856 759 LOGICAL , INTENT(in ) :: before 857 INTEGER , INTENT(in ) :: nb , ndir 858 ! 859 LOGICAL :: western_side, eastern_side,northern_side,southern_side 760 ! 860 761 !!---------------------------------------------------------------------- 861 762 ! … … 863 764 ptab(i1:i2,j1:j2) = sshn(i1:i2,j1:j2) 864 765 ELSE 865 western_side = (nb == 1).AND.(ndir == 1) 866 eastern_side = (nb == 1).AND.(ndir == 2) 867 southern_side = (nb == 2).AND.(ndir == 1) 868 northern_side = (nb == 2).AND.(ndir == 2) 869 !! clem ghost 870 IF(western_side) hbdy_w(1:nbghostcells,j1:j2) = ptab(i1:i2,j1:j2) * tmask(i1:i2,j1:j2,1) 871 IF(eastern_side) hbdy_e(1:nbghostcells,j1:j2) = ptab(i1:i2,j1:j2) * tmask(i1:i2,j1:j2,1) 872 IF(southern_side) hbdy_s(i1:i2,1:nbghostcells) = ptab(i1:i2,j1:j2) * tmask(i1:i2,j1:j2,1) 873 IF(northern_side) hbdy_n(i1:i2,1:nbghostcells) = ptab(i1:i2,j1:j2) * tmask(i1:i2,j1:j2,1) 766 hbdy(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) * tmask(i1:i2,j1:j2,1) 874 767 ENDIF 875 768 ! 876 769 END SUBROUTINE interpsshn 877 770 878 SUBROUTINE interpun( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before , nb, ndir)771 SUBROUTINE interpun( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before ) 879 772 !!---------------------------------------------------------------------- 880 773 !! *** ROUTINE interpun *** … … 884 777 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) :: ptab 885 778 LOGICAL, INTENT(in) :: before 886 INTEGER, INTENT(in) :: nb , ndir887 779 !! 888 780 INTEGER :: ji,jj,jk 889 REAL(wp) :: zrhoy 781 REAL(wp) :: zrhoy, zhtot 890 782 ! vertical interpolation: 891 783 REAL(wp), DIMENSION(k1:k2) :: tabin, h_in 892 784 REAL(wp), DIMENSION(1:jpk) :: h_out 893 INTEGER :: N_in, N_out , iref785 INTEGER :: N_in, N_out 894 786 REAL(wp) :: h_diff 895 LOGICAL :: western_side, eastern_side896 787 !!--------------------------------------------- 897 788 ! … … 902 793 ptab(ji,jj,jk,1) = (e2u(ji,jj) * e3u_n(ji,jj,jk) * un(ji,jj,jk)*umask(ji,jj,jk)) 903 794 # if defined key_vertical 904 ptab(ji,jj,jk,2) = (umask(ji,jj,jk) * e2u(ji,jj) * e3u_n(ji,jj,jk)) 795 ! Interpolate thicknesses (masked for subsequent extrapolation) 796 ptab(ji,jj,jk,2) = umask(ji,jj,jk) * e2u(ji,jj) * e3u_n(ji,jj,jk) 905 797 # endif 906 798 END DO 907 799 END DO 908 800 END DO 801 # if defined key_vertical 802 ! Extrapolate thicknesses in partial bottom cells: 803 ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 804 IF (ln_zps) THEN 805 DO jj=j1,j2 806 DO ji=i1,i2 807 jk = mbku(ji,jj) 808 ptab(ji,jj,jk,2) = 0._wp 809 END DO 810 END DO 811 END IF 812 ! Save ssh at last level: 813 ptab(i1:i2,j1:j2,k2,2) = 0._wp 814 IF (.NOT.ln_linssh) THEN 815 ! This vertical sum below should be replaced by the sea-level at U-points (optimization): 816 DO jk=1,jpk 817 ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) + e3u_n(i1:i2,j1:j2,jk) * umask(i1:i2,j1:j2,jk) 818 END DO 819 ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hu_0(i1:i2,j1:j2) 820 END IF 821 # endif 822 ! 909 823 ELSE 910 824 zrhoy = Agrif_rhoy() 911 825 # if defined key_vertical 912 826 ! VERTICAL REFINEMENT BEGIN 913 western_side = (nb == 1).AND.(ndir == 1) 914 eastern_side = (nb == 1).AND.(ndir == 2)827 828 IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 915 829 916 830 DO ji=i1,i2 917 iref = ji918 IF (western_side) iref = MAX(2,ji)919 IF (eastern_side) iref = MIN(nlci-2,ji)920 831 DO jj=j1,j2 921 N_in = 0 922 DO jk=k1,k2 923 IF (ptab(ji,jj,jk,2) == 0) EXIT 924 N_in = N_in + 1 925 tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2) 926 h_in(N_in) = ptab(ji,jj,jk,2)/(e2u(ji,jj)*zrhoy) 832 ua(ji,jj,:) = 0._wp 833 N_in = mbku_parent(ji,jj) 834 zhtot = 0._wp 835 DO jk=1,N_in 836 IF (jk==N_in) THEN 837 h_in(jk) = hu0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot 838 ELSE 839 h_in(jk) = ptab(ji,jj,jk,2)/(e2u(ji,jj)*zrhoy) 840 ENDIF 841 zhtot = zhtot + h_in(jk) 842 tabin(jk) = ptab(ji,jj,jk,1)/(e2u(ji,jj)*zrhoy*h_in(jk)) 927 843 ENDDO 928 929 IF (N_in == 0) THEN 930 ua(ji,jj,:) = 0._wp 931 CYCLE 932 ENDIF 933 844 934 845 N_out = 0 935 846 DO jk=1,jpk 936 if (umask( iref,jj,jk) == 0) EXIT847 if (umask(ji,jj,jk) == 0) EXIT 937 848 N_out = N_out + 1 938 h_out(N_out) = e3u_a( iref,jj,jk)849 h_out(N_out) = e3u_a(ji,jj,jk) 939 850 ENDDO 940 941 IF (N_out == 0) THEN 942 ua(ji,jj,:) = 0._wp 943 CYCLE 851 IF (N_in*N_out > 0) THEN 852 CALL reconstructandremap(tabin(1:N_in),h_in(1:N_in),ua(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1) 944 853 ENDIF 945 946 IF (N_in * N_out > 0) THEN947 h_diff = sum(h_out(1:N_out))-sum(h_in(1:N_in))948 ! Should be able to remove the next IF/ELSEIF statement once scale factors are dealt with properly949 if (h_diff < -1.e4) then950 print *,'CHECK YOUR BATHY ...', h_diff, sum(h_out(1:N_out)), sum(h_in(1:N_in))951 ! stop952 endif953 ENDIF954 call reconstructandremap(tabin(1:N_in),h_in(1:N_in),ua(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out)955 854 ENDDO 956 855 ENDDO … … 968 867 END SUBROUTINE interpun 969 868 970 SUBROUTINE interpvn( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before , nb, ndir)869 SUBROUTINE interpvn( ptab, i1, i2, j1, j2, k1, k2, m1, m2, before ) 971 870 !!---------------------------------------------------------------------- 972 871 !! *** ROUTINE interpvn *** … … 976 875 REAL(wp), DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) :: ptab 977 876 LOGICAL, INTENT(in) :: before 978 INTEGER, INTENT(in) :: nb , ndir979 877 ! 980 878 INTEGER :: ji,jj,jk … … 983 881 REAL(wp), DIMENSION(k1:k2) :: tabin, h_in 984 882 REAL(wp), DIMENSION(1:jpk) :: h_out 985 INTEGER :: N_in, N_out, jref 986 REAL(wp) :: h_diff 987 LOGICAL :: northern_side,southern_side 883 INTEGER :: N_in, N_out 884 REAL(wp) :: h_diff, zhtot 988 885 !!--------------------------------------------- 989 886 ! … … 994 891 ptab(ji,jj,jk,1) = (e1v(ji,jj) * e3v_n(ji,jj,jk) * vn(ji,jj,jk)*vmask(ji,jj,jk)) 995 892 # if defined key_vertical 893 ! Interpolate thicknesses (masked for subsequent extrapolation) 996 894 ptab(ji,jj,jk,2) = vmask(ji,jj,jk) * e1v(ji,jj) * e3v_n(ji,jj,jk) 997 895 # endif … … 999 897 END DO 1000 898 END DO 899 # if defined key_vertical 900 ! Extrapolate thicknesses in partial bottom cells: 901 ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 902 IF (ln_zps) THEN 903 DO jj=j1,j2 904 DO ji=i1,i2 905 jk = mbkv(ji,jj) 906 ptab(ji,jj,jk,2) = 0._wp 907 END DO 908 END DO 909 END IF 910 ! Save ssh at last level: 911 ptab(i1:i2,j1:j2,k2,2) = 0._wp 912 IF (.NOT.ln_linssh) THEN 913 ! This vertical sum below should be replaced by the sea-level at V-points (optimization): 914 DO jk=1,jpk 915 ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) + e3v_n(i1:i2,j1:j2,jk) * vmask(i1:i2,j1:j2,jk) 916 END DO 917 ptab(i1:i2,j1:j2,k2,2) = ptab(i1:i2,j1:j2,k2,2) - hv_0(i1:i2,j1:j2) 918 END IF 919 # endif 1001 920 ELSE 1002 921 zrhox = Agrif_rhox() 1003 922 # if defined key_vertical 1004 923 1005 southern_side = (nb == 2).AND.(ndir == 1) 1006 northern_side = (nb == 2).AND.(ndir == 2) 924 IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 1007 925 1008 926 DO jj=j1,j2 1009 jref = jj1010 IF (southern_side) jref = MAX(2,jj)1011 IF (northern_side) jref = MIN(nlcj-2,jj)1012 927 DO ji=i1,i2 1013 N_in = 0 1014 DO jk=k1,k2 1015 if (ptab(ji,jj,jk,2) == 0) EXIT 1016 N_in = N_in + 1 1017 tabin(jk) = ptab(ji,jj,jk,1)/ptab(ji,jj,jk,2) 1018 h_in(N_in) = ptab(ji,jj,jk,2)/(e1v(ji,jj)*zrhox) 1019 END DO 1020 IF (N_in == 0) THEN 1021 va(ji,jj,:) = 0._wp 1022 CYCLE 1023 ENDIF 928 va(ji,jj,:) = 0._wp 929 N_in = mbkv_parent(ji,jj) 930 zhtot = 0._wp 931 DO jk=1,N_in 932 IF (jk==N_in) THEN 933 h_in(jk) = hv0_parent(ji,jj) + ptab(ji,jj,k2,2) - zhtot 934 ELSE 935 h_in(jk) = ptab(ji,jj,jk,2)/(e1v(ji,jj)*zrhox) 936 ENDIF 937 zhtot = zhtot + h_in(jk) 938 tabin(jk) = ptab(ji,jj,jk,1)/(e1v(ji,jj)*zrhox*h_in(jk)) 939 ENDDO 1024 940 1025 941 N_out = 0 1026 942 DO jk=1,jpk 1027 if (vmask(ji,j ref,jk) == 0) EXIT943 if (vmask(ji,jj,jk) == 0) EXIT 1028 944 N_out = N_out + 1 1029 h_out(N_out) = e3v_a(ji,jref,jk) 1030 END DO 1031 IF (N_out == 0) THEN 1032 va(ji,jj,:) = 0._wp 1033 CYCLE 945 h_out(N_out) = e3v_a(ji,jj,jk) 946 END DO 947 IF (N_in*N_out > 0) THEN 948 call reconstructandremap(tabin(1:N_in),h_in(1:N_in),va(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out,1) 1034 949 ENDIF 1035 call reconstructandremap(tabin(1:N_in),h_in(1:N_in),va(ji,jj,1:N_out),h_out(1:N_out),N_in,N_out)1036 950 END DO 1037 951 END DO … … 1045 959 END SUBROUTINE interpvn 1046 960 1047 SUBROUTINE interpunb( ptab, i1, i2, j1, j2, before , nb, ndir)961 SUBROUTINE interpunb( ptab, i1, i2, j1, j2, before) 1048 962 !!---------------------------------------------------------------------- 1049 963 !! *** ROUTINE interpunb *** … … 1052 966 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1053 967 LOGICAL , INTENT(in ) :: before 1054 INTEGER , INTENT(in ) :: nb , ndir1055 968 ! 1056 969 INTEGER :: ji, jj 1057 970 REAL(wp) :: zrhoy, zrhot, zt0, zt1, ztcoeff 1058 LOGICAL :: western_side, eastern_side,northern_side,southern_side1059 971 !!---------------------------------------------------------------------- 1060 972 ! … … 1062 974 ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * hu_n(i1:i2,j1:j2) * un_b(i1:i2,j1:j2) 1063 975 ELSE 1064 western_side = (nb == 1).AND.(ndir == 1)1065 eastern_side = (nb == 1).AND.(ndir == 2)1066 southern_side = (nb == 2).AND.(ndir == 1)1067 northern_side = (nb == 2).AND.(ndir == 2)1068 976 zrhoy = Agrif_Rhoy() 1069 977 zrhot = Agrif_rhot() … … 1071 979 zt0 = REAL(Agrif_NbStepint() , wp) / zrhot 1072 980 zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot 1073 ! Polynomial interpolation coefficients: 1074 IF( bdy_tinterp == 1 ) THEN 1075 ztcoeff = zrhot * ( zt1**2._wp * ( zt1 - 1._wp) & 1076 & - zt0**2._wp * ( zt0 - 1._wp) ) 1077 ELSEIF( bdy_tinterp == 2 ) THEN 1078 ztcoeff = zrhot * ( zt1 * ( zt1 - 1._wp)**2._wp & 1079 & - zt0 * ( zt0 - 1._wp)**2._wp ) 1080 ELSE 1081 ztcoeff = 1 1082 ENDIF 1083 ! 1084 IF(western_side) ubdy_w(1:nbghostcells,j1:j2) = ubdy_w(1:nbghostcells,j1:j2) + ztcoeff * ptab(i1:i2,j1:j2) 1085 IF(eastern_side) ubdy_e(1:nbghostcells,j1:j2) = ubdy_e(1:nbghostcells,j1:j2) + ztcoeff * ptab(i1:i2,j1:j2) 1086 IF(southern_side) ubdy_s(i1:i2,1:nbghostcells) = ubdy_s(i1:i2,1:nbghostcells) + ztcoeff * ptab(i1:i2,j1:j2) 1087 IF(northern_side) ubdy_n(i1:i2,1:nbghostcells) = ubdy_n(i1:i2,1:nbghostcells) + ztcoeff * ptab(i1:i2,j1:j2) 1088 ! 1089 IF( bdy_tinterp == 0 .OR. bdy_tinterp == 2) THEN 1090 IF(western_side) ubdy_w(1:nbghostcells,j1:j2) = ubdy_w(1:nbghostcells,j1:j2) / (zrhoy*e2u(i1:i2,j1:j2)) * umask(i1:i2,j1:j2,1) 1091 IF(eastern_side) ubdy_e(1:nbghostcells,j1:j2) = ubdy_e(1:nbghostcells,j1:j2) / (zrhoy*e2u(i1:i2,j1:j2)) * umask(i1:i2,j1:j2,1) 1092 IF(southern_side) ubdy_s(i1:i2,1:nbghostcells) = ubdy_s(i1:i2,1:nbghostcells) / (zrhoy*e2u(i1:i2,j1:j2)) * umask(i1:i2,j1:j2,1) 1093 IF(northern_side) ubdy_n(i1:i2,1:nbghostcells) = ubdy_n(i1:i2,1:nbghostcells) / (zrhoy*e2u(i1:i2,j1:j2)) * umask(i1:i2,j1:j2,1) 1094 ENDIF 1095 ENDIF 981 ! 982 DO ji = i1, i2 983 DO jj = j1, j2 984 IF ( utint_stage(ji,jj) < (bdy_tinterp + 1) ) THEN 985 IF ( utint_stage(ji,jj) == 1 ) THEN 986 ztcoeff = zrhot * ( zt1**2._wp * ( zt1 - 1._wp) & 987 & - zt0**2._wp * ( zt0 - 1._wp) ) 988 ELSEIF( utint_stage(ji,jj) == 2 ) THEN 989 ztcoeff = zrhot * ( zt1 * ( zt1 - 1._wp)**2._wp & 990 & - zt0 * ( zt0 - 1._wp)**2._wp ) 991 ELSEIF( utint_stage(ji,jj) == 0 ) THEN 992 ztcoeff = 1._wp 993 ELSE 994 ztcoeff = 0._wp 995 ENDIF 996 ! 997 ubdy(ji,jj) = ubdy(ji,jj) + ztcoeff * ptab(ji,jj) 998 ! 999 IF (( utint_stage(ji,jj) == 2 ).OR.( utint_stage(ji,jj) == 0 )) THEN 1000 ubdy(ji,jj) = ubdy(ji,jj) / (zrhoy*e2u(ji,jj)) * umask(ji,jj,1) 1001 ENDIF 1002 ! 1003 utint_stage(ji,jj) = utint_stage(ji,jj) + 1 1004 ENDIF 1005 END DO 1006 END DO 1007 END IF 1096 1008 ! 1097 1009 END SUBROUTINE interpunb 1098 1010 1099 1011 1100 SUBROUTINE interpvnb( ptab, i1, i2, j1, j2, before , nb, ndir)1012 SUBROUTINE interpvnb( ptab, i1, i2, j1, j2, before ) 1101 1013 !!---------------------------------------------------------------------- 1102 1014 !! *** ROUTINE interpvnb *** … … 1105 1017 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1106 1018 LOGICAL , INTENT(in ) :: before 1107 INTEGER , INTENT(in ) :: nb , ndir 1108 ! 1109 INTEGER :: ji,jj 1019 ! 1020 INTEGER :: ji, jj 1110 1021 REAL(wp) :: zrhox, zrhot, zt0, zt1, ztcoeff 1111 LOGICAL :: western_side, eastern_side,northern_side,southern_side1112 1022 !!---------------------------------------------------------------------- 1113 1023 ! … … 1115 1025 ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * hv_n(i1:i2,j1:j2) * vn_b(i1:i2,j1:j2) 1116 1026 ELSE 1117 western_side = (nb == 1).AND.(ndir == 1)1118 eastern_side = (nb == 1).AND.(ndir == 2)1119 southern_side = (nb == 2).AND.(ndir == 1)1120 northern_side = (nb == 2).AND.(ndir == 2)1121 1027 zrhox = Agrif_Rhox() 1122 1028 zrhot = Agrif_rhot() 1123 1029 ! Time indexes bounds for integration 1124 1030 zt0 = REAL(Agrif_NbStepint() , wp) / zrhot 1125 zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot 1126 IF( bdy_tinterp == 1 ) THEN 1127 ztcoeff = zrhot * ( zt1**2._wp * ( zt1 - 1._wp) & 1128 & - zt0**2._wp * ( zt0 - 1._wp) ) 1129 ELSEIF( bdy_tinterp == 2 ) THEN 1130 ztcoeff = zrhot * ( zt1 * ( zt1 - 1._wp)**2._wp & 1131 & - zt0 * ( zt0 - 1._wp)**2._wp ) 1132 ELSE 1133 ztcoeff = 1 1134 ENDIF 1135 !! clem ghost 1136 IF(western_side) vbdy_w(1:nbghostcells,j1:j2) = vbdy_w(1:nbghostcells,j1:j2) + ztcoeff * ptab(i1:i2,j1:j2) 1137 IF(eastern_side) vbdy_e(1:nbghostcells,j1:j2) = vbdy_e(1:nbghostcells,j1:j2) + ztcoeff * ptab(i1:i2,j1:j2) 1138 IF(southern_side) vbdy_s(i1:i2,1:nbghostcells) = vbdy_s(i1:i2,1:nbghostcells) + ztcoeff * ptab(i1:i2,j1:j2) 1139 IF(northern_side) vbdy_n(i1:i2,1:nbghostcells) = vbdy_n(i1:i2,1:nbghostcells) + ztcoeff * ptab(i1:i2,j1:j2) 1140 ! 1141 IF( bdy_tinterp == 0 .OR. bdy_tinterp == 2) THEN 1142 IF(western_side) vbdy_w(1:nbghostcells,j1:j2) = vbdy_w(1:nbghostcells,j1:j2) / (zrhox*e1v(i1:i2,j1:j2)) * vmask(i1:i2,j1:j2,1) 1143 IF(eastern_side) vbdy_e(1:nbghostcells,j1:j2) = vbdy_e(1:nbghostcells,j1:j2) / (zrhox*e1v(i1:i2,j1:j2)) * vmask(i1:i2,j1:j2,1) 1144 IF(southern_side) vbdy_s(i1:i2,1:nbghostcells) = vbdy_s(i1:i2,1:nbghostcells) / (zrhox*e1v(i1:i2,j1:j2)) * vmask(i1:i2,j1:j2,1) 1145 IF(northern_side) vbdy_n(i1:i2,1:nbghostcells) = vbdy_n(i1:i2,1:nbghostcells) / (zrhox*e1v(i1:i2,j1:j2)) * vmask(i1:i2,j1:j2,1) 1146 ENDIF 1031 zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot 1032 ! 1033 DO ji = i1, i2 1034 DO jj = j1, j2 1035 IF ( vtint_stage(ji,jj) < (bdy_tinterp + 1) ) THEN 1036 IF ( vtint_stage(ji,jj) == 1 ) THEN 1037 ztcoeff = zrhot * ( zt1**2._wp * ( zt1 - 1._wp) & 1038 & - zt0**2._wp * ( zt0 - 1._wp) ) 1039 ELSEIF( vtint_stage(ji,jj) == 2 ) THEN 1040 ztcoeff = zrhot * ( zt1 * ( zt1 - 1._wp)**2._wp & 1041 & - zt0 * ( zt0 - 1._wp)**2._wp ) 1042 ELSEIF( vtint_stage(ji,jj) == 0 ) THEN 1043 ztcoeff = 1._wp 1044 ELSE 1045 ztcoeff = 0._wp 1046 ENDIF 1047 ! 1048 vbdy(ji,jj) = vbdy(ji,jj) + ztcoeff * ptab(ji,jj) 1049 ! 1050 IF (( vtint_stage(ji,jj) == 2 ).OR.( vtint_stage(ji,jj) == 0 )) THEN 1051 vbdy(ji,jj) = vbdy(ji,jj) / (zrhox*e1v(ji,jj)) * vmask(ji,jj,1) 1052 ENDIF 1053 ! 1054 vtint_stage(ji,jj) = vtint_stage(ji,jj) + 1 1055 ENDIF 1056 END DO 1057 END DO 1147 1058 ENDIF 1148 1059 ! … … 1150 1061 1151 1062 1152 SUBROUTINE interpub2b( ptab, i1, i2, j1, j2, before , nb, ndir)1063 SUBROUTINE interpub2b( ptab, i1, i2, j1, j2, before ) 1153 1064 !!---------------------------------------------------------------------- 1154 1065 !! *** ROUTINE interpub2b *** … … 1157 1068 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1158 1069 LOGICAL , INTENT(in ) :: before 1159 INTEGER , INTENT(in ) :: nb , ndir1160 1070 ! 1161 1071 INTEGER :: ji,jj 1162 REAL(wp) :: zrhot, zt0, zt1,zat 1163 LOGICAL :: western_side, eastern_side,northern_side,southern_side 1072 REAL(wp) :: zrhot, zt0, zt1, zat 1164 1073 !!---------------------------------------------------------------------- 1165 1074 IF( before ) THEN … … 1170 1079 ENDIF 1171 1080 ELSE 1172 western_side = (nb == 1).AND.(ndir == 1)1173 eastern_side = (nb == 1).AND.(ndir == 2)1174 southern_side = (nb == 2).AND.(ndir == 1)1175 northern_side = (nb == 2).AND.(ndir == 2)1176 zrhot = Agrif_rhot()1177 ! Time indexes bounds for integration1178 zt0 = REAL(Agrif_NbStepint() , wp) / zrhot1179 zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot1180 ! Polynomial interpolation coefficients:1181 zat = zrhot * ( zt1**2._wp * (-2._wp*zt1 + 3._wp) &1182 & - zt0**2._wp * (-2._wp*zt0 + 3._wp) )1183 !! clem ghost1184 IF(western_side ) ubdy_w(1:nbghostcells,j1:j2) = zat * ptab(i1:i2,j1:j2)1185 IF(eastern_side ) ubdy_e(1:nbghostcells,j1:j2) = zat * ptab(i1:i2,j1:j2)1186 IF(southern_side) ubdy_s(i1:i2,1:nbghostcells) = zat * ptab(i1:i2,j1:j2)1187 IF(northern_side) ubdy_n(i1:i2,1:nbghostcells) = zat * ptab(i1:i2,j1:j2)1188 ENDIF1189 !1190 END SUBROUTINE interpub2b1191 1192 1193 SUBROUTINE interpvb2b( ptab, i1, i2, j1, j2, before, nb, ndir )1194 !!----------------------------------------------------------------------1195 !! *** ROUTINE interpvb2b ***1196 !!----------------------------------------------------------------------1197 INTEGER , INTENT(in ) :: i1, i2, j1, j21198 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab1199 LOGICAL , INTENT(in ) :: before1200 INTEGER , INTENT(in ) :: nb , ndir1201 !1202 INTEGER :: ji,jj1203 REAL(wp) :: zrhot, zt0, zt1,zat1204 LOGICAL :: western_side, eastern_side,northern_side,southern_side1205 !!----------------------------------------------------------------------1206 !1207 IF( before ) THEN1208 IF ( ln_bt_fw ) THEN1209 ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2)1210 ELSE1211 ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2)1212 ENDIF1213 ELSE1214 western_side = (nb == 1).AND.(ndir == 1)1215 eastern_side = (nb == 1).AND.(ndir == 2)1216 southern_side = (nb == 2).AND.(ndir == 1)1217 northern_side = (nb == 2).AND.(ndir == 2)1218 1081 zrhot = Agrif_rhot() 1219 1082 ! Time indexes bounds for integration … … 1224 1087 & - zt0**2._wp * (-2._wp*zt0 + 3._wp) ) 1225 1088 ! 1226 IF(western_side ) vbdy_w(1:nbghostcells,j1:j2) = zat * ptab(i1:i2,j1:j2) 1227 IF(eastern_side ) vbdy_e(1:nbghostcells,j1:j2) = zat * ptab(i1:i2,j1:j2) 1228 IF(southern_side) vbdy_s(i1:i2,1:nbghostcells) = zat * ptab(i1:i2,j1:j2) 1229 IF(northern_side) vbdy_n(i1:i2,1:nbghostcells) = zat * ptab(i1:i2,j1:j2) 1089 ubdy(i1:i2,j1:j2) = zat * ptab(i1:i2,j1:j2) 1090 ! 1091 ! Update interpolation stage: 1092 utint_stage(i1:i2,j1:j2) = 1 1093 ENDIF 1094 ! 1095 END SUBROUTINE interpub2b 1096 1097 1098 SUBROUTINE interpvb2b( ptab, i1, i2, j1, j2, before ) 1099 !!---------------------------------------------------------------------- 1100 !! *** ROUTINE interpvb2b *** 1101 !!---------------------------------------------------------------------- 1102 INTEGER , INTENT(in ) :: i1, i2, j1, j2 1103 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1104 LOGICAL , INTENT(in ) :: before 1105 ! 1106 INTEGER :: ji,jj 1107 REAL(wp) :: zrhot, zt0, zt1, zat 1108 !!---------------------------------------------------------------------- 1109 ! 1110 IF( before ) THEN 1111 IF ( ln_bt_fw ) THEN 1112 ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2) 1113 ELSE 1114 ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2) 1115 ENDIF 1116 ELSE 1117 zrhot = Agrif_rhot() 1118 ! Time indexes bounds for integration 1119 zt0 = REAL(Agrif_NbStepint() , wp) / zrhot 1120 zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot 1121 ! Polynomial interpolation coefficients: 1122 zat = zrhot * ( zt1**2._wp * (-2._wp*zt1 + 3._wp) & 1123 & - zt0**2._wp * (-2._wp*zt0 + 3._wp) ) 1124 ! 1125 vbdy(i1:i2,j1:j2) = zat * ptab(i1:i2,j1:j2) 1126 ! 1127 ! update interpolation stage: 1128 vtint_stage(i1:i2,j1:j2) = 1 1230 1129 ENDIF 1231 1130 ! … … 1233 1132 1234 1133 1235 SUBROUTINE interpe3t( ptab, i1, i2, j1, j2, k1, k2, before , nb, ndir)1134 SUBROUTINE interpe3t( ptab, i1, i2, j1, j2, k1, k2, before ) 1236 1135 !!---------------------------------------------------------------------- 1237 1136 !! *** ROUTINE interpe3t *** … … 1240 1139 REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab 1241 1140 LOGICAL , INTENT(in ) :: before 1242 INTEGER , INTENT(in ) :: nb , ndir1243 1141 ! 1244 1142 INTEGER :: ji, jj, jk 1245 LOGICAL :: western_side, eastern_side, northern_side, southern_side1246 1143 !!---------------------------------------------------------------------- 1247 1144 ! … … 1249 1146 ptab(i1:i2,j1:j2,k1:k2) = tmask(i1:i2,j1:j2,k1:k2) * e3t_0(i1:i2,j1:j2,k1:k2) 1250 1147 ELSE 1251 western_side = (nb == 1).AND.(ndir == 1)1252 eastern_side = (nb == 1).AND.(ndir == 2)1253 southern_side = (nb == 2).AND.(ndir == 1)1254 northern_side = (nb == 2).AND.(ndir == 2)1255 1148 ! 1256 1149 DO jk = k1, k2 1257 1150 DO jj = j1, j2 1258 1151 DO ji = i1, i2 1259 !1260 1152 IF( ABS( ptab(ji,jj,jk) - tmask(ji,jj,jk) * e3t_0(ji,jj,jk) ) > 1.D-2) THEN 1261 IF (western_side.AND.(ptab(i1+nbghostcells-1,jj,jk)>0._wp)) THEN 1262 WRITE(numout,*) 'ERROR bathymetry merge at the western border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk 1263 WRITE(numout,*) ptab(ji,jj,jk), e3t_0(ji,jj,jk) 1264 kindic_agr = kindic_agr + 1 1265 ELSEIF (eastern_side.AND.(ptab(i2-nbghostcells+1,jj,jk)>0._wp)) THEN 1266 WRITE(numout,*) 'ERROR bathymetry merge at the eastern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk 1267 WRITE(numout,*) ptab(ji,jj,jk), e3t_0(ji,jj,jk) 1268 kindic_agr = kindic_agr + 1 1269 ELSEIF (southern_side.AND.(ptab(ji,j1+nbghostcells-1,jk)>0._wp)) THEN 1270 WRITE(numout,*) 'ERROR bathymetry merge at the southern border ji,jj,jk', ji+nimpp-1,jj+njmpp-1,jk 1271 WRITE(numout,*) ptab(ji,jj,jk), e3t_0(ji,jj,jk) 1272 kindic_agr = kindic_agr + 1 1273 ELSEIF (northern_side.AND.(ptab(ji,j2-nbghostcells+1,jk)>0._wp)) THEN 1274 WRITE(numout,*) 'ERROR bathymetry merge at the northen border ji,jj,jk', ji+nimpp-1,jj+njmpp-1,jk 1275 WRITE(numout,*) ptab(ji,jj,jk), e3t_0(ji,jj,jk) 1276 kindic_agr = kindic_agr + 1 1277 ENDIF 1153 WRITE(numout,*) ' Agrif error for e3t_0: parent , child, i, j, k ', & 1154 & ptab(ji,jj,jk), tmask(ji,jj,jk) * e3t_0(ji,jj,jk), & 1155 & ji+nimpp-1, jj+njmpp-1, jk 1156 kindic_agr = kindic_agr + 1 1278 1157 ENDIF 1279 1158 END DO … … 1284 1163 ! 1285 1164 END SUBROUTINE interpe3t 1286 1287 1288 SUBROUTINE interpumsk( ptab, i1, i2, j1, j2, k1, k2, before, nb, ndir )1289 !!----------------------------------------------------------------------1290 !! *** ROUTINE interpumsk ***1291 !!----------------------------------------------------------------------1292 INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k21293 REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab1294 LOGICAL , INTENT(in ) :: before1295 INTEGER , INTENT(in ) :: nb , ndir1296 !1297 INTEGER :: ji, jj, jk1298 LOGICAL :: western_side, eastern_side1299 !!----------------------------------------------------------------------1300 !1301 IF( before ) THEN1302 ptab(i1:i2,j1:j2,k1:k2) = umask(i1:i2,j1:j2,k1:k2)1303 ELSE1304 western_side = (nb == 1).AND.(ndir == 1)1305 eastern_side = (nb == 1).AND.(ndir == 2)1306 DO jk = k1, k21307 DO jj = j1, j21308 DO ji = i1, i21309 ! Velocity mask at boundary edge points:1310 IF (ABS(ptab(ji,jj,jk) - umask(ji,jj,jk)) > 1.D-2) THEN1311 IF (western_side) THEN1312 WRITE(numout,*) 'ERROR with umask at the western border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk1313 WRITE(numout,*) ' masks: parent, child ', ptab(ji,jj,jk), umask(ji,jj,jk)1314 kindic_agr = kindic_agr + 11315 ELSEIF (eastern_side) THEN1316 WRITE(numout,*) 'ERROR with umask at the eastern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk1317 WRITE(numout,*) ' masks: parent, child ', ptab(ji,jj,jk), umask(ji,jj,jk)1318 kindic_agr = kindic_agr + 11319 ENDIF1320 ENDIF1321 END DO1322 END DO1323 END DO1324 !1325 ENDIF1326 !1327 END SUBROUTINE interpumsk1328 1329 1330 SUBROUTINE interpvmsk( ptab, i1, i2, j1, j2, k1, k2, before, nb, ndir )1331 !!----------------------------------------------------------------------1332 !! *** ROUTINE interpvmsk ***1333 !!----------------------------------------------------------------------1334 INTEGER , INTENT(in ) :: i1,i2,j1,j2,k1,k21335 REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2), INTENT(inout) :: ptab1336 LOGICAL , INTENT(in ) :: before1337 INTEGER , INTENT(in ) :: nb , ndir1338 !1339 INTEGER :: ji, jj, jk1340 LOGICAL :: northern_side, southern_side1341 !!----------------------------------------------------------------------1342 !1343 IF( before ) THEN1344 ptab(i1:i2,j1:j2,k1:k2) = vmask(i1:i2,j1:j2,k1:k2)1345 ELSE1346 southern_side = (nb == 2).AND.(ndir == 1)1347 northern_side = (nb == 2).AND.(ndir == 2)1348 DO jk = k1, k21349 DO jj = j1, j21350 DO ji = i1, i21351 ! Velocity mask at boundary edge points:1352 IF (ABS(ptab(ji,jj,jk) - vmask(ji,jj,jk)) > 1.D-2) THEN1353 IF (southern_side) THEN1354 WRITE(numout,*) 'ERROR with vmask at the southern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk1355 WRITE(numout,*) ' masks: parent, child ', ptab(ji,jj,jk), vmask(ji,jj,jk)1356 kindic_agr = kindic_agr + 11357 ELSEIF (northern_side) THEN1358 WRITE(numout,*) 'ERROR with vmask at the northern border ji,jj,jk ', ji+nimpp-1,jj+njmpp-1,jk1359 WRITE(numout,*) ' masks: parent, child ', ptab(ji,jj,jk), vmask(ji,jj,jk)1360 kindic_agr = kindic_agr + 11361 ENDIF1362 ENDIF1363 END DO1364 END DO1365 END DO1366 !1367 ENDIF1368 !1369 END SUBROUTINE interpvmsk1370 1165 1371 1166 … … 1377 1172 REAL(wp),DIMENSION(i1:i2,j1:j2,k1:k2,m1:m2), INTENT(inout) :: ptab 1378 1173 LOGICAL , INTENT(in ) :: before 1379 REAL(wp), DIMENSION(k1:k2) :: tabin, h_in 1380 REAL(wp), DIMENSION(1:jpk) :: h_out 1381 INTEGER :: N_in, N_out, ji, jj, jk 1174 ! 1175 INTEGER :: ji, jj, jk 1176 INTEGER :: N_in, N_out 1177 REAL(wp), DIMENSION(k1:k2) :: tabin, z_in 1178 REAL(wp), DIMENSION(1:jpk) :: z_out 1382 1179 !!---------------------------------------------------------------------- 1383 1180 ! … … 1390 1187 END DO 1391 1188 END DO 1392 #ifdef key_vertical 1189 1190 # if defined key_vertical 1191 ! Interpolate thicknesses 1192 ! Warning: these are masked, hence extrapolated prior interpolation. 1393 1193 DO jk=k1,k2 1394 1194 DO jj=j1,j2 1395 1195 DO ji=i1,i2 1396 ptab(ji,jj,jk,2) = wmask(ji,jj,jk) * e3w_n(ji,jj,jk)1196 ptab(ji,jj,jk,2) = tmask(ji,jj,jk) * e3t_n(ji,jj,jk) 1397 1197 END DO 1398 1198 END DO 1399 1199 END DO 1400 #endif 1200 1201 ! Extrapolate thicknesses in partial bottom cells: 1202 ! Set them to Agrif_SpecialValue (0.). Correct bottom thicknesses are retrieved later on 1203 IF (ln_zps) THEN 1204 DO jj=j1,j2 1205 DO ji=i1,i2 1206 jk = mbkt(ji,jj) 1207 ptab(ji,jj,jk,2) = 0._wp 1208 END DO 1209 END DO 1210 END IF 1211 1212 ! Save ssh at last level: 1213 IF (.NOT.ln_linssh) THEN 1214 ptab(i1:i2,j1:j2,k2,2) = sshn(i1:i2,j1:j2)*tmask(i1:i2,j1:j2,1) 1215 ELSE 1216 ptab(i1:i2,j1:j2,k2,2) = 0._wp 1217 END IF 1218 # endif 1401 1219 ELSE 1402 1220 #ifdef key_vertical 1403 avm_k(i1:i2,j1:j2,1:jpk) = 0. 1404 DO jj=j1,j2 1405 DO ji=i1,i2 1406 N_in = 0 1407 DO jk=k1,k2 !k2 = jpk of parent grid 1408 IF (ptab(ji,jj,jk,2) == 0) EXIT 1409 N_in = N_in + 1 1410 tabin(jk) = ptab(ji,jj,jk,1) 1411 h_in(N_in) = ptab(ji,jj,jk,2) 1221 IF (ln_linssh) ptab(i1:i2,j1:j2,k2,2) = 0._wp 1222 avm_k(i1:i2,j1:j2,k1:k2) = 0._wp 1223 1224 DO jj = j1, j2 1225 DO ji =i1, i2 1226 N_in = mbkt_parent(ji,jj) 1227 IF ( tmask(ji,jj,1) == 0._wp) N_in = 0 1228 z_in(N_in+1) = ht0_parent(ji,jj) + ptab(ji,jj,k2,2) 1229 DO jk = N_in, 1, -1 ! Parent vertical grid 1230 z_in(jk) = z_in(jk+1) - ptab(ji,jj,jk,2) 1231 tabin(jk) = ptab(ji,jj,jk,1) 1412 1232 END DO 1413 1233 N_out = 0 1414 DO jk=1,jpk ! jpk of child grid 1415 IF (wmask(ji,jj,jk) == 0) EXIT 1234 z_out(1) = 0._wp 1235 DO jk = 2, jpk ! Child vertical grid 1236 IF (tmask(ji,jj,jk) == 0._wp) EXIT 1416 1237 N_out = N_out + 1 1417 h_out(jk) = e3t_n(ji,jj,jk)1238 z_out(jk) = z_out(jk-1) + e3t_n(ji,jj,jk-1) 1418 1239 ENDDO 1419 IF (N_in > 0) THEN1420 CALL re constructandremap(tabin(1:N_in),h_in,avm_k(ji,jj,1:N_out),h_out,N_in,N_out)1240 IF (N_in*N_out > 0) THEN 1241 CALL remap_linear(tabin(1:N_in),z_in(1:N_in),avm_k(ji,jj,1:N_out),z_out(1:N_out),N_in,N_out,1) 1421 1242 ENDIF 1422 1243 ENDDO … … 1428 1249 ! 1429 1250 END SUBROUTINE interpavm 1251 1252 # if defined key_vertical 1253 SUBROUTINE interpmbkt( ptab, i1, i2, j1, j2, before ) 1254 !!---------------------------------------------------------------------- 1255 !! *** ROUTINE interpsshn *** 1256 !!---------------------------------------------------------------------- 1257 INTEGER , INTENT(in ) :: i1, i2, j1, j2 1258 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1259 LOGICAL , INTENT(in ) :: before 1260 ! 1261 !!---------------------------------------------------------------------- 1262 ! 1263 IF( before) THEN 1264 ptab(i1:i2,j1:j2) = REAL(mbkt(i1:i2,j1:j2),wp) 1265 ELSE 1266 mbkt_parent(i1:i2,j1:j2) = NINT(ptab(i1:i2,j1:j2)) 1267 ENDIF 1268 ! 1269 END SUBROUTINE interpmbkt 1270 1271 SUBROUTINE interpht0( ptab, i1, i2, j1, j2, before ) 1272 !!---------------------------------------------------------------------- 1273 !! *** ROUTINE interpsshn *** 1274 !!---------------------------------------------------------------------- 1275 INTEGER , INTENT(in ) :: i1, i2, j1, j2 1276 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1277 LOGICAL , INTENT(in ) :: before 1278 ! 1279 !!---------------------------------------------------------------------- 1280 ! 1281 IF( before) THEN 1282 ptab(i1:i2,j1:j2) = ht_0(i1:i2,j1:j2) 1283 ELSE 1284 ht0_parent(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) 1285 ENDIF 1286 ! 1287 END SUBROUTINE interpht0 1288 #endif 1430 1289 1431 1290 #else
Note: See TracChangeset
for help on using the changeset viewer.