- Timestamp:
- 2017-06-06T15:55:44+02:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_r7881_HPC09_ZDF/NEMOGCM/NEMO/OPA_SRC/SBC/sbcisf.F90
r7816 r8143 5 5 !! shelf 6 6 !!====================================================================== 7 !! History : 3.2 8 !! X.X 9 !! 3.4 7 !! History : 3.2 ! 2011-02 (C.Harris ) Original code isf cav 8 !! X.X ! 2006-02 (C. Wang ) Original code bg03 9 !! 3.4 ! 2013-03 (P. Mathiot) Merging + parametrization 10 10 !!---------------------------------------------------------------------- 11 11 12 12 !!---------------------------------------------------------------------- 13 !! sbc_isf 13 !! sbc_isf : update sbc under ice shelf 14 14 !!---------------------------------------------------------------------- 15 USE oce 16 USE dom_oce 17 USE phycst 18 USE eosbn2 19 USE sbc_oce 20 USE zdf bfr !15 USE oce ! ocean dynamics and tracers 16 USE dom_oce ! ocean space and time domain 17 USE phycst ! physical constants 18 USE eosbn2 ! equation of state 19 USE sbc_oce ! surface boundary condition: ocean fields 20 USE zdfdrg ! vertical physics: top/bottom drag coef. 21 21 ! 22 USE in_out_manager 23 USE iom 24 USE fldread 25 USE lbclnk 26 USE wrk_nemo 27 USE timing 28 USE lib_fortran 22 USE in_out_manager ! I/O manager 23 USE iom ! I/O manager library 24 USE fldread ! read input field at current time step 25 USE lbclnk ! 26 USE wrk_nemo ! Memory allocation 27 USE timing ! Timing 28 USE lib_fortran ! glob_sum 29 29 30 30 IMPLICIT NONE … … 77 77 CONTAINS 78 78 79 SUBROUTINE sbc_isf( kt)79 SUBROUTINE sbc_isf( kt ) 80 80 !!--------------------------------------------------------------------- 81 81 !! *** ROUTINE sbc_isf *** … … 94 94 INTEGER :: ji, jj, jk ! loop index 95 95 INTEGER :: ikt, ikb ! loop index 96 REAL(wp), DIMENSION (:,:), POINTER :: zt_frz, zdep! freezing temperature (zt_frz) at depth (zdep)96 REAL(wp), DIMENSION(jpi,jpj) :: zt_frz, zdep ! freezing temperature (zt_frz) at depth (zdep) 97 97 REAL(wp), DIMENSION(:,:,:), POINTER :: zfwfisf3d, zqhcisf3d, zqlatisf3d 98 98 REAL(wp), DIMENSION(:,: ), POINTER :: zqhcisf2d … … 100 100 ! 101 101 IF( MOD( kt-1, nn_fsbc) == 0 ) THEN 102 ! allocation103 CALL wrk_alloc( jpi,jpj, zt_frz, zdep )104 102 105 103 ! compute salt and heat flux … … 204 202 CALL wrk_dealloc( jpi,jpj, zqhcisf2d ) 205 203 END IF 206 ! deallocation207 CALL wrk_dealloc( jpi,jpj, zt_frz, zdep )208 204 ! 209 205 END IF … … 254 250 END FUNCTION 255 251 252 256 253 SUBROUTINE sbc_isf_init 257 254 !!--------------------------------------------------------------------- … … 289 286 290 287 IF ( lwp ) WRITE(numout,*) 291 IF ( lwp ) WRITE(numout,*) 'sbc_isf: heat flux of the ice shelf' 292 IF ( lwp ) WRITE(numout,*) '~~~~~~~~~' 293 IF ( lwp ) WRITE(numout,*) 'sbcisf :' 294 IF ( lwp ) WRITE(numout,*) '~~~~~~~~' 288 IF ( lwp ) WRITE(numout,*) 'sbc_isf_init : heat flux of the ice shelf' 289 IF ( lwp ) WRITE(numout,*) '~~~~~~~~~~~' 295 290 IF ( lwp ) WRITE(numout,*) ' nn_isf = ', nn_isf 296 291 IF ( lwp ) WRITE(numout,*) ' nn_isfblk = ', nn_isfblk … … 299 294 IF ( lwp ) WRITE(numout,*) ' rn_gammat0 = ', rn_gammat0 300 295 IF ( lwp ) WRITE(numout,*) ' rn_gammas0 = ', rn_gammas0 301 IF ( lwp ) WRITE(numout,*) ' rn_ tfri2 = ', rn_tfri2296 IF ( lwp ) WRITE(numout,*) ' rn_Cd0 = ', r_Cdmin_top 302 297 ! 303 298 ! Allocate public variable … … 305 300 ! 306 301 ! initialisation 307 qisf (:,:) = 0._wp ;fwfisf (:,:) = 0._wp308 risf_tsc(:,:,:) = 0._wp ;fwfisf_b(:,:) = 0._wp302 qisf (:,:) = 0._wp ; fwfisf (:,:) = 0._wp 303 risf_tsc(:,:,:) = 0._wp ; fwfisf_b(:,:) = 0._wp 309 304 ! 310 305 ! define isf tbl tickness, top and bottom indice … … 312 307 CASE ( 1 ) 313 308 rhisf_tbl(:,:) = rn_hisf_tbl 314 misfkt (:,:)= mikt(:,:) ! same indice for bg03 et cav => used in isfdiv309 misfkt (:,:) = mikt(:,:) ! same indice for bg03 et cav => used in isfdiv 315 310 316 311 CASE ( 2 , 3 ) … … 346 341 DO jj = 1, jpj 347 342 ik = 2 343 !!gm potential bug: use gdepw_0 not _n 348 344 DO WHILE ( ik <= mbkt(ji,jj) .AND. gdepw_n(ji,jj,ik) < rzisf_tbl(ji,jj) ) ; ik = ik + 1 ; END DO 349 345 misfkt(ji,jj) = ik-1 … … 354 350 ! as in nn_isf == 1 355 351 rhisf_tbl(:,:) = rn_hisf_tbl 356 misfkt (:,:)= mikt(:,:) ! same indice for bg03 et cav => used in isfdiv352 misfkt (:,:) = mikt(:,:) ! same indice for bg03 et cav => used in isfdiv 357 353 358 354 ! load variable used in fldread (use for temporal interpolation of isf fwf forcing) … … 377 373 ! determine the deepest level influenced by the boundary layer 378 374 DO jk = ikt+1, mbkt(ji,jj) 379 IF ( (SUM(e3t_n(ji,jj,ikt:jk-1)) < rhisf_tbl(ji,jj)) .AND. (tmask(ji,jj,jk) == 1) )ikb = jk375 IF( (SUM(e3t_n(ji,jj,ikt:jk-1)) < rhisf_tbl(ji,jj)) .AND. (tmask(ji,jj,jk) == 1) ) ikb = jk 380 376 END DO 381 377 rhisf_tbl(ji,jj) = MIN(rhisf_tbl(ji,jj), SUM(e3t_n(ji,jj,ikt:ikb))) ! limit the tbl to water thickness. … … 390 386 END SUBROUTINE sbc_isf_init 391 387 388 392 389 SUBROUTINE sbc_isf_bg03(kt) 393 390 !!--------------------------------------------------------------------- … … 402 399 !! interaction for climate models", Ocean Modelling 5(2003) 157-170. 403 400 !! (hereafter BG) 404 !! History : 405 !! 06-02 (C. Wang) Original code 401 !! History : 06-02 (C. Wang) Original code 406 402 !!---------------------------------------------------------------------- 407 403 INTEGER, INTENT ( in ) :: kt … … 415 411 !!---------------------------------------------------------------------- 416 412 417 IF ( nn_timing == 1 )CALL timing_start('sbc_isf_bg03')413 IF( nn_timing == 1 ) CALL timing_start('sbc_isf_bg03') 418 414 ! 419 415 DO ji = 1, jpi … … 441 437 !add to salinity trend 442 438 ELSE 443 qisf(ji,jj) = 0._wp ;fwfisf(ji,jj) = 0._wp439 qisf(ji,jj) = 0._wp ; fwfisf(ji,jj) = 0._wp 444 440 END IF 445 441 END DO 446 442 END DO 447 443 ! 448 IF( nn_timing == 1 ) CALL timing_stop('sbc_isf_bg03')444 IF( nn_timing == 1 ) CALL timing_stop('sbc_isf_bg03') 449 445 ! 450 446 END SUBROUTINE sbc_isf_bg03 447 451 448 452 449 SUBROUTINE sbc_isf_cav( kt ) … … 463 460 !! emp, emps : update freshwater flux below ice shelf 464 461 !!--------------------------------------------------------------------- 465 INTEGER, INTENT(in) :: kt! ocean time step462 INTEGER, INTENT(in) :: kt ! ocean time step 466 463 ! 467 464 INTEGER :: ji, jj ! dummy loop indices 468 465 INTEGER :: nit 466 LOGICAL :: lit 469 467 REAL(wp) :: zlamb1, zlamb2, zlamb3 470 468 REAL(wp) :: zeps1,zeps2,zeps3,zeps4,zeps6,zeps7 … … 472 470 REAL(wp) :: zeps = 1.e-20_wp 473 471 REAL(wp) :: zerr 474 REAL(wp), DIMENSION(:,:), POINTER :: zfrz 475 REAL(wp), DIMENSION(:,:), POINTER :: zgammat, zgammas 476 REAL(wp), DIMENSION(:,:), POINTER :: zfwflx, zhtflx, zhtflx_b 477 LOGICAL :: lit 472 REAL(wp), DIMENSION(jpi,jpj) :: zfrz 473 REAL(wp), DIMENSION(jpi,jpj) :: zgammat, zgammas 474 REAL(wp), DIMENSION(jpi,jpj) :: zfwflx, zhtflx, zhtflx_b 478 475 !!--------------------------------------------------------------------- 479 476 ! coeficient for linearisation of potential tfreez … … 484 481 IF( nn_timing == 1 ) CALL timing_start('sbc_isf_cav') 485 482 ! 486 CALL wrk_alloc( jpi,jpj, zfrz , zgammat, zgammas )487 CALL wrk_alloc( jpi,jpj, zfwflx, zhtflx , zhtflx_b )488 489 483 ! initialisation 490 484 zgammat(:,:) = rn_gammat0 ; zgammas (:,:) = rn_gammas0 … … 578 572 CALL iom_put('isfgammas', zgammas) 579 573 ! 580 CALL wrk_dealloc( jpi,jpj, zfrz , zgammat, zgammas )581 CALL wrk_dealloc( jpi,jpj, zfwflx, zhtflx , zhtflx_b )582 !583 574 IF( nn_timing == 1 ) CALL timing_stop('sbc_isf_cav') 584 575 ! … … 600 591 INTEGER :: ikt 601 592 INTEGER :: ji, jj ! loop index 602 REAL(wp), DIMENSION(:,:), POINTER :: zustar ! U, V at T point and friction velocity603 593 REAL(wp) :: zdku, zdkv ! U, V shear 604 594 REAL(wp) :: zPr, zSc, zRc ! Prandtl, Scmidth and Richardson number … … 614 604 REAL(wp), PARAMETER :: znu = 1.95e-6_wp ! kinamatic viscosity of sea water (m2.s-1) 615 605 REAL(wp), DIMENSION(2) :: zts, zab 606 REAL(wp), DIMENSION(jpi,jpj) :: zustar ! U, V at T point and friction velocity 616 607 !!--------------------------------------------------------------------- 617 CALL wrk_alloc( jpi,jpj, zustar )618 608 ! 619 609 SELECT CASE ( nn_gammablk ) … … 626 616 !! Jenkins et al., 2010, JPO, p2298-2312 627 617 !! Adopted by Asay-Davis et al. (2015) 628 629 !! compute ustar (eq. 24) 630 zustar(:,:) = SQRT( rn_tfri2 * (utbl(:,:) * utbl(:,:) + vtbl(:,:) * vtbl(:,:) + rn_tfeb2) ) 618 !!gm I don't understand the u* expression in those papers... (see for example zdfglf module) 619 !! for me ustar= Cd0 * |U| not (Cd0)^1/2 * |U| .... which is what you can find in Jenkins et al. 620 621 !! compute ustar (eq. 24) !! NB: here r_Cdmin_top = rn_Cd0 read in namdrg_top namelist) 622 zustar(:,:) = SQRT( r_Cdmin_top * (utbl(:,:) * utbl(:,:) + vtbl(:,:) * vtbl(:,:) + r_ke0_top) ) 631 623 632 624 !! Compute gammats … … 638 630 !! as MOL depends of flux and flux depends of MOL, best will be iteration (TO DO) 639 631 !! compute ustar 640 zustar(:,:) = SQRT( r n_tfri2 * (utbl(:,:) * utbl(:,:) + vtbl(:,:) * vtbl(:,:) + rn_tfeb2) )632 zustar(:,:) = SQRT( r_Cdmin_top * (utbl(:,:) * utbl(:,:) + vtbl(:,:) * vtbl(:,:) + r_ke0_top) ) 641 633 642 634 !! compute Pr and Sc number (can be improved) … … 649 641 650 642 !! compute gamma 651 DO ji =2,jpi652 DO jj =2,jpj643 DO ji = 2, jpi 644 DO jj = 2, jpj 653 645 ikt = mikt(ji,jj) 654 646 655 IF (zustar(ji,jj) == 0._wp) THEN ! only for kt = 1 I think647 IF( zustar(ji,jj) == 0._wp ) THEN ! only for kt = 1 I think 656 648 pgt = rn_gammat0 657 649 pgs = rn_gammas0 658 650 ELSE 659 651 !! compute Rc number (as done in zdfric.F90) 652 !!gm better to do it like in the new zdfric.F90 i.e. avm weighted Ri computation 653 !!gm moreover, use Max(rn2,0) to take care of static instabilities.... 660 654 zcoef = 0.5_wp / e3w_n(ji,jj,ikt) 661 655 ! ! shear of horizontal velocity … … 703 697 CALL lbc_lnk(pgs(:,:),'T',1.) 704 698 END SELECT 705 CALL wrk_dealloc( jpi,jpj, zustar )706 699 ! 707 700 END SUBROUTINE sbc_isf_gammats 708 701 702 709 703 SUBROUTINE sbc_isf_tbl( pvarin, pvarout, cd_ptin ) 710 704 !!---------------------------------------------------------------------- … … 714 708 !! 715 709 !!---------------------------------------------------------------------- 716 REAL(wp), DIMENSION(:,:,:), INTENT( in ) :: pvarin 717 REAL(wp), DIMENSION(:,:) , INTENT( out ) :: pvarout 718 CHARACTER(len=1), INTENT( in ) :: cd_ptin ! point of variable in/out 719 ! 720 REAL(wp) :: ze3, zhk 710 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pvarin 711 REAL(wp), DIMENSION(:,:) , INTENT( out) :: pvarout 712 CHARACTER(len=1), INTENT(in ) :: cd_ptin ! point of variable in/out 713 ! 714 INTEGER :: ji, jj, jk ! loop index 715 INTEGER :: ikt, ikb ! top and bottom index of the tbl 716 REAL(wp) :: ze3, zhk 721 717 REAL(wp), DIMENSION(:,:), POINTER :: zhisf_tbl ! thickness of the tbl 722 723 INTEGER :: ji, jj, jk ! loop index724 INTEGER :: ikt, ikb ! top and bottom index of the tbl725 718 !!---------------------------------------------------------------------- 726 719 ! allocation … … 736 729 ikt = miku(ji,jj) ; ikb = miku(ji,jj) 737 730 ! thickness of boundary layer at least the top level thickness 738 zhisf_tbl(ji,jj) = MAX( rhisf_tbl_0(ji,jj), e3u_n(ji,jj,ikt))731 zhisf_tbl(ji,jj) = MAX( rhisf_tbl_0(ji,jj) , e3u_n(ji,jj,ikt) ) 739 732 740 733 ! determine the deepest level influenced by the boundary layer … … 755 748 END DO 756 749 END DO 757 DO jj = 2,jpj 758 DO ji = 2,jpi 750 DO jj = 2, jpj 751 DO ji = 2, jpi 752 !!gm a wet-point only average should be used here !!! 759 753 pvarout(ji,jj) = 0.5_wp * (pvarout(ji,jj) + pvarout(ji-1,jj)) 760 754 END DO … … 786 780 END DO 787 781 END DO 788 DO jj = 2,jpj 789 DO ji = 2,jpi 782 DO jj = 2, jpj 783 DO ji = 2, jpi 784 !!gm a wet-point only average should be used here !!! 790 785 pvarout(ji,jj) = 0.5_wp * (pvarout(ji,jj) + pvarout(ji,jj-1)) 791 786 END DO … … 882 877 ! 883 878 END SUBROUTINE sbc_isf_div 879 884 880 !!====================================================================== 885 881 END MODULE sbcisf
Note: See TracChangeset
for help on using the changeset viewer.
