Ignore:
Timestamp:
Mar 20, 2014, 10:57:19 AM (10 years ago)
Author:
Laurent Fairhead
Message:

Merged trunk changes r1920:1997 into testing branch

Location:
LMDZ5/branches/testing
Files:
5 deleted
7 edited
1 copied

Legend:

Unmodified
Added
Removed
  • LMDZ5/branches/testing

  • LMDZ5/branches/testing/libf/phy1d/1DUTILS.h

    r1921 r1999  
    44c
    55c
    6       SUBROUTINE conf_unicol( tapedef )
     6      SUBROUTINE conf_unicol
    77c
    88#ifdef CPP_IOIPSL
     
    1515c-----------------------------------------------------------------------
    1616c     Auteurs :   A. Lahellec  .
    17 c
    18 c     Arguments :
    19 c
    20 c     tapedef   :
    21 
    22        INTEGER tapedef
    2317c
    2418c   Declarations :
     
    367361c   Variables locales pour NetCDF:
    368362c   ------------------------------
    369       INTEGER nid, nvarid
    370       INTEGER idim_s
    371       INTEGER ierr, ierr_file
    372363      INTEGER iq
    373364      INTEGER length
     
    378369      character*80 abort_message
    379370      LOGICAL found
    380 c
    381       INTEGER nb
    382371
    383372      modname = 'dyn1deta0 : '
     
    508497c   ----------
    509498      CHARACTER*(*) fichnom
    510       REAL time
    511499cAl1 plev tronque pour .nc mais plev(klev+1):=0
    512500      real :: plev(klon,klev),play (klon,klev),phi(klon,klev)
     
    520508c   Variables locales pour NetCDF:
    521509c   ------------------------------
    522       INTEGER nid, nvarid
    523       INTEGER idim_s
    524       INTEGER ierr, ierr_file
     510      INTEGER nid
     511      INTEGER ierr
    525512      INTEGER iq,l
    526513      INTEGER length
     
    535522      DATA nb / 0 /
    536523
    537       REAL zan0,zjulian,hours
    538       INTEGER yyears0,jjour0, mmois0
    539       character*30 unites
    540 
    541 cDbg
    542524      CALL open_restartphy(fichnom)
    543525      print*,'redm1 ',fichnom,klon,klev,nqtot
     
    550532      ierr = NF_OPEN(fichnom, NF_WRITE, nid)
    551533      IF (ierr .NE. NF_NOERR) THEN
    552          PRINT*, "Pb. d ouverture "//fichnom
    553          CALL abort
     534         abort_message="Pb. d ouverture "//fichnom
     535         CALL abort_gcm('Modele 1D',abort_message,1)
    554536      ENDIF
    555537
     
    661643!   traitement des point normaux
    662644         DO j=2,jm-1
    663             ig=2+(j-2)*(im-1)
     645            ig=2+(j-2)*(im-1)
    664646            CALL SCOPY(im-1,pfi(ig,ifield),1,pdyn(1,j,ifield),1)
    665             pdyn(im,j,ifield)=pdyn(1,j,ifield)
     647            pdyn(im,j,ifield)=pdyn(1,j,ifield)
    666648         ENDDO
    667649      ENDDO
     
    683665!         ierr    = severity of situation ( = 0 normal )
    684666 
    685       character*20 modname
     667      character(len=*) modname
    686668      integer ierr
    687       character*80 message
     669      character(len=*) message
    688670 
    689671      write(*,*) 'in abort_gcm'
     
    764746      RETURN
    765747      END
    766       subroutine wrgradsfi(if,nl,field,name,titlevar)
    767       implicit none
    768  
    769 !   Declarations
    770  
    771 #include "gradsdef.h"
    772  
    773 !   arguments
    774       integer if,nl
    775       real field(imx*jmx*lmx)
    776       character*10 name,file
    777       character*10 titlevar
    778  
    779 !   local
    780  
    781       integer im,jm,lm,i,j,l,iv,iii,iji,iif,ijf
    782  
    783       logical writectl
    784  
    785  
    786       writectl=.false.
    787  
    788 !     print*,if,iid(if),jid(if),ifd(if),jfd(if)
    789       iii=iid(if)
    790       iji=jid(if)
    791       iif=ifd(if)
    792       ijf=jfd(if)
    793       im=iif-iii+1
    794       jm=ijf-iji+1
    795       lm=lmd(if)
    796 
    797 
    798 !     print*,'im,jm,lm,name,firsttime(if)'
    799 !     print*,im,jm,lm,name,firsttime(if)
    800  
    801       if(firsttime(if)) then
    802          if(name.eq.var(1,if)) then
    803             firsttime(if)=.false.
    804             ivar(if)=1
    805          print*,'fin de l initialiation de l ecriture du fichier'
    806          print*,file
    807            print*,'fichier no: ',if
    808            print*,'unit ',unit(if)
    809            print*,'nvar  ',nvar(if)
    810            print*,'vars ',(var(iv,if),iv=1,nvar(if))
    811          else
    812             ivar(if)=ivar(if)+1
    813             nvar(if)=ivar(if)
    814             var(ivar(if),if)=name
    815             tvar(ivar(if),if)=trim(titlevar)
    816             nld(ivar(if),if)=nl
    817             print*,'initialisation ecriture de ',var(ivar(if),if)
    818             print*,'if ivar(if) nld ',if,ivar(if),nld(ivar(if),if)
    819          endif
    820          writectl=.true.
    821          itime(if)=1
    822       else
    823          ivar(if)=mod(ivar(if),nvar(if))+1
    824          if (ivar(if).eq.nvar(if)) then
    825             writectl=.true.
    826             itime(if)=itime(if)+1
    827          endif
    828  
    829          if(var(ivar(if),if).ne.name) then
    830            print*,'Il faut stoker la meme succession de champs a chaque'
    831            print*,'pas de temps'
    832            print*,'fichier no: ',if
    833            print*,'unit ',unit(if)
    834            print*,'nvar  ',nvar(if)
    835            print*,'vars ',(var(iv,if),iv=1,nvar(if))
    836  
    837            stop
    838          endif
    839       endif
    840  
    841 !     print*,'ivar(if),nvar(if),var(ivar(if),if),writectl'
    842 !     print*,ivar(if),nvar(if),var(ivar(if),if),writectl
    843       do l=1,nl
    844          irec(if)=irec(if)+1
    845 !        print*,'Ecrit rec=',irec(if),iii,iif,iji,ijf,
    846 !    s (l-1)*imd(if)*jmd(if)+(iji-1)*imd(if)+iii
    847 !    s ,(l-1)*imd(if)*jmd(if)+(ijf-1)*imd(if)+iif
    848          write(unit(if)+1,rec=irec(if))
    849      s   ((field((l-1)*imd(if)*jmd(if)+(j-1)*imd(if)+i)
    850      s   ,i=iii,iif),j=iji,ijf)
    851       enddo
    852       if (writectl) then
    853  
    854       file=fichier(if)
    855 !   WARNING! on reecrase le fichier .ctl a chaque ecriture
    856       open(unit(if),file=trim(file)//'.ctl',
    857      &         form='formatted',status='unknown')
    858       write(unit(if),'(a5,1x,a40)')
    859      &       'DSET ','^'//trim(file)//'.dat'
    860  
    861       write(unit(if),'(a12)') 'UNDEF 1.0E30'
    862       write(unit(if),'(a5,1x,a40)') 'TITLE ',title(if)
    863       call formcoord(unit(if),im,xd(iii,if),1.,.false.,'XDEF')
    864       call formcoord(unit(if),jm,yd(iji,if),1.,.true.,'YDEF')
    865       call formcoord(unit(if),lm,zd(1,if),1.,.false.,'ZDEF')
    866       write(unit(if),'(a4,i10,a30)')
    867      &       'TDEF ',itime(if),' LINEAR 07AUG1998 30MN '
    868       write(unit(if),'(a4,2x,i5)') 'VARS',nvar(if)
    869       do iv=1,nvar(if)
    870 !        print*,'if,var(iv,if),nld(iv,if),nld(iv,if)-1/nld(iv,if)'
    871 !        print*,if,var(iv,if),nld(iv,if),nld(iv,if)-1/nld(iv,if)
    872          write(unit(if),1000) var(iv,if),nld(iv,if)-1/nld(iv,if)
    873      &     ,99,tvar(iv,if)
    874       enddo
    875       write(unit(if),'(a7)') 'ENDVARS'
    876 !
    877 1000  format(a5,3x,i4,i3,1x,a39)
    878  
    879       close(unit(if))
    880  
    881       endif ! writectl
    882  
    883       return
    884  
    885       END
    886  
    887       subroutine inigrads(if,im
    888      s  ,x,fx,xmin,xmax,jm,y,ymin,ymax,fy,lm,z,fz
    889      s  ,dt,file,titlel)
    890  
    891  
    892       implicit none
    893  
    894       integer if,im,jm,lm,i,j,l
    895       real x(im),y(jm),z(lm),fx,fy,fz,dt
    896       real xmin,xmax,ymin,ymax
    897       integer nf
    898  
    899       character file*10,titlel*40
    900  
    901 #include "gradsdef.h"
    902  
    903       data unit/24,32,34,36,38,40,42,44,46,48/
    904       data nf/0/
    905  
    906       if (if.le.nf) stop'verifier les appels a inigrads'
    907  
    908       print*,'Entree dans inigrads'
    909  
    910       nf=if
    911       title(if)=titlel
    912       ivar(if)=0
    913  
    914       fichier(if)=trim(file)
    915  
    916       firsttime(if)=.true.
    917       dtime(if)=dt
    918  
    919       iid(if)=1
    920       ifd(if)=im
    921       imd(if)=im
    922       do i=1,im
    923          xd(i,if)=x(i)*fx
    924          if(xd(i,if).lt.xmin) iid(if)=i+1
    925          if(xd(i,if).le.xmax) ifd(if)=i
    926       enddo
    927       print*,'On stoke du point ',iid(if),'  a ',ifd(if),' en x'
    928  
    929       jid(if)=1
    930       jfd(if)=jm
    931       jmd(if)=jm
    932       do j=1,jm
    933          yd(j,if)=y(j)*fy
    934          if(yd(j,if).gt.ymax) jid(if)=j+1
    935          if(yd(j,if).ge.ymin) jfd(if)=j
    936       enddo
    937       print*,'On stoke du point ',jid(if),'  a ',jfd(if),' en y'
    938 
    939       print*,'Open de dat'
    940       print*,'file=',file
    941       print*,'fichier(if)=',fichier(if)
    942  
    943       print*,4*(ifd(if)-iid(if))*(jfd(if)-jid(if))
    944       print*,trim(file)//'.dat'
    945  
    946       OPEN (unit(if)+1,FILE=trim(file)//'.dat',
    947      s   FORM='UNFORMATTED',
    948      s   ACCESS='DIRECT'
    949      s  ,RECL=4*(ifd(if)-iid(if)+1)*(jfd(if)-jid(if)+1))
    950  
    951       print*,'Open de dat ok'
    952  
    953       lmd(if)=lm
    954       do l=1,lm
    955          zd(l,if)=z(l)*fz
    956       enddo
    957  
    958       irec(if)=0
    959 !CR
    960 !      print*,if,imd(if),jmd(if),lmd(if)
    961 !      print*,'if,imd(if),jmd(if),lmd(if)'
    962  
    963       return
    964       end
     748 
    965749      SUBROUTINE gr_dyn_fi(nfield,im,jm,ngrid,pdyn,pfi)
    966750      IMPLICIT NONE
     
    992776      DO ifield=1,nfield
    993777         DO j=2,jm-1
    994             ig=2+(j-2)*(im-1)
     778            ig=2+(j-2)*(im-1)
    995779            CALL SCOPY(im-1,pdyn(1,j,ifield),1,pfi(ig,ifield),1)
    996780         ENDDO
     
    1151935
    1152936!======================================================================
    1153        SUBROUTINE read_tsurf1d(knon,knindex,sst_out)
     937       SUBROUTINE read_tsurf1d(knon,sst_out)
    1154938
    1155939! This subroutine specifies the surface temperature to be used in 1D simulations
     
    1158942
    1159943      INTEGER, INTENT(IN)                  :: knon     ! nomber of points on compressed grid
    1160       INTEGER, DIMENSION(klon), INTENT(IN) :: knindex  ! grid point number for compressed grid
    1161944      REAL, DIMENSION(klon), INTENT(OUT)   :: sst_out  ! tsurf used to force the single-column model
    1162945
     
    12201003
    12211004       SUBROUTINE advect_va(llm,omega,d_t_va,d_q_va,d_u_va,d_v_va,
    1222      !                q,temp,u,v,
    1223      !            play,plev)
     1005     s                q,temp,u,v,play)
    12241006!itlmd
    12251007!----------------------------------------------------------------------
     
    12371019        real  q(llm,3),temp(llm)
    12381020        real  u(llm),v(llm)
    1239         real  play(llm),plev(llm+1)
     1021        real  play(llm)
    12401022! interne
    12411023        integer l
     
    13231105        real dph(llm),dqdp(llm),dtdp(llm)
    13241106! interne
    1325         integer l,k
    1326         real alpha,omdn,omup
     1107        integer k
     1108        real omdn,omup
    13271109
    13281110!        dudp=0.
     
    14031185      character*80 fich_toga
    14041186
    1405       integer no,l,k,ip
    1406       real riy,rim,rid,rih,bid
     1187      integer k,ip
     1188      real bid
    14071189
    14081190      integer iy,im,id,ih
     
    14221204
    14231205       do k = 1, nlev_toga
    1424          read(21,230) plev_toga(k,ip), t_toga(k,ip), q_toga(k,ip) 
     1206         read(21,230) plev_toga(k,ip), t_toga(k,ip), q_toga(k,ip)
    14251207     :       ,u_toga(k,ip), v_toga(k,ip), w_toga(k,ip)
    14261208     :       ,ht_toga(k,ip), vt_toga(k,ip), hq_toga(k,ip), vq_toga(k,ip)
     
    14431225
    14441226  223 format(4i3,6f8.2)
    1445   226 format(f7.1,1x,10f8.2)
    1446   227 format(f7.1,1x,1p,4e11.3)
    14471227  230 format(6f9.3,4e11.3)
    14481228
     
    14621242      character*80 fich_sandu
    14631243
    1464       integer no,l,k,ip
    1465       real riy,rim,rid,rih,bid
    1466 
     1244      integer ip
    14671245      integer iy,im,id,ih
    14681246
    1469        real plev_min
    1470 
    1471        plev_min = 55000.  ! pas de tendance de vap. d eau au-dessus de 55 hPa
     1247      real plev_min
     1248
     1249      print*,'nlev_sandu',nlev_sandu
     1250      plev_min = 55000.  ! pas de tendance de vap. d eau au-dessus de 55 hPa
    14721251
    14731252      open(21,file=trim(fich_sandu),form='formatted')
     
    14821261
    14831262  223 format(4i3,f8.2)
    1484   226 format(f7.1,1x,10f8.2)
    1485   227 format(f7.1,1x,1p,4e11.3)
    1486   230 format(6f9.3,4e11.3)
    14871263
    14881264          return
     
    15041280      character*80 fich_astex
    15051281
    1506       integer no,l,k,ip
    1507       real riy,rim,rid,rih,bid
    1508 
     1282      integer ip
    15091283      integer iy,im,id,ih
    15101284
    15111285       real plev_min
    15121286
     1287      print*,'nlev_astex',nlev_astex
    15131288       plev_min = 55000.  ! pas de tendance de vap. d eau au-dessus de 55 hPa
    15141289
     
    15281303
    15291304  223 format(4i3,e13.2,f7.2,f7.3,f7.2,f7.3,f7.2)
    1530   226 format(f7.1,1x,10f8.2)
    1531   227 format(f7.1,1x,1p,4e11.3)
    1532   230 format(6f9.3,4e11.3)
    15331305
    15341306          return
     
    15511323      character*80 :: fich_twpice
    15521324      real*8 time(ntime)
    1553       real*8 lat, lon, alt, phis       
     1325      real*8 lat, lon, alt, phis
    15541326      real*8 lev(nlevel)
    15551327      real*8 plev(nlevel,ntime)
     
    15611333      real*8 T_adv_h(nlevel,ntime)
    15621334      real*8 T_adv_v(nlevel,ntime), q_adv_h(nlevel,ntime)
    1563       real*8 q_adv_v(nlevel,ntime)     
     1335      real*8 q_adv_v(nlevel,ntime)
    15641336      real*8 s(nlevel,ntime), s_adv_h(nlevel,ntime)
    15651337      real*8 s_adv_v(nlevel,ntime)
     
    19761748         integer ierr
    19771749
    1978          integer i
    19791750         integer timevar,levvar
    19801751         integer timelen,levlen
     
    20501821       real omega_mod(llm),o3mmr_mod(llm)
    20511822
    2052        integer l,k,k1,k2,kp
    2053        real aa,frac,frac1,frac2,fact
     1823       integer l,k,k1,k2
     1824       real frac,frac1,frac2,fact
    20541825
    20551826       do l = 1, llm
     
    21681939       real o3mmr_mod(llm),ql_mod(llm),qt_mod(llm)
    21691940
    2170        integer l,k,k1,k2,kp
    2171        real aa,frac,frac1,frac2,fact
     1941       integer l,k,k1,k2
     1942       real frac,frac1,frac2,fact
    21721943
    21731944       do l = 1, llm
     
    24442215        real ts_prof
    24452216! local:
    2446         integer it_sandu1, it_sandu2,k
     2217        integer it_sandu1, it_sandu2
    24472218        real timeit,time_sandu1,time_sandu2,frac
    24482219! Check that initial day of the simulation consistent with SANDU period:
     
    25112282      character*80 fich_armcu
    25122283
    2513       integer no,l,k,ip
    2514       real riy,rim,rid,rih,bid
     2284      integer ip
    25152285
    25162286      integer iy,im,id,ih,in
     2287
     2288      print*,'nlev_armcu',nlev_armcu
    25172289
    25182290      open(21,file=trim(fich_armcu),form='formatted')
     
    25292301
    25302302  223 format(5i3,5f8.3)
    2531   226 format(f7.1,1x,10f8.2)
    2532   227 format(f7.1,1x,1p,4e11.3)
    2533   230 format(6f9.3,4e11.3)
    25342303
    25352304          return
     
    25712340       real hq_mod(llm),vq_mod(llm)
    25722341 
    2573        integer l,k,k1,k2,kp
    2574        real aa,frac,frac1,frac2,fact
     2342       integer l,k,k1,k2
     2343       real frac,frac1,frac2,fact
    25752344 
    25762345       do l = 1, llm
     
    26842453        real div_prof,ts_prof,ug_prof,vg_prof,ufa_prof,vfa_prof
    26852454! local:
    2686         integer it_astex1, it_astex2,k
     2455        integer it_astex1, it_astex2
    26872456        real timeit,time_astex1,time_astex2,frac
    26882457
     
    29682737
    29692738!=====================================================================
    2970       subroutine readprofiles(nlev_max,kmax,height,
     2739      subroutine readprofiles(nlev_max,kmax,ntrac,height,
    29712740     .           thlprof,qtprof,uprof,
    29722741     .           vprof,e12prof,ugprof,vgprof,
    29732742     .           wfls,dqtdxls,dqtdyls,dqtdtls,
    2974      .           thlpcar)
     2743     .           thlpcar,tracer,nt1,nt2)
    29752744      implicit none
    29762745
    2977         integer nlev_max,kmax,kmax2
     2746        integer nlev_max,kmax,kmax2,ntrac
    29782747        logical :: llesread = .true.
    29792748
     
    29822751     .       ugprof(nlev_max),vgprof(nlev_max),wfls(nlev_max),
    29832752     .       dqtdxls(nlev_max),dqtdyls(nlev_max),dqtdtls(nlev_max),
    2984      .           thlpcar(nlev_max)
     2753     .           thlpcar(nlev_max),tracer(nlev_max,ntrac)
    29852754
    29862755        integer, parameter :: ilesfile=1
    2987         integer :: ierr,irad,imax,jtot,k
    2988         logical :: lmoist,lcoriol,ltimedep
    2989         real :: xsize,ysize
    2990         real :: ustin,wsvsurf,timerad
    2991         character(80) :: chmess
     2756        integer :: ierr,k,itrac,nt1,nt2
    29922757
    29932758        if(.not.(llesread)) return
     
    30162781        close(ilesfile)
    30172782
     2783       open(ilesfile,file='trac.inp.001',status='old',iostat=ierr)
     2784        if (ierr /= 0) then
     2785            print*,'WARNING : trac.inp does not exist'
     2786        else
     2787        read (ilesfile,*) kmax2,nt1,nt2
     2788        if (nt2>ntrac) then
     2789          stop'Augmenter le nombre de traceurs dans traceur.def'
     2790        endif
     2791        if (kmax .ne. kmax2) then
     2792          print *, 'fichiers prof.inp et lscale.inp incompatibles :'
     2793          print *, 'nbre de niveaux : ',kmax,' et ',kmax2
     2794          stop 'lecture profiles'
     2795        endif
     2796        do k=1,kmax
     2797          read (ilesfile,*) height(k),(tracer(k,itrac),itrac=nt1,nt2)
     2798        end do
     2799        close(ilesfile)
     2800        endif
     2801
    30182802        return
    30192803        end
     
    30242808      implicit none
    30252809
    3026         integer nlev_max,kmax,kmax2
     2810        integer nlev_max,kmax
    30272811        logical :: llesread = .true.
    30282812
     
    30332817
    30342818        integer, parameter :: ilesfile=1
    3035         integer :: ierr,irad,imax,jtot,k
    3036         logical :: lmoist,lcoriol,ltimedep
    3037         real :: xsize,ysize
    3038         real :: ustin,wsvsurf,timerad
    3039         character(80) :: chmess
     2819        integer :: k,ierr
    30402820
    30412821        if(.not.(llesread)) return
     
    30602840      implicit none
    30612841
    3062         integer nlev_max,kmax,kmax2
     2842        integer nlev_max,kmax
    30632843        logical :: llesread = .true.
    30642844
     
    30692849
    30702850        integer, parameter :: ilesfile=1
    3071         integer :: ierr,irad,imax,jtot,k
    3072         logical :: lmoist,lcoriol,ltimedep
    3073         real :: xsize,ysize
    3074         real :: ustin,wsvsurf,timerad
    3075         character(80) :: chmess
     2851        integer :: ierr,k
    30762852
    30772853        if(.not.(llesread)) return
     
    30982874      implicit none
    30992875
    3100         integer nlev_max,kmax,kmax2
     2876        integer nlev_max,kmax
    31012877        logical :: llesread = .true.
    31022878
     
    31082884        integer, parameter :: ilesfile=1
    31092885        integer, parameter :: ifile=2
    3110         integer :: ierr,irad,imax,jtot,k
    3111         logical :: lmoist,lcoriol,ltimedep
    3112         real :: xsize,ysize
    3113         real :: ustin,wsvsurf,timerad
    3114         character(80) :: chmess
     2886        integer :: ierr,jtot,k
    31152887
    31162888        if(.not.(llesread)) return
     
    31632935
    31642936      integer ntime,nlevel
    3165       integer l,k
    31662937      character*80 :: fich_amma
    3167       real*8 time(ntime)
    3168       real*8 zz(nlevel)
     2938      real*8 zz(nlevel)
    31692939
    31702940      real*8 temp(nlevel),pp(nlevel)
     
    31732943      real*8 dw(nlevel,ntime)
    31742944      real*8 dt(nlevel,ntime)
    3175       real*8 dq(nlevel,ntime)   
     2945      real*8 dq(nlevel,ntime)
    31762946      real*8 flat(ntime),sens(ntime)
    31772947
     
    35033273
    35043274      integer ntime,nlevel
    3505       integer l,k
    35063275      character*80 :: fich_fire
    3507       real*8 time(ntime)
    3508       real*8 zz(nlevel)
     3276      real*8 zz(nlevel)
    35093277
    35103278      real*8 thl(nlevel)
     
    35133281      real*8 ug(nlevel,ntime),vg(nlevel,ntime),wls(nlevel,ntime)
    35143282      real*8 dqtdx(nlevel,ntime),dqtdy(nlevel,ntime)
    3515       real*8 dqtdt(nlevel,ntime),thl_rad(nlevel,ntime) 
     3283      real*8 dqtdt(nlevel,ntime),thl_rad(nlevel,ntime)
    35163284
    35173285      integer nid, ierr
  • LMDZ5/branches/testing/libf/phy1d/1D_decl_cases.h

    r1910 r1999  
    1616        real    sec_print
    1717!!
    18         integer nn
    19         integer it_toga1, it_toga2
    20         real time_toga1,time_toga2
    21 
    2218        real ts_toga(nt_toga)
    2319        real plev_toga(nlev_toga,nt_toga),w_toga(nlev_toga,nt_toga)
     
    3430        real hq_prof(nlev_toga),vq_prof(nlev_toga)
    3531
    36         real plev_mod(llm),w_mod(llm), t_mod(llm),q_mod(llm)
     32        real w_mod(llm), t_mod(llm),q_mod(llm)
    3733        real u_mod(llm),v_mod(llm), ht_mod(llm),vt_mod(llm)
    3834        real hq_mod(llm),vq_mod(llm),qv_mod(llm),ql_mod(llm),qt_mod(llm)
     
    8783        character*80 :: fich_amma
    8884! Option du cas AMMA ou on impose la discretisation verticale (Ap,Bp)
    89         logical  :: fixe_disvert=.true.
    9085        integer nlev_amma, nt_amma
    9186!       parameter (nlev_amma=29, nt_amma=48)  ! Fleur, juillet 2012
     
    10499!profils initiaux:
    105100        real plev_amma(nlev_amma)
    106         real tv_amma(nlev_amma),rho_amma(nlev_amma)
    107         real thv_amma(nlev_amma)
    108101       
    109102        real z_amma(nlev_amma)
     
    111104        real u_amma(nlev_amma)
    112105        real v_amma(nlev_amma)
    113 
    114         real thvsurf_amma,tvsurf_amma,rhosurf_amma,thsurf
    115106
    116107        real th_ammai(nlev_amma),q_ammai(nlev_amma)
     
    130121
    131122!champs interpoles
    132         real plev_profamma(nlev_amma),vitw_profamma(nlev_amma)
     123        real vitw_profamma(nlev_amma)
    133124        real ht_profamma(nlev_amma)
    134125        real hq_profamma(nlev_amma)
     
    148139        integer year_ini_fire, day_ini_fire, mth_ini_fire
    149140        real heure_ini_fire
    150         real day_ju_ini_fire   ! Julian day of fire first day
    151141        parameter (year_ini_fire=1987)
    152142        parameter (mth_ini_fire=7)
     
    154144        parameter (heure_ini_fire=0.) !0h en secondes
    155145
    156 !profils initiaux:
    157         real z_fire(nlev_fire)
    158         real thl_fire(nlev_fire),qt_fire(nlev_fire)
    159         real u_fire(nlev_fire), v_fire(nlev_fire)
    160         real tke_fire(nlev_fire)
    161        
    162 !forcings
    163         real ugeo_fire(nlev_fire),vgeo_fire(nlev_fire)
    164         real wls_fire(nlev_fire),dqtdx_fire(nlev_fire)
    165         real dqtdy_fire(nlev_fire)
    166         real dqtdt_fire(nlev_fire),thl_rad_fire(nlev_fire)
    167          
    168         real ugeo_mod(llm),vgeo_mod(llm),wls_mod(llm)
    169         real dqtdx_mod(llm),dqtdy_mod(llm),dqtdt_mod(llm)
    170         real thl_rad_mod(llm)
    171146!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    172147! Declarations specifiques au cas GCSSold
     
    180155        real  Ts_gcssold
    181156        real  dtime_frcg
    182 cAl1     logical :: imp_fcg_gcssold
    183 c        logical :: ts_fcg_gcssold
    184 c        logical :: Tp_fcg_gcssold
    185157        logical :: Turb_fcg_gcssold
    186         common /turb_forcing/ dtime_frcg,
    187      $      Turb_fcg_gcssold, hthturb_gcssold, hqturb_gcssold
     158
     159        common /turb_forcing/
     160     s  dtime_frcg,hthturb_gcssold, hqturb_gcssold,Turb_fcg_gcssold
    188161!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    189162! Declarations specifiques au cas Arm_cu
     
    206179        real adv_qt_armcu(nt_armcu)
    207180        real theta_mod(llm),rv_mod(llm),play_mod(llm)
    208         real d_t_dyn_ls(llm),d_q_dyn_ls(llm)
    209181! profc comme "profil armcu"
    210         real h_profc,play_profc,t_profc,th_profc,plev_profc
    211         real u_profc,v_profc,qv_profc,rv_profc
    212182       
    213183! forcages interpoles dans le temps
     
    229199        logical  :: trouve_700=.true.
    230200        parameter (dt_sandu=6.*3600.)   ! forcages donnes ttes les 6 heures par ifa_sandu.txt
    231 !       parameter (tau_sandu=3600.)  ! temps de relaxation u,v,thetal,qt vers profil init et au dessus 700hPa
     201        parameter (tau_sandu=3600.)  ! temps de relaxation u,v,thetal,qt vers profil init et au dessus 700hPa
    232202!!
    233         integer it_sandu1, it_sandu2
    234         real time_sandu1,time_sandu2
    235 
    236203        real ts_sandu(nt_sandu)
    237204! profs comme "profil sandu"
     
    242209        real omega_profs(nlev_sandu),o3mmr_profs(nlev_sandu)
    243210
     211        real, dimension(llm) :: relax_u,relax_v,relax_thl
     212        real, dimension(llm,2) :: relax_q
     213
    244214        real thl_mod(llm),omega_mod(llm),o3mmr_mod(llm),tke_mod(llm)
    245 ! pour relaxer u,v,thl et qt vers les profils initiaux au dessus de 700hPa
    246         real relax_u(llm),relax_v(llm),relax_thl(llm),relax_q(llm,2)
    247215!vertical advection computation
    248216        real d_t_z(llm), d_q_z(llm)
     
    260228        parameter (mth_ini_astex=6)
    261229        parameter (day_ini_astex=13)  ! 165 = 13 juin 1992
    262         real dt_astex, tau_astex
     230        real dt_astex
    263231        parameter (dt_astex=3600.)    ! forcages donnes ttes les heures par ifa_astex.txt
    264         integer it_astex1, it_astex2
    265         real time_astex1,time_astex2
    266232        real ts_astex(nt_astex),div_astex(nt_astex),ug_astex(nt_astex)
    267233        real vg_astex(nt_astex),ufa_astex(nt_astex),vfa_astex(nt_astex)
  • LMDZ5/branches/testing/libf/phy1d/1D_read_forc_cases.h

    r1910 r1999  
    33! forcing_radconv = .T. : Pure radiative-convective equilibrium:
    44!----------------------------------------------------------------------
     5
     6
     7      nq1=0
     8      nq2=0
    59
    610      if (forcing_les .or. forcing_radconv
     
    2428!----------------------------------------------------------------------
    2529
    26       call readprofiles(nlev_max,kmax,height,
     30      call readprofiles(nlev_max,kmax,nqtot,height,
    2731     .           tttprof,qtprof,uprof,vprof,
    2832     .           e12prof,ugprof,vgprof,
    2933     .           wfls,dqtdxls,dqtdyls,dqtdtls,
    30      .           thlpcar)
     34     .           thlpcar,qprof,nq1,nq2)
    3135      endif
    3236
     
    6266        ug(l)   = ugprof(kmax)-frac*( ugprof(kmax)- ugprof(kmax-1))
    6367        vg(l)   = vgprof(kmax)-frac*( vgprof(kmax)- vgprof(kmax-1))
     68        IF (nq2>0) q(l,nq1:nq2)=qprof(kmax,nq1:nq2)
     69     s               -frac*(qprof(kmax,nq1:nq2)-qprof(kmax-1,nq1:nq2))
    6470        omega(l)=   wfls(kmax)-frac*(   wfls(kmax)-   wfls(kmax-1))
    6571
     
    8591            ug(l)   = ugprof(k)-frac*( ugprof(k)- ugprof(k-1))
    8692            vg(l)   = vgprof(k)-frac*( vgprof(k)- vgprof(k-1))
     93            IF (nq2>0) q(l,nq1:nq2)=qprof(k,nq1:nq2)
     94     s                   -frac*(qprof(k,nq1:nq2)-qprof(k-1,nq1:nq2))
    8795            omega(l)=   wfls(k)-frac*(   wfls(k)-   wfls(k-1))
    8896            dq_dyn(l,1)=dqtdtls(k)-frac*(dqtdtls(k)-dqtdtls(k-1))
     
    104112            vg(l)   = vgprof(1)
    105113            omega(l)=   wfls(1)
     114            IF (nq2>0) q(l,nq1:nq2)=qprof(1,nq1:nq2)
    106115            dq_dyn(l,1)  =dqtdtls(1)
    107116            dt_cooling(l)=thlpcar(1)
  • LMDZ5/branches/testing/libf/phy1d/1Dconv.h

    r1910 r1999  
    2222      REAL hplaym(100) !pression en hPa milieux des couches Meso-NH
    2323
    24       integer i,j,k,ii,ll,in
    25       REAL tsol,qsol
     24      integer i,j,k,ll,in
    2625
    2726      CHARACTER*80 file_forctl,file_fordat
    2827
    29       COMMON/com1_phys_gcss/klev,play,JM,coef1,coef2
    30       COMMON/com2_phys_gcss/nblvlm,playm,hplaym
     28      COMMON/com1_phys_gcss/play,coef1,coef2,JM,klev
     29      COMMON/com2_phys_gcss/playm,hplaym,nblvlm
    3130
    3231c======================================================================
     
    170169             hqTurbbef(i)=hqTurbaft(i)
    171170          enddo
    172                tsbef = tsaft
     171          tsbef = tsaft
    173172          timebef=pasprev*dt
    174173          timeaft=timebef+dt
     
    213212               print*,'hqTurb_mes ',(hqTurb_mes(i),i=1,nblvlm)
    214213                  endif
    215                IF (ts_fcg) print*,'ts_subr', ts_subr
     214          IF (ts_fcg) print*,'ts_subr', ts_subr
    216215c*** on interpole les champs meso_NH sur les niveaux de pression***
    217216c*** gcm . on obtient le nouveau champ after                    ***
     
    263262               hqTurb(ll)=hqTurbaft(ll)
    264263          enddo
    265                ts_subr = tsaft
     264          ts_subr = tsaft
    266265      else ! temps.ge.pasmax
    267266c*** on interpole sur les pas de temps de 10mn du gcm a partir   ***
     
    282281             endif ! Turb_fcg
    283282         enddo
    284               ts_subr = ((timeaft-time)*tsbef + (time-timebef)*tsaft)/dt
     283         ts_subr = ((timeaft-time)*tsbef + (time-timebef)*tsaft)/dt
    285284       endif ! temps.ge.pasmax
    286285c
     
    440439                endif ! Turb_fcg
    441440          enddo
    442                ts_subr = tsaft
     441          ts_subr = tsaft
    443442          close(99)
    444443          close(98)
     
    505504      REAL hplaym(100)!pression en hecto-Pa des milieux de couche Meso-NH
    506505
    507       COMMON/com1_phys_gcss/klev,play,JM,coef1,coef2
    508       COMMON/com2_phys_gcss/nblvlm,playm,hplaym
    509 
    510       integer i,k,klevgcm
     506      COMMON/com1_phys_gcss/play,coef1,coef2,JM,klev
     507      COMMON/com2_phys_gcss/playm,hplaym,nblvlm
     508
     509      integer k,klevgcm
    511510      real playgcm(klevgcm) ! pression en milieu de couche du gcm
    512511      real psolgcm
     
    577576      REAL playm(100)  !pression en Pa milieu de chaque couche Meso-NH
    578577      REAL hplaym(100) !pression en hPa des milieux de couche Meso-NH
    579       COMMON/com2_phys_gcss/nblvlm,playm,hplaym
    580 
    581       INTEGER i,lu,k,mlz,mlzh,j
     578      COMMON/com2_phys_gcss/playm,hplaym,nblvlm
     579
     580      INTEGER i,lu,mlz,mlzh
    582581 
    583582      character*80 file_forctl
     
    644643      real ts
    645644c
    646       INTEGER i, k
     645      INTEGER k
    647646c
    648647      LOGICAL imp_fcg,ts_fcg,Turb_fcg
     
    725724      REAL hplaym(100)!pression en hPa des milieux de couche Meso-NH
    726725
    727       COMMON/com1_phys_gcss/klev,play,JM,coef1,coef2
    728       COMMON/com2_phys_gcss/nblvlm,playm,hplaym
     726      COMMON/com1_phys_gcss/play,coef1,coef2,JM,klev
     727      COMMON/com2_phys_gcss/playm,hplaym,nblvlm
    729728
    730729      REAL psol
    731730      REAL val
    732       INTEGER k, mlz, mlzh
     731      INTEGER k, mlz
    733732
    734733
  • LMDZ5/branches/testing/libf/phy1d/compar1d.h

    r1910 r1999  
    2727      logical :: ok_old_disvert
    2828
    29       common/com_par1d/forcing_type,nat_surf,tsurf,rugos,               &
     29      common/com_par1d/
     30     & nat_surf,tsurf,rugos,                                            &
    3031     & qsol,qsurf,psurf,zsurf,albedo,time,time_ini,xlat,xlon,airefi,    &
    3132     & wtsurf,wqsurf,restart_runoff,xagesno,qsolinp,zpicinp,            &
     33     & forcing_type,
    3234     & restart,ok_old_disvert
    3335
  • LMDZ5/branches/testing/libf/phy1d/lmdz1d.F

    r1921 r1999  
    8181      integer :: an
    8282 
    83 !
    84       real :: paire    = 1.     ! aire de la maille
    85 !**      common /flux_arp/fsens,flat,ok_flux_surf
    86 
    8783!---------------------------------------------------------------------
    8884!  Declarations related to forcing and initial profiles
     
    9086
    9187        integer :: kmax = llm
    92         integer nlev_max,llm700
    93         parameter (nlev_max = 1000)
    94         real timestep, frac, timeit
     88        integer llm700,nq1,nq2
     89        INTEGER, PARAMETER :: nlev_max=1000, nqmx=1000
     90        real timestep, frac
    9591        real height(nlev_max),tttprof(nlev_max),qtprof(nlev_max),
    9692     .              uprof(nlev_max),vprof(nlev_max),e12prof(nlev_max),
    9793     .              ugprof(nlev_max),vgprof(nlev_max),wfls(nlev_max),
    9894     .              dqtdxls(nlev_max),dqtdyls(nlev_max),
    99      .              dqtdtls(nlev_max),thlpcar(nlev_max)
    100 
    101         real    :: fff
     95     .              dqtdtls(nlev_max),thlpcar(nlev_max),
     96     .              qprof(nlev_max,nqmx)
     97
    10298c        integer :: forcing_type
    10399        logical :: forcing_les     = .false.
     
    143139!---------------------------------------------------------------------
    144140
    145       integer :: iq
    146141      real :: phi(llm)
    147142      real :: teta(llm),tetal(llm),temp(llm),u(llm),v(llm),w(llm)
     
    151146      real :: sfdt, cfdt
    152147      real :: du_phys(llm),dv_phys(llm),dt_phys(llm)
    153       real :: du_dyn(llm),dv_dyn(llm),dt_dyn(llm)
    154       real :: dt_cooling(llm),d_t_cool(llm),d_th_adv(llm)
    155       real :: dq_cooling(llm),d_q_cool(llm)
    156       real :: tmpvar(llm)
     148      real :: dt_dyn(llm)
     149      real :: dt_cooling(llm),d_th_adv(llm)
    157150      real :: alpha
     151      real :: ttt
    158152
    159153      REAL, ALLOCATABLE, DIMENSION(:,:):: q
     
    161155      REAL, ALLOCATABLE, DIMENSION(:,:):: dq_dyn
    162156      REAL, ALLOCATABLE, DIMENSION(:,:):: d_q_adv
     157!     REAL, ALLOCATABLE, DIMENSION(:):: d_th_adv
    163158
    164159!---------------------------------------------------------------------
     
    204199!  Fichiers et d'autres variables
    205200!---------------------------------------------------------------------
    206       real ttt,bow,q1
    207       integer :: ierr,k,l,i,it=1,mxcalc
     201      integer :: k,l,i,it=1,mxcalc
    208202      integer jjmp1
    209203      parameter (jjmp1=jjm+1-1/jjm)
     
    230224!---------------------------------------------------------------------
    231225cAl1
    232         call conf_unicol(99)
     226        call conf_unicol
    233227cAl1 moves this gcssold var from common fcg_gcssold to
    234228        Turb_fcg_gcssold = xTurb_fcg_gcssold
     
    357351c      Le numero du jour est dans "day". L heure est traitee separement.
    358352c      La date complete est dans "daytime" (l'unite est le jour).
    359       fnday=nday
     353      if (nday>0) then
     354         fnday=nday
     355      else
     356         fnday=-nday/float(day_step)
     357      endif
     358
    360359c Special case for arm_cu which lasts less than one day : 53100s !! (MPL 20111026)
    361360      IF(forcing_type .EQ. 61) fnday=53100./86400.
     
    369368      itau_phy = 0
    370369      call ymds2ju(annee_ref,mois,day_ref,heure,day)
    371       day_ini = day
    372       day_end = day_ini + nday
     370      day_ini = int(day)
     371      day_end = day_ini + fnday
    373372
    374373      IF (forcing_type .eq.2) THEN
     
    422421      call infotrac_init
    423422
     423      if (nqtot>nqmx) STOP'Augmenter nqmx dans lmdz1d.F'
    424424      allocate(q(llm,nqtot)) ; q(:,:)=0.
    425425      allocate(dq(llm,nqtot))
    426426      allocate(dq_dyn(llm,nqtot))
    427427      allocate(d_q_adv(llm,nqtot))
     428!     allocate(d_th_adv(llm))
    428429
    429430c
     
    463464!! mpl et jyg le 22/08/2012 :
    464465!!  pour que les cas a flux de surface imposes marchent
    465       IF(.NOT.ok_flux_surf) THEN
     466      IF(.NOT.ok_flux_surf.or.max(abs(wtsurf),abs(wqsurf))>0.) THEN
    466467       fsens=-wtsurf*rcpd*rho(1)
    467468       flat=-wqsurf*rlvtt*rho(1)
    468469       print *,'Flux: ok_flux wtsurf wqsurf',ok_flux_surf,wtsurf,wqsurf
    469470      ENDIF
     471      print*,'Flux sol ',fsens,flat
    470472!!      ok_flux_surf=.false.
    471473!!      fsens=-wtsurf*rcpd*rho(1)
     
    851853!
    852854        du_age(1:mxcalc)= -2.*sfdt/timestep*
    853      :          (sfdt*(u(1:mxcalc)-ug(1:mxcalc)) -
    854      :           cfdt*(v(1:mxcalc)-vg(1:mxcalc))  )
     855     s          (sfdt*(u(1:mxcalc)-ug(1:mxcalc)) -
     856     s           cfdt*(v(1:mxcalc)-vg(1:mxcalc))  )
    855857!!     : fcoriolis*(v(1:mxcalc)-vg(1:mxcalc))
    856858!
    857859       dv_age(1:mxcalc)= -2.*sfdt/timestep*
    858      :          (cfdt*(u(1:mxcalc)-ug(1:mxcalc)) +
    859      :           sfdt*(v(1:mxcalc)-vg(1:mxcalc))  )
     860     s          (cfdt*(u(1:mxcalc)-ug(1:mxcalc)) +
     861     s           sfdt*(v(1:mxcalc)-vg(1:mxcalc))  )
    860862!!     : -fcoriolis*(u(1:mxcalc)-ug(1:mxcalc))
    861863!
     
    870872!! Increment state variables
    871873!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     874! pour les cas sandu et astex, on reclacule u,v,q,temp et teta dans 1D_nudge_sandu_astex.h
     875! au dessus de 700hpa, on relaxe vers les profils initiaux
     876      if (forcing_sandu .OR. forcing_astex) then
     877#include "1D_nudge_sandu_astex.h"
     878      else
    872879        u(1:mxcalc)=u(1:mxcalc) + timestep*(
    873      :              du_phys(1:mxcalc)
    874      :             +du_age(1:mxcalc) )           
     880     s              du_phys(1:mxcalc)
     881     s             +du_age(1:mxcalc) )           
    875882        v(1:mxcalc)=v(1:mxcalc) + timestep*(
    876                   dv_phys(1:mxcalc)
    877      :             +dv_age(1:mxcalc) )
     883     s              dv_phys(1:mxcalc)
     884     s             +dv_age(1:mxcalc) )
    878885        q(1:mxcalc,:)=q(1:mxcalc,:)+timestep*(
    879      :                dq(1:mxcalc,:)
    880      :               +d_q_adv(1:mxcalc,:) )
     886     s                dq(1:mxcalc,:)
     887     s               +d_q_adv(1:mxcalc,:) )
    881888
    882889        if (prt_level.ge.1) then
     
    893900     .             +d_th_adv(1:mxcalc)
    894901     .             +dt_cooling(1:mxcalc))  ! Taux de chauffage ou refroid.
     902
     903      endif  ! forcing_sandu or forcing_astex
    895904
    896905        teta=temp*(pzero/play)**rkappa
Note: See TracChangeset for help on using the changeset viewer.