Changeset 11574 for NEMO/branches/2019/dev_r11233_AGRIF-05_jchanut_vert_coord_interp/src/NST/agrif_oce_interp.F90
- Timestamp:
- 2019-09-19T12:08:31+02:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11233_AGRIF-05_jchanut_vert_coord_interp/src/NST/agrif_oce_interp.F90
r10068 r11574 37 37 PRIVATE 38 38 39 PUBLIC Agrif_dyn, Agrif_ssh, Agrif_dyn_ts, Agrif_ ssh_ts, Agrif_dta_ts39 PUBLIC Agrif_dyn, Agrif_ssh, Agrif_dyn_ts, Agrif_dyn_ts_flux, Agrif_ssh_ts, Agrif_dta_ts 40 40 PUBLIC Agrif_tra, Agrif_avm 41 41 PUBLIC interpun , interpvn … … 43 43 PUBLIC interpunb, interpvnb , interpub2b, interpvb2b 44 44 PUBLIC interpe3t, interpumsk, interpvmsk 45 46 INTEGER :: bdy_tinterp = 047 45 48 46 # include "vectopt_loop_substitute.h90" … … 78 76 ! 79 77 INTEGER :: ji, jj, jk ! dummy loop indices 80 INTEGER :: j1, j2, i1, i281 78 INTEGER :: ibdy1, jbdy1, ibdy2, jbdy2 82 79 REAL(wp), DIMENSION(jpi,jpj) :: zub, zvb … … 93 90 Agrif_UseSpecialValue = .FALSE. 94 91 ! 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 92 ! --- 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 93 ibdy1 = 2 94 ibdy2 = 1+nbghostcells 95 ! 96 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport 97 DO ji = mi0(ibdy1), mi1(ibdy2) 98 ua_b(ji,:) = 0._wp 99 110 100 DO jk = 1, jpkm1 111 101 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 102 ua_b(ji,jj) = ua_b(ji,jj) + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 103 END DO 104 END DO 105 116 106 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 107 ua_b(ji,jj) = ua_b(ji,jj) * r1_hu_a(ji,jj) 108 END DO 109 END DO 110 ENDIF 111 ! 112 DO ji = mi0(ibdy1), mi1(ibdy2) 113 zub(ji,:) = 0._wp ! Correct transport 130 114 DO jk = 1, jpkm1 131 115 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)116 zub(ji,jj) = zub(ji,jj) & 117 & + e3u_a(ji,jj,jk) * ua(ji,jj,jk)*umask(ji,jj,jk) 134 118 END DO 135 119 END DO 136 120 DO jj=1,jpj 137 zub( ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu_a(ibdy1:ibdy2,jj)121 zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj) 138 122 END DO 139 123 140 124 DO jk = 1, jpkm1 141 125 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 126 ua(ji,jj,jk) = ( ua(ji,jj,jk) + ua_b(ji,jj)-zub(ji,jj)) * umask(ji,jj,jk) 127 END DO 128 END DO 129 END DO 146 130 147 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 148 zvb(ibdy1:ibdy2,:) = 0._wp 131 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 132 DO ji = mi0(ibdy1), mi1(ibdy2) 133 zvb(ji,:) = 0._wp 149 134 DO jk = 1, jpkm1 150 135 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) 136 zvb(ji,jj) = zvb(ji,jj) + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 153 137 END DO 154 138 END DO 155 139 DO jj = 1, jpj 156 zvb( ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv_a(ibdy1:ibdy2,jj)140 zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj) 157 141 END DO 158 142 DO jk = 1, jpkm1 159 143 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 144 va(ji,jj,jk) = ( va(ji,jj,jk) + va_b(ji,jj)-zvb(ji,jj))*vmask(ji,jj,jk) 145 END DO 146 END DO 147 END DO 172 148 ENDIF 173 149 174 150 ! --- 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._wp151 ibdy1 = jpiglo-1-nbghostcells 152 ibdy2 = jpiglo-2 153 ! 154 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport 155 DO ji = mi0(ibdy1), mi1(ibdy2) 156 ua_b(ji,:) = 0._wp 181 157 DO jk = 1, jpkm1 182 158 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)159 ua_b(ji,jj) = ua_b(ji,jj) & 160 & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 185 161 END DO 186 162 END DO 187 163 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 164 ua_b(ji,jj) = ua_b(ji,jj) * r1_hu_a(ji,jj) 165 END DO 166 END DO 167 ENDIF 168 ! 169 DO ji = mi0(ibdy1), mi1(ibdy2) 170 zub(ji,:) = 0._wp ! Correct transport 201 171 DO jk = 1, jpkm1 202 172 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)173 zub(ji,jj) = zub(ji,jj) & 174 & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 205 175 END DO 206 176 END DO 207 177 DO jj=1,jpj 208 zub( ibdy1:ibdy2,jj) = zub(ibdy1:ibdy2,jj) * r1_hu_a(ibdy1:ibdy2,jj)178 zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj) 209 179 END DO 210 180 211 181 DO jk = 1, jpkm1 212 182 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 183 ua(ji,jj,jk) = ( ua(ji,jj,jk) & 184 & + ua_b(ji,jj)-zub(ji,jj))*umask(ji,jj,jk) 185 END DO 186 END DO 187 END DO 217 188 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 189 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 190 ibdy1 = jpiglo-nbghostcells 191 ibdy2 = jpiglo-1 192 DO ji = mi0(ibdy1), mi1(ibdy2) 193 zvb(ji,:) = 0._wp 222 194 DO jk = 1, jpkm1 223 195 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)196 zvb(ji,jj) = zvb(ji,jj) & 197 & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 226 198 END DO 227 199 END DO 228 200 DO jj = 1, jpj 229 zvb( ibdy1:ibdy2,jj) = zvb(ibdy1:ibdy2,jj) * r1_hv_a(ibdy1:ibdy2,jj)201 zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj) 230 202 END DO 231 203 DO jk = 1, jpkm1 232 204 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 205 va(ji,jj,jk) = ( va(ji,jj,jk) & 206 & + va_b(ji,jj)-zvb(ji,jj)) * vmask(ji,jj,jk) 207 END DO 208 END DO 209 END DO 245 210 ENDIF 246 211 247 212 ! --- 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._wp213 jbdy1 = 2 214 jbdy2 = 1+nbghostcells 215 ! 216 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport 217 DO jj = mj0(jbdy1), mj1(jbdy2) 218 va_b(:,jj) = 0._wp 254 219 DO jk = 1, jpkm1 255 220 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)221 va_b(ji,jj) = va_b(ji,jj) & 222 & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 258 223 END DO 259 224 END DO 260 225 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 226 va_b(ji,jj) = va_b(ji,jj) * r1_hv_a(ji,jj) 227 END DO 228 END DO 229 ENDIF 230 ! 231 DO jj = mj0(jbdy1), mj1(jbdy2) 232 zvb(:,jj) = 0._wp ! Correct transport 274 233 DO jk=1,jpkm1 275 234 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)235 zvb(ji,jj) = zvb(ji,jj) & 236 & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 278 237 END DO 279 238 END DO 280 239 DO ji = 1, jpi 281 zvb(ji,j bdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv_a(ji,jbdy1:jbdy2)240 zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj) 282 241 END DO 283 242 284 243 DO jk = 1, jpkm1 285 244 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 245 va(ji,jj,jk) = ( va(ji,jj,jk) & 246 & + va_b(ji,jj) - zvb(ji,jj) ) * vmask(ji,jj,jk) 247 END DO 248 END DO 249 END DO 290 250 291 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 292 zub(:,jbdy1:jbdy2) = 0._wp 251 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 252 DO jj = mj0(jbdy1), mj1(jbdy2) 253 zub(:,jj) = 0._wp 293 254 DO jk = 1, jpkm1 294 255 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)256 zub(ji,jj) = zub(ji,jj) & 257 & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 297 258 END DO 298 259 END DO 299 260 DO ji = 1, jpi 300 zub(ji,j bdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu_a(ji,jbdy1:jbdy2)261 zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj) 301 262 END DO 302 263 303 264 DO jk = 1, jpkm1 304 265 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 266 ua(ji,jj,jk) = ( ua(ji,jj,jk) & 267 & + ua_b(ji,jj) - zub(ji,jj) ) * umask(ji,jj,jk) 268 END DO 269 END DO 270 END DO 317 271 ENDIF 318 272 319 273 ! --- 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._wp274 jbdy1 = jpjglo-1-nbghostcells 275 jbdy2 = jpjglo-2 276 ! 277 IF( .NOT.ln_dynspg_ts ) THEN ! Store transport 278 DO jj = mj0(jbdy1), mj1(jbdy2) 279 va_b(:,jj) = 0._wp 326 280 DO jk = 1, jpkm1 327 281 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)282 va_b(ji,jj) = va_b(ji,jj) & 283 & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 330 284 END DO 331 285 END DO 332 286 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 287 va_b(ji,jj) = va_b(ji,jj) * r1_hv_a(ji,jj) 288 END DO 289 END DO 290 ENDIF 291 ! 292 DO jj = mj0(jbdy1), mj1(jbdy2) 293 zvb(:,jj) = 0._wp ! Correct transport 346 294 DO jk=1,jpkm1 347 295 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)296 zvb(ji,jj) = zvb(ji,jj) & 297 & + e3v_a(ji,jj,jk) * va(ji,jj,jk) * vmask(ji,jj,jk) 350 298 END DO 351 299 END DO 352 300 DO ji = 1, jpi 353 zvb(ji,j bdy1:jbdy2) = zvb(ji,jbdy1:jbdy2) * r1_hv_a(ji,jbdy1:jbdy2)301 zvb(ji,jj) = zvb(ji,jj) * r1_hv_a(ji,jj) 354 302 END DO 355 303 356 304 DO jk = 1, jpkm1 357 305 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 306 va(ji,jj,jk) = ( va(ji,jj,jk) & 307 & + va_b(ji,jj) - zvb(ji,jj) ) * vmask(ji,jj,jk) 308 END DO 309 END DO 310 END DO 362 311 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 312 IF( ln_dynspg_ts ) THEN ! Set tangential velocities to time splitting estimate 313 jbdy1 = jpjglo-nbghostcells 314 jbdy2 = jpjglo-1 315 DO jj = mj0(jbdy1), mj1(jbdy2) 316 zub(:,jj) = 0._wp 367 317 DO jk = 1, jpkm1 368 318 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)319 zub(ji,jj) = zub(ji,jj) & 320 & + e3u_a(ji,jj,jk) * ua(ji,jj,jk) * umask(ji,jj,jk) 371 321 END DO 372 322 END DO 373 323 DO ji = 1, jpi 374 zub(ji,j bdy1:jbdy2) = zub(ji,jbdy1:jbdy2) * r1_hu_a(ji,jbdy1:jbdy2)324 zub(ji,jj) = zub(ji,jj) * r1_hu_a(ji,jj) 375 325 END DO 376 326 377 327 DO jk = 1, jpkm1 378 328 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 329 ua(ji,jj,jk) = ( ua(ji,jj,jk) & 330 & + ua_b(ji,jj) - zub(ji,jj) ) * umask(ji,jj,jk) 331 END DO 332 END DO 333 END DO 391 334 ENDIF 392 335 ! … … 401 344 !! 402 345 INTEGER :: ji, jj 346 INTEGER :: istart, iend, jstart, jend 403 347 !!---------------------------------------------------------------------- 404 348 ! 405 349 IF( Agrif_Root() ) RETURN 406 350 ! 407 IF((nbondi == -1).OR.(nbondi == 2)) THEN 351 !--- West ---! 352 istart = 2 353 iend = nbghostcells+1 354 DO ji = mi0(istart), mi1(iend) 408 355 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 356 va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 357 ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 358 END DO 359 END DO 360 ! 361 !--- East ---! 362 istart = jpiglo-nbghostcells 363 iend = jpiglo-1 364 DO ji = mi0(istart), mi1(iend) 419 365 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 366 va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 367 END DO 368 END DO 369 istart = jpiglo-nbghostcells-1 370 iend = jpiglo-2 371 DO ji = mi0(istart), mi1(iend) 372 DO jj=1,jpj 373 ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 374 END DO 375 END DO 376 ! 377 !--- South ---! 378 jstart = 2 379 jend = nbghostcells+1 380 DO jj = mj0(jstart), mj1(jend) 430 381 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 382 ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 383 va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 384 END DO 385 END DO 386 ! 387 !--- North ---! 388 jstart = jpjglo-nbghostcells 389 jend = jpjglo-1 390 DO jj = mj0(jstart), mj1(jend) 441 391 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 392 ua_e(ji,jj) = ubdy(ji,jj) * hur_e(ji,jj) 393 END DO 394 END DO 395 jstart = jpjglo-nbghostcells-1 396 jend = jpjglo-2 397 DO jj = mj0(jstart), mj1(jend) 398 DO ji=1,jpi 399 va_e(ji,jj) = vbdy(ji,jj) * hvr_e(ji,jj) 400 END DO 401 END DO 450 402 ! 451 403 END SUBROUTINE Agrif_dyn_ts 452 404 405 SUBROUTINE Agrif_dyn_ts_flux( jn, zu, zv ) 406 !!---------------------------------------------------------------------- 407 !! *** ROUTINE Agrif_dyn_ts_flux *** 408 !!---------------------------------------------------------------------- 409 INTEGER, INTENT(in) :: jn 410 REAL(wp), DIMENSION(jpi,jpj), INTENT(inout) :: zu, zv 411 !! 412 INTEGER :: ji, jj 413 INTEGER :: istart, iend, jstart, jend 414 !!---------------------------------------------------------------------- 415 ! 416 IF( Agrif_Root() ) RETURN 417 ! 418 !--- West ---! 419 istart = 2 420 iend = nbghostcells+1 421 DO ji = mi0(istart), mi1(iend) 422 DO jj=1,jpj 423 zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 424 zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 425 END DO 426 END DO 427 ! 428 !--- East ---! 429 istart = jpiglo-nbghostcells 430 iend = jpiglo-1 431 DO ji = mi0(istart), mi1(iend) 432 DO jj=1,jpj 433 zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 434 END DO 435 END DO 436 istart = jpiglo-nbghostcells-1 437 iend = jpiglo-2 438 DO ji = mi0(istart), mi1(iend) 439 DO jj=1,jpj 440 zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 441 END DO 442 END DO 443 ! 444 !--- South ---! 445 jstart = 2 446 jend = nbghostcells+1 447 DO jj = mj0(jstart), mj1(jend) 448 DO ji=1,jpi 449 zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 450 zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 451 END DO 452 END DO 453 ! 454 !--- North ---! 455 jstart = jpjglo-nbghostcells 456 jend = jpjglo-1 457 DO jj = mj0(jstart), mj1(jend) 458 DO ji=1,jpi 459 zu(ji,jj) = ubdy(ji,jj) * e2u(ji,jj) 460 END DO 461 END DO 462 jstart = jpjglo-nbghostcells-1 463 jend = jpjglo-2 464 DO jj = mj0(jstart), mj1(jend) 465 DO ji=1,jpi 466 zv(ji,jj) = vbdy(ji,jj) * e1v(ji,jj) 467 END DO 468 END DO 469 ! 470 END SUBROUTINE Agrif_dyn_ts_flux 453 471 454 472 SUBROUTINE Agrif_dta_ts( kt ) … … 470 488 ! 471 489 ! Interpolate barotropic fluxes 472 Agrif_SpecialValue =0._wp490 Agrif_SpecialValue = 0._wp 473 491 Agrif_UseSpecialValue = ln_spc_dyn 492 ! 493 ! Set bdy time interpolation stage to 0 (latter incremented locally do deal with corners) 494 utint_stage(:,:) = 0 495 vtint_stage(:,:) = 0 474 496 ! 475 497 IF( ll_int_cons ) THEN ! Conservative interpolation 476 498 ! order matters here !!!!!! 477 499 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 ) 479 bdy_tinterp = 1500 CALL Agrif_Bc_variable( vb2b_interp_id, calledweight=1._wp, procname=interpvb2b ) 501 ! 480 502 CALL Agrif_Bc_variable( unb_id , calledweight=1._wp, procname=interpunb ) ! After 481 503 CALL Agrif_Bc_variable( vnb_id , calledweight=1._wp, procname=interpvnb ) 482 bdy_tinterp = 2504 ! 483 505 CALL Agrif_Bc_variable( unb_id , calledweight=0._wp, procname=interpunb ) ! Before 484 506 CALL Agrif_Bc_variable( vnb_id , calledweight=0._wp, procname=interpvnb ) 485 507 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 508 ! 509 ubdy(:,:) = 0._wp ; vbdy(:,:) = 0._wp 491 510 CALL Agrif_Bc_variable( unb_id, procname=interpunb ) 492 511 CALL Agrif_Bc_variable( vnb_id, procname=interpvnb ) … … 503 522 INTEGER, INTENT(in) :: kt 504 523 ! 505 INTEGER :: ji, jj, indx, indy 524 INTEGER :: ji, jj 525 INTEGER :: istart, iend, jstart, jend 506 526 !!---------------------------------------------------------------------- 507 527 ! … … 516 536 ! 517 537 ! --- West --- ! 518 IF((nbondi == -1).OR.(nbondi == 2)) THEN 519 indx = 1+nbghostcells 538 istart = 2 539 iend = 1 + nbghostcells 540 DO ji = mi0(istart), mi1(iend) 520 541 DO jj = 1, jpj 521 DO ji = 2, indx 522 ssha(ji,jj) = hbdy_w(ji-1,jj) 523 ENDDO 542 ssha(ji,jj) = hbdy(ji,jj) 524 543 ENDDO 525 END IF544 ENDDO 526 545 ! 527 546 ! --- East --- ! 528 IF((nbondi == 1).OR.(nbondi == 2)) THEN 529 indx = nlci-nbghostcells 547 istart = jpiglo - nbghostcells 548 iend = jpiglo - 1 549 DO ji = mi0(istart), mi1(iend) 530 550 DO jj = 1, jpj 531 DO ji = indx, nlci-1 532 ssha(ji,jj) = hbdy_e(ji-indx+1,jj) 533 ENDDO 551 ssha(ji,jj) = hbdy(ji,jj) 534 552 ENDDO 535 END IF553 ENDDO 536 554 ! 537 555 ! --- 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 556 jstart = 2 557 jend = 1 + nbghostcells 558 DO jj = mj0(jstart), mj1(jend) 559 DO ji = 1, jpi 560 ssha(ji,jj) = hbdy(ji,jj) 544 561 ENDDO 545 END IF562 ENDDO 546 563 ! 547 564 ! --- 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 565 jstart = jpjglo - nbghostcells 566 jend = jpjglo - 1 567 DO jj = mj0(jstart), mj1(jend) 568 DO ji = 1, jpi 569 ssha(ji,jj) = hbdy(ji,jj) 554 570 ENDDO 555 END IF571 ENDDO 556 572 ! 557 573 END SUBROUTINE Agrif_ssh … … 564 580 INTEGER, INTENT(in) :: jn 565 581 !! 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)582 INTEGER :: ji, jj 583 INTEGER :: istart, iend, jstart, jend 584 !!---------------------------------------------------------------------- 569 585 ! 570 586 IF( Agrif_Root() ) RETURN 571 587 ! 572 588 ! --- West --- ! 573 IF((nbondi == -1).OR.(nbondi == 2)) THEN 574 indx = 1+nbghostcells 589 istart = 2 590 iend = 1+nbghostcells 591 DO ji = mi0(istart), mi1(iend) 575 592 DO jj = 1, jpj 576 DO ji = 2, indx 577 ssha_e(ji,jj) = hbdy_w(ji-1,jj) 578 ENDDO 593 ssha_e(ji,jj) = hbdy(ji,jj) 579 594 ENDDO 580 END IF595 ENDDO 581 596 ! 582 597 ! --- East --- ! 583 IF((nbondi == 1).OR.(nbondi == 2)) THEN 584 indx = nlci-nbghostcells 598 istart = jpiglo - nbghostcells 599 iend = jpiglo - 1 600 DO ji = mi0(istart), mi1(iend) 585 601 DO jj = 1, jpj 586 DO ji = indx, nlci-1 587 ssha_e(ji,jj) = hbdy_e(ji-indx+1,jj) 588 ENDDO 602 ssha_e(ji,jj) = hbdy(ji,jj) 589 603 ENDDO 590 END IF604 ENDDO 591 605 ! 592 606 ! --- 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 607 jstart = 2 608 jend = 1+nbghostcells 609 DO jj = mj0(jstart), mj1(jend) 610 DO ji = 1, jpi 611 ssha_e(ji,jj) = hbdy(ji,jj) 599 612 ENDDO 600 END IF613 ENDDO 601 614 ! 602 615 ! --- 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 616 jstart = jpjglo - nbghostcells 617 jend = jpjglo - 1 618 DO jj = mj0(jstart), mj1(jend) 619 DO ji = 1, jpi 620 ssha_e(ji,jj) = hbdy(ji,jj) 609 621 ENDDO 610 END IF622 ENDDO 611 623 ! 612 624 END SUBROUTINE Agrif_ssh_ts … … 634 646 635 647 636 SUBROUTINE interptsn( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before , nb, ndir)648 SUBROUTINE interptsn( ptab, i1, i2, j1, j2, k1, k2, n1, n2, before ) 637 649 !!---------------------------------------------------------------------- 638 650 !! *** ROUTINE interptsn *** … … 641 653 INTEGER , INTENT(in ) :: i1, i2, j1, j2, k1, k2, n1, n2 642 654 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 655 ! 656 INTEGER :: ji, jj, jk, jn ! dummy loop indices 657 INTEGER :: N_in, N_out 649 658 ! vertical interpolation: 650 659 REAL(wp), DIMENSION(i1:i2,j1:j2,1:jpk,n1:n2) :: ptab_child … … 652 661 REAL(wp), DIMENSION(k1:k2) :: h_in 653 662 REAL(wp), DIMENSION(1:jpk) :: h_out 654 REAL(wp) :: h_diff655 663 656 664 IF( before ) THEN … … 676 684 ELSE 677 685 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 686 # if defined key_vertical 682 687 DO jj=j1,j2 683 688 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 689 N_in = 0 691 690 DO jk=k1,k2 !k2 = jpk of parent grid … … 697 696 N_out = 0 698 697 DO jk=1,jpk ! jpk of child grid 699 IF (tmask( iref,jref,jk) == 0) EXIT698 IF (tmask(ji,jj,jk) == 0) EXIT 700 699 N_out = N_out + 1 701 h_out(jk) = e3t_n( iref,jref,jk)700 h_out(jk) = e3t_n(ji,jj,jk) 702 701 ENDDO 703 702 IF (N_in > 0) THEN … … 716 715 END DO 717 716 718 IF ( .NOT.lk_agrif_clp ) THEN719 !720 imin = i1 ; imax = i2721 jmin = j1 ; jmax = j2722 !723 ! Remove CORNERS724 IF((nbondj == -1).OR.(nbondj == 2)) jmin = 2 + nbghostcells725 IF((nbondj == +1).OR.(nbondj == 2)) jmax = nlcj - nbghostcells - 1726 IF((nbondi == -1).OR.(nbondi == 2)) imin = 2 + nbghostcells727 IF((nbondi == +1).OR.(nbondi == 2)) imax = nlci - nbghostcells - 1728 !729 IF( eastern_side ) THEN730 zrho = Agrif_Rhox()731 z1 = ( zrho - 1._wp ) * 0.5_wp732 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 - z7736 !737 ibdy = nlci-nbghostcells738 DO jn = 1, jpts739 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, jpkm1741 DO jj = jmin,jmax742 IF( umask(ibdy-1,jj,jk) == 0._wp ) THEN743 tsa(ibdy,jj,jk,jn) = tsa(ibdy+1,jj,jk,jn) * tmask(ibdy,jj,jk)744 ELSE745 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 ) THEN747 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 ENDIF750 ENDIF751 END DO752 END DO753 ! 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 DO756 ENDIF757 !758 IF( northern_side ) THEN759 zrho = Agrif_Rhoy()760 z1 = ( zrho - 1._wp ) * 0.5_wp761 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 - z7765 !766 jbdy = nlcj-nbghostcells767 DO jn = 1, jpts768 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, jpkm1770 DO ji = imin,imax771 IF( vmask(ji,jbdy-1,jk) == 0._wp ) THEN772 tsa(ji,jbdy,jk,jn) = tsa(ji,jbdy+1,jk,jn) * tmask(ji,jbdy,jk)773 ELSE774 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 ) THEN776 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 ENDIF779 ENDIF780 END DO781 END DO782 ! 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 DO785 ENDIF786 !787 IF( western_side ) THEN788 zrho = Agrif_Rhox()789 z1 = ( zrho - 1._wp ) * 0.5_wp790 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 - z7794 !795 ibdy = 1+nbghostcells796 DO jn = 1, jpts797 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, jpkm1799 DO jj = jmin,jmax800 IF( umask(ibdy,jj,jk) == 0._wp ) THEN801 tsa(ibdy,jj,jk,jn) = tsa(ibdy-1,jj,jk,jn) * tmask(ibdy,jj,jk)802 ELSE803 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 ) THEN805 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 ENDIF808 ENDIF809 END DO810 END DO811 ! 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 DO814 ENDIF815 !816 IF( southern_side ) THEN817 zrho = Agrif_Rhoy()818 z1 = ( zrho - 1._wp ) * 0.5_wp819 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 - z7823 !824 jbdy=1+nbghostcells825 DO jn = 1, jpts826 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, jpkm1828 DO ji = imin,imax829 IF( vmask(ji,jbdy,jk) == 0._wp ) THEN830 tsa(ji,jbdy,jk,jn)=tsa(ji,jbdy-1,jk,jn) * tmask(ji,jbdy,jk)831 ELSE832 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 ) THEN834 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 ENDIF837 ENDIF838 END DO839 END DO840 ! 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 DO843 ENDIF844 !845 ENDIF846 717 ENDIF 847 718 ! 848 719 END SUBROUTINE interptsn 849 720 850 SUBROUTINE interpsshn( ptab, i1, i2, j1, j2, before , nb, ndir)721 SUBROUTINE interpsshn( ptab, i1, i2, j1, j2, before ) 851 722 !!---------------------------------------------------------------------- 852 723 !! *** ROUTINE interpsshn *** … … 855 726 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 856 727 LOGICAL , INTENT(in ) :: before 857 INTEGER , INTENT(in ) :: nb , ndir 858 ! 859 LOGICAL :: western_side, eastern_side,northern_side,southern_side 728 ! 860 729 !!---------------------------------------------------------------------- 861 730 ! … … 863 732 ptab(i1:i2,j1:j2) = sshn(i1:i2,j1:j2) 864 733 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) 734 hbdy(i1:i2,j1:j2) = ptab(i1:i2,j1:j2) * tmask(i1:i2,j1:j2,1) 874 735 ENDIF 875 736 ! … … 1045 906 END SUBROUTINE interpvn 1046 907 1047 SUBROUTINE interpunb( ptab, i1, i2, j1, j2, before , nb, ndir)908 SUBROUTINE interpunb( ptab, i1, i2, j1, j2, before) 1048 909 !!---------------------------------------------------------------------- 1049 910 !! *** ROUTINE interpunb *** … … 1052 913 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1053 914 LOGICAL , INTENT(in ) :: before 1054 INTEGER , INTENT(in ) :: nb , ndir1055 915 ! 1056 916 INTEGER :: ji, jj 1057 917 REAL(wp) :: zrhoy, zrhot, zt0, zt1, ztcoeff 1058 LOGICAL :: western_side, eastern_side,northern_side,southern_side1059 918 !!---------------------------------------------------------------------- 1060 919 ! … … 1062 921 ptab(i1:i2,j1:j2) = e2u(i1:i2,j1:j2) * hu_n(i1:i2,j1:j2) * un_b(i1:i2,j1:j2) 1063 922 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 923 zrhoy = Agrif_Rhoy() 1069 924 zrhot = Agrif_rhot() … … 1071 926 zt0 = REAL(Agrif_NbStepint() , wp) / zrhot 1072 927 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 928 ! 929 DO ji = i1, i2 930 DO jj = j1, j2 931 IF ( utint_stage(ji,jj) == 1 ) THEN 932 ztcoeff = zrhot * ( zt1**2._wp * ( zt1 - 1._wp) & 933 & - zt0**2._wp * ( zt0 - 1._wp) ) 934 ELSEIF( utint_stage(ji,jj) == 2 ) THEN 935 ztcoeff = zrhot * ( zt1 * ( zt1 - 1._wp)**2._wp & 936 & - zt0 * ( zt0 - 1._wp)**2._wp ) 937 ELSEIF( utint_stage(ji,jj) == 0 ) THEN 938 ztcoeff = 1._wp 939 ELSE 940 ztcoeff = 0._wp 941 ENDIF 942 ! 943 ubdy(ji,jj) = ubdy(ji,jj) + ztcoeff * ptab(ji,jj) 944 ! 945 IF (( utint_stage(ji,jj) == 2 ).OR.( utint_stage(ji,jj) == 0 )) THEN 946 ubdy(ji,jj) = ubdy(ji,jj) / (zrhoy*e2u(ji,jj)) * umask(ji,jj,1) 947 utint_stage(ji,jj) = 3 948 ELSE 949 utint_stage(ji,jj) = utint_stage(ji,jj) + 1 950 ENDIF 951 END DO 952 END DO 953 END IF 1096 954 ! 1097 955 END SUBROUTINE interpunb 1098 956 1099 957 1100 SUBROUTINE interpvnb( ptab, i1, i2, j1, j2, before , nb, ndir)958 SUBROUTINE interpvnb( ptab, i1, i2, j1, j2, before ) 1101 959 !!---------------------------------------------------------------------- 1102 960 !! *** ROUTINE interpvnb *** … … 1105 963 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1106 964 LOGICAL , INTENT(in ) :: before 1107 INTEGER , INTENT(in ) :: nb , ndir 1108 ! 1109 INTEGER :: ji,jj 965 ! 966 INTEGER :: ji, jj 1110 967 REAL(wp) :: zrhox, zrhot, zt0, zt1, ztcoeff 1111 LOGICAL :: western_side, eastern_side,northern_side,southern_side1112 968 !!---------------------------------------------------------------------- 1113 969 ! … … 1115 971 ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * hv_n(i1:i2,j1:j2) * vn_b(i1:i2,j1:j2) 1116 972 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 973 zrhox = Agrif_Rhox() 1122 974 zrhot = Agrif_rhot() 1123 975 ! Time indexes bounds for integration 1124 976 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 977 zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot 978 ! 979 DO ji = i1, i2 980 DO jj = j1, j2 981 IF ( vtint_stage(ji,jj) == 1 ) THEN 982 ztcoeff = zrhot * ( zt1**2._wp * ( zt1 - 1._wp) & 983 & - zt0**2._wp * ( zt0 - 1._wp) ) 984 ELSEIF( vtint_stage(ji,jj) == 2 ) THEN 985 ztcoeff = zrhot * ( zt1 * ( zt1 - 1._wp)**2._wp & 986 & - zt0 * ( zt0 - 1._wp)**2._wp ) 987 ELSEIF( vtint_stage(ji,jj) == 0 ) THEN 988 ztcoeff = 1._wp 989 ELSE 990 ztcoeff = 0._wp 991 ENDIF 992 ! 993 vbdy(ji,jj) = vbdy(ji,jj) + ztcoeff * ptab(ji,jj) 994 ! 995 IF (( vtint_stage(ji,jj) == 2 ).OR.( vtint_stage(ji,jj) == 0 )) THEN 996 vbdy(ji,jj) = vbdy(ji,jj) / (zrhox*e1v(ji,jj)) * vmask(ji,jj,1) 997 vtint_stage(ji,jj) = 3 998 ELSE 999 vtint_stage(ji,jj) = vtint_stage(ji,jj) + 1 1000 ENDIF 1001 END DO 1002 END DO 1147 1003 ENDIF 1148 1004 ! … … 1150 1006 1151 1007 1152 SUBROUTINE interpub2b( ptab, i1, i2, j1, j2, before , nb, ndir)1008 SUBROUTINE interpub2b( ptab, i1, i2, j1, j2, before ) 1153 1009 !!---------------------------------------------------------------------- 1154 1010 !! *** ROUTINE interpub2b *** … … 1157 1013 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1158 1014 LOGICAL , INTENT(in ) :: before 1159 INTEGER , INTENT(in ) :: nb , ndir1160 1015 ! 1161 1016 INTEGER :: ji,jj 1162 REAL(wp) :: zrhot, zt0, zt1,zat 1163 LOGICAL :: western_side, eastern_side,northern_side,southern_side 1017 REAL(wp) :: zrhot, zt0, zt1, zat 1164 1018 !!---------------------------------------------------------------------- 1165 1019 IF( before ) THEN … … 1170 1024 ENDIF 1171 1025 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 1026 zrhot = Agrif_rhot() 1219 1027 ! Time indexes bounds for integration … … 1224 1032 & - zt0**2._wp * (-2._wp*zt0 + 3._wp) ) 1225 1033 ! 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) 1034 ubdy(i1:i2,j1:j2) = zat * ptab(i1:i2,j1:j2) 1035 ! 1036 ! Update interpolation stage: 1037 utint_stage(i1:i2,j1:j2) = 1 1038 ENDIF 1039 ! 1040 END SUBROUTINE interpub2b 1041 1042 1043 SUBROUTINE interpvb2b( ptab, i1, i2, j1, j2, before ) 1044 !!---------------------------------------------------------------------- 1045 !! *** ROUTINE interpvb2b *** 1046 !!---------------------------------------------------------------------- 1047 INTEGER , INTENT(in ) :: i1, i2, j1, j2 1048 REAL(wp), DIMENSION(i1:i2,j1:j2), INTENT(inout) :: ptab 1049 LOGICAL , INTENT(in ) :: before 1050 ! 1051 INTEGER :: ji,jj 1052 REAL(wp) :: zrhot, zt0, zt1, zat 1053 !!---------------------------------------------------------------------- 1054 ! 1055 IF( before ) THEN 1056 IF ( ln_bt_fw ) THEN 1057 ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vb2_b(i1:i2,j1:j2) 1058 ELSE 1059 ptab(i1:i2,j1:j2) = e1v(i1:i2,j1:j2) * vn_adv(i1:i2,j1:j2) 1060 ENDIF 1061 ELSE 1062 zrhot = Agrif_rhot() 1063 ! Time indexes bounds for integration 1064 zt0 = REAL(Agrif_NbStepint() , wp) / zrhot 1065 zt1 = REAL(Agrif_NbStepint()+1, wp) / zrhot 1066 ! Polynomial interpolation coefficients: 1067 zat = zrhot * ( zt1**2._wp * (-2._wp*zt1 + 3._wp) & 1068 & - zt0**2._wp * (-2._wp*zt0 + 3._wp) ) 1069 ! 1070 vbdy(i1:i2,j1:j2) = zat * ptab(i1:i2,j1:j2) 1071 ! 1072 ! update interpolation stage: 1073 vtint_stage(i1:i2,j1:j2) = 1 1230 1074 ENDIF 1231 1075 !
Note: See TracChangeset
for help on using the changeset viewer.
