Ignore:
Timestamp:
Apr 9, 2009, 12:11:35 PM (15 years ago)
Author:
Laurent Fairhead
Message:

Réintegration dans le tronc des modifications issues de la branche LMDZ-dev
comprises entre la révision 1074 et 1145
Validation: une simulation de 1 jour en séquentiel sur PC donne les mêmes
résultats entre la trunk et la dev
LF

Location:
LMDZ4/trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • LMDZ4/trunk

  • LMDZ4/trunk/libf/dyn3d/iniacademic.F

    r524 r1146  
    44c
    55c
    6       SUBROUTINE iniacademic(nq,vcov,ucov,teta,q,masse,ps,phis,time_0)
     6      SUBROUTINE iniacademic(vcov,ucov,teta,q,masse,ps,phis,time_0)
     7
     8      USE filtreg_mod
     9      USE infotrac, ONLY : nqtot
    710
    811c%W%    %G%
    … …  
    4245#include "temps.h"
    4346#include "control.h"
     47#include "iniprint.h"
    4448
    4549c   Arguments:
    4650c   ----------
    4751
    48       integer nq
    4952      real time_0
    5053
    … …  
    5255      REAL vcov(ip1jm,llm),ucov(ip1jmp1,llm) ! vents covariants
    5356      REAL teta(ip1jmp1,llm)                 ! temperature potentielle
    54       REAL q(ip1jmp1,llm,nq)               ! champs advectes
     57      REAL q(ip1jmp1,llm,nqtot)               ! champs advectes
    5558      REAL ps(ip1jmp1)                       ! pression  au sol
    5659      REAL masse(ip1jmp1,llm)                ! masse d'air
     60      REAL phis(ip1jmp1)                     ! geopotentiel au sol
     61
     62c   Local:
     63c   ------
     64
    5765      REAL p (ip1jmp1,llmp1  )               ! pression aux interfac.des couches
    5866      REAL pks(ip1jmp1)                      ! exner au  sol
    5967      REAL pk(ip1jmp1,llm)                   ! exner au milieu des couches
    6068      REAL pkf(ip1jmp1,llm)                  ! exner filt.au milieu des couches
    61       REAL phis(ip1jmp1)                     ! geopotentiel au sol
    6269      REAL phi(ip1jmp1,llm)                  ! geopotentiel
    63 
    64 
    65 
    66 
    67 
    68 c   Local:
    69 c   ------
    70 
    7170      REAL ddsin,tetarappelj,tetarappell,zsig
    7271      real tetajl(jjp1,llm)
    … …  
    7978
    8079c-----------------------------------------------------------------------
     80! 1. Initializations for Earth-like case
     81! --------------------------------------
     82      if (planet_type=="earth") then
     83c
     84        time_0=0.
    8185
    82 c
    83       time_0=0.
     86        im         = iim
     87        jm         = jjm
     88        day_ini    = 0
     89        omeg       = 4.*asin(1.)/86400.
     90        rad    = 6371229.
     91        g      = 9.8
     92        daysec = 86400.
     93        dtvr    = daysec/FLOAT(day_step)
     94        zdtvr=dtvr
     95        kappa  = 0.2857143
     96        cpp    = 1004.70885
     97        preff     = 101325.
     98        pa        =  50000.
     99        etot0      = 0.
     100        ptot0      = 0.
     101        ztot0      = 0.
     102        stot0      = 0.
     103        ang0       = 0.
    84104
    85       im         = iim
    86       jm         = jjm
    87       day_ini    = 0
    88       omeg       = 4.*asin(1.)/86400.
    89       rad    = 6371229.
    90       g      = 9.8
    91       daysec = 86400.
    92       dtvr    = daysec/FLOAT(day_step)
    93       zdtvr=dtvr
    94       kappa  = 0.2857143
    95       cpp    = 1004.70885
    96       preff     = 101325.
    97       pa        =  50 000.
    98       etot0      = 0.
    99       ptot0      = 0.
    100       ztot0      = 0.
    101       stot0      = 0.
    102       ang0       = 0.
    103       pa         = 0.
     105        CALL iniconst
     106        CALL inigeom
     107        CALL inifilr
    104108
    105       CALL inicons0
    106       CALL inigeom
    107       CALL inifilr
    108 
    109       ps=0.
    110       phis=0.
     109        ps=0.
     110        phis=0.
    111111c---------------------------------------------------------------------
    112112
    113       taurappel=10.*daysec
     113        taurappel=10.*daysec
    114114
    115115c---------------------------------------------------------------------
    … …  
    117117c   --------------------------------------
    118118
    119       DO l=1,llm
    120        zsig=ap(l)/preff+bp(l)
    121        if (zsig.gt.0.3) then
    122          lsup=l
    123          tetarappell=1./8.*(-log(zsig)-.5)
    124          DO j=1,jjp1
    125             ddsin=sin(rlatu(j))-sin(pi/20.)
    126             tetajl(j,l)=300.*(1+1./18.*(1.-3.*ddsin*ddsin)+tetarappell)
    127          ENDDO
    128         else
     119        DO l=1,llm
     120         zsig=ap(l)/preff+bp(l)
     121         if (zsig.gt.0.3) then
     122           lsup=l
     123           tetarappell=1./8.*(-log(zsig)-.5)
     124           DO j=1,jjp1
     125             ddsin=sin(rlatu(j))-sin(pi/20.)
     126             tetajl(j,l)=300.*(1+1./18.*(1.-3.*ddsin*ddsin)+tetarappell)
     127           ENDDO
     128          else
    129129c   Choix isotherme au-dessus de 300 mbar
    130          do j=1,jjp1
    131             tetajl(j,l)=tetajl(j,lsup)*(0.3/zsig)**kappa
    132          enddo
    133         endif
    134       ENDDO
     130           do j=1,jjp1
     131             tetajl(j,l)=tetajl(j,lsup)*(0.3/zsig)**kappa
     132           enddo
     133          endif ! of if (zsig.gt.0.3)
     134        ENDDO ! of DO l=1,llm
    135135
    136       do l=1,llm
    137          do j=1,jjp1
    138             do i=1,iip1
    139                ij=(j-1)*iip1+i
    140                tetarappel(ij,l)=tetajl(j,l)
    141             enddo
    142          enddo
    143       enddo
     136        do l=1,llm
     137           do j=1,jjp1
     138              do i=1,iip1
     139                 ij=(j-1)*iip1+i
     140                 tetarappel(ij,l)=tetajl(j,l)
     141              enddo
     142           enddo
     143        enddo
    144144
    145 c     call dump2d(jjp1,llm,tetajl,'TEQ   ')
     145c       call dump2d(jjp1,llm,tetajl,'TEQ   ')
    146146
    147       ps=1.e5
    148       phis=0.
    149       CALL pression ( ip1jmp1, ap, bp, ps, p       )
    150       CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
    151       CALL massdair(p,masse)
     147        ps=1.e5
     148        phis=0.
     149        CALL pression ( ip1jmp1, ap, bp, ps, p       )
     150        CALL exner_hyb( ip1jmp1, ps, p,alpha,beta, pks, pk, pkf )
     151        CALL massdair(p,masse)
    152152
    153153c  intialisation du vent et de la temperature
    154       teta(:,:)=tetarappel(:,:)
    155       CALL geopot(ip1jmp1,teta,pk,pks,phis,phi)
    156       call ugeostr(phi,ucov)
    157       vcov=0.
    158       q(:,:,1   )=1.e-10
    159       q(:,:,2   )=1.e-15
    160       q(:,:,3:nq)=0.
     154        teta(:,:)=tetarappel(:,:)
     155        CALL geopot(ip1jmp1,teta,pk,pks,phis,phi)
     156        call ugeostr(phi,ucov)
     157        vcov=0.
     158        q(:,:,1   )=1.e-10
     159        q(:,:,2   )=1.e-15
     160        q(:,:,3:nqtot)=0.
    161161
    162162
    163 c   perturbation al\351atoire sur la temp\351rature
    164       idum  = -1
    165       zz = ran1(idum)
    166       idum  = 0
    167       do l=1,llm
    168          do ij=iip2,ip1jm
    169             teta(ij,l)=teta(ij,l)*(1.+0.005*ran1(idum))
    170          enddo
    171       enddo
     163c   perturbation aleatoire sur la temperature
     164        idum  = -1
     165        zz = ran1(idum)
     166        idum  = 0
     167        do l=1,llm
     168           do ij=iip2,ip1jm
     169              teta(ij,l)=teta(ij,l)*(1.+0.005*ran1(idum))
     170           enddo
     171        enddo
    172172
    173       do l=1,llm
    174          do ij=1,ip1jmp1,iip1
    175             teta(ij+iim,l)=teta(ij,l)
    176          enddo
    177       enddo
     173        do l=1,llm
     174           do ij=1,ip1jmp1,iip1
     175              teta(ij+iim,l)=teta(ij,l)
     176           enddo
     177        enddo
    178178
    179179
    … …  
    185185
    186186c   initialisation d'un traceur sur une colonne
    187       j=jjp1*3/4
    188       i=iip1/2
    189       ij=(j-1)*iip1+i
    190       q(ij,:,3)=1.
    191 
     187        j=jjp1*3/4
     188        i=iip1/2
     189        ij=(j-1)*iip1+i
     190        q(ij,:,3)=1.
     191     
     192      else
     193        write(lunout,*)"iniacademic: planet types other than earth",
     194     &                 " not implemented (yet)."
     195        stop
     196      endif ! of if (planet_type=="earth")
    192197      return
    193198      END
Note: See TracChangeset for help on using the changeset viewer.