Ignore:
Timestamp:
2020-07-30T15:40:36+02:00 (8 weeks ago)
Author:
mathiot
Message:

ticket #1900: add subroutine icb_utl_zavg, fix issue in icb_utl_getkb, add subroutine test_icb_utl_getkb to test awkward cases (will be removed after the review), simplify icbdyn.F90 icbthm.F90 according to the new routine created.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2020/tickets_icb_1900/src/OCE/ICB/icbutl.F90

    r13357 r13359  
    4141   PUBLIC   icb_utl_copy          ! routine called in icbstp module 
    4242   PUBLIC   icb_utl_getkb         ! routine called in icbdyn and icbthm modules 
     43   PUBLIC   test_icb_utl_getkb    ! routine called in icbdyn and icbthm modules 
     44   PUBLIC   icb_utl_zavg          ! routine called in icbdyn and icbthm modules 
    4345   PUBLIC   icb_utl_interp        ! routine called in icbdyn, icbthm modules 
    4446   PUBLIC   icb_utl_bilin_h       ! routine called in icbdyn module 
     
    482484      !! 
    483485      !!---------------------------------------------------------------------- 
    484       INTEGER, INTENT(out)               :: kb 
     486      INTEGER,                INTENT(out):: kb 
    485487      REAL(wp), DIMENSION(:), INTENT(in) :: pe3 
    486488      REAL(wp),               INTENT(in) :: pD 
     
    488490      INTEGER  :: jk 
    489491      REAL(wp) :: zdepw 
    490       !! 
    491       zdepw = 0.0 
    492       kb = 1 
     492      !!---------------------------------------------------------------------- 
     493      !! 
     494      zdepw = pe3(1) ; kb = 2 
    493495      DO WHILE ( zdepw <  pD) 
    494496         zdepw = zdepw + pe3(kb) 
    495497         kb = kb + 1 
    496498      END DO 
    497       kb = kb - 1 
     499      kb = MIN(kb - 1,jpk) 
     500   END SUBROUTINE 
     501 
     502   SUBROUTINE icb_utl_zavg(pzavg, pdat, pe3, pD, kb ) 
     503      !!---------------------------------------------------------------------- 
     504      !!                ***  ROUTINE icb_utl_getkb         *** 
     505      !! 
     506      !! ** Purpose :   compute the vertical average of ocean properties affected by icb 
     507      !! 
     508      !!---------------------------------------------------------------------- 
     509      INTEGER,                INTENT(in ) :: kb        ! deepest level affected by icb 
     510      REAL(wp), DIMENSION(:), INTENT(in ) :: pe3, pdat ! vertical profile 
     511      REAL(wp),               INTENT(in ) :: pD        ! draft 
     512      REAL(wp),               INTENT(out) :: pzavg     ! z average 
     513      !!---------------------------------------------------------------------- 
     514      INTEGER  :: jk 
     515      REAL(wp) :: zdep 
     516      !!---------------------------------------------------------------------- 
     517      pzavg = 0.0 ; zdep = 0.0 
     518      DO jk = 1,kb-1 
     519         pzavg = pzavg + pe3(jk)*pdat(jk) 
     520         zdep  = zdep  + pe3(jk) 
     521      END DO 
     522      ! if kb is limited by mbkt  => bottom value is used between bathy and icb tail 
     523      ! if kb not limited by mbkt => ocean value over mask is used (ie 0.0 for u, v) 
     524      pzavg = ( pzavg + (pD - zdep)*pdat(kb)) / pD 
    498525   END SUBROUTINE 
    499526 
     
    855882   END FUNCTION icb_utl_heat 
    856883 
     884   SUBROUTINE test_icb_utl_getkb 
     885      INTEGER :: ikb 
     886      REAL(wp) :: zD, zout 
     887      REAL(wp), DIMENSION(jpk) :: ze3, zin 
     888      WRITE(numout,*) 'Test icb_utl_getkb : ' 
     889      zD = 0.0 ; ze3= 20.0 
     890      WRITE(numout,*) 'INPUT : zD = ',zD,' ze3 = ',ze3(1) 
     891      CALL icb_utl_getkb(ikb, ze3, zD) 
     892      WRITE(numout,*) 'OUTPUT : kb = ',ikb 
     893 
     894      zD = 8000000.0 ; ze3= 20.0 
     895      WRITE(numout,*) 'INPUT : zD = ',zD,' ze3 = ',ze3(1) 
     896      CALL icb_utl_getkb(ikb, ze3, zD) 
     897      WRITE(numout,*) 'OUTPUT : kb = ',ikb 
     898 
     899      zD = 80.0 ; ze3= 20.0 
     900      WRITE(numout,*) 'INPUT : zD = ',zD,' ze3 = ',ze3(1) 
     901      CALL icb_utl_getkb(ikb, ze3, zD) 
     902      WRITE(numout,*) 'OUTPUT : kb = ',ikb 
     903 
     904      zD = 85.0 ; ze3= 20.0 
     905      WRITE(numout,*) 'INPUT : zD = ',zD,' ze3 = ',ze3(1) 
     906      CALL icb_utl_getkb(ikb, ze3, zD) 
     907      WRITE(numout,*) 'OUTPUT : kb = ',ikb 
     908 
     909      zD = 75.0 ; ze3= 20.0 
     910      WRITE(numout,*) 'INPUT : zD = ',zD,' ze3 = ',ze3(1) 
     911      CALL icb_utl_getkb(ikb, ze3, zD) 
     912      WRITE(numout,*) 'OUTPUT : kb = ',ikb 
     913 
     914      WRITE(numout,*) '==================================' 
     915      WRITE(numout,*) 'Test icb_utl_zavg' 
     916      zD = 0.0 ; ze3= 20.0 ; zin=1.0 
     917      CALL icb_utl_getkb(ikb, ze3, zD) 
     918      CALL icb_utl_zavg(zout, zin, ze3, zD, ikb) 
     919      WRITE(numout,*) 'INPUT  : zD = ',zD,' ze3 = ',ze3(1),' zin = ', zin, ' ikb = ',ikb 
     920      WRITE(numout,*) 'OUTPUT : zout = ',zout 
     921 
     922      zD = 50.0 ; ze3= 20.0 ; zin=1.0; zin(3:jpk) = 0.0 
     923      CALL icb_utl_getkb(ikb, ze3, zD) 
     924      CALL icb_utl_zavg(zout, zin, ze3, zD, ikb) 
     925      WRITE(numout,*) 'INPUT  : zD = ',zD,' ze3 = ',ze3(1),' zin = ', zin, ' ikb = ',ikb 
     926      WRITE(numout,*) 'OUTPUT : zout = ',zout 
     927      CALL FLUSH(numout) 
     928 
     929      zD = 80.0 ; ze3= 20.0 ; zin=1.0; zin(3:jpk) = 0.0 
     930      CALL icb_utl_getkb(ikb, ze3, zD) 
     931      CALL icb_utl_zavg(zout, zin, ze3, zD, ikb) 
     932      WRITE(numout,*) 'INPUT  : zD = ',zD,' ze3 = ',ze3(1),' zin = ', zin, ' ikb = ',ikb 
     933      WRITE(numout,*) 'OUTPUT : zout = ',zout 
     934 
     935      zD = 80 ; ze3= 20.0 ; zin=1.0 ; zin(3:jpk) = 0.0 
     936      CALL icb_utl_getkb(ikb, ze3, zD) 
     937      ikb = 2 
     938      CALL icb_utl_zavg(zout, zin, ze3, zD, ikb) 
     939      WRITE(numout,*) 'INPUT  : zD = ',zD,' ze3 = ',ze3(1),' zin = ', zin, ' ikb = ',ikb 
     940      WRITE(numout,*) 'OUTPUT : zout = ',zout 
     941 
     942      CALL FLUSH(numout) 
     943 
     944   END SUBROUTINE test_icb_utl_getkb 
     945 
    857946   !!====================================================================== 
    858947END MODULE icbutl 
Note: See TracChangeset for help on using the changeset viewer.