Changeset 7698 for trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynhpg.F90
- Timestamp:
- 2017-02-18T10:02:03+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/DYN/dynhpg.F90
r7646 r7698 84 84 !!---------------------------------------------------------------------- 85 85 INTEGER, INTENT(in) :: kt ! ocean time-step index 86 INTEGER :: jk, jj, ji 86 87 REAL(wp), POINTER, DIMENSION(:,:,:) :: ztrdu, ztrdv 87 88 !!---------------------------------------------------------------------- … … 91 92 IF( l_trddyn ) THEN ! Temporary saving of ua and va trends (l_trddyn) 92 93 CALL wrk_alloc( jpi,jpj,jpk, ztrdu, ztrdv ) 93 ztrdu(:,:,:) = ua(:,:,:) 94 ztrdv(:,:,:) = va(:,:,:) 94 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji) 95 DO jk = 1, jpk 96 DO jj = 1, jpj 97 DO ji = 1, jpi 98 ztrdu(ji,jj,jk) = ua(ji,jj,jk) 99 ztrdv(ji,jj,jk) = va(ji,jj,jk) 100 END DO 101 END DO 102 END DO 95 103 ENDIF 96 104 ! … … 105 113 ! 106 114 IF( l_trddyn ) THEN ! save the hydrostatic pressure gradient trends for momentum trend diagnostics 107 ztrdu(:,:,:) = ua(:,:,:) - ztrdu(:,:,:) 108 ztrdv(:,:,:) = va(:,:,:) - ztrdv(:,:,:) 115 !$OMP PARALLEL DO schedule(static) private(jk, jj, ji) 116 DO jk = 1, jpk 117 DO jj = 1, jpj 118 DO ji = 1, jpi 119 ztrdu(ji,jj,jk) = ua(ji,jj,jk) - ztrdu(ji,jj,jk) 120 ztrdv(ji,jj,jk) = va(ji,jj,jk) - ztrdv(ji,jj,jk) 121 END DO 122 END DO 123 END DO 109 124 CALL trd_dyn( ztrdu, ztrdv, jpdyn_hpg, kt ) 110 125 CALL wrk_dealloc( jpi,jpj,jpk, ztrdu, ztrdv ) … … 198 213 ! 199 214 ! initialisation of ice shelf load 200 IF ( .NOT. ln_isfcav ) riceload(:,:)=0.0 215 IF ( .NOT. ln_isfcav ) THEN 216 !$OMP PARALLEL DO schedule(static) private(jj, ji) 217 DO jj = 1, jpj 218 DO ji = 1, jpi 219 riceload(ji,jj)=0.0 220 END DO 221 END DO 222 END IF 201 223 IF ( ln_isfcav ) THEN 202 224 CALL wrk_alloc( jpi,jpj, 2, ztstop) … … 212 234 213 235 ! assume water displaced by the ice shelf is at T=-1.9 and S=34.4 (rude) 214 ztstop(:,:,1)=-1.9_wp ; ztstop(:,:,2)=34.4_wp 236 !$OMP PARALLEL DO schedule(static) private(jj, ji) 237 DO jj = 1, jpj 238 DO ji = 1, jpi 239 ztstop(ji,jj,1)=-1.9_wp 240 ztstop(ji,jj,2)=34.4_wp 241 END DO 242 END DO 215 243 216 244 ! compute density of the water displaced by the ice shelf … … 226 254 ! divided by 2 later 227 255 ziceload = 0._wp 256 !$OMP PARALLEL 257 !$OMP DO schedule(static) private(jj,ji,ikt,jk) 228 258 DO jj = 1, jpj 229 259 DO ji = 1, jpi … … 238 268 END DO 239 269 END DO 240 riceload(:,:)=ziceload(:,:) ! need to be saved for diaar5 270 !$OMP DO schedule(static) private(jj, ji) 271 DO jj = 1, jpj 272 DO ji = 1, jpi 273 riceload(ji,jj)=ziceload(ji,jj) ! need to be saved for diaar5 274 END DO 275 END DO 276 !$OMP END PARALLEL 241 277 242 278 CALL wrk_dealloc( jpi,jpj, 2, ztstop) … … 282 318 283 319 ! Surface value 320 !$OMP PARALLEL 321 !$OMP DO schedule(static) private(ji,jj, zcoef1) 284 322 DO jj = 2, jpjm1 285 323 DO ji = fs_2, fs_jpim1 ! vector opt. … … 297 335 ! interior value (2=<jk=<jpkm1) 298 336 DO jk = 2, jpkm1 337 !$OMP DO schedule(static) private(ji,jj, zcoef1) 299 338 DO jj = 2, jpjm1 300 339 DO ji = fs_2, fs_jpim1 ! vector opt. … … 313 352 END DO 314 353 END DO 315 END DO 354 !$OMP END DO NOWAIT 355 END DO 356 !$OMP END PARALLEL 316 357 ! 317 358 CALL wrk_dealloc( jpi,jpj,jpk, zhpi, zhpj ) … … 351 392 352 393 ! Surface value (also valid in partial step case) 394 !$OMP PARALLEL 395 !$OMP DO schedule(static) private(ji,jj,zcoef1) 353 396 DO jj = 2, jpjm1 354 397 DO ji = fs_2, fs_jpim1 ! vector opt. … … 365 408 ! interior value (2=<jk=<jpkm1) 366 409 DO jk = 2, jpkm1 410 !$OMP DO schedule(static) private(ji,jj, zcoef1) 367 411 DO jj = 2, jpjm1 368 412 DO ji = fs_2, fs_jpim1 ! vector opt. … … 384 428 385 429 ! partial steps correction at the last level (use gru & grv computed in zpshde.F90) 430 !$OMP DO schedule(static) private(ji,jj,iku,ikv,zcoef2,zcoef3) 386 431 DO jj = 2, jpjm1 387 432 DO ji = 2, jpim1 … … 404 449 END DO 405 450 END DO 451 !$OMP END PARALLEL 406 452 ! 407 453 CALL wrk_dealloc( jpi,jpj,jpk, zhpi, zhpj )
Note: See TracChangeset
for help on using the changeset viewer.