Changeset 3198 for LMDZ6/trunk/libf


Ignore:
Timestamp:
Feb 12, 2018, 1:24:03 AM (7 years ago)
Author:
fhourdin
Message:

Retour vers l'insensibilite au decoupage en sous domaine.
Les routines gwd_rando incluait le calcul de niveaux de reference
sur la base d'un profile pris au milieu du domaine (en klon/2).
Rempace par un test en presnivs.

Une autre intercation entre routines concernant la tke a fait apparaitre
que la tke n'était pas passee correctement au niveau klev+1 au moment
du regroupement des mailles sous les sous surface.

Ces changements garantissent la convergence numerique si
addtkeoro=0
iflag_pbl<12
et
ok_gwd_rando=n
La convergence n'est pas garantie pour les dernieres versions des physiq.def
mais les differences devraient etre mineures.

Location:
LMDZ6/trunk/libf/phylmd
Files:
5 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/phylmd/acama_gwd_rando_m.F90

    r2665 r3198  
    2020    use dimphy, only: klon, klev
    2121    use assert_m, only: assert
     22    USE ioipsl_getin_p_mod, ONLY : getin_p
     23    USE vertical_layers_mod, ONLY : presnivs
     24
    2225    include "YOMCST.h"
    2326    include "clesphys.h"
    … …  
    111114    REAL BV(KLON, KLEV + 1) ! Brunt Vaisala freq. (BVF) at 1/2 levels
    112115    REAL BVSEC ! Security to avoid negative BVF
     116
     117    REAL, DIMENSION(klev+1) ::HREF
     118    LOGICAL, SAVE :: gwd_reproductibilite_mpiomp=.true.
     119    LOGICAL, SAVE :: firstcall = .TRUE.
     120  !$OMP THREADPRIVATE(firstcall,gwd_reproductibilite_mpiomp)
     121
     122    CHARACTER (LEN=20) :: modname='flott_gwd_rando'
     123    CHARACTER (LEN=80) :: abort_message
     124
     125
     126
     127  IF (firstcall) THEN
     128    ! Cle introduite pour resoudre un probleme de non reproductibilite
     129    ! Le but est de pouvoir tester de revenir a la version precedenete
     130    ! A eliminer rapidement
     131    CALL getin_p('gwd_reproductibilite_mpiomp',gwd_reproductibilite_mpiomp)
     132    IF (NW+4*(NA-1)+NA>=KLEV) THEN
     133       abort_message = 'NW+3*NA>=KLEV Probleme pour generation des ondes'
     134       CALL abort_physic (modname,abort_message,1)
     135    ENDIF
     136    firstcall=.false.
     137!    CALL iophys_ini
     138  ENDIF
    113139
    114140    !-----------------------------------------------------------------
    … …  
    205231    ! Launching altitude
    206232
     233    IF (gwd_reproductibilite_mpiomp) THEN
     234       ! Reprend la formule qui calcule PH en fonction de PP=play
     235       DO LL = 2, KLEV
     236          HREF(LL) = EXP((LOG(presnivs(LL)) + LOG(presnivs(LL - 1))) / 2.)
     237       end DO
     238       HREF(KLEV + 1) = 0.
     239       HREF(1) = 2. * presnivs(1) - HREF(2)
     240    ELSE
     241       HREF(1:KLEV)=PH(KLON/2,1:KLEV)
     242    ENDIF
     243
    207244    LAUNCH=0
    208245    LTROP =0
    209246    DO LL = 1, KLEV
    210        IF (PH(KLON / 2, LL) / PH(KLON / 2, 1) > XLAUNCH) LAUNCH = LL
     247       IF (HREF(LL) / HREF(1) > XLAUNCH) LAUNCH = LL
    211248    ENDDO
    212249    DO LL = 1, KLEV
    213        IF (PH(KLON / 2, LL) / PH(KLON / 2, 1) > XTROP) LTROP = LL
     250       IF (HREF(LL) / HREF(1) > XTROP) LTROP = LL
    214251    ENDDO
     252    !LAUNCH=22 ; LTROP=33
     253!   print*,'LAUNCH=',LAUNCH,'LTROP=',LTROP
     254
    215255
    216256!   PRINT *,'LAUNCH IN ACAMARA:',LAUNCH
    … …  
    293333
    294334    JW = 0
    295     DO JP = 1, NP
    296        DO JK = 1, NK
    297           DO JO = 1, NO
    298              JW = JW + 1
     335    DO JW = 1, NW
    299336             ! Angle
    300337             DO II = 1, KLON
    … …  
    340377                ! RUW0(JW, II) = RUWFRT
    341378             ENDDO
    342           end DO
    343        end DO
    344379    end DO
    345380

  • LMDZ6/trunk/libf/phylmd/flott_gwd_rando_m.F90

    r2665 r3198  
    1818      use dimphy, only: klon, klev
    1919      use assert_m, only: assert
     20      USE ioipsl_getin_p_mod, ONLY : getin_p
     21      USE vertical_layers_mod, ONLY : presnivs
     22
    2023      include "YOMCST.h"
    2124      include "clesphys.h"
    … …  
    103106    REAL PH(KLON, KLEV + 1) ! Pressure at 1/2 levels
    104107    REAL PSEC ! Security to avoid division by 0 pressure
    105     REAL PHM1(KLON, KLEV + 1) ! 1/Press at 1/2 levels
    106108    REAL BV(KLON, KLEV + 1) ! Brunt Vaisala freq. (BVF) at 1/2 levels
    107109    REAL BVSEC ! Security to avoid negative BVF
     110    REAL RAN_NUM_1,RAN_NUM_2,RAN_NUM_3
     111
     112    REAL, DIMENSION(klev+1) ::HREF
     113
     114    LOGICAL, SAVE :: gwd_reproductibilite_mpiomp=.true.
     115    LOGICAL, SAVE :: firstcall = .TRUE.
     116  !$OMP THREADPRIVATE(firstcall,gwd_reproductibilite_mpiomp)
     117
     118    CHARACTER (LEN=20) :: modname='flott_gwd_rando'
     119    CHARACTER (LEN=80) :: abort_message
     120
     121
     122
     123  IF (firstcall) THEN
     124    ! Cle introduite pour resoudre un probleme de non reproductibilite
     125    ! Le but est de pouvoir tester de revenir a la version precedenete
     126    ! A eliminer rapidement
     127    CALL getin_p('gwd_reproductibilite_mpiomp',gwd_reproductibilite_mpiomp)
     128    IF (NW+3*NA>=KLEV) THEN
     129       abort_message = 'NW+3*NA>=KLEV Probleme pour generation des ondes'
     130       CALL abort_physic (modname,abort_message,1)
     131    ENDIF
     132    firstcall=.false.
     133  ENDIF
     134
    108135
    109136    !-----------------------------------------------------------------
    … …  
    156183    ZOISEC = 1.E-6 ! Security FOR 0 INTRINSIC FREQ
    157184
     185IF (1==0) THEN
    158186    !ONLINE
    159187        call assert(klon == (/size(pp, 1), size(tt, 1), size(uu, 1), &
    … …  
    167195          "FLOTT_GWD_RANDO klev")
    168196    !END ONLINE
     197ENDIF
    169198
    170199    IF(DELTAT < DTIME)THEN
    … …  
    183212    DO LL = 2, KLEV
    184213       PH(:, LL) = EXP((LOG(PP(:, LL)) + LOG(PP(:, LL - 1))) / 2.)
    185        PHM1(:, LL) = 1. / PH(:, LL)
    186     end DO
    187 
     214    end DO
    188215    PH(:, KLEV + 1) = 0.
    189     PHM1(:, KLEV + 1) = 1. / PSEC
    190216    PH(:, 1) = 2. * PP(:, 1) - PH(:, 2)
    191217
    192218    ! Launching altitude
     219
     220    !Pour revenir a la version non reproductible en changeant le nombre de process
     221    IF (gwd_reproductibilite_mpiomp) THEN
     222       ! Reprend la formule qui calcule PH en fonction de PP=play
     223       DO LL = 2, KLEV
     224          HREF(LL) = EXP((LOG(presnivs(LL)) + LOG(presnivs(LL - 1))) / 2.)
     225       end DO
     226       HREF(KLEV + 1) = 0.
     227       HREF(1) = 2. * presnivs(1) - HREF(2)
     228    ELSE
     229       HREF(1:KLEV)=PH(KLON/2,1:KLEV)
     230    ENDIF
    193231
    194232    LAUNCH=0
    195233    LTROP =0
    196234    DO LL = 1, KLEV
    197        IF (PH(KLON / 2, LL) / PH(KLON / 2, 1) > XLAUNCH) LAUNCH = LL
     235       IF (HREF(LL) / HREF(1) > XLAUNCH) LAUNCH = LL
    198236    ENDDO
    199237    DO LL = 1, KLEV
    200        IF (PH(KLON / 2, LL) / PH(KLON / 2, 1) > XTROP) LTROP = LL
     238       IF (HREF(LL) / HREF(1) > XTROP) LTROP = LL
    201239    ENDDO
     240    !LAUNCH=22 ; LTROP=33
     241!   print*,'LAUNCH=',LAUNCH,'LTROP=',LTROP
    202242
    203243    ! Log pressure vert. coordinate
    … …  
    245285    ! waves characteristics in an almost stochastic way
    246286
    247     JW = 0
    248     DO JP = 1, NP
    249        DO JK = 1, NK
    250           DO JO = 1, NO
    251              JW = JW + 1
     287    DO JW = 1, NW
    252288             ! Angle
    253289             DO II = 1, KLON
    254290                ! Angle (0 or PI so far)
    255                 ZP(JW, II) = (SIGN(1., 0.5 - MOD(TT(II, JW) * 10., 1.)) + 1.) &
     291                RAN_NUM_1=MOD(TT(II, JW) * 10., 1.)
     292                RAN_NUM_2= MOD(TT(II, JW) * 100., 1.)
     293                ZP(JW, II) = (SIGN(1., 0.5 - RAN_NUM_1) + 1.) &
    256294                     * RPI / 2.
    257295                ! Horizontal wavenumber amplitude
    258                 ZK(JW, II) = KMIN + (KMAX - KMIN) * MOD(TT(II, JW) * 100., 1.)
     296                ZK(JW, II) = KMIN + (KMAX - KMIN) *RAN_NUM_2
    259297                ! Horizontal phase speed
    260298                CPHA = 0.
    261299                DO JJ = 1, NA
     300                    RAN_NUM_3=MOD(TT(II, JW+3*JJ)**2, 1.)
    262301                    CPHA = CPHA + &
    263          CMAX*2.*(MOD(TT(II, JW+3*JJ)**2, 1.)-0.5)*SQRT(3.)/SQRT(NA*1.)
     302                    CMAX*2.*(RAN_NUM_3 -0.5)*SQRT(3.)/SQRT(NA*1.)
    264303                END DO
    265304                IF (CPHA.LT.0.)  THEN
    … …  
    276315                RUW0(JW, II) = RUWMAX
    277316             ENDDO
    278           end DO
    279        end DO
    280     end DO
     317    ENDDO
    281318
    282319    ! 4. COMPUTE THE FLUXES
    … …  
    417454    ENDDO
    418455
     456
    419457  END SUBROUTINE FLOTT_GWD_RANDO
    420458

  • LMDZ6/trunk/libf/phylmd/pbl_surface_mod.F90

    r3179 r3198  
    12941294       ENDDO
    12951295!!! jyg le 07/02/2012 et le 10/04/2013
    1296         DO k = 1, klev
     1296        DO k = 1, klev+1
    12971297          DO j = 1, knon
    12981298             i = ni(j)
    … …  
    13001300!!             ytke(j,k)   = tke(i,k,nsrf)
    13011301             ytke(j,k)   = tke_x(i,k,nsrf)
     1302          ENDDO
     1303        ENDDO
    13021304!>jyg
     1305        DO k = 1, klev
     1306          DO j = 1, knon
     1307             i = ni(j)
    13031308!FC
    13041309             y_treedrg(j,k) =  treedrg(i,k,nsrf)
    … …  
    24082413       IF (iflag_split .eq.0) THEN
    24092414        wake_dltke(:,:,nsrf) = 0.
    2410         DO k = 1, klev
     2415        DO k = 1, klev+1
    24112416           DO j = 1, knon
    24122417              i = ni(j)
    … …  
    24212426
    24222427       ELSE  ! (iflag_split .eq.0)
    2423         DO k = 1, klev
     2428        DO k = 1, klev+1
    24242429          DO j = 1, knon
    24252430            i = ni(j)

  • LMDZ6/trunk/libf/phylmd/physiq_mod.F90

    r3180 r3198  
    43904390
    43914391
    4392     CALL tend_to_tke(pdtphys,paprs,exner,t_seri,u_seri,v_seri,dtadd,duadd,dvadd,pbl_tke)
     4392    CALL tend_to_tke(pdtphys,paprs,exner,t_seri,u_seri,v_seri,dtadd,duadd,dvadd,pctsrf,pbl_tke)
    43934393
    43944394

  • LMDZ6/trunk/libf/phylmd/tend_to_tke.F90

    r3188 r3198  
    3232!**************************************************************************************
    3333
    34  SUBROUTINE tend_to_tke(dt,plev,exner,temp,windu,windv,dt_a,du_a,dv_a,tke)
     34 SUBROUTINE tend_to_tke(dt,plev,exner,temp,windu,windv,dt_a,du_a,dv_a,pctsrf,tke)
    3535
    3636 USE dimphy, ONLY: klon, klev
    … …  
    5353  REAL du_a(klon,klev)      ! Zonal wind speed tendency [m/s], grid-cell average or for a one subsurface
    5454  REAL dv_a(klon,klev)      ! Meridional wind speed tendency [m/s], grid-cell average or for a one subsurface
     55  REAL pctsrf(klon,nbsrf+1)       ! Turbulent Kinetic energy [m2/s2], grid-cell average or for a subsurface
    5556
    5657! Inputs/Outputs
    … …  
    119120
    120121
    121  DO isrf=1,nbsrf
     122 DO isrf=1,nsrf
    122123    DO k=1,klev
    123        tke(:,k,isrf)= tke(:,k,isrf)+tendu(:,k)+tendv(:,k)+tendt(:,k)
    124        tke(:,k,isrf)= max(tke(:,k,isrf),1.e-10)
     124       DO i=1,klon
     125          IF (pctsrf(i,isrf)>0.) THEN
     126            tke(i,k,isrf)= tke(i,k,isrf)+tendu(i,k)+tendv(i,k)+tendt(i,k)
     127            tke(i,k,isrf)= max(tke(i,k,isrf),1.e-10)
     128          ENDIF
     129       ENDDO
    125130    ENDDO
    126131 ENDDO
    127 
    128 ! dtke_t(:,:)=tendt(:,:)
    129 ! dtke_u(:,:)=tendu(:,:)
    130 ! dtke_v(:,:)=tendv(:,:)
    131132
    132133
Note: See TracChangeset for help on using the changeset viewer.