Changeset 15422
- Timestamp:
- 2021-10-21T11:19:25+02:00 (3 years ago)
- Location:
- NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE
- Files:
-
- 2 added
- 15 edited
-
BDY/bdyini.F90 (modified) (1 diff)
-
DIA/diaharm.F90 (modified) (14 diffs)
-
DIA/diaharm_fast.F90 (added)
-
DOM/dtatsd.F90 (modified) (9 diffs)
-
DYN/dynnxt.F90 (modified) (7 diffs)
-
DYN/dynspg.F90 (modified) (3 diffs)
-
SBC/fldread.F90 (modified) (1 diff)
-
SBC/sbctide.F90 (modified) (3 diffs)
-
SBC/tide_FES14.h90 (added)
-
SBC/tide_mod.F90 (modified) (3 diffs)
-
SBC/tideini.F90 (modified) (4 diffs)
-
USR/usrdef_istate.F90 (modified) (1 diff)
-
USR/usrdef_sbc.F90 (modified) (3 diffs)
-
par_oce.F90 (modified) (1 diff)
-
step.F90 (modified) (1 diff)
-
step_oce.F90 (modified) (1 diff)
-
stpctl.F90 (modified) (6 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/BDY/bdyini.F90
r14075 r15422 167 167 ! Check and write out namelist parameters 168 168 ! ----------------------------------------- 169 IF( jperio /= 0 ) CALL ctl_stop( 'bdy_segs: Cyclic or symmetric,', &170 & ' and general open boundary condition are not compatible' )169 ! IF( jperio /= 0 ) CALL ctl_stop( 'bdy_segs: Cyclic or symmetric,', & 170 ! & ' and general open boundary condition are not compatible' ) 171 171 172 172 IF(lwp) WRITE(numout,*) 'Number of open boundary sets : ', nb_bdy -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/DIA/diaharm.F90
r14075 r15422 5 5 !!====================================================================== 6 6 !! History : 3.1 ! 2007 (O. Le Galloudec, J. Chanut) Original code 7 !! 8 !! NB: 2017-12 : add 3D harmonic analysis of velocities 9 !! integration of Maria Luneva's development 10 !! 'key_3Ddiaharm 7 11 !!---------------------------------------------------------------------- 8 12 USE oce ! ocean dynamics and tracers variables … … 13 17 USE sbctide ! Tidal forcing or not 14 18 ! 19 # if defined key_3Ddiaharm 20 USE zdf_oce 21 #endif 22 ! 15 23 USE in_out_manager ! I/O units 16 24 USE iom ! I/0 library … … 33 41 INTEGER :: nb_ana ! Number of harmonics to analyse 34 42 35 INTEGER , ALLOCATABLE, DIMENSION(:) :: name 36 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ana_temp 37 REAL(wp), ALLOCATABLE, DIMENSION(:) :: ana_freq, ut, vt, ft 38 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: out_eta, out_u, out_v 43 INTEGER , ALLOCATABLE, DIMENSION(:) :: name 44 REAL(wp), ALLOCATABLE, DIMENSION(:) :: ana_freq, ut, vt, ft 45 # if defined key_3Ddiaharm 46 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:,:) :: ana_temp 47 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: out_eta, out_u, out_v, out_w, out_dzi 48 # else 49 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ana_temp 50 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: out_eta, out_u, out_v 51 # endif 39 52 40 53 INTEGER :: ninco, nsparse … … 76 89 WRITE(numout,*) 77 90 WRITE(numout,*) 'dia_harm_init: Tidal harmonic analysis initialization' 91 # if defined key_3Ddiaharm 92 WRITE(numout,*) ' - 3D harmonic analysis of currents activated (key_3Ddiaharm)' 93 #endif 78 94 WRITE(numout,*) '~~~~~~~~~~~~~ ' 79 95 ENDIF … … 155 171 ! Initialize temporary arrays: 156 172 ! ---------------------------- 157 ALLOCATE( ana_temp(jpi,jpj,2*nb_ana,3) ) 158 ana_temp(:,:,:,:) = 0._wp 173 # if defined key_3Ddiaharm 174 ALLOCATE( ana_temp( jpi, jpj, 2*nb_ana, 5, jpk ) ) 175 ana_temp(:,:,:,:,:) = 0._wp 176 # else 177 ALLOCATE( ana_temp( jpi, jpj, 2*nb_ana, 3 ) ) 178 ana_temp(:,:,:,: ) = 0._wp 179 #endif 159 180 160 181 ENDIF … … 175 196 ! 176 197 INTEGER :: ji, jj, jh, jc, nhc 198 # if defined key_3Ddiaharm 199 INTEGER :: jk 200 # endif 177 201 REAL(wp) :: ztime, ztemp 178 202 !!-------------------------------------------------------------------- … … 190 214 & +(1.-MOD(jc,2))* ft(jh) *SIN(ana_freq(jh)*ztime + vt(jh) + ut(jh))) 191 215 ! 216 ! ssh, ub, vb are stored at the last level of 5d array 192 217 DO jj = 2, jpjm1 193 218 DO ji = 2, jpim1 194 ana_temp(ji,jj,nhc,1) = ana_temp(ji,jj,nhc,1) + ztemp * sshn(ji,jj) * ssmask (ji,jj) ! elevation 195 ana_temp(ji,jj,nhc,2) = ana_temp(ji,jj,nhc,2) + ztemp * un_b(ji,jj) * ssumask(ji,jj) ! u-vel 196 ana_temp(ji,jj,nhc,3) = ana_temp(ji,jj,nhc,3) + ztemp * vn_b(ji,jj) * ssvmask(ji,jj) ! v-vel 219 ! Elevation and currents 220 # if defined key_3Ddiaharm 221 ana_temp(ji,jj,nhc,1,jpk) = ana_temp(ji,jj,nhc,1,jpk) + ztemp*sshn(ji,jj)*ssmask (ji,jj) 222 ana_temp(ji,jj,nhc,2,jpk) = ana_temp(ji,jj,nhc,2,jpk) + ztemp*un_b(ji,jj)*ssumask(ji,jj) 223 ana_temp(ji,jj,nhc,3,jpk) = ana_temp(ji,jj,nhc,3,jpk) + ztemp*vn_b(ji,jj)*ssvmask(ji,jj) 224 225 ana_temp(ji,jj,nhc,5,jpk) = ana_temp(ji,jj,nhc,5,jpk) & 226 & + ztemp*bfrva(ji,jj)*vn(ji,jj,mbkv(ji,jj))*ssvmask(ji,jj) 227 ana_temp(ji,jj,nhc,4,jpk) = ana_temp(ji,jj,nhc,4,jpk) & 228 & + ztemp*bfrua(ji,jj)*un(ji,jj,mbku(ji,jj))*ssumask(ji,jj) 229 # else 230 ana_temp(ji,jj,nhc,1) = ana_temp(ji,jj,nhc,1) + ztemp*sshn(ji,jj)*ssmask (ji,jj) 231 ana_temp(ji,jj,nhc,2) = ana_temp(ji,jj,nhc,2) + ztemp*un_b(ji,jj)*ssumask(ji,jj) 232 ana_temp(ji,jj,nhc,3) = ana_temp(ji,jj,nhc,3) + ztemp*vn_b(ji,jj)*ssvmask(ji,jj) 233 # endif 197 234 END DO 198 235 END DO 236 ! 237 # if defined key_3Ddiaharm 238 ! 3d velocity and density: 239 DO jk=1,jpk-1 240 DO jj = 1,jpj 241 DO ji = 1,jpi 242 ! density and velocity 243 ana_temp(ji,jj,nhc,1,jk) = ana_temp(ji,jj,nhc,1,jk) + ztemp*rhd(ji,jj,jk) 244 ana_temp(ji,jj,nhc,2,jk) = ana_temp(ji,jj,nhc,2,jk) + ztemp*(un(ji,jj,jk)-un_b(ji,jj)) & 245 & *umask(ji,jj,jk) 246 ana_temp(ji,jj,nhc,3,jk) = ana_temp(ji,jj,nhc,3,jk) + ztemp*(vn(ji,jj,jk)-vn_b(ji,jj)) & 247 & *vmask(ji,jj,jk) 248 ana_temp(ji,jj,nhc,4,jk) = ana_temp(ji,jj,nhc,4,jk) + ztemp*wn(ji,jj,jk) 249 250 ana_temp(ji,jj,nhc,5,jk) = ana_temp(ji,jj,nhc,5,jk) - 0.5*grav*ztemp*(rhd(ji,jj,jk)+rhd(ji,jj,jk+1))/max(rn2(ji,jj,jk),1.e-8_wp) 251 END DO 252 END DO 253 ENDDO 254 # endif 199 255 END DO 200 256 END DO … … 218 274 !!-------------------------------------------------------------------- 219 275 INTEGER :: ji, jj, jh, jc, jn, nhan 276 # if defined key_3Ddiaharm 277 INTEGER :: jk 278 # endif 220 279 INTEGER :: ksp, kun, keq 221 280 REAL(wp) :: ztime, ztime_ini, ztime_end, z1_han … … 226 285 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 227 286 228 ALLOCATE( out_eta(jpi,jpj,2*nb_ana), out_u(jpi,jpj,2*nb_ana), out_v(jpi,jpj,2*nb_ana) )229 230 287 ztime_ini = nit000_han*rdt ! Initial time in seconds at the beginning of analysis 231 288 ztime_end = nitend_han*rdt ! Final time in seconds at the end of analysis … … 233 290 z1_han = 1._wp / REAL(nhan-1) 234 291 292 # if defined key_3Ddiaharm 293 ALLOCATE( out_eta(jpi,jpj,jpk,2*nb_ana), & 294 & out_u (jpi,jpj,jpk,2*nb_ana), & 295 & out_v (jpi,jpj,jpk,2*nb_ana), & 296 & out_w (jpi,jpj,jpk,2*nb_ana), & 297 & out_dzi(jpi,jpj,jpk,2*nb_ana) ) 298 # else 299 ALLOCATE( out_eta(jpi,jpj,2*nb_ana), & 300 & out_u (jpi,jpj,2*nb_ana), & 301 & out_v (jpi,jpj,2*nb_ana) ) 302 # endif 303 304 IF(lwp) WRITE(numout,*) 'ANA F OLD', ft 305 IF(lwp) WRITE(numout,*) 'ANA U OLD', ut 306 IF(lwp) WRITE(numout,*) 'ANA V OLD', vt 307 235 308 ninco = 2*nb_ana 236 309 … … 260 333 CALL SUR_DETERMINE_INIT 261 334 262 ! Elevation: 335 ! Density and Elevation: 336 # if defined key_3Ddiaharm 337 DO jk=1,jpk 338 # endif 263 339 DO jj = 2, jpjm1 264 340 DO ji = 2, jpim1 265 341 266 342 ! Fill input array 343 # if defined key_3Ddiaharm 344 ztmp4(1:nb_ana*2) = ana_temp(ji,jj,1:nb_ana*2,1,jk) 345 # else 267 346 ztmp4(1:nb_ana*2) = ana_temp(ji,jj,1:nb_ana*2,1) 347 # endif 268 348 CALL SUR_DETERMINE 269 349 270 350 ! Fill output array 271 351 DO jh = 1, nb_ana 272 out_eta(ji,jj,jh ) = ztmp7((jh-1)*2+1) * ssmask(ji,jj) 273 out_eta(ji,jj,jh+nb_ana) = -ztmp7((jh-1)*2+2) * ssmask(ji,jj) 352 # if defined key_3Ddiaharm 353 out_eta(ji,jj,jk,jh ) = ztmp7((jh-1)*2+1) * ssmask(ji,jj) 354 out_eta(ji,jj,jk,jh+nb_ana) = -ztmp7((jh-1)*2+2) * ssmask(ji,jj) 355 # else 356 out_eta(ji,jj, jh ) = ztmp7((jh-1)*2+1) * ssmask(ji,jj) 357 out_eta(ji,jj, jh+nb_ana) = -ztmp7((jh-1)*2+2) * ssmask(ji,jj) 358 # endif 274 359 END DO 275 360 END DO … … 281 366 282 367 ! Fill input array 368 # if defined key_3Ddiaharm 369 ztmp4(1:nb_ana*2) = ana_temp(ji,jj,1:nb_ana*2,2,jk) 370 # else 283 371 ztmp4(1:nb_ana*2) = ana_temp(ji,jj,1:nb_ana*2,2) 372 # endif 284 373 CALL SUR_DETERMINE 285 374 286 375 ! Fill output array 287 376 DO jh = 1, nb_ana 288 out_u(ji,jj, jh) = ztmp7((jh-1)*2+1) * ssumask(ji,jj) 289 out_u(ji,jj,nb_ana+jh) = -ztmp7((jh-1)*2+2) * ssumask(ji,jj) 377 # if defined key_3Ddiaharm 378 out_u(ji,jj,jk, jh) = ztmp7((jh-1)*2+1) * ssumask(ji,jj) 379 out_u(ji,jj,jk,nb_ana+jh) = -ztmp7((jh-1)*2+2) * ssumask(ji,jj) 380 # else 381 out_u(ji,jj, jh) = ztmp7((jh-1)*2+1) * ssumask(ji,jj) 382 out_u(ji,jj, nb_ana+jh) = -ztmp7((jh-1)*2+2) * ssumask(ji,jj) 383 # endif 384 290 385 END DO 291 386 … … 298 393 299 394 ! Fill input array 395 # if defined key_3Ddiaharm 396 ztmp4(1:nb_ana*2) = ana_temp(ji,jj,1:nb_ana*2,3,jk) 397 # else 300 398 ztmp4(1:nb_ana*2) = ana_temp(ji,jj,1:nb_ana*2,3) 399 # endif 301 400 CALL SUR_DETERMINE 302 401 303 402 ! Fill output array 304 403 DO jh = 1, nb_ana 305 out_v(ji,jj, jh) = ztmp7((jh-1)*2+1) * ssvmask(ji,jj) 306 out_v(ji,jj,nb_ana+jh) = -ztmp7((jh-1)*2+2) * ssvmask(ji,jj) 404 # if defined key_3Ddiaharm 405 out_v(ji,jj,jk, jh) = ztmp7((jh-1)*2+1) * ssvmask(ji,jj) 406 out_v(ji,jj,jk,nb_ana+jh) = -ztmp7((jh-1)*2+2) * ssvmask(ji,jj) 407 # else 408 out_v(ji,jj, jh) = ztmp7((jh-1)*2+1) * ssvmask(ji,jj) 409 out_v(ji,jj, nb_ana+jh) = -ztmp7((jh-1)*2+2) * ssvmask(ji,jj) 410 # endif 307 411 END DO 308 412 309 413 END DO 310 414 END DO 415 416 # if defined key_3Ddiaharm 417 ! w- velocity 418 DO jj = 1, jpj 419 DO ji = 1, jpi 420 ! Fill input array 421 kun=0 422 DO jh = 1,nb_ana 423 DO jc = 1,2 424 kun = kun + 1 425 ztmp4(kun)=ana_temp(ji,jj,kun,4,jk) 426 END DO 427 END DO 428 429 CALL SUR_DETERMINE(jj+1) 430 431 ! Fill output array 432 DO jh = 1, nb_ana 433 ana_amp(ji,jj,jh,1)=ztmp7((jh-1)*2+1) 434 ana_amp(ji,jj,jh,2)=ztmp7((jh-1)*2+2) 435 END DO 436 437 END DO 438 END DO 439 440 DO jj = 1, jpj 441 DO ji = 1, jpi 442 DO jh = 1, nb_ana 443 X1=ana_amp(ji,jj,jh,1) 444 X2=-ana_amp(ji,jj,jh,2) 445 out_w(ji,jj,jk, jh)=X1 * tmask_i(ji,jj) 446 out_w(ji,jj,jk,nb_ana+jh)=X2 * tmask_i(ji,jj) 447 END DO 448 END DO 449 END DO 450 451 ! dzi- isopycnal displacements 452 DO jj = 1, jpj 453 DO ji = 1, jpi 454 ! Fill input array 455 kun=0 456 DO jh = 1,nb_ana 457 DO jc = 1,2 458 kun = kun + 1 459 ztmp4(kun)=ana_temp(ji,jj,kun,5,jk) 460 END DO 461 END DO 462 463 CALL SUR_DETERMINE(jj+1) 464 465 ! Fill output array 466 DO jh = 1, nb_ana 467 ana_amp(ji,jj,jh,1)=ztmp7((jh-1)*2+1) 468 ana_amp(ji,jj,jh,2)=ztmp7((jh-1)*2+2) 469 END DO 470 471 END DO 472 END DO 473 474 DO jj = 1, jpj 475 DO ji = 1, jpi 476 DO jh = 1, nb_ana 477 X1=ana_amp(ji,jj,jh,1) 478 X2=-ana_amp(ji,jj,jh,2) 479 out_dzi(ji,jj,jk, jh)=X1 * tmask_i(ji,jj) 480 out_dzi(ji,jj,jk,nb_ana+jh)=X2 * tmask_i(ji,jj) 481 END DO 482 END DO 483 END DO 484 485 ENDDO ! jk 486 # endif 311 487 ! 312 488 ! clem: we could avoid this call if all the loops were from 1:jpi and 1:jpj … … 328 504 !!-------------------------------------------------------------------- 329 505 INTEGER :: jh 330 !!---------------------------------------------------------------------- 506 507 # if defined key_3Ddiaharm 508 CHARACTER(LEN=lc) :: cdfile_name_W ! name of the file created (W-points) 509 INTEGER :: jk 510 REAL(WP), ALLOCATABLE, DIMENSION (:,:,:) :: z3real, z3im 511 REAL(WP), ALLOCATABLE, DIMENSION (:,:) :: z2real, z2im 512 # endif 513 !!---------------------------------------------------------------------- 514 515 #if defined key_dimgout 516 cdfile_name_T = TRIM(cexper)//'_Tidal_harmonics_gridT.dimgproc' 517 cdfile_name_U = TRIM(cexper)//'_Tidal_harmonics_gridU.dimgproc' 518 cdfile_name_V = TRIM(cexper)//'_Tidal_harmonics_gridV.dimgproc' 519 # if defined key_3Ddiaharm 520 cdfile_name_W = TRIM(cexper)//'_Tidal_harmonics_gridW.dimgproc' 521 # endif 522 #endif 331 523 332 524 IF(lwp) WRITE(numout,*) ' ' 333 525 IF(lwp) WRITE(numout,*) 'dia_wri_harm : Write harmonic analysis results' 526 #if defined key_dimgout 527 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~ Output files: ', TRIM(cdfile_name_T) 528 IF(lwp) WRITE(numout,*) ' ', TRIM(cdfile_name_U) 529 IF(lwp) WRITE(numout,*) ' ', TRIM(cdfile_name_V) 530 # if defined key_3Ddiaharm 531 IF(lwp) WRITE(numout,*) ' ', TRIM(cdfile_name_W) 532 # endif 533 #endif 334 534 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 335 535 336 ! A) Elevation 536 # if defined key_3Ddiaharm 537 ALLOCATE(z3real(jpi,jpj,jpk),z3im(jpi,jpj,jpk),z2real(jpi,jpj),z2im(jpi,jpj)) 538 # endif 539 540 ! A) density and elevation 337 541 !///////////// 542 #if defined key_dimgout 543 cltext='density amplitude and phase; elevation is level=jpk ' 544 CALL dia_wri_dimg(TRIM(cdfile_name_T), TRIM(cltext), out_eta, 2*nb_ana, '2') 545 #else 546 # if defined key_3Ddiaharm 547 z3real(:,:,:) = 0._wp; z3im(:,:,:) = 0._wp 548 # endif 338 549 DO jh = 1, nb_ana 550 # if defined key_3Ddiaharm 551 DO jk=1,jpkm1 552 z3real(:,:,jk)=out_eta(:,:,jk,jh) 553 z3im (:,:,jk)=out_eta(:,:,jk,jh+nb_ana) 554 ENDDO 555 z2real(:,:)=out_eta(:,:,jpk,jh); z2im(:,:)=out_eta(:,:,jpk,jh+nb_ana) 556 CALL iom_put( TRIM(tname(jh))//'x_ro', z3real(:,:,:) ) 557 CALL iom_put( TRIM(tname(jh))//'y_ro', z3im (:,:,:) ) 558 CALL iom_put( TRIM(tname(jh))//'x' , z2real(:,: ) ) 559 CALL iom_put( TRIM(tname(jh))//'y' , z2im (:,: ) ) 560 # else 561 WRITE(numout,*) "OUTPUT ORI: ", TRIM(tname(jh))//'x', ' & ', TRIM(tname(jh))//'y', MAXVAL(out_eta(:,:,jh)) 339 562 CALL iom_put( TRIM(tname(jh))//'x', out_eta(:,:,jh) ) 340 563 CALL iom_put( TRIM(tname(jh))//'y', out_eta(:,:,jh+nb_ana) ) 341 END DO 342 343 ! B) ubar 564 # endif 565 END DO 566 #endif 567 568 ! B) u 344 569 !///////// 570 #if defined key_dimgout 571 cltext='3d u amplitude and phase; ubar is the last level' 572 CALL dia_wri_dimg(TRIM(cdfile_name_U), TRIM(cltext), out_u, 2*nb_ana, '2') 573 #else 574 # if defined key_3Ddiaharm 575 z3real(:,:,:) = 0._wp; z3im(:,:,:) = 0._wp 576 # endif 345 577 DO jh = 1, nb_ana 346 CALL iom_put( TRIM(tname(jh))//'x_u', out_u(:,:,jh) ) 347 CALL iom_put( TRIM(tname(jh))//'y_u', out_u(:,:,jh+nb_ana) ) 348 END DO 349 350 ! C) vbar 578 # if defined key_3Ddiaharm 579 DO jk=1,jpkm1 580 z3real(:,:,jk)=out_u(:,:,jk,jh) 581 z3im (:,:,jk)=out_u(:,:,jk,jh+nb_ana) 582 ENDDO 583 z2real(:,:)=out_u(:,:,jpk,jh); z2im(:,:)=out_u(:,:,jpk,jh+nb_ana) 584 CALL iom_put( TRIM(tname(jh))//'x_u3d', z3real(:,:,:) ) 585 CALL iom_put( TRIM(tname(jh))//'y_u3d', z3im (:,:,:) ) 586 CALL iom_put( TRIM(tname(jh))//'x_u2d', z2real(:,:) ) 587 CALL iom_put( TRIM(tname(jh))//'y_u2d', z2im (:,:) ) 588 z2real(:,:)=out_w(:,:,jpk,jh); z2im(:,:)=out_w(:,:,jpk,jh+nb_ana) 589 CALL iom_put( TRIM(tname(jh))//'x_tabx', z2real(:,:) ) 590 CALL iom_put( TRIM(tname(jh))//'y_tabx', z2im (:,:) ) 591 # else 592 CALL iom_put( TRIM(tname(jh))//'x_u2d', out_u(:,:,jh) ) 593 CALL iom_put( TRIM(tname(jh))//'y_u2d', out_u(:,:,nb_ana+jh) ) 594 # endif 595 END DO 596 #endif 597 598 ! C) v 351 599 !///////// 600 #if defined key_dimgout 601 cltext='3d v amplitude and phase; vbar is the last level' 602 CALL dia_wri_dimg(TRIM(cdfile_name_V), TRIM(cltext), out_v, 2*nb_ana, '2') 603 #else 604 # if defined key_3Ddiaharm 605 z3real(:,:,:) = 0._wp; z3im(:,:,:) = 0._wp 606 # endif 352 607 DO jh = 1, nb_ana 353 CALL iom_put( TRIM(tname(jh))//'x_v', out_v(:,:,jh ) ) 354 CALL iom_put( TRIM(tname(jh))//'y_v', out_v(:,:,jh+nb_ana) ) 355 END DO 356 ! 608 # if defined key_3Ddiaharm 609 DO jk=1,jpkm1 610 z3real(:,:,jk)=out_v(:,:,jk,jh) 611 z3im (:,:,jk)=out_v(:,:,jk,jh+nb_ana) 612 ENDDO 613 z2real(:,:)=out_v(:,:,jpk,jh); z2im(:,:)=out_v(:,:,jpk,jh+nb_ana) 614 CALL iom_put( TRIM(tname(jh))//'x_v3d', z3real(:,:,:) ) 615 CALL iom_put( TRIM(tname(jh))//'y_v3d', z3im (:,:,:) ) 616 CALL iom_put( TRIM(tname(jh))//'x_v2d' , z2real(:,:) ) 617 CALL iom_put( TRIM(tname(jh))//'y_v2d' , z2im (:,:) ) 618 z2real(:,:)=out_dzi(:,:,jpk,jh); z2im(:,:)=out_dzi(:,:,jpk,jh+nb_ana) 619 CALL iom_put( TRIM(tname(jh))//'x_taby', z2real(:,:) ) 620 CALL iom_put( TRIM(tname(jh))//'y_taby', z2im (:,:) ) 621 # else 622 CALL iom_put( TRIM(tname(jh))//'x_v2d', out_v(:,:,jh ) ) 623 CALL iom_put( TRIM(tname(jh))//'y_v2d', out_v(:,:,jh+nb_ana) ) 624 # endif 625 END DO 626 627 #endif 628 ! D) w 629 # if defined key_3Ddiaharm 630 # if defined key_dimgout 631 cltext='3d w amplitude and phase; vort_baro is the last level' 632 CALL dia_wri_dimg(TRIM(cdfile_name_W), TRIM(cltext), out_w, 2*nb_ana, '2') 633 # else 634 DO jh = 1, nb_ana 635 DO jk=1,jpkm1 636 z3real(:,:,jk)=out_w(:,:,jk,jh) 637 z3im(:,:,jk)=out_w(:,:,jk,jh+nb_ana) 638 ENDDO 639 CALL iom_put( TRIM(tname(jh))//'x_w3d', z3real(:,:,:) ) 640 CALL iom_put( TRIM(tname(jh))//'y_w3d', z3im(:,:,:) ) 641 END DO 642 # endif 643 644 ! E) dzi + tau_bot 645 # if defined key_dimgout 646 cltext='dzi=g*ro/N2 amplitude and phase' 647 CALL dia_wri_dimg(TRIM(cdfile_name_W), TRIM(cltext), out_w, 2*nb_ana, '2') 648 # else 649 DO jh = 1, nb_ana 650 DO jk=1,jpkm1 651 z3real(:,:,jk)=out_dzi(:,:,jk,jh) 652 z3im(:,:,jk)=out_dzi(:,:,jk,jh+nb_ana) 653 ENDDO 654 CALL iom_put( TRIM(tname(jh))//'x_dzi', z3real(:,:,:) ) 655 CALL iom_put( TRIM(tname(jh))//'y_dzi', z3im(:,:,:) ) 656 END DO 657 # endif 658 # endif 659 660 ! 661 # if defined key_3Ddiaharm 662 DEALLOCATE(z3real, z3im, z2real,z2im) 663 # endif 664 357 665 END SUBROUTINE dia_wri_harm 358 666 -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/DOM/dtatsd.F90
r14075 r15422 21 21 ! 22 22 USE in_out_manager ! I/O manager 23 USE iom ! IOM library 23 24 USE lib_mpp ! MPP library 24 25 … … 31 32 ! !!* namtsd namelist : Temperature & Salinity Data * 32 33 LOGICAL , PUBLIC :: ln_tsd_init !: T & S data flag 34 LOGICAL , PUBLIC :: ln_tsd_interp !: vertical interpolation flag 33 35 LOGICAL , PUBLIC :: ln_tsd_dmp !: internal damping toward input data flag 34 36 35 37 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf_tsd ! structure of input SST (file informations, fields read) 38 INTEGER :: jpk_init, inum_dta 39 INTEGER :: id, linum ! local integers 40 INTEGER :: zdim(4) 36 41 37 42 !!---------------------------------------------------------------------- … … 53 58 LOGICAL, INTENT(in), OPTIONAL :: ld_tradmp ! force the initialization when tradp is used 54 59 ! 55 INTEGER :: ios, ierr0, ierr1, ierr2, ierr3 ! local integers60 INTEGER :: ios, ierr0, ierr1, ierr2, ierr3, ierr4, ierr5 ! local integers 56 61 !! 57 62 CHARACTER(len=100) :: cn_dir ! Root directory for location of ssr files 58 TYPE(FLD_N), DIMENSION( jpts) :: slf_i ! array of namelist informations on the fields to read59 TYPE(FLD_N) :: sn_tem, sn_sal 63 TYPE(FLD_N), DIMENSION(jpts+2):: slf_i ! array of namelist informations on the fields to read 64 TYPE(FLD_N) :: sn_tem, sn_sal, sn_dep, sn_msk 60 65 !! 61 66 NAMELIST/namtsd/ ln_tsd_init, ln_tsd_dmp, cn_dir, sn_tem, sn_sal … … 63 68 ! 64 69 ! Initialisation 65 ierr0 = 0 ; ierr1 = 0 ; ierr2 = 0 ; ierr3 = 0 70 ierr0 = 0 ; ierr1 = 0 ; ierr2 = 0 ; ierr3 = 0 ; ierr4 = 0 ; ierr5 = 0 66 71 ! 67 72 REWIND( numnam_ref ) ! Namelist namtsd in reference namelist : … … 80 85 WRITE(numout,*) '~~~~~~~~~~~~ ' 81 86 WRITE(numout,*) ' Namelist namtsd' 82 WRITE(numout,*) ' Initialisation of ocean T & S with T &S input data ln_tsd_init = ', ln_tsd_init 83 WRITE(numout,*) ' damping of ocean T & S toward T &S input data ln_tsd_dmp = ', ln_tsd_dmp 87 WRITE(numout,*) ' Initialisation of ocean T & S with T & S input data ln_tsd_init = ', ln_tsd_init 88 WRITE(numout,*) ' Interpolation of initial conditions in the vertical ln_tsd_interp = ', ln_tsd_interp 89 WRITE(numout,*) ' damping of ocean T & S toward T & S input data ln_tsd_dmp = ', ln_tsd_dmp 84 90 WRITE(numout,*) 85 91 IF( .NOT.ln_tsd_init .AND. .NOT.ln_tsd_dmp ) THEN … … 94 100 ln_tsd_init = .FALSE. 95 101 ENDIF 102 IF( ln_tsd_interp .AND. ln_tsd_dmp ) THEN 103 CALL ctl_stop( 'dta_tsd_init: Tracer damping and vertical interpolation not yet configured' ) ; RETURN 104 ENDIF 105 IF( ln_tsd_interp .AND. LEN(TRIM(sn_msk%wname)) > 0 ) THEN 106 CALL ctl_stop( 'dta_tsd_init: Using vertical interpolation and weights files not recommended' ) ; RETURN 107 ENDIF 96 108 ! 97 109 ! ! allocate the arrays (if necessary) 98 110 IF( ln_tsd_init .OR. ln_tsd_dmp ) THEN 99 111 ! 100 ALLOCATE( sf_tsd(jpts), STAT=ierr0 ) 112 IF( ln_tsd_interp ) THEN 113 ALLOCATE( sf_tsd(jpts+2), STAT=ierr0 ) ! to carry the addtional depth information 114 ELSE 115 ALLOCATE( sf_tsd(jpts ), STAT=ierr0 ) 116 ENDIF 101 117 IF( ierr0 > 0 ) THEN 102 118 CALL ctl_stop( 'dta_tsd_init: unable to allocate sf_tsd structure' ) ; RETURN 103 119 ENDIF 104 120 ! 105 ALLOCATE( sf_tsd(jp_tem)%fnow(jpi,jpj,jpk) , STAT=ierr0 ) 106 IF( sn_tem%ln_tint ) ALLOCATE( sf_tsd(jp_tem)%fdta(jpi,jpj,jpk,2) , STAT=ierr1 ) 107 ALLOCATE( sf_tsd(jp_sal)%fnow(jpi,jpj,jpk) , STAT=ierr2 ) 108 IF( sn_sal%ln_tint ) ALLOCATE( sf_tsd(jp_sal)%fdta(jpi,jpj,jpk,2) , STAT=ierr3 ) 109 ! 110 IF( ierr0 + ierr1 + ierr2 + ierr3 > 0 ) THEN 121 IF( ln_tsd_interp ) THEN 122 CALL iom_open ( trim(cn_dir) // trim(sn_dep%clname), inum_dta ) 123 id = iom_varid( inum_dta, sn_dep%clvar, zdim ) 124 jpk_init = zdim(3) 125 IF(lwp) WRITE(numout,*) 'Dimension of vertical coordinate in ICs: ', jpk_init 126 CALL iom_close( inum_dta ) ! Close the input file 127 ! 128 ALLOCATE( sf_tsd(jp_tem)%fnow(jpi,jpj,jpk_init ) , STAT=ierr0 ) 129 IF( sn_tem%ln_tint ) ALLOCATE( sf_tsd(jp_tem)%fdta(jpi,jpj,jpk_init,2) , STAT=ierr1 ) 130 ALLOCATE( sf_tsd(jp_sal)%fnow(jpi,jpj,jpk_init ) , STAT=ierr2 ) 131 IF( sn_sal%ln_tint ) ALLOCATE( sf_tsd(jp_sal)%fdta(jpi,jpj,jpk_init,2) , STAT=ierr3 ) 132 ALLOCATE( sf_tsd(jp_dep)%fnow(jpi,jpj,jpk_init ) , STAT=ierr4 ) 133 ALLOCATE( sf_tsd(jp_msk)%fnow(jpi,jpj,jpk_init ) , STAT=ierr5 ) 134 ELSE 135 ALLOCATE( sf_tsd(jp_tem)%fnow(jpi,jpj,jpk) , STAT=ierr0 ) 136 IF( sn_tem%ln_tint ) ALLOCATE( sf_tsd(jp_tem)%fdta(jpi,jpj,jpk,2) , STAT=ierr1 ) 137 ALLOCATE( sf_tsd(jp_sal)%fnow(jpi,jpj,jpk) , STAT=ierr2 ) 138 IF( sn_sal%ln_tint ) ALLOCATE( sf_tsd(jp_sal)%fdta(jpi,jpj,jpk,2) , STAT=ierr3 ) 139 ENDIF ! ln_tsd_interp 140 ! 141 IF( ierr0 + ierr1 + ierr2 + ierr3 + ierr4 + ierr5 > 0 ) THEN 111 142 CALL ctl_stop( 'dta_tsd : unable to allocate T & S data arrays' ) ; RETURN 112 143 ENDIF 113 144 ! ! fill sf_tsd with sn_tem & sn_sal and control print 114 145 slf_i(jp_tem) = sn_tem ; slf_i(jp_sal) = sn_sal 146 IF( ln_tsd_interp ) slf_i(jp_dep) = sn_dep ; slf_i(jp_msk) = sn_msk 115 147 CALL fld_fill( sf_tsd, slf_i, cn_dir, 'dta_tsd', 'Temperature & Salinity data', 'namtsd', no_print ) 116 148 ! … … 138 170 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT( out) :: ptsd ! T & S data 139 171 ! 140 INTEGER :: ji, jj, jk, jl, jk k ! dummy loop indicies141 INTEGER :: ik, il0, il1, ii0, ii1, ij0, ij1 ! local integers172 INTEGER :: ji, jj, jk, jl, jk_init ! dummy loop indices 173 INTEGER :: ik, il0, il1, ii0, ii1, ij0, ij1 ! local integers 142 174 REAL(wp):: zl, zi ! local scalars 143 REAL(wp), DIMENSION(jpk) :: ztp, zsp ! 1D workspace144 175 !!---------------------------------------------------------------------- 145 176 ! … … 176 207 !!gm end 177 208 ! 178 ptsd(:,:,:,jp_tem) = sf_tsd(jp_tem)%fnow(:,:,:) ! NO mask 179 ptsd(:,:,:,jp_sal) = sf_tsd(jp_sal)%fnow(:,:,:) 180 ! 181 IF( ln_sco ) THEN !== s- or mixed s-zps-coordinate ==! 182 ! 183 IF( kt == nit000 .AND. lwp )THEN 184 WRITE(numout,*) 185 WRITE(numout,*) 'dta_tsd: interpolates T & S data onto the s- or mixed s-z-coordinate mesh' 186 ENDIF 187 ! 188 DO jj = 1, jpj ! vertical interpolation of T & S 189 DO ji = 1, jpi 190 DO jk = 1, jpk ! determines the intepolated T-S profiles at each (i,j) points 209 IF( kt == nit000 .AND. lwp )THEN 210 WRITE(numout,*) 211 WRITE(numout,*) 'dta_tsd: interpolates T & S data onto current mesh' 212 ENDIF 213 ! 214 IF( ln_tsd_interp ) THEN ! probably should use pointers in the following to make more readable 215 ! 216 DO jk = 1, jpk ! determines the intepolated T-S profiles at each (i,j) points 217 DO jj= 1, jpj 218 DO ji= 1, jpi 191 219 zl = gdept_0(ji,jj,jk) 192 IF( zl < gdept_1d(1 ) ) THEN ! above the first level of data 193 ztp(jk) = ptsd(ji,jj,1 ,jp_tem) 194 zsp(jk) = ptsd(ji,jj,1 ,jp_sal) 195 ELSEIF( zl > gdept_1d(jpk) ) THEN ! below the last level of data 196 ztp(jk) = ptsd(ji,jj,jpkm1,jp_tem) 197 zsp(jk) = ptsd(ji,jj,jpkm1,jp_sal) 198 ELSE ! inbetween : vertical interpolation between jkk & jkk+1 199 DO jkk = 1, jpkm1 ! when gdept(jkk) < zl < gdept(jkk+1) 200 IF( (zl-gdept_1d(jkk)) * (zl-gdept_1d(jkk+1)) <= 0._wp ) THEN 201 zi = ( zl - gdept_1d(jkk) ) / (gdept_1d(jkk+1)-gdept_1d(jkk)) 202 ztp(jk) = ptsd(ji,jj,jkk,jp_tem) + ( ptsd(ji,jj,jkk+1,jp_tem) - ptsd(ji,jj,jkk,jp_tem) ) * zi 203 zsp(jk) = ptsd(ji,jj,jkk,jp_sal) + ( ptsd(ji,jj,jkk+1,jp_sal) - ptsd(ji,jj,jkk,jp_sal) ) * zi 220 IF( zl < sf_tsd(jp_dep)%fnow(ji,jj,1) ) THEN ! above the first level of data 221 ptsd(ji,jj,jk,jp_tem) = sf_tsd(jp_tem)%fnow(ji,jj,1) 222 ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,1) 223 ELSEIF( zl > sf_tsd(jp_dep)%fnow(ji,jj,jpk_init) ) THEN ! below the last level of data 224 ptsd(ji,jj,jk,jp_tem) = sf_tsd(jp_tem)%fnow(ji,jj,jpk_init) 225 ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,jpk_init) 226 ELSE ! inbetween : vertical interpolation between jk_init & jk_init+1 227 DO jk_init = 1, jpk_init-1 ! when gdept(jk_init) < zl < gdept(jk_init+1) 228 IF( sf_tsd(jp_msk)%fnow(ji,jj,jk_init+1) == 0 ) THEN ! if there is no data fill down 229 sf_tsd(jp_tem)%fnow(ji,jj,jk_init+1) = sf_tsd(jp_tem)%fnow(ji,jj,jk_init) 230 sf_tsd(jp_sal)%fnow(ji,jj,jk_init+1) = sf_tsd(jp_sal)%fnow(ji,jj,jk_init) 231 ENDIF 232 IF((zl-sf_tsd(jp_dep)%fnow(ji,jj,jk_init)) * (zl-sf_tsd(jp_dep)%fnow(ji,jj,jk_init+1)) <= 0._wp ) THEN 233 zi = ( zl - sf_tsd(jp_dep)%fnow(ji,jj,jk_init) ) / & 234 & (sf_tsd(jp_dep)%fnow(ji,jj,jk_init+1)-sf_tsd(jp_dep)%fnow(ji,jj,jk_init)) 235 ptsd(ji,jj,jk,jp_tem) = sf_tsd(jp_tem)%fnow(ji,jj,jk_init) + & 236 & (sf_tsd(jp_tem)%fnow(ji,jj,jk_init+1)-sf_tsd(jp_tem)%fnow(ji,jj,jk_init)) * zi 237 ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,jk_init) + & 238 & (sf_tsd(jp_sal)%fnow(ji,jj,jk_init+1)-sf_tsd(jp_sal)%fnow(ji,jj,jk_init)) * zi 204 239 ENDIF 205 240 END DO 206 241 ENDIF 207 242 END DO 208 DO jk = 1, jpkm1209 ptsd(ji,jj,jk,jp_tem) = ztp(jk) * tmask(ji,jj,jk) ! mask required for mixed zps-s-coord210 ptsd(ji,jj,jk,jp_sal) = zsp(jk) * tmask(ji,jj,jk)211 END DO212 ptsd(ji,jj,jpk,jp_tem) = 0._wp213 ptsd(ji,jj,jpk,jp_sal) = 0._wp214 243 END DO 215 244 END DO 216 245 ! 246 ptsd(:,:,:,jp_tem) = ptsd(:,:,:,jp_tem) *tmask(:,:,:) 247 ptsd(:,:,:,jp_sal) = ptsd(:,:,:,jp_sal) *tmask(:,:,:) 217 248 ELSE !== z- or zps- coordinate ==! 218 249 ! 219 ptsd(:,:,:,jp_tem) = ptsd(:,:,:,jp_tem) * tmask(:,:,:) ! Mask220 ptsd(:,:,:,jp_sal) = ptsd(:,:,:,jp_sal) * tmask(:,:,:)250 ptsd(:,:,:,jp_tem) = sf_tsd(jp_tem)%fnow(:,:,:) * tmask(:,:,:) ! Mask 251 ptsd(:,:,:,jp_sal) = sf_tsd(jp_sal)%fnow(:,:,:) * tmask(:,:,:) 221 252 ! 222 253 IF( ln_zps ) THEN ! zps-coordinate (partial steps) interpolation at the last ocean level … … 248 279 DEALLOCATE( sf_tsd(jp_sal)%fnow ) ! S arrays in the structure 249 280 IF( sf_tsd(jp_sal)%ln_tint ) DEALLOCATE( sf_tsd(jp_sal)%fdta ) 281 IF( ln_tsd_interp ) DEALLOCATE( sf_tsd(jp_dep)%fnow ) ! T arrays in the structure 282 IF( ln_tsd_interp ) DEALLOCATE( sf_tsd(jp_msk)%fnow ) ! T arrays in the structure 250 283 DEALLOCATE( sf_tsd ) ! the structure itself 251 284 ENDIF -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/DYN/dynnxt.F90
r14075 r15422 33 33 USE dynadv ! dynamics: vector invariant versus flux form 34 34 USE dynspg_ts ! surface pressure gradient: split-explicit scheme 35 USE dynspg 35 36 USE domvvl ! variable volume 36 37 USE bdy_oce , ONLY : ln_bdy … … 58 59 PUBLIC dyn_nxt ! routine called by step.F90 59 60 61 !! Substitution 62 # include "vectopt_loop_substitute.h90" 60 63 !!---------------------------------------------------------------------- 61 64 !! NEMO/OCE 4.0 , NEMO Consortium (2018) … … 182 185 vn(:,:,jk) = va(:,:,jk) 183 186 END DO 187 ! limit velocities 188 IF (ln_ulimit) THEN 189 call dyn_limit_velocity (kt) 190 ENDIF 184 191 IF( .NOT.ln_linssh ) THEN ! e3._b <-- e3._n 185 192 !!gm BUG ???? I don't understand why it is not : e3._n <-- e3._a … … 214 221 END DO 215 222 END DO 223 ! limit velocities 224 IF (ln_ulimit) THEN 225 call dyn_limit_velocity (kt) 226 ENDIF 216 227 ! ! ================! 217 228 ELSE ! Variable volume ! … … 263 274 END DO 264 275 END DO 276 ! limit velocities 277 IF (ln_ulimit) THEN 278 call dyn_limit_velocity (kt) 279 ENDIF 265 280 ! 266 281 ELSE ! Asselin filter applied on thickness weighted velocity … … 290 305 END DO 291 306 END DO 307 ! limit velocities 308 IF (ln_ulimit) THEN 309 call dyn_limit_velocity (kt) 310 ENDIF 292 311 e3u_b(:,:,1:jpkm1) = ze3u_f(:,:,1:jpkm1) ! e3u_b <-- filtered scale factor 293 312 e3v_b(:,:,1:jpkm1) = ze3v_f(:,:,1:jpkm1) … … 401 420 END SUBROUTINE dyn_nxt 402 421 422 SUBROUTINE dyn_limit_velocity (kt) 423 ! limits maximum values of un and vn by volume courant number 424 INTEGER, INTENT( in ) :: kt ! ocean time-step index 425 ! 426 INTEGER :: ji, jj, jk ! dummy loop indices 427 REAL(wp) :: zzu,zplim,zmlim,isp,ism,zcn,ze3e1,zzcn,zcnn,idivp,idivm 428 429 ! limit fluxes 430 zcn =cn_ulimit !0.9 ! maximum velocity inverse courant number 431 zcnn = cnn_ulimit !0.54 ! how much to reduce cn by in divergen flow 432 433 DO jk = 1, jpkm1 434 DO jj = 2, jpjm1 435 DO ji = fs_2, fs_jpim1 ! vect. opt. 436 ! U direction 437 zzu = un(ji,jj,jk) 438 ze3e1 = e3u_n(ji ,jj,jk) * e2u(ji ,jj) 439 ! ips is 1 if flow is positive othersize zero 440 isp = 0.5 * (sign(1.0_wp,zzu) + 1.0_wp ) 441 ism = -0.5 * (sign(1.0_wp,zzu) - 1.0_wp ) 442 ! idev is 1 if divergent flow otherwise zero 443 idivp = -isp * 0.5 * (sign(1.0_wp, un(ji-1,jj,jk)) - 1.0_wp ) 444 idivm = ism * 0.5 * (sign(1.0_wp, un(ji+1,jj,jk)) + 1.0_wp ) 445 zzcn = (idivp+idivm)*(zcnn-1.0_wp)+1.0_wp 446 zzcn = zzcn * zcn 447 zplim = zzcn * (e3t_n(ji ,jj,jk) * e1t(ji ,jj) * e2t(ji ,jj)) / & 448 (2.0*rdt * ze3e1)*umask(ji,jj,jk) 449 zmlim = -zzcn * (e3t_n(ji+1,jj,jk) * e1t(ji+1,jj) * e2t(ji+1,jj)) / & 450 (2.0*rdt * ze3e1)*umask(ji,jj,jk) 451 ! limit currents 452 un(ji,jj,jk) = min ( zzu,zplim) * isp + max(zzu,zmlim) *ism 453 454 ! V direction 455 zzu = vn(ji,jj,jk) 456 ze3e1 = e3v_n(ji ,jj,jk) * e1v(ji ,jj) 457 isp = 0.5 * (sign(1.0_wp,zzu) + 1.0_wp ) 458 ism = -0.5 * (sign(1.0_wp,zzu) - 1.0_wp ) 459 ! idev is 1 if divergent flow otherwise zero 460 idivp = -isp * 0.5 * (sign(1.0_wp, vn(ji,jj-1,jk)) - 1.0_wp ) 461 idivm = ism * 0.5 * (sign(1.0_wp, vn(ji,jj+1,jk)) + 1.0_wp ) 462 zzcn = (idivp+idivm)*(zcnn-1.0_wp)+1.0_wp 463 zzcn = zzcn * zcn 464 zplim = zzcn * (e3t_n(ji,jj ,jk) * e1t(ji,jj ) * e2t(ji,jj )) / & 465 (2.0*rdt * ze3e1)*vmask(ji,jj,jk) 466 zmlim = -zzcn * (e3t_n(ji,jj+1,jk) * e1t(ji,jj+1) * e2t(ji,jj+1)) / & 467 (2.0*rdt * ze3e1)*vmask(ji,jj,jk) 468 ! limit currents 469 vn(ji,jj,jk) = min ( zzu,zplim) * isp + max(zzu,zmlim) *ism 470 ENDDO 471 ENDDO 472 ENDDO 473 CALL lbc_lnk_multi( 'dynnxt', un(:,:,:), 'U', -1., vn(:,:,:), 'V', -1. ) 474 475 END SUBROUTINE dyn_limit_velocity 476 403 477 !!========================================================================= 404 478 END MODULE dynnxt -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/DYN/dynspg.F90
r14075 r15422 39 39 INTEGER :: nspg = 0 ! type of surface pressure gradient scheme defined from lk_dynspg_... 40 40 41 LOGICAL, PUBLIC :: ln_ulimit 42 REAL(wp), PUBLIC :: cn_ulimit,cnn_ulimit 43 41 44 ! ! Parameter to control the surface pressure gradient scheme 42 45 INTEGER, PARAMETER :: np_TS = 1 ! split-explicit time stepping (Time-Splitting) … … 191 194 NAMELIST/namdyn_spg/ ln_dynspg_exp , ln_dynspg_ts, & 192 195 & ln_bt_fw, ln_bt_av , ln_bt_auto , & 193 & nn_baro , rn_bt_cmax, nn_bt_flt, rn_bt_alpha 196 & nn_baro , rn_bt_cmax, nn_bt_flt, rn_bt_alpha, & 197 & ln_ulimit, cn_ulimit, cnn_ulimit 194 198 !!---------------------------------------------------------------------- 195 199 ! … … 213 217 WRITE(numout,*) ' Explicit free surface ln_dynspg_exp = ', ln_dynspg_exp 214 218 WRITE(numout,*) ' Free surface with time splitting ln_dynspg_ts = ', ln_dynspg_ts 219 WRITE(numout,*) ' Limit velocities ln_ulimit = ', ln_ulimit 220 WRITE(numout,*) ' Limit velocities max inverse Courant number = ', cn_ulimit 221 WRITE(numout,*) ' Limit velocities multiplier for divergant flow = ', cnn_ulimit 215 222 ENDIF 216 223 ! ! Control of surface pressure gradient scheme options -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/SBC/fldread.F90
r14075 r15422 901 901 IF( zdepth(jk) <= pdta_read_z(jb,1,1) ) THEN ! above the first level of external data 902 902 pdta(jb,1,jk) = pdta_read(jb,1,1) 903 ELSEIF( zdepth(jk) > pdta_read_z(jb,1,ipkmax) ) THEN ! below the last level of external data /= FillValue904 pdta(jb,1,jk) = pdta_read(jb,1,ipkmax)903 ELSEIF( zdepth(jk) > pdta_read_z(jb,1,ipkmax) ) THEN ! below the last level of external data 904 pdta(jb,1,jk) = pdta_read(jb,1,ipkmax) 905 905 ELSE ! inbetween: vertical interpolation between jkb & jkb+1 906 906 DO jkb = 1, ipkmax-1 -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/SBC/sbctide.F90
r14075 r15422 16 16 USE ioipsl ! NetCDF IPSL library 17 17 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 18 USE bdytides 18 19 19 20 IMPLICIT NONE … … 108 109 109 110 DO jk = 1, nb_harmo 110 zcons = 0.7_wp * Wave(ntide(jk))%equitide * ftide(jk) 111 ! love number now provides in tide namelist 112 zcons = dn_love_number * Wave(ntide(jk))%equitide * ftide(jk) 111 113 DO ji = 1, jpi 112 114 DO jj = 1, jpj … … 119 121 IF ( Wave(ntide(jk))%nutide == 1 ) THEN ; zcs = zcons * SIN( 2._wp*zlat ) 120 122 ELSEIF( Wave(ntide(jk))%nutide == 2 ) THEN ; zcs = zcons * COS( zlat )**2 123 ! Add tide potential for long period tides 124 ELSEIF( Wave(ntide(jk))%nutide == 0 ) THEN ; zcs = zcons * (0.5_wp-1.5_wp*SIN(zlat)**2._wp) 121 125 ELSE ; zcs = 0._wp 122 126 ENDIF -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/SBC/tide_mod.F90
r14075 r15422 16 16 PUBLIC tide_init_Wave ! called by tideini and diaharm modules 17 17 18 INTEGER, PUBLIC, PARAMETER :: jpmax_harmo = 19!: maximum number of harmonic18 INTEGER, PUBLIC, PARAMETER :: jpmax_harmo = 31 !: maximum number of harmonic 19 19 20 20 TYPE, PUBLIC :: tide … … 41 41 42 42 SUBROUTINE tide_init_Wave 43 # if defined key_FES14_tides 44 # include "tide_FES14.h90" 45 # else 43 46 # include "tide.h90" 47 # endif 44 48 END SUBROUTINE tide_init_Wave 45 49 … … 331 335 zf = zf * zf1 * zf1 332 336 ! 337 CASE( 20 ) !== formule 20, compound waves ( 78 x 78 x 78 x 78 ) 338 zf1 = nodal_factort(78) 339 zf = zf1 * zf1 * zf1 * zf1 340 ! 333 341 CASE( 73 ) !== formule 73 334 342 zs = sin(sh_I) -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/SBC/tideini.F90
r14075 r15422 33 33 INTEGER , PUBLIC :: kt_tide !: 34 34 REAL(wp), PUBLIC :: rdttideramp !: 35 REAL(wp), PUBLIC :: dn_love_number !: 35 36 REAL(wp), PUBLIC :: rn_scal_load !: 36 37 CHARACTER(lc), PUBLIC :: cn_tide_load !: … … 54 55 ! 55 56 NAMELIST/nam_tide/ln_tide, ln_tide_pot, ln_scal_load, ln_read_load, cn_tide_load, & 56 & ln_tide_ramp, rn_scal_load, rdttideramp, clname57 & ln_tide_ramp, rn_scal_load, rdttideramp, dn_love_number, clname 57 58 !!---------------------------------------------------------------------- 58 59 ! … … 78 79 WRITE(numout,*) ' Read load potential from file ln_read_load = ', ln_read_load 79 80 WRITE(numout,*) ' Apply ramp on tides at startup ln_tide_ramp = ', ln_tide_ramp 81 WRITE(numout,*) ' dn_love_number = ', dn_love_number 80 82 WRITE(numout,*) ' Fraction of SSH used in scal. approx. rn_scal_load = ', rn_scal_load 81 83 WRITE(numout,*) ' Duration (days) of ramp rdttideramp = ', rdttideramp … … 99 101 END DO 100 102 ! 103 IF (ln_tide .and.lwp) WRITE(numout,*) 'nb_harmo = ', nb_harmo 104 101 105 ! Ensure that tidal components have been set in namelist_cfg 102 106 IF( nb_harmo == 0 ) CALL ctl_stop( 'tide_init : No tidal components set in nam_tide' ) -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/USR/usrdef_istate.F90
r14075 r15422 62 62 DO jj = 1, jpj 63 63 DO ji = 1, jpi 64 pts(ji,jj,jk,jp_tem) = ( ( 16. - 12. * TANH( (pdept(ji,jj,jk) - 400) / 700 ) ) & 65 & * (-TANH( (500. - pdept(ji,jj,jk)) / 150. ) + 1.) / 2. & 66 & + ( 15. * ( 1. - TANH( (pdept(ji,jj,jk)-50.) / 1500.) ) & 67 & - 1.4 * TANH((pdept(ji,jj,jk)-100.) / 100.) & 68 & + 7. * (1500. - pdept(ji,jj,jk) ) / 1500.) & 69 & * (-TANH( (pdept(ji,jj,jk) - 500.) / 150.) + 1.) / 2. ) * ptmask(ji,jj,jk) 70 71 pts(ji,jj,jk,jp_sal) = ( ( 36.25 - 1.13 * TANH( (pdept(ji,jj,jk) - 305) / 460 ) ) & 72 & * (-TANH((500. - pdept(ji,jj,jk)) / 150.) + 1.) / 2 & 73 & + ( 35.55 + 1.25 * (5000. - pdept(ji,jj,jk)) / 5000. & 74 & - 1.62 * TANH( (pdept(ji,jj,jk) - 60. ) / 650. ) & 75 & + 0.2 * TANH( (pdept(ji,jj,jk) - 35. ) / 100. ) & 76 & + 0.2 * TANH( (pdept(ji,jj,jk) - 1000.) / 5000.) ) & 77 & * (-TANH( (pdept(ji,jj,jk) - 500.) / 150.) + 1.) / 2 ) * ptmask(ji,jj,jk) 64 pts(ji,jj,jk,jp_tem) = 20._wp * ptmask(ji,jj,jk) 65 pts(ji,jj,jk,jp_sal) = 36.25_wp * ptmask(ji,jj,jk) 78 66 END DO 79 67 END DO -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/USR/usrdef_sbc.F90
r14075 r15422 1 1 MODULE usrdef_sbc 2 2 !!====================================================================== 3 !! *** MODULE usrdef_sbc ***4 !! 5 !! === GYRE configuration ===3 !! *** MODULE usrdef_sbc *** 4 !! 5 !! === AMM7_SURGE configuration === 6 6 !! 7 7 !! User defined : surface forcing of a user configuration 8 8 !!====================================================================== 9 9 !! History : 4.0 ! 2016-03 (S. Flavoni, G. Madec) user defined interface 10 !! 4.0 ! 2017-12 (C. O'Neill) add necessary options for surge work - either no fluxes 11 !! (for tide-only run) or wind and pressure only 10 12 !!---------------------------------------------------------------------- 11 13 12 14 !!---------------------------------------------------------------------- 13 !! usrdef_sbc : user defined surface bounday conditions in GYRE case15 !! usrdef_sbc : user defined surface bounday conditions in LOCK_EXCHANGE case 14 16 !!---------------------------------------------------------------------- 15 USE oce ! ocean dynamics and tracers 16 USE dom_oce ! ocean space and time domain 17 USE sbc_oce ! Surface boundary condition: ocean fields 18 USE phycst ! physical constants 17 USE oce ! ocean dynamics and tracers 18 USE dom_oce ! ocean space and time domain 19 USE sbc_oce ! Surface boundary condition: ocean fields 20 USE sbc_ice ! Surface boundary condition: ocean fields 21 USE fldread ! read input fields 22 USE phycst ! physical constants 23 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 19 24 ! 20 USE in_out_manager ! I/O manager 21 USE lib_mpp ! distribued memory computing library 22 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 23 USE lib_fortran ! 25 USE in_out_manager ! I/O manager 26 USE iom 27 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 28 USE lib_mpp ! distribued memory computing library 29 !USE wrk_nemo ! work arrays 30 USE timing ! Timing 31 USE prtctl ! Print control 24 32 25 33 IMPLICIT NONE 26 34 PRIVATE 27 35 28 PUBLIC usrdef_sbc_oce ! routine called in sbcmod module 29 PUBLIC usrdef_sbc_ice_tau ! routine called by icestp.F90 for ice dynamics 30 PUBLIC usrdef_sbc_ice_flx ! routine called by icestp.F90 for ice thermo 36 PUBLIC usrdef_sbc_oce ! routine called in sbcmod module 37 PUBLIC usrdef_sbc_ice_tau ! routine called by sbcice_lim.F90 for ice dynamics 38 PUBLIC usrdef_sbc_ice_flx ! routine called by sbcice_lim.F90 for ice thermo 39 ! !!* Namelist namsbc_usr 40 REAL(wp) :: rn_vfac ! multiplication factor for ice/ocean velocity in the calculation of wind stress (clem) 41 REAL(wp) :: rn_charn_const 42 LOGICAL :: ln_use_sbc ! Surface fluxes on or not 31 43 32 44 !! * Substitutions … … 43 55 !! *** ROUTINE usrdef_sbc *** 44 56 !! 45 !! ** Purpose : provide at each time-step the GYREsurface boundary57 !! ** Purpose : provide at each time-step the surface boundary 46 58 !! condition, i.e. the momentum, heat and freshwater fluxes. 47 59 !! 48 !! ** Method : a nalytical seasonal cycle for GYRE configuration.60 !! ** Method : all 0 fields, for AMM7_SURGE case 49 61 !! CAUTION : never mask the surface stress field ! 50 62 !! 51 !! ** Action : - set the ocean surface boundary condition, i.e.63 !! ** Action : - if tide-only case - set to ZERO all the ocean surface boundary condition, i.e. 52 64 !! utau, vtau, taum, wndm, qns, qsr, emp, sfx 53 !! 54 !! Reference : Hazeleger, W., and S. Drijfhout, JPO, 30, 677-695, 2000. 65 !! - if tide+surge case - read in wind and air pressure !! 55 66 !!---------------------------------------------------------------------- 56 67 INTEGER, INTENT(in) :: kt ! ocean time step 57 !! 58 INTEGER :: ji, jj ! dummy loop indices 59 INTEGER :: zyear0 ! initial year 60 INTEGER :: zmonth0 ! initial month 61 INTEGER :: zday0 ! initial day 62 INTEGER :: zday_year0 ! initial day since january 1st 63 REAL(wp) :: ztau , ztau_sais ! wind intensity and of the seasonal cycle 64 REAL(wp) :: ztime ! time in hour 65 REAL(wp) :: ztimemax , ztimemin ! 21th June, and 21th decem. if date0 = 1st january 66 REAL(wp) :: ztimemax1, ztimemin1 ! 21th June, and 21th decem. if date0 = 1st january 67 REAL(wp) :: ztimemax2, ztimemin2 ! 21th June, and 21th decem. if date0 = 1st january 68 REAL(wp) :: ztaun ! intensity 69 REAL(wp) :: zemp_s, zemp_n, zemp_sais, ztstar 70 REAL(wp) :: zcos_sais1, zcos_sais2, ztrp, zconv, t_star 71 REAL(wp) :: zsumemp, zsurf 72 REAL(wp) :: zrhoa = 1.22 ! Air density kg/m3 73 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient 74 REAL(wp) :: ztx, zty, zmod, zcoef ! temporary variables 75 REAL(wp) :: zyydd ! number of days in one year 68 69 INTEGER :: ios ! Local integer output status for namelist read 70 ! 71 CHARACTER(len=100) :: cn_dir ! Root directory for location of flux files 72 TYPE(FLD_N) :: sn_wndi, sn_wndj ! informations about the fields to be read 73 74 NAMELIST/namsbc_usr/ ln_use_sbc, cn_dir , rn_vfac, & 75 & sn_wndi, sn_wndj, rn_charn_const 76 76 !!--------------------------------------------------------------------- 77 zyydd = REAL(nyear_len(1),wp) 77 ! 78 IF( kt == nit000 ) THEN 79 80 81 REWIND( numnam_cfg ) ! Namelist namsbc_usr in configuration namelist 82 READ ( numnam_cfg, namsbc_usr, IOSTAT = ios, ERR = 902 ) 83 902 IF( ios /= 0 ) CALL ctl_nam ( ios , 'namsbc_surge in configuration namelist' ) 78 84 79 ! ---------------------------- ! 80 ! heat and freshwater fluxes ! 81 ! ---------------------------- ! 82 !same temperature, E-P as in HAZELEGER 2000 85 IF(lwm) WRITE( numond, namsbc_usr ) 86 IF(lwp) WRITE(numout,*)' usr_sbc : AMM7_SURGE tide only case: NO surface forcing' 87 IF(lwp) WRITE(numout,*)' ~~~~~~~~~~~ utau = vtau = taum = wndm = qns = qsr = emp = sfx = 0' 83 88 84 zyear0 = ndate0 / 10000 ! initial year 85 zmonth0 = ( ndate0 - zyear0 * 10000 ) / 100 ! initial month 86 zday0 = ndate0 - zyear0 * 10000 - zmonth0 * 100 ! initial day betwen 1 and 30 87 zday_year0 = ( zmonth0 - 1 ) * 30.+zday0 ! initial day betwen 1 and 360 88 89 ! current day (in hours) since january the 1st of the current year 90 ztime = REAL( kt ) * rdt / (rmmss * rhhmm) & ! total incrementation (in hours) 91 & - (nyear - 1) * rjjhh * zyydd ! minus years since beginning of experiment (in hours) 92 93 ztimemax1 = ((5.*30.)+21.)* 24. ! 21th june at 24h in hours 94 ztimemin1 = ztimemax1 + rjjhh * zyydd / 2 ! 21th december in hours 95 ztimemax2 = ((6.*30.)+21.)* 24. ! 21th july at 24h in hours 96 ztimemin2 = ztimemax2 - rjjhh * zyydd / 2 ! 21th january in hours 97 ! ! NB: rjjhh * zyydd / 4 = one seasonal cycle in hours 98 99 ! amplitudes 100 zemp_S = 0.7 ! intensity of COS in the South 101 zemp_N = 0.8 ! intensity of COS in the North 102 zemp_sais = 0.1 103 zTstar = 28.3 ! intemsity from 28.3 a -5 deg 104 105 ! 1/2 period between 21th June and 21th December and between 21th July and 21th January 106 zcos_sais1 = COS( (ztime - ztimemax1) / (ztimemin1 - ztimemax1) * rpi ) 107 zcos_sais2 = COS( (ztime - ztimemax2) / (ztimemax2 - ztimemin2) * rpi ) 108 109 ztrp= - 40.e0 ! retroaction term on heat fluxes (W/m2/K) 110 zconv = 3.16e-5 ! convertion factor: 1 m/yr => 3.16e-5 mm/s 111 DO jj = 1, jpj 112 DO ji = 1, jpi 113 ! domain from 15 deg to 50 deg between 27 and 28 degC at 15N, -3 114 ! and 13 degC at 50N 53.5 + or - 11 = 1/4 period : 115 ! 64.5 in summer, 42.5 in winter 116 t_star = zTstar * ( 1. + 1. / 50. * zcos_sais2 ) & 117 & * COS( rpi * (gphit(ji,jj) - 5.) & 118 & / ( 53.5 * ( 1 + 11 / 53.5 * zcos_sais2 ) * 2.) ) 119 ! 23.5 deg : tropics 120 qsr (ji,jj) = 230 * COS( 3.1415 * ( gphit(ji,jj) - 23.5 * zcos_sais1 ) / ( 0.9 * 180 ) ) 121 qns (ji,jj) = ztrp * ( tsb(ji,jj,1,jp_tem) - t_star ) - qsr(ji,jj) 122 IF( gphit(ji,jj) >= 14.845 .AND. 37.2 >= gphit(ji,jj) ) THEN ! zero at 37.8 deg, max at 24.6 deg 123 emp (ji,jj) = zemp_S * zconv & 124 & * SIN( rpi / 2 * (gphit(ji,jj) - 37.2) / (24.6 - 37.2) ) & 125 & * ( 1 - zemp_sais / zemp_S * zcos_sais1) 126 ELSE 127 emp (ji,jj) = - zemp_N * zconv & 128 & * SIN( rpi / 2 * (gphit(ji,jj) - 37.2) / (46.8 - 37.2) ) & 129 & * ( 1 - zemp_sais / zemp_N * zcos_sais1 ) 130 ENDIF 131 END DO 132 END DO 133 134 zsumemp = GLOB_SUM( 'usrdef_sbc', emp (:,:) ) 135 zsurf = GLOB_SUM( 'usrdef_sbc', tmask(:,:,1) ) 136 zsumemp = zsumemp / zsurf ! Default GYRE configuration 137 138 ! freshwater (mass flux) and update of qns with heat content of emp 139 emp (:,:) = emp(:,:) - zsumemp * tmask(:,:,1) ! freshwater flux (=0 in domain average) 140 sfx (:,:) = 0.0_wp ! no salt flux 141 qns (:,:) = qns(:,:) - emp(:,:) * sst_m(:,:) * rcp ! evap and precip are at SST 142 143 144 ! ---------------------------- ! 145 ! momentum fluxes ! 146 ! ---------------------------- ! 147 ! same wind as in Wico 148 !test date0 : ndate0 = 010203 149 zyear0 = ndate0 / 10000 150 zmonth0 = ( ndate0 - zyear0 * 10000 ) / 100 151 zday0 = ndate0 - zyear0 * 10000 - zmonth0 * 100 152 !Calculates nday_year, day since january 1st 153 zday_year0 = (zmonth0-1)*30.+zday0 154 155 !accumulates days of previous months of this year 156 ! day (in hours) since january the 1st 157 ztime = FLOAT( kt ) * rdt / (rmmss * rhhmm) & ! incrementation in hour 158 & - (nyear - 1) * rjjhh * zyydd ! - nber of hours the precedent years 159 ztimemax = ((5.*30.)+21.)* 24. ! 21th june in hours 160 ztimemin = ztimemax + rjjhh * zyydd / 2 ! 21th december in hours 161 ! ! NB: rjjhh * zyydd / 4 = 1 seasonal cycle in hours 162 163 ! mean intensity at 0.105 ; srqt(2) because projected with 45deg angle 164 ztau = 0.105 / SQRT( 2. ) 165 ! seasonal oscillation intensity 166 ztau_sais = 0.015 167 ztaun = ztau - ztau_sais * COS( (ztime - ztimemax) / (ztimemin - ztimemax) * rpi ) 168 DO jj = 1, jpj 169 DO ji = 1, jpi 170 ! domain from 15deg to 50deg and 1/2 period along 14deg 171 ! so 5/4 of half period with seasonal cycle 172 utau(ji,jj) = - ztaun * SIN( rpi * (gphiu(ji,jj) - 15.) / (29.-15.) ) 173 vtau(ji,jj) = ztaun * SIN( rpi * (gphiv(ji,jj) - 15.) / (29.-15.) ) 174 END DO 175 END DO 176 177 ! module of wind stress and wind speed at T-point 178 zcoef = 1. / ( zrhoa * zcdrag ) 179 DO jj = 2, jpjm1 180 DO ji = fs_2, fs_jpim1 ! vect. opt. 181 ztx = utau(ji-1,jj ) + utau(ji,jj) 182 zty = vtau(ji ,jj-1) + vtau(ji,jj) 183 zmod = 0.5 * SQRT( ztx * ztx + zty * zty ) 184 taum(ji,jj) = zmod 185 wndm(ji,jj) = SQRT( zmod * zcoef ) 186 END DO 187 END DO 188 CALL lbc_lnk_multi( 'usrdef_sbc', taum(:,:), 'T', 1. , wndm(:,:), 'T', 1. ) 189 190 ! ---------------------------------- ! 191 ! control print at first time-step ! 192 ! ---------------------------------- ! 193 IF( kt == nit000 .AND. lwp ) THEN 194 WRITE(numout,*) 195 WRITE(numout,*)'usrdef_sbc_oce : analytical surface fluxes for GYRE configuration' 196 WRITE(numout,*)'~~~~~~~~~~~ ' 197 WRITE(numout,*)' nyear = ', nyear 198 WRITE(numout,*)' nmonth = ', nmonth 199 WRITE(numout,*)' nday = ', nday 200 WRITE(numout,*)' nday_year = ', nday_year 201 WRITE(numout,*)' ztime = ', ztime 202 WRITE(numout,*)' ztimemax = ', ztimemax 203 WRITE(numout,*)' ztimemin = ', ztimemin 204 WRITE(numout,*)' ztimemax1 = ', ztimemax1 205 WRITE(numout,*)' ztimemin1 = ', ztimemin1 206 WRITE(numout,*)' ztimemax2 = ', ztimemax2 207 WRITE(numout,*)' ztimemin2 = ', ztimemin2 208 WRITE(numout,*)' zyear0 = ', zyear0 209 WRITE(numout,*)' zmonth0 = ', zmonth0 210 WRITE(numout,*)' zday0 = ', zday0 211 WRITE(numout,*)' zday_year0 = ', zday_year0 212 WRITE(numout,*)' zyydd = ', zyydd 213 WRITE(numout,*)' zemp_S = ', zemp_S 214 WRITE(numout,*)' zemp_N = ', zemp_N 215 WRITE(numout,*)' zemp_sais = ', zemp_sais 216 WRITE(numout,*)' zTstar = ', zTstar 217 WRITE(numout,*)' zsumemp = ', zsumemp 218 WRITE(numout,*)' zsurf = ', zsurf 219 WRITE(numout,*)' ztrp = ', ztrp 220 WRITE(numout,*)' zconv = ', zconv 221 WRITE(numout,*)' ndastp = ', ndastp 222 WRITE(numout,*)' adatrj = ', adatrj 89 utau(:,:) = 0._wp 90 vtau(:,:) = 0._wp 91 taum(:,:) = 0._wp 92 wndm(:,:) = 0._wp 93 ! 94 emp (:,:) = 0._wp 95 sfx (:,:) = 0._wp 96 qns (:,:) = 0._wp 97 qsr (:,:) = 0._wp 98 ! 223 99 ENDIF 224 100 ! … … 231 107 232 108 233 SUBROUTINE usrdef_sbc_ice_flx( kt , phs, phi)109 SUBROUTINE usrdef_sbc_ice_flx( kt ) 234 110 INTEGER, INTENT(in) :: kt ! ocean time step 235 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: phs ! snow thickness236 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: phi ! ice thickness237 111 END SUBROUTINE usrdef_sbc_ice_flx 238 112 -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/par_oce.F90
r14075 r15422 64 64 INTEGER, PUBLIC, PARAMETER :: jp_tem = 1 !: indice for temperature 65 65 INTEGER, PUBLIC, PARAMETER :: jp_sal = 2 !: indice for salinity 66 INTEGER, PUBLIC, PARAMETER :: jp_dep = 3 !: indice for depth 67 INTEGER, PUBLIC, PARAMETER :: jp_msk = 4 !: indice for mask 66 68 67 69 !!---------------------------------------------------------------------- -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/step.F90
r14075 r15422 206 206 CALL dia_ar5 ( kstp ) ! ar5 diag 207 207 IF( ln_diaptr ) CALL dia_ptr ! Poleward adv/ldf TRansports diagnostics 208 IF( lk_diaharm_fast ) & 209 & CALL dia_harm_fast( kstp ) ! Tidal harmonic analysis - restart and faster version 208 210 IF( ln_diaharm ) CALL dia_harm( kstp ) ! Tidal harmonic analysis 209 211 CALL dia_wri ( kstp ) ! ocean model: outputs -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/step_oce.F90
r14075 r15422 79 79 USE diahth ! thermocline depth (dia_hth routine) 80 80 USE diahsb ! heat, salt and volume budgets (dia_hsb routine) 81 USE diaharm_fast ! harmonic analysis of tides (harm_ana routine) 81 82 USE diaharm 82 83 USE diacfl -
NEMO/branches/UKMO/r14075_India_uncoupled/src/OCE/stpctl.F90
r14075 r15422 126 126 zmax(1) = MAXVAL( ABS( sshn(:,:) ) ) ! ssh max 127 127 ENDIF 128 zmax(2) = MAXVAL( ABS( un(:,:,:) ) ) ! velocity max (zonal only)128 zmax(2) = MAXVAL( un(:,:,:)*un(:,:,:) + vn(:,:,:)*vn(:,:,:) ) ! velocity max 129 129 llmsk(:,:,:) = tmask(:,:,:) == 1._wp 130 130 IF( COUNT( llmsk(:,:,:) ) > 0 ) THEN ! avoid huge values sent back for land processors... … … 167 167 ! !== error handling ==! 168 168 IF( zmax(1) > 20._wp .OR. & ! too large sea surface height ( > 20 m ) 169 & zmax(2) > 10._wp .OR. & ! too large velocity ( > 10 m/s)169 & zmax(2) > 100._wp .OR. & ! too large velocity ( > 10 m/s) 170 170 & zmax(3) >= 0._wp .OR. & ! negative or zero sea surface salinity 171 171 & zmax(4) >= 100._wp .OR. & ! too large sea surface salinity ( > 100 ) … … 177 177 IF( lwm .AND. kt /= nitend ) istatus = NF90_CLOSE(idrun) 178 178 CALL mpp_maxloc( 'stpctl', ABS(sshn) , ssmask(:,:) , zzz, ih(1:2) ) ; ih(3) = 0 179 CALL mpp_maxloc( 'stpctl', ABS(un), umask (:,:,:), zzz, iu )179 CALL mpp_maxloc( 'stpctl', un*un + vn*vn , umask (:,:,:), zzz, iu ) 180 180 CALL mpp_minloc( 'stpctl', tsn(:,:,:,jp_sal), tmask (:,:,:), zzz, is1 ) 181 181 CALL mpp_maxloc( 'stpctl', tsn(:,:,:,jp_sal), tmask (:,:,:), zzz, is2 ) … … 192 192 ELSE 193 193 ih(1:2)= MAXLOC( ABS( sshn(:,:) ) ) + (/ nimpp - 1, njmpp - 1 /) ; ih(3) = 0 194 iu(:) = MAXLOC( ABS( un (:,:,:) )) + (/ nimpp - 1, njmpp - 1, 0 /)194 iu(:) = MAXLOC( un(:,:,:)*un(:,:,:) + vn(:,:,:)*vn(:,:,:) ) + (/ nimpp - 1, njmpp - 1, 0 /) 195 195 is1(:) = MINLOC( tsn(:,:,:,jp_sal), mask = tmask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /) 196 196 is2(:) = MAXLOC( tsn(:,:,:,jp_sal), mask = tmask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /) … … 200 200 WRITE(ctmp1,*) ' stp_ctl: |ssh| > 20 m or |U| > 10 m/s or S <= 0 or S >= 100 or NaN encounter in the tests' 201 201 CALL wrt_line(ctmp2, kt, ' |ssh| max ', zmax(1), ih , iareasum(1), iareamin(1), iareamax(1) ) 202 CALL wrt_line(ctmp3, kt, ' |U|max ', zmax(2), iu , iareasum(2), iareamin(2), iareamax(2) )202 CALL wrt_line(ctmp3, kt, ' Vel2 max ', zmax(2), iu , iareasum(2), iareamin(2), iareamax(2) ) 203 203 CALL wrt_line(ctmp4, kt, ' Sal min ', - zmax(3), is1, iareasum(3), iareamin(3), iareamax(3) ) 204 204 CALL wrt_line(ctmp5, kt, ' Sal max ', zmax(4), is2, iareasum(4), iareamin(4), iareamax(4) ) … … 226 226 ENDIF 227 227 ! 228 9500 FORMAT(' it :', i8, ' |ssh|_max: ', D23.16, ' |U|_max: ', D23.16,' S_min: ', D23.16,' S_max: ', D23.16)228 9500 FORMAT(' it :', i8, ' |ssh|_max: ', D23.16, ' Vel2_max: ', D23.16,' S_min: ', D23.16,' S_max: ', D23.16) 229 229 ! 230 230 END SUBROUTINE stp_ctl
Note: See TracChangeset
for help on using the changeset viewer.
