Changeset 12453
- Timestamp:
- 2020-02-25T16:29:34+01:00 (3 years ago)
- Location:
- NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE
- Files:
-
- 2 added
- 14 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE/BDY/bdyini.F90
r11715 r12453 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/r12083_India_uncoupled/src/OCE/DIA/diaharm.F90
r11715 r12453 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 = 1,jpj 193 218 DO ji = 1,jpi 194 ! Elevation 195 ana_temp(ji,jj,nhc,1) = ana_temp(ji,jj,nhc,1) + ztemp*sshn(ji,jj)*ssmask (ji,jj) 196 ana_temp(ji,jj,nhc,2) = ana_temp(ji,jj,nhc,2) + ztemp*un_b(ji,jj)*ssumask(ji,jj) 197 ana_temp(ji,jj,nhc,3) = ana_temp(ji,jj,nhc,3) + ztemp*vn_b(ji,jj)*ssvmask(ji,jj) 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 198 234 END DO 199 235 END DO 200 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 201 255 END DO 202 256 END DO … … 221 275 !!-------------------------------------------------------------------- 222 276 INTEGER :: ji, jj, jh, jc, jn, nhan, jl 277 # if defined key_3Ddiaharm 278 INTEGER :: jk 279 # endif 223 280 INTEGER :: ksp, kun, keq 224 281 REAL(wp) :: ztime, ztime_ini, ztime_end … … 234 291 ztime_end = nitend_han*rdt ! Final time in seconds at the end of analysis 235 292 nhan = (nitend_han-nit000_han+1)/nstep_han ! Number of dumps used for analysis 293 294 # if defined key_3Ddiaharm 295 ALLOCATE( out_eta(jpi,jpj,jpk,2*nb_ana), & 296 & out_u (jpi,jpj,jpk,2*nb_ana), & 297 & out_v (jpi,jpj,jpk,2*nb_ana), & 298 & out_w (jpi,jpj,jpk,2*nb_ana), & 299 & out_dzi(jpi,jpj,jpk,2*nb_ana) ) 300 # else 301 ALLOCATE( out_eta(jpi,jpj,2*nb_ana), & 302 & out_u (jpi,jpj,2*nb_ana), & 303 & out_v (jpi,jpj,2*nb_ana) ) 304 # endif 305 306 IF(lwp) WRITE(numout,*) 'ANA F OLD', ft 307 IF(lwp) WRITE(numout,*) 'ANA U OLD', ut 308 IF(lwp) WRITE(numout,*) 'ANA V OLD', vt 236 309 237 310 ninco = 2*nb_ana … … 257 330 nsparse = ksp 258 331 259 ! Elevation: 332 ! Density and Elevation: 333 # if defined key_3Ddiaharm 334 DO jk=1,jpk 335 # endif 260 336 DO jj = 1, jpj 261 337 DO ji = 1, jpi … … 265 341 DO jc = 1, 2 266 342 kun = kun + 1 343 # if defined key_3Ddiaharm 344 ztmp4(kun)=ana_temp(ji,jj,kun,1,jk) 345 # else 267 346 ztmp4(kun)=ana_temp(ji,jj,kun,1) 347 # endif 268 348 END DO 269 349 END DO … … 278 358 END DO 279 359 END DO 280 281 ALLOCATE( out_eta(jpi,jpj,2*nb_ana), &282 & out_u (jpi,jpj,2*nb_ana), &283 & out_v (jpi,jpj,2*nb_ana) )284 360 285 361 DO jj = 1, jpj … … 288 364 X1 = ana_amp(ji,jj,jh,1) 289 365 X2 =-ana_amp(ji,jj,jh,2) 290 out_eta(ji,jj,jh ) = X1 * tmask_i(ji,jj) 291 out_eta(ji,jj,jh+nb_ana) = X2 * tmask_i(ji,jj) 292 END DO 293 END DO 294 END DO 295 296 ! ubar: 366 # if defined key_3Ddiaharm 367 out_eta(ji,jj,jk,jh ) = X1 * tmask_i(ji,jj) 368 out_eta(ji,jj,jk,jh+nb_ana) = X2 * tmask_i(ji,jj) 369 # else 370 out_eta(ji,jj ,jh ) = X1 * tmask_i(ji,jj) 371 out_eta(ji,jj ,jh+nb_ana) = X2 * tmask_i(ji,jj) 372 # endif 373 END DO 374 END DO 375 END DO 376 377 ! u-component of velocity 297 378 DO jj = 1, jpj 298 379 DO ji = 1, jpi … … 302 383 DO jc = 1,2 303 384 kun = kun + 1 385 # if defined key_3Ddiaharm 386 ztmp4(kun)=ana_temp(ji,jj,kun,2,jk) 387 # else 304 388 ztmp4(kun)=ana_temp(ji,jj,kun,2) 389 # endif 305 390 END DO 306 391 END DO … … 322 407 X1= ana_amp(ji,jj,jh,1) 323 408 X2=-ana_amp(ji,jj,jh,2) 324 out_u(ji,jj, jh) = X1 * ssumask(ji,jj) 325 out_u(ji,jj,nb_ana+jh) = X2 * ssumask(ji,jj) 409 # if defined key_3Ddiaharm 410 out_u(ji,jj,jk, jh) = X1 * ssumask(ji,jj) 411 out_u(ji,jj,jk,nb_ana+jh) = X2 * ssumask(ji,jj) 412 # else 413 out_u(ji,jj, jh) = X1 * ssumask(ji,jj) 414 out_u(ji,jj, nb_ana+jh) = X2 * ssumask(ji,jj) 415 # endif 326 416 ENDDO 327 417 ENDDO 328 418 ENDDO 329 419 330 ! v bar:420 ! v- velocity 331 421 DO jj = 1, jpj 332 422 DO ji = 1, jpi … … 336 426 DO jc = 1,2 337 427 kun = kun + 1 428 # if defined key_3Ddiaharm 429 ztmp4(kun)=ana_temp(ji,jj,kun,3,jk) 430 # else 338 431 ztmp4(kun)=ana_temp(ji,jj,kun,3) 432 # endif 339 433 END DO 340 434 END DO … … 356 450 X1=ana_amp(ji,jj,jh,1) 357 451 X2=-ana_amp(ji,jj,jh,2) 358 out_v(ji,jj, jh)=X1 * ssvmask(ji,jj) 359 out_v(ji,jj,nb_ana+jh)=X2 * ssvmask(ji,jj) 360 END DO 361 END DO 362 END DO 452 # if defined key_3Ddiaharm 453 out_v(ji,jj,jk, jh)=X1 * ssvmask(ji,jj) 454 out_v(ji,jj,jk,nb_ana+jh)=X2 * ssvmask(ji,jj) 455 # else 456 out_v(ji,jj, jh)=X1 * ssvmask(ji,jj) 457 out_v(ji,jj, nb_ana+jh)=X2 * ssvmask(ji,jj) 458 # endif 459 END DO 460 END DO 461 END DO 462 463 # if defined key_3Ddiaharm 464 ! w- velocity 465 DO jj = 1, jpj 466 DO ji = 1, jpi 467 ! Fill input array 468 kun=0 469 DO jh = 1,nb_ana 470 DO jc = 1,2 471 kun = kun + 1 472 ztmp4(kun)=ana_temp(ji,jj,kun,4,jk) 473 END DO 474 END DO 475 476 CALL SUR_DETERMINE(jj+1) 477 478 ! Fill output array 479 DO jh = 1, nb_ana 480 ana_amp(ji,jj,jh,1)=ztmp7((jh-1)*2+1) 481 ana_amp(ji,jj,jh,2)=ztmp7((jh-1)*2+2) 482 END DO 483 484 END DO 485 END DO 486 487 DO jj = 1, jpj 488 DO ji = 1, jpi 489 DO jh = 1, nb_ana 490 X1=ana_amp(ji,jj,jh,1) 491 X2=-ana_amp(ji,jj,jh,2) 492 out_w(ji,jj,jk, jh)=X1 * tmask_i(ji,jj) 493 out_w(ji,jj,jk,nb_ana+jh)=X2 * tmask_i(ji,jj) 494 END DO 495 END DO 496 END DO 497 498 ! dzi- isopycnal displacements 499 DO jj = 1, jpj 500 DO ji = 1, jpi 501 ! Fill input array 502 kun=0 503 DO jh = 1,nb_ana 504 DO jc = 1,2 505 kun = kun + 1 506 ztmp4(kun)=ana_temp(ji,jj,kun,5,jk) 507 END DO 508 END DO 509 510 CALL SUR_DETERMINE(jj+1) 511 512 ! Fill output array 513 DO jh = 1, nb_ana 514 ana_amp(ji,jj,jh,1)=ztmp7((jh-1)*2+1) 515 ana_amp(ji,jj,jh,2)=ztmp7((jh-1)*2+2) 516 END DO 517 518 END DO 519 END DO 520 521 DO jj = 1, jpj 522 DO ji = 1, jpi 523 DO jh = 1, nb_ana 524 X1=ana_amp(ji,jj,jh,1) 525 X2=-ana_amp(ji,jj,jh,2) 526 out_dzi(ji,jj,jk, jh)=X1 * tmask_i(ji,jj) 527 out_dzi(ji,jj,jk,nb_ana+jh)=X2 * tmask_i(ji,jj) 528 END DO 529 END DO 530 END DO 531 532 ENDDO ! jk 533 # endif 363 534 ! 364 535 CALL dia_wri_harm ! Write results in files … … 379 550 cdfile_name_V ! name of the file created (V-points) 380 551 INTEGER :: jh 381 !!---------------------------------------------------------------------- 552 553 # if defined key_3Ddiaharm 554 CHARACTER(LEN=lc) :: cdfile_name_W ! name of the file created (W-points) 555 INTEGER :: jk 556 REAL(WP), ALLOCATABLE, DIMENSION (:,:,:) :: z3real, z3im 557 REAL(WP), ALLOCATABLE, DIMENSION (:,:) :: z2real, z2im 558 # endif 559 !!---------------------------------------------------------------------- 560 561 #if defined key_dimgout 562 cdfile_name_T = TRIM(cexper)//'_Tidal_harmonics_gridT.dimgproc' 563 cdfile_name_U = TRIM(cexper)//'_Tidal_harmonics_gridU.dimgproc' 564 cdfile_name_V = TRIM(cexper)//'_Tidal_harmonics_gridV.dimgproc' 565 # if defined key_3Ddiaharm 566 cdfile_name_W = TRIM(cexper)//'_Tidal_harmonics_gridW.dimgproc' 567 # endif 568 #endif 382 569 383 570 IF(lwp) WRITE(numout,*) ' ' 384 571 IF(lwp) WRITE(numout,*) 'dia_wri_harm : Write harmonic analysis results' 572 #if defined key_dimgout 573 IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~ Output files: ', TRIM(cdfile_name_T) 574 IF(lwp) WRITE(numout,*) ' ', TRIM(cdfile_name_U) 575 IF(lwp) WRITE(numout,*) ' ', TRIM(cdfile_name_V) 576 # if defined key_3Ddiaharm 577 IF(lwp) WRITE(numout,*) ' ', TRIM(cdfile_name_W) 578 # endif 579 #endif 385 580 IF(lwp) WRITE(numout,*) ' ' 386 581 387 ! A) Elevation 582 # if defined key_3Ddiaharm 583 ALLOCATE(z3real(jpi,jpj,jpk),z3im(jpi,jpj,jpk),z2real(jpi,jpj),z2im(jpi,jpj)) 584 # endif 585 586 ! A) density and elevation 388 587 !///////////// 389 588 ! 589 #if defined key_dimgout 590 cltext='density amplitude and phase; elevation is level=jpk ' 591 CALL dia_wri_dimg(TRIM(cdfile_name_T), TRIM(cltext), out_eta, 2*nb_ana, '2') 592 #else 593 # if defined key_3Ddiaharm 594 z3real(:,:,:) = 0._wp; z3im(:,:,:) = 0._wp 595 # endif 390 596 DO jh = 1, nb_ana 597 # if defined key_3Ddiaharm 598 DO jk=1,jpkm1 599 z3real(:,:,jk)=out_eta(:,:,jk,jh) 600 z3im (:,:,jk)=out_eta(:,:,jk,jh+nb_ana) 601 ENDDO 602 z2real(:,:)=out_eta(:,:,jpk,jh); z2im(:,:)=out_eta(:,:,jpk,jh+nb_ana) 603 CALL iom_put( TRIM(tname(jh))//'x_ro', z3real(:,:,:) ) 604 CALL iom_put( TRIM(tname(jh))//'y_ro', z3im (:,:,:) ) 605 CALL iom_put( TRIM(tname(jh))//'x' , z2real(:,: ) ) 606 CALL iom_put( TRIM(tname(jh))//'y' , z2im (:,: ) ) 607 # else 608 WRITE(numout,*) "OUTPUT ORI: ", TRIM(tname(jh))//'x', ' & ', TRIM(tname(jh))//'y', MAXVAL(out_eta(:,:,jh)) 391 609 CALL iom_put( TRIM(tname(jh))//'x', out_eta(:,:,jh) ) 392 610 CALL iom_put( TRIM(tname(jh))//'y', out_eta(:,:,nb_ana+jh) ) 393 END DO 394 395 ! B) ubar 611 # endif 612 END DO 613 #endif 614 615 ! B) u 396 616 !///////// 397 617 ! 618 #if defined key_dimgout 619 cltext='3d u amplitude and phase; ubar is the last level' 620 CALL dia_wri_dimg(TRIM(cdfile_name_U), TRIM(cltext), out_u, 2*nb_ana, '2') 621 #else 622 # if defined key_3Ddiaharm 623 z3real(:,:,:) = 0._wp; z3im(:,:,:) = 0._wp 624 # endif 398 625 DO jh = 1, nb_ana 399 CALL iom_put( TRIM(tname(jh))//'x_u', out_u(:,:,jh) ) 400 CALL iom_put( TRIM(tname(jh))//'y_u', out_u(:,:,nb_ana+jh) ) 401 END DO 402 403 ! C) vbar 626 # if defined key_3Ddiaharm 627 DO jk=1,jpkm1 628 z3real(:,:,jk)=out_u(:,:,jk,jh) 629 z3im (:,:,jk)=out_u(:,:,jk,jh+nb_ana) 630 ENDDO 631 z2real(:,:)=out_u(:,:,jpk,jh); z2im(:,:)=out_u(:,:,jpk,jh+nb_ana) 632 CALL iom_put( TRIM(tname(jh))//'x_u3d', z3real(:,:,:) ) 633 CALL iom_put( TRIM(tname(jh))//'y_u3d', z3im (:,:,:) ) 634 CALL iom_put( TRIM(tname(jh))//'x_u2d', z2real(:,:) ) 635 CALL iom_put( TRIM(tname(jh))//'y_u2d', z2im (:,:) ) 636 z2real(:,:)=out_w(:,:,jpk,jh); z2im(:,:)=out_w(:,:,jpk,jh+nb_ana) 637 CALL iom_put( TRIM(tname(jh))//'x_tabx', z2real(:,:) ) 638 CALL iom_put( TRIM(tname(jh))//'y_tabx', z2im (:,:) ) 639 # else 640 CALL iom_put( TRIM(tname(jh))//'x_u2d', out_u(:,:,jh) ) 641 CALL iom_put( TRIM(tname(jh))//'y_u2d', out_u(:,:,nb_ana+jh) ) 642 # endif 643 END DO 644 #endif 645 646 ! C) v 404 647 !///////// 405 648 ! 649 #if defined key_dimgout 650 cltext='3d v amplitude and phase; vbar is the last level' 651 CALL dia_wri_dimg(TRIM(cdfile_name_V), TRIM(cltext), out_v, 2*nb_ana, '2') 652 #else 653 # if defined key_3Ddiaharm 654 z3real(:,:,:) = 0._wp; z3im(:,:,:) = 0._wp 655 # endif 406 656 DO jh = 1, nb_ana 407 CALL iom_put( TRIM(tname(jh))//'x_v', out_v(:,:,jh ) ) 408 CALL iom_put( TRIM(tname(jh))//'y_v', out_v(:,:,jh+nb_ana) ) 409 END DO 410 ! 657 # if defined key_3Ddiaharm 658 DO jk=1,jpkm1 659 z3real(:,:,jk)=out_v(:,:,jk,jh) 660 z3im (:,:,jk)=out_v(:,:,jk,jh+nb_ana) 661 ENDDO 662 z2real(:,:)=out_v(:,:,jpk,jh); z2im(:,:)=out_v(:,:,jpk,jh+nb_ana) 663 CALL iom_put( TRIM(tname(jh))//'x_v3d', z3real(:,:,:) ) 664 CALL iom_put( TRIM(tname(jh))//'y_v3d', z3im (:,:,:) ) 665 CALL iom_put( TRIM(tname(jh))//'x_v2d' , z2real(:,:) ) 666 CALL iom_put( TRIM(tname(jh))//'y_v2d' , z2im (:,:) ) 667 z2real(:,:)=out_dzi(:,:,jpk,jh); z2im(:,:)=out_dzi(:,:,jpk,jh+nb_ana) 668 CALL iom_put( TRIM(tname(jh))//'x_taby', z2real(:,:) ) 669 CALL iom_put( TRIM(tname(jh))//'y_taby', z2im (:,:) ) 670 # else 671 CALL iom_put( TRIM(tname(jh))//'x_v2d', out_v(:,:,jh ) ) 672 CALL iom_put( TRIM(tname(jh))//'y_v2d', out_v(:,:,jh+nb_ana) ) 673 # endif 674 END DO 675 676 #endif 677 ! D) w 678 # if defined key_3Ddiaharm 679 # if defined key_dimgout 680 cltext='3d w amplitude and phase; vort_baro is the last level' 681 CALL dia_wri_dimg(TRIM(cdfile_name_W), TRIM(cltext), out_w, 2*nb_ana, '2') 682 # else 683 DO jh = 1, nb_ana 684 DO jk=1,jpkm1 685 z3real(:,:,jk)=out_w(:,:,jk,jh) 686 z3im(:,:,jk)=out_w(:,:,jk,jh+nb_ana) 687 ENDDO 688 CALL iom_put( TRIM(tname(jh))//'x_w3d', z3real(:,:,:) ) 689 CALL iom_put( TRIM(tname(jh))//'y_w3d', z3im(:,:,:) ) 690 END DO 691 # endif 692 693 ! E) dzi + tau_bot 694 # if defined key_dimgout 695 cltext='dzi=g*ro/N2 amplitude and phase' 696 CALL dia_wri_dimg(TRIM(cdfile_name_W), TRIM(cltext), out_w, 2*nb_ana, '2') 697 # else 698 DO jh = 1, nb_ana 699 DO jk=1,jpkm1 700 z3real(:,:,jk)=out_dzi(:,:,jk,jh) 701 z3im(:,:,jk)=out_dzi(:,:,jk,jh+nb_ana) 702 ENDDO 703 CALL iom_put( TRIM(tname(jh))//'x_dzi', z3real(:,:,:) ) 704 CALL iom_put( TRIM(tname(jh))//'y_dzi', z3im(:,:,:) ) 705 END DO 706 # endif 707 # endif 708 709 ! 710 # if defined key_3Ddiaharm 711 DEALLOCATE(z3real, z3im, z2real,z2im) 712 # endif 713 411 714 END SUBROUTINE dia_wri_harm 412 715 -
NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE/DOM/dtatsd.F90
r11715 r12453 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 NAMELIST/namtsd/ ln_tsd_init, ln_tsd_ dmp, cn_dir, sn_tem, sn_sal66 NAMELIST/namtsd/ ln_tsd_init, ln_tsd_interp, ln_tsd_dmp, cn_dir, sn_tem, sn_sal, sn_dep, sn_msk 62 67 !!---------------------------------------------------------------------- 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 ! 142 IF( ierr0 + ierr1 + ierr2 + ierr3 + ierr4 + ierr5 > 0 ) THEN 111 143 CALL ctl_stop( 'dta_tsd : unable to allocate T & S data arrays' ) ; RETURN 112 144 ENDIF 113 145 ! ! fill sf_tsd with sn_tem & sn_sal and control print 114 146 slf_i(jp_tem) = sn_tem ; slf_i(jp_sal) = sn_sal 147 IF( ln_tsd_interp ) slf_i(jp_dep) = sn_dep ; slf_i(jp_msk) = sn_msk 115 148 CALL fld_fill( sf_tsd, slf_i, cn_dir, 'dta_tsd', 'Temperature & Salinity data', 'namtsd', no_print ) 116 149 ! … … 138 171 REAL(wp), DIMENSION(jpi,jpj,jpk,jpts), INTENT( out) :: ptsd ! T & S data 139 172 ! 140 INTEGER :: ji, jj, jk, jl, jk k ! dummy loop indicies141 INTEGER :: ik, il0, il1, ii0, ii1, ij0, ij1 ! local integers173 INTEGER :: ji, jj, jk, jl, jk_init ! dummy loop indices 174 INTEGER :: ik, il0, il1, ii0, ii1, ij0, ij1 ! local integers 142 175 REAL(wp):: zl, zi ! local scalars 143 REAL(wp), DIMENSION(jpk) :: ztp, zsp ! 1D workspace144 176 !!---------------------------------------------------------------------- 145 177 ! … … 176 208 !!gm end 177 209 ! 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 210 IF( kt == nit000 .AND. lwp )THEN 211 WRITE(numout,*) 212 WRITE(numout,*) 'dta_tsd: interpolates T & S data onto current mesh' 213 ENDIF 214 ! 215 IF( ln_tsd_interp ) THEN ! probably should use pointers in the following to make more readable 216 ! 217 DO jk = 1, jpk ! determines the intepolated T-S profiles at each (i,j) points 218 DO jj= 1, jpj 219 DO ji= 1, jpi 191 220 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 221 IF( zl < sf_tsd(jp_dep)%fnow(ji,jj,1) ) THEN ! above the first level of data 222 ptsd(ji,jj,jk,jp_tem) = sf_tsd(jp_tem)%fnow(ji,jj,1) 223 ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,1) 224 ELSEIF( zl > sf_tsd(jp_dep)%fnow(ji,jj,jpk_init) ) THEN ! below the last level of data 225 ptsd(ji,jj,jk,jp_tem) = sf_tsd(jp_tem)%fnow(ji,jj,jpk_init) 226 ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,jpk_init) 227 ELSE ! inbetween : vertical interpolation between jk_init & jk_init+1 228 DO jk_init = 1, jpk_init-1 ! when gdept(jk_init) < zl < gdept(jk_init+1) 229 IF( sf_tsd(jp_msk)%fnow(ji,jj,jk_init+1) == 0 ) THEN ! if there is no data fill down 230 sf_tsd(jp_tem)%fnow(ji,jj,jk_init+1) = sf_tsd(jp_tem)%fnow(ji,jj,jk_init) 231 sf_tsd(jp_sal)%fnow(ji,jj,jk_init+1) = sf_tsd(jp_sal)%fnow(ji,jj,jk_init) 232 ENDIF 233 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 234 zi = ( zl - sf_tsd(jp_dep)%fnow(ji,jj,jk_init) ) / & 235 & (sf_tsd(jp_dep)%fnow(ji,jj,jk_init+1)-sf_tsd(jp_dep)%fnow(ji,jj,jk_init)) 236 ptsd(ji,jj,jk,jp_tem) = sf_tsd(jp_tem)%fnow(ji,jj,jk_init) + & 237 & (sf_tsd(jp_tem)%fnow(ji,jj,jk_init+1)-sf_tsd(jp_tem)%fnow(ji,jj,jk_init)) * zi 238 ptsd(ji,jj,jk,jp_sal) = sf_tsd(jp_sal)%fnow(ji,jj,jk_init) + & 239 & (sf_tsd(jp_sal)%fnow(ji,jj,jk_init+1)-sf_tsd(jp_sal)%fnow(ji,jj,jk_init)) * zi 204 240 ENDIF 205 241 END DO 206 242 ENDIF 207 243 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 244 END DO 215 245 END DO 216 246 ! 247 ptsd(:,:,:,jp_tem) = ptsd(:,:,:,jp_tem) *tmask(:,:,:) 248 ptsd(:,:,:,jp_sal) = ptsd(:,:,:,jp_sal) *tmask(:,:,:) 217 249 ELSE !== z- or zps- coordinate ==! 218 250 ! 219 ptsd(:,:,:,jp_tem) = ptsd(:,:,:,jp_tem) * tmask(:,:,:) ! Mask220 ptsd(:,:,:,jp_sal) = ptsd(:,:,:,jp_sal)* tmask(:,:,:)251 ptsd(:,:,:,jp_tem) = sf_tsd(jp_tem)%fnow(:,:,:) * tmask(:,:,:) ! Mask 252 ptsd(:,:,:,jp_sal) = sf_tsd(jp_sal)%fnow(:,:,:) * tmask(:,:,:) 221 253 ! 222 254 IF( ln_zps ) THEN ! zps-coordinate (partial steps) interpolation at the last ocean level … … 248 280 DEALLOCATE( sf_tsd(jp_sal)%fnow ) ! S arrays in the structure 249 281 IF( sf_tsd(jp_sal)%ln_tint ) DEALLOCATE( sf_tsd(jp_sal)%fdta ) 282 IF( ln_tsd_interp ) DEALLOCATE( sf_tsd(jp_dep)%fnow ) ! T arrays in the structure 283 IF( ln_tsd_interp ) DEALLOCATE( sf_tsd(jp_msk)%fnow ) ! T arrays in the structure 250 284 DEALLOCATE( sf_tsd ) ! the structure itself 251 285 ENDIF -
NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE/DYN/dynnxt.F90
r11715 r12453 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 … … 178 179 vn(:,:,jk) = va(:,:,jk) 179 180 END DO 181 ! limit velocities 182 IF (ln_ulimit) THEN 183 call dyn_limit_velocity (kt) 184 ENDIF 180 185 IF( .NOT.ln_linssh ) THEN ! e3._b <-- e3._n 181 186 !!gm BUG ???? I don't understand why it is not : e3._n <-- e3._a … … 210 215 END DO 211 216 END DO 217 ! limit velocities 218 IF (ln_ulimit) THEN 219 call dyn_limit_velocity (kt) 220 ENDIF 212 221 ! ! ================! 213 222 ELSE ! Variable volume ! … … 275 284 END DO 276 285 END DO 286 ! limit velocities 287 IF (ln_ulimit) THEN 288 call dyn_limit_velocity (kt) 289 ENDIF 277 290 ! 278 291 ELSE ! Asselin filter applied on thickness weighted velocity … … 302 315 END DO 303 316 END DO 317 ! limit velocities 318 IF (ln_ulimit) THEN 319 call dyn_limit_velocity (kt) 320 ENDIF 304 321 e3u_b(:,:,1:jpkm1) = ze3u_f(:,:,1:jpkm1) ! e3u_b <-- filtered scale factor 305 322 e3v_b(:,:,1:jpkm1) = ze3v_f(:,:,1:jpkm1) … … 377 394 END SUBROUTINE dyn_nxt 378 395 396 SUBROUTINE dyn_limit_velocity (kt) 397 ! limits maxming vlaues of un and vn by volume courant number 398 INTEGER, INTENT( in ) :: kt ! ocean time-step index 399 ! 400 INTEGER :: ji, jj, jk ! dummy loop indices 401 REAL(wp) :: zzu,zplim,zmlim,isp,ism,zcn,ze3e1,zzcn,zcnn,idivp,idivm 402 403 ! limit fluxes 404 zcn =cn_ulimit !0.9 ! maximum velocity inverse courant number 405 zcnn = cnn_ulimit !0.54 ! how much to reduce cn by in divergen flow 406 407 DO jk = 1, jpkm1 408 DO jj = 1, jpjm1 409 DO ji = 1, jpim1 410 ! U direction 411 zzu = un(ji,jj,jk) 412 ze3e1 = e3u_n(ji ,jj,jk) * e2u(ji ,jj) 413 ! ips is 1 if flow is positive othersize zero 414 isp = 0.5 * (sign(1.0_wp,zzu) + 1.0_wp ) 415 ism = -0.5 * (sign(1.0_wp,zzu) - 1.0_wp ) 416 ! idev is 1 if divergent flow otherwise zero 417 idivp = -isp * 0.5 * (sign(1.0_wp, un(ji-1,jj,jk)) - 1.0_wp ) 418 idivm = ism * 0.5 * (sign(1.0_wp, un(ji+1,jj,jk)) + 1.0_wp ) 419 zzcn = (idivp+idivm)*(zcnn-1.0_wp)+1.0_wp 420 zzcn = zzcn * zcn 421 zplim = zzcn * (e3t_n(ji ,jj,jk) * e1t(ji ,jj) * e2t(ji ,jj)) / & 422 (2.0*rdt * ze3e1)*umask(ji,jj,jk) 423 zmlim = -zzcn * (e3t_n(ji+1,jj,jk) * e1t(ji+1,jj) * e2t(ji+1,jj)) / & 424 (2.0*rdt * ze3e1)*umask(ji,jj,jk) 425 ! limit currents 426 un(ji,jj,jk) = min ( zzu,zplim) * isp + max(zzu,zmlim) *ism 427 ! V direction 428 zzu = vn(ji,jj,jk) 429 ze3e1 = e3v_n(ji ,jj,jk) * e1v(ji ,jj) 430 isp = 0.5 * (sign(1.0_wp,zzu) + 1.0_wp ) 431 ism = -0.5 * (sign(1.0_wp,zzu) - 1.0_wp ) 432 ! idev is 1 if divergent flow otherwise zero 433 idivp = -isp * 0.5 * (sign(1.0_wp, vn(ji,jj-1,jk)) - 1.0_wp ) 434 idivm = ism * 0.5 * (sign(1.0_wp, vn(ji,jj+1,jk)) + 1.0_wp ) 435 zzcn = (idivp+idivm)*(zcnn-1.0_wp)+1.0_wp 436 zzcn = zzcn * zcn 437 zplim = zzcn * (e3t_n(ji,jj ,jk) * e1t(ji,jj ) * e2t(ji,jj )) / & 438 (2.0*rdt * ze3e1)*vmask(ji,jj,jk) 439 zmlim = -zzcn * (e3t_n(ji,jj+1,jk) * e1t(ji,jj+1) * e2t(ji,jj+1)) / & 440 (2.0*rdt * ze3e1)*vmask(ji,jj,jk) 441 vn(ji,jj,jk) = min ( zzu,zplim) * isp + max(zzu,zmlim) *ism 442 ENDDO 443 ENDDO 444 ENDDO 445 446 END SUBROUTINE dyn_limit_velocity 447 379 448 !!========================================================================= 380 449 END MODULE dynnxt -
NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE/DYN/dynspg.F90
r11715 r12453 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/r12083_India_uncoupled/src/OCE/SBC/sbctide.F90
r11715 r12453 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/r12083_India_uncoupled/src/OCE/SBC/tide_mod.F90
r11715 r12453 16 16 PUBLIC tide_init_Wave ! called by tideini and diaharm modules 17 17 18 # if defined key_FES14_tides 19 INTEGER, PUBLIC, PARAMETER :: jpmax_harmo = 34 !: maximum number of harmonic 20 # else 18 21 INTEGER, PUBLIC, PARAMETER :: jpmax_harmo = 19 !: maximum number of harmonic 22 # endif 19 23 20 24 TYPE, PUBLIC :: tide … … 41 45 42 46 SUBROUTINE tide_init_Wave 47 # if defined key_FES14_tides 48 # include "tide_FES14.h90" 49 # else 43 50 # include "tide.h90" 51 # endif 44 52 END SUBROUTINE tide_init_Wave 45 53 … … 331 339 zf = zf * zf1 * zf1 332 340 ! 341 CASE( 20 ) !== formule 20, compound waves ( 78 x 78 x 78 x 78 ) 342 zf1 = nodal_factort(78) 343 zf = zf1 * zf1 * zf1 * zf1 344 ! 333 345 CASE( 73 ) !== formule 73 334 346 zs = sin(sh_I) -
NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE/SBC/tideini.F90
r11715 r12453 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/r12083_India_uncoupled/src/OCE/USR/usrdef_istate.F90
r11715 r12453 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/r12083_India_uncoupled/src/OCE/USR/usrdef_sbc.F90
r11715 r12453 3 3 !! *** MODULE usrdef_sbc *** 4 4 !! 5 !! === GYREconfiguration ===5 !! === WAD_TEST_CASES configuration === 6 6 !! 7 7 !! User defined : surface forcing of a user configuration … … 11 11 12 12 !!---------------------------------------------------------------------- 13 !! usrdef_sbc : user defined surface bounday conditions in GYREcase13 !! usrdef_sbc : user defined surface bounday conditions in WAD_TEST_CASES case 14 14 !!---------------------------------------------------------------------- 15 15 USE oce ! ocean dynamics and tracers … … 21 21 USE lib_mpp ! distribued memory computing library 22 22 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 23 USE lib_fortran ! 23 USE lib_fortran ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 24 24 25 25 IMPLICIT NONE … … 43 43 !! *** ROUTINE usrdef_sbc *** 44 44 !! 45 !! ** Purpose : provide at each time-step the GYREsurface boundary45 !! ** Purpose : provide at each time-step the surface boundary 46 46 !! condition, i.e. the momentum, heat and freshwater fluxes. 47 47 !! 48 !! ** Method : a nalytical seasonal cycle for GYRE configuration.48 !! ** Method : all 0 fields, for WAD_TEST_CASES case 49 49 !! CAUTION : never mask the surface stress field ! 50 50 !! 51 !! ** Action : - set t he ocean surface boundary condition, i.e.51 !! ** Action : - set to ZERO all the ocean surface boundary condition, i.e. 52 52 !! utau, vtau, taum, wndm, qns, qsr, emp, sfx 53 53 !! 54 !! Reference : Hazeleger, W., and S. Drijfhout, JPO, 30, 677-695, 2000.55 54 !!---------------------------------------------------------------------- 56 55 INTEGER, INTENT(in) :: kt ! ocean time step 57 !!58 INTEGER :: ji, jj ! dummy loop indices59 INTEGER :: zyear0 ! initial year60 INTEGER :: zmonth0 ! initial month61 INTEGER :: zday0 ! initial day62 INTEGER :: zday_year0 ! initial day since january 1st63 REAL(wp) :: ztau , ztau_sais ! wind intensity and of the seasonal cycle64 REAL(wp) :: ztime ! time in hour65 REAL(wp) :: ztimemax , ztimemin ! 21th June, and 21th decem. if date0 = 1st january66 REAL(wp) :: ztimemax1, ztimemin1 ! 21th June, and 21th decem. if date0 = 1st january67 REAL(wp) :: ztimemax2, ztimemin2 ! 21th June, and 21th decem. if date0 = 1st january68 REAL(wp) :: ztaun ! intensity69 REAL(wp) :: zemp_s, zemp_n, zemp_sais, ztstar70 REAL(wp) :: zcos_sais1, zcos_sais2, ztrp, zconv, t_star71 REAL(wp) :: zsumemp, zsurf72 REAL(wp) :: zrhoa = 1.22 ! Air density kg/m373 REAL(wp) :: zcdrag = 1.5e-3 ! drag coefficient74 REAL(wp) :: ztx, zty, zmod, zcoef ! temporary variables75 REAL(wp) :: zyydd ! number of days in one year76 56 !!--------------------------------------------------------------------- 77 zyydd = REAL(nyear_len(1),wp) 78 79 ! ---------------------------- ! 80 ! heat and freshwater fluxes ! 81 ! ---------------------------- ! 82 !same temperature, E-P as in HAZELEGER 2000 83 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 57 ! 58 IF( kt == nit000 ) THEN 59 ! 60 IF(lwp) WRITE(numout,*)' usr_sbc : WAD_TEST_CASES case: NO surface forcing' 61 IF(lwp) WRITE(numout,*)' ~~~~~~~~~~~ utau = vtau = taum = wndm = qns = qsr = emp = sfx = 0' 62 ! 63 utau(:,:) = 0._wp 64 vtau(:,:) = 0._wp 65 taum(:,:) = 0._wp 66 wndm(:,:) = 0._wp 67 ! 68 emp (:,:) = 0._wp 69 sfx (:,:) = 0._wp 70 qns (:,:) = 0._wp 71 qsr (:,:) = 0._wp 72 ! 223 73 ENDIF 224 74 ! -
NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE/par_oce.F90
r11715 r12453 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/r12083_India_uncoupled/src/OCE/step.F90
r11715 r12453 209 209 CALL dia_ar5 ( kstp ) ! ar5 diag 210 210 IF( ln_diaharm ) CALL dia_harm( kstp ) ! Tidal harmonic analysis 211 IF( lk_diaharm_fast ) & 212 & CALL dia_harm_fast( kstp ) ! Tidal harmonic analysis - restart and faster version 211 213 CALL dia_wri ( kstp ) ! ocean model: outputs 212 214 ! -
NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE/step_oce.F90
r11715 r12453 80 80 USE diahsb ! heat, salt and volume budgets (dia_hsb routine) 81 81 USE diaharm 82 USE diaharm_fast ! harmonic analysis of tides (harm_ana routine) 82 83 USE diacfl 83 84 USE diaobs ! Observation operator -
NEMO/branches/UKMO/r12083_India_uncoupled/src/OCE/stpctl.F90
r11715 r12453 115 115 zmax(1) = MAXVAL( ABS( sshn(:,:) ) ) ! ssh max 116 116 ENDIF 117 zmax(2) = MAXVAL( ABS( un(:,:,:) ) ) ! velocity max (zonal only)117 zmax(2) = MAXVAL( un(:,:,:)*un(:,:,:) + vn(:,:,:)*vn(:,:,:) ) ! velocity max 118 118 zmax(3) = MAXVAL( -tsn(:,:,:,jp_sal) , mask = tmask(:,:,:) == 1._wp ) ! minus salinity max 119 119 zmax(4) = MAXVAL( tsn(:,:,:,jp_sal) , mask = tmask(:,:,:) == 1._wp ) ! salinity max … … 149 149 IF( ( ln_ctl .OR. lsomeoce ) .AND. ( & ! have use mpp_max (because ln_ctl=.T.) or contains some ocean points 150 150 & zmax(1) > 20._wp .OR. & ! too large sea surface height ( > 20 m ) 151 & zmax(2) > 10._wp .OR. & ! too large velocity ( > 10 m/s)151 & zmax(2) > 100._wp .OR. & ! too large velocity ( > 10 m/s) 152 152 & zmax(3) >= 0._wp .OR. & ! negative or zero sea surface salinity 153 153 & zmax(4) >= 100._wp .OR. & ! too large sea surface salinity ( > 100 ) … … 156 156 IF( lk_mpp .AND. ln_ctl ) THEN 157 157 CALL mpp_maxloc( 'stpctl', ABS(sshn) , ssmask(:,:) , zzz, ih ) 158 CALL mpp_maxloc( 'stpctl', ABS(un), umask (:,:,:), zzz, iu )158 CALL mpp_maxloc( 'stpctl', un*un + vn*vn , umask (:,:,:), zzz, iu ) 159 159 CALL mpp_minloc( 'stpctl', tsn(:,:,:,jp_sal), tmask (:,:,:), zzz, is1 ) 160 160 CALL mpp_maxloc( 'stpctl', tsn(:,:,:,jp_sal), tmask (:,:,:), zzz, is2 ) 161 161 ELSE 162 162 ih(:) = MAXLOC( ABS( sshn(:,:) ) ) + (/ nimpp - 1, njmpp - 1 /) 163 iu(:) = MAXLOC( ABS( un (:,:,:) )) + (/ nimpp - 1, njmpp - 1, 0 /)163 iu(:) = MAXLOC( un(:,:,:)*un(:,:,:) + vn(:,:,:)*vn(:,:,:) ) + (/ nimpp - 1, njmpp - 1, 0 /) 164 164 is1(:) = MINLOC( tsn(:,:,:,jp_sal), mask = tmask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /) 165 165 is2(:) = MAXLOC( tsn(:,:,:,jp_sal), mask = tmask(:,:,:) == 1._wp ) + (/ nimpp - 1, njmpp - 1, 0 /) … … 187 187 ! 188 188 9100 FORMAT (' kt=',i8,' |ssh| max: ',1pg11.4,', at i j : ',2i5) 189 9200 FORMAT (' kt=',i8,' |U|max: ',1pg11.4,', at i j k: ',3i5)189 9200 FORMAT (' kt=',i8,' Vel2 max: ',1pg11.4,', at i j k: ',3i5) 190 190 9300 FORMAT (' kt=',i8,' S min: ',1pg11.4,', at i j k: ',3i5) 191 191 9400 FORMAT (' kt=',i8,' S max: ',1pg11.4,', at i j k: ',3i5) 192 9500 FORMAT(' it :', i8, ' |ssh|_max: ', D23.16, ' |U|_max: ', D23.16,' S_min: ', D23.16,' S_max: ', D23.16)192 9500 FORMAT(' it :', i8, ' |ssh|_max: ', D23.16, ' Vel2_max: ', D23.16,' S_min: ', D23.16,' S_max: ', D23.16) 193 193 ! 194 194 END SUBROUTINE stp_ctl
Note: See TracChangeset
for help on using the changeset viewer.