Ignore:
Timestamp:
Oct 27, 2019, 5:48:03 PM (5 years ago)
Author:
fhourdin
Message:

Plein de petites corrections pour le format commun.
Semble fonctionner correctement pour ARMCU/REF et RICO/REF

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/phylmd/dyn1d/scm.F90

    r3541 r3592  
    1414       zgam, zmax0, zmea, zpic, zsig, &
    1515       zstd, zthe, zval, ale_bl, ale_bl_trig, alp_bl, ql_ancien, qs_ancien, &
    16        prlw_ancien, prsw_ancien, prw_ancien
     16       prlw_ancien, prsw_ancien, prw_ancien, &
     17       u10m,v10m,ale_wake,ale_bl_stat
     18
    1719 
    1820   USE dimphy
    … …  
    518520
    519521#include "1D_read_forc_cases.h"
     522   print*,'A d_t_adv ',d_t_adv(1:20)*86400
    520523
    521524      if (forcing_GCM2SCM) then
    … …  
    714717        v_ancien(1,:)=v(:)
    715718 
     719u10m=0.
     720v10m=0.
     721ale_wake=0.
     722ale_bl_stat=0.
     723
    716724!------------------------------------------------------------------------
    717725! Make file containing restart for the physics (startphy.nc)
    … …  
    842850
    843851#include "1D_interp_cases.h"
    844 
    845       if (forcing_GCM2SCM) then
    846         write (*,*) 'forcing_GCM2SCM not yet implemented'
    847         stop 'in time loop'
    848       endif ! forcing_GCM2SCM
     852   ! Vertical advection
     853!  call lstendH(llm,nqtot,omega,dt_dyn,dq_dyn,q,temp,u,v,play)
     854!  print*,'B d_t_adv ',d_t_adv(1:20)*86400
     855!  print*,'B dt_dyn ',dt_dyn(1:20)*86400
     856!  print*,'B dt omega ',omega
     857   teta=temp*(pzero/play)**rkappa
     858   do l=2,llm-1
     859     dt_dyn(l)=-(omega(l)*(teta(l+1)-teta(l-1))/(play(l+1)-play(l-1)))/(pzero/play(l))**rkappa
     860     dq_dyn(l,1)=-omega(l)*(q(l+1,1)-q(l-1,1))/(play(l+1)-play(l-1))
     861   enddo
     862   d_t_adv(:)=d_t_adv(:)+dt_dyn(:)
     863   d_q_adv(:,1)=d_q_adv(:,1)+dq_dyn(:,1)
    849864
    850865!---------------------------------------------------------------------
    … …  
    876891     &   presnivs(l),u(l),v(l),temp(l),q(l,1),q(l,2),omega2(l),l=1,llm)
    877892       endif
     893
     894       CALL iophys_ecrit('dtadv',klev,'dtadv','K/day',86400*d_t_adv)
     895       CALL iophys_ecrit('dtdyn',klev,'dtdyn','K/day',86400*dt_dyn)
    878896
    879897!---------------------------------------------------------------------
    … …  
    911929
    912930       fcoriolis=2.*sin(rpi*xlat/180.)*romega
    913        if (forcing_radconv .or. forcing_fire) then
    914          fcoriolis=0.0
    915          dt_cooling=0.0
    916          d_t_adv=0.0
    917          d_q_adv=0.0
    918        endif
    919 !      print*, 'calcul de fcoriolis ', fcoriolis
    920 
    921        if (forcing_toga .or. forcing_GCSSold .or. forcing_twpice            &
    922      &    .or.forcing_amma .or. forcing_type.eq.101) then
    923          fcoriolis=0.0 ; ug=0. ; vg=0.
    924        endif
    925 
    926        if(forcing_rico) then
    927           dt_cooling=0.
    928        endif
    929 
     931
     932!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     933!! CONSERVE EN ATTENDANT QUE LE CAS EN QUESTION FONCTIONNE EN STD !!
     934!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     935!       if (forcing_radconv .or. forcing_fire) then
     936!         fcoriolis=0.0
     937!         dt_cooling=0.0
     938!         d_t_adv=0.0
     939!         d_q_adv=0.0
     940!       endif
     941!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     942
     943!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     944!! CONSERVE EN ATTENDANT QUE LE CAS EN QUESTION FONCTIONNE EN STD !!
     945!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     946!      if (forcing_toga .or. forcing_GCSSold .or. forcing_twpice            &
     947!    &    .or.forcing_amma .or. forcing_type.eq.101) then
     948!        fcoriolis=0.0 ; ug=0. ; vg=0.
     949!      endif
     950!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     951
     952!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     953!! CONSERVE EN ATTENDANT QUE LE CAS EN QUESTION FONCTIONNE EN STD !!
     954!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     955!      if(forcing_rico) then
     956!         dt_cooling=0.
     957!      endif
     958!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     959
     960!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     961!! CONSERVE EN ATTENDANT QUE LE CAS EN QUESTION FONCTIONNE EN STD !!
     962!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    930963!CRio:Attention modif sp??cifique cas de Caroline
    931       if (forcing_type==-1) then
    932          fcoriolis=0.
    933 !Nudging
    934        
     964!     if (forcing_type==-1) then
     965!        fcoriolis=0.
     966!       
    935967!on calcule dt_cooling
    936         do l=1,llm
    937         if (play(l).ge.20000.) then
    938             dt_cooling(l)=-1.5/86400.
    939         elseif ((play(l).ge.10000.).and.((play(l).lt.20000.))) then
    940             dt_cooling(l)=-1.5/86400.*(play(l)-10000.)/(10000.)-1./86400.*(20000.-play(l))/10000.*(temp(l)-200.)
    941         else
    942             dt_cooling(l)=-1.*(temp(l)-200.)/86400.
    943         endif
    944         enddo
    945 
    946       endif     
    947 !RC
    948       if (forcing_sandu) then
    949          ug(1:llm)=u_mod(1:llm)
    950          vg(1:llm)=v_mod(1:llm)
    951       endif
     968!       do l=1,llm
     969!       if (play(l).ge.20000.) then
     970!           dt_cooling(l)=-1.5/86400.
     971!       elseif ((play(l).ge.10000.).and.((play(l).lt.20000.))) then
     972!           dt_cooling(l)=-1.5/86400.*(play(l)-10000.)/(10000.)-1./86400.*(20000.-play(l))/10000.*(temp(l)-200.)
     973!       else
     974!           dt_cooling(l)=-1.*(temp(l)-200.)/86400.
     975!       endif
     976!       enddo
     977!
     978!     endif     
     979!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     980
     981!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     982!! CONSERVE EN ATTENDANT QUE LE CAS EN QUESTION FONCTIONNE EN STD !!
     983!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     984!     if (forcing_sandu) then
     985!        ug(1:llm)=u_mod(1:llm)
     986!        vg(1:llm)=v_mod(1:llm)
     987!     endif
     988!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    952989
    953990      IF (prt_level >= 5) print*, 'fcoriolis, xlat,mxcalc ', &
    … …  
    9921029      endif
    9931030!
    994        if (forcing_fire) THEN
    995 
    996 !let ww=if ( alt le 1100 ) then alt*-0.00001 else 0
    997 !let wt=if ( alt le 1100 ) then min( -3.75e-5 , -7.5e-8*alt)  else 0
    998 !let wq=if ( alt le 1100 ) then max( 1.5e-8 , 3e-11*alt)  else 0
    999            d_t_adv=0.
    1000            d_q_adv=0.
    1001            teta=temp*(pzero/play)**rkappa
    1002            d_t_adv=0.
    1003            d_q_adv=0.
    1004            do l=2,llm-1
    1005               if (zlay(l)<=1100) then
    1006                   wwww=-0.00001*zlay(l)
    1007                   d_t_adv(l)=-wwww*(teta(l)-teta(l+1))/(zlay(l)-zlay(l+1)) /(pzero/play(l))**rkappa
    1008                   d_q_adv(l,1:2)=-wwww*(q(l,1:2)-q(l+1,1:2))/(zlay(l)-zlay(l+1))
    1009                   d_t_adv(l)=d_t_adv(l)+min(-3.75e-5 , -7.5e-8*zlay(l))
    1010                   d_q_adv(l,1)=d_q_adv(l,1)+max( 1.5e-8 , 3e-11*zlay(l))
    1011               endif
    1012            enddo
    1013 
    1014         endif
     1031!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1032!! CONSERVE EN ATTENDANT QUE LE CAS EN QUESTION FONCTIONNE EN STD !!
     1033!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1034!       if (forcing_fire) THEN
     1035!               print*,'Enlever cette section rapidement'
     1036!               stop
     1037!               
     1038!
     1039!!let ww=if ( alt le 1100 ) then alt*-0.00001 else 0
     1040!!let wt=if ( alt le 1100 ) then min( -3.75e-5 , -7.5e-8*alt)  else 0
     1041!!let wq=if ( alt le 1100 ) then max( 1.5e-8 , 3e-11*alt)  else 0
     1042!           d_t_adv=0.
     1043!           d_q_adv=0.
     1044!           teta=temp*(pzero/play)**rkappa
     1045!           d_t_adv=0.
     1046!           d_q_adv=0.
     1047!           do l=2,llm-1
     1048!              if (zlay(l)<=1100) then
     1049!                  wwww=-0.00001*zlay(l)
     1050!                  d_t_adv(l)=-wwww*(teta(l)-teta(l+1))/(zlay(l)-zlay(l+1)) /(pzero/play(l))**rkappa
     1051!                  d_q_adv(l,1:2)=-wwww*(q(l,1:2)-q(l+1,1:2))/(zlay(l)-zlay(l+1))
     1052!                  d_t_adv(l)=d_t_adv(l)+min(-3.75e-5 , -7.5e-8*zlay(l))
     1053!                  d_q_adv(l,1)=d_q_adv(l,1)+max( 1.5e-8 , 3e-11*zlay(l))
     1054!              endif
     1055!           enddo
     1056!
     1057!        endif
     1058!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    10151059
    10161060!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    … …  
    10261070    IF (flag_inhib_forcing == 0) then ! if tendency of forcings should be added
    10271071
     1072!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1073!! CONSERVE EN ATTENDANT QUE LE CAS EN QUESTION FONCTIONNE EN STD !!
     1074!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    10281075! pour les cas sandu et astex, on reclacule u,v,q,temp et teta dans 1D_nudge_sandu_astex.h
    10291076! au dessus de 700hpa, on relaxe vers les profils initiaux
    1030       if (forcing_sandu .OR. forcing_astex) then
    1031 #include "1D_nudge_sandu_astex.h"
    1032       else
     1077!     if (forcing_sandu .OR. forcing_astex) then
     1078!#include "1D_nudge_sandu_astex.h"
     1079!      else
     1080!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    10331081        u(1:mxcalc)=u(1:mxcalc) + timestep*(                                &
    10341082     &              du_phys(1:mxcalc)                                       &
    … …  
    10621110     &             +dt_cooling(1:mxcalc))  ! Taux de chauffage ou refroid.
    10631111
    1064       endif  ! forcing_sandu or forcing_astex
     1112       print*,'MXCALC d_t_adv ',mxcalc,d_t_adv(1:20)*86400
     1113
     1114!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1115!! CONSERVE EN ATTENDANT QUE LE CAS EN QUESTION FONCTIONNE EN STD !!
     1116!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     1117!     endif  ! forcing_sandu or forcing_astex
     1118!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    10651119
    10661120        teta=temp*(pzero/play)**rkappa
Note: See TracChangeset for help on using the changeset viewer.