Ignore:
Timestamp:
Mar 19, 2012, 11:27:18 AM (13 years ago)
Author:
emillour
Message:

Generic GCM:
Some cleanup and bug fixing:

  • "cloudfrac" was not well written to restartfi (wrong size).
  • missing save attribute for "reffrad" in physiq.F90.
  • cleanup recomputation of surface pressure in newstart and change loop order in interp_horiz (which "fixes" an odd behaviour which fills some arrays with zeros, but only when using some versions of ifort!)

EM

Location:
trunk/LMDZ.GENERIC/libf/dyn3d
Files:
3 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.GENERIC/libf/dyn3d/interp_horiz.F

    r135 r588  
    1919c  """""""""
    2020       
    21        INTEGER imo, jmo ! dimensions ancienne grille (input)
    22        INTEGER imn,jmn  ! dimensions nouvelle grille (input)
     21       INTEGER,INTENT(IN) :: imo, jmo ! dimensions ancienne grille (input)
     22       INTEGER,INTENT(IN) :: imn,jmn  ! dimensions nouvelle grille (input)
    2323
    24        REAL rlonuo(imo+1)     !  Latitude et
    25        REAL rlatvo(jmo)       !  longitude des
    26        REAL rlonun(imn+1)     !  bord des
    27        REAL rlatvn(jmn)     !  cases "scalaires" (input)
     24       REAL,INTENT(IN) :: rlonuo(imo+1)     !  Latitude et
     25       REAL,INTENT(IN) :: rlatvo(jmo)       !  longitude des
     26       REAL,INTENT(IN) :: rlonun(imn+1)     !  bord des
     27       REAL,INTENT(IN) :: rlatvn(jmn)     !  cases "scalaires" (input)
    2828
    29        INTEGER lm ! dimension verticale (input)
    30        REAL varo (imo+1, jmo+1,lm) ! var dans l'ancienne grille (input)
    31        REAL varn (imn+1,jmn+1,lm) ! var dans la nouvelle grille (output)
     29       INTEGER,INTENT(IN) :: lm ! dimension verticale (input)
     30       REAL,INTENT(IN) :: varo (imo+1, jmo+1,lm) ! var dans l'ancienne grille (input)
     31       REAL,INTENT(OUT) :: varn (imn+1,jmn+1,lm) ! var dans la nouvelle grille (output)
    3232
    3333c Autres variables
    … …  
    3939       INTEGER ii,jj,l
    4040
    41        REAL airen ((imnmx2+1)*(jmnmx2+1)) ! aire dans la nouvelle grille
     41       REAL,SAVE :: airen ((imnmx2+1)*(jmnmx2+1)) ! aire dans la nouvelle grille
    4242       REAL airentotn   ! aire totale pole nord dans la nouvelle grille
    4343       REAL airentots   ! aire totale pole sud dans la nouvelle grille
    … …  
    4848c + des pouiemes (cas ou une maille est a cheval sur 2 ou 4 mailles)
    4949c  Il y a un test dans iniinterp_h pour s'assurer que ktotal < kmax
    50        INTEGER kmax, k, ktotal
     50       INTEGER kmax, k
     51       integer,save :: ktotal
    5152       parameter (kmax = 360*179 + 200000)
    5253c      parameter (kmax = 360*179 + 40000)
    5354
    54        INTEGER iik(kmax), jjk(kmax),jk(kmax),ik(kmax)
    55        REAL intersec(kmax)
    56        REAL R
     55       INTEGER,SAVE :: iik(kmax), jjk(kmax),jk(kmax),ik(kmax)
     56       REAL,SAVE :: intersec(kmax)
     57       REAL r
    5758       REAL totn, tots
    58        integer prev_sumdim
    59        save prev_sumdim
    60        data prev_sumdim /0/
    61        
     59       integer,save :: prev_sumdim=0
    6260
    63        logical firsttest, aire_ok
    64        save firsttest
    65        data firsttest /.true./
    66        data aire_ok /.true./
     61       logical,save :: firsttest=.true. , aire_ok=.true.
    6762
    68        integer imoS,jmoS,imnS,jmnS
    69        save imoS,jmoS,imnS,jmnS
    70        save ktotal,iik,jjk,jk,ik,intersec,airen
    71        REAL pi
     63       integer,save :: imoS,jmoS,imnS,jmnS
    7264
    7365c Test dimensions imnmx2 jmnmx2
    … …  
    115107      end if
    116108
    117       do l=1,lm
    118        do jj =1 , jmn+1
    119         do ii=1, imn+1
    120           varn(ii,jj,l) =0.
    121         end do
    122        end do
    123       end do
     109! initialize varn() to zero
     110      varn(1:imn+1,1:jmn+1,1:lm)=0.
    124111       
    125112c Interpolation horizontale
    … …  
    128115c de l'ancienne et la  nouvelle grille
    129116c
    130      
    131       do k=1,ktotal
    132         do l=1,lm
     117! Ehouarn 2012: for some strange reason, with ifort v12.x,
     118!               when the order of the loop below is changed
     119!               values of varn(:,:,l=2...) are then sometimes remain zero!   
     120      do l=1,lm
     121        do k=1,ktotal
    133122         varn(iik(k),jjk(k),l) = varn(iik(k),jjk(k),l)
    134123     &   + varo(ik(k), jk(k),l)*intersec(k)/airen(iik(k)
    … …  
    161150         end do
    162151
    163       ENDDO
     152      ENDDO ! of DO l=1, lm
    164153           
    165154
    … …  
    167156c  TEST  TEST  TEST  TEST  TEST  TEST  TEST  TEST  TEST  TEST
    168157      if (firsttest) then
    169       pi=2.*asin(1.)
    170158      firsttest = .false.
    171159      write (*,*) 'INTERP. HORIZ. : TEST SUR LES AIRES:'
    … …  
    192180      end do
    193181      if (aire_ok) write(*,*) 'INTERP. HORIZ. : AIRES OK'
    194       endif
     182      endif ! of if (firsttest)
    195183
    196184c FIN TEST  FIN TEST  FIN TEST  FIN TEST  FIN TEST  FIN TEST  FIN TEST
    … …  
    198186
    199187
    200         return
    201188        end

  • trunk/LMDZ.GENERIC/libf/dyn3d/lect_start_archive.F

    r253 r588  
    4646      INTEGER   imold,jmold,lmold,nsoilold,nqold
    4747
    48 c et autres:
    49 c----------
    50       INTEGER lnblnk
    51       EXTERNAL lnblnk
    5248
    5349c Variables pour les lectures des fichiers "ini"
    … …  
    109105c     REAL year_day,periheli,aphelie,peri_day
    110106c     REAL obliquit,z0,emin_turb,lmixmin
    111 c     REAL emissiv,emisice(2),albedice(2),tauvis
     107c     REAL emissiv,emisice(2),albedice(2)
    112108c     REAL iceradius(2) , dtemisice(2)
    113109
    … …  
    994990      co2icetotal = 0.
    995991      if (igcm_co2_ice.ne.0) then
    996       DO j=1,jjp1
     992        ! recast surface CO2 ice on new grid
     993        call interp_horiz(qsurfold(1,1,igcm_co2_ice),
     994     &                  qsurfs(1,1,igcm_co2_ice),
     995     &                  imold,jmold,iim,jjm,1,
     996     &                  rlonuold,rlatvold,rlonu,rlatv)
     997       DO j=1,jjp1
    997998         DO i=1,iim
    998999           !co2icetotal = co2icetotal + co2iceS(i,j)*aire(i,j)
    9991000           co2icetotal=co2icetotal+qsurfS(i,j,igcm_co2_ice)*aire(i,j)
    10001001         END DO
    1001       END DO
     1002       END DO
     1003      else
     1004        write(*,*) "Warning: No co2_ice surface tracer"
    10021005      endif
    10031006
    … …  
    10091012      write(*,*)'Ancienne grille: masse de la glace CO2:',co2icetotalold
    10101013      write(*,*)'Nouvelle grille: masse de la glace CO2:',co2icetotal
     1014      if (co2icetotalold.ne.0.) then
    10111015      write(*,*)'Ratio new ice./old ice =',co2icetotal/co2icetotalold
     1016      endif
    10121017      write(*,*)
    10131018

  • trunk/LMDZ.GENERIC/libf/dyn3d/newstart.F

    r535 r588  
    307307      pa=tab_cntrl(17)  ! reference pressure at which coord is purely pressure
    308308      !NB: in start_archive files tab_cntrl(17)=tab_cntrl(18)=0
    309       if (preff.eq.0) then
    310         preff=610
    311         pa=20
    312       endif
    313309      write(*,*) "Newstart: preff=",preff," pa=",pa
    314310      yes=' '
    … …  
    16711667            cloudfrac(ig,l)=0.5
    16721668         ENDDO
     1669! Ehouarn, march 2012: also add some initialisation for hice
     1670         hice(ig)=0.0
    16731671      ENDDO
    16741672     
    … …  
    16821680!      ENDDO
    16831681
    1684 
    16851682c=======================================================================
    16861683c   Correct pressure on the new grid (menu 0)
    16871684c=======================================================================
     1685
    16881686
    16891687      if ((choix_1.eq.0).and.(.not.flagps0)) then
    16901688        r = 1000.*8.31/mugaz
    16911689
    1692         phi0=0.0
    16931690        do j=1,jjp1
    16941691          do i=1,iip1
    1695              phi0 = phi0+phis(i,j)*aire(i,j)
    1696           end do
    1697         end do
    1698         phi0=phi0/airetot
    1699 
    1700         do j=1,jjp1
    1701           do i=1,iip1
    1702              ps(i,j) = ps(i,j) *
    1703 !     .            exp((phisold_newgrid(i,j)-phis(i,j)) /
    1704      .            exp((phi0-phis(i,j)) /
     1692             ps(i,j) = ps(i,j) *
     1693     .            exp((phisold_newgrid(i,j)-phis(i,j)) /
    17051694     .                                  (t(i,j,1) * r))
    17061695          end do
    17071696        end do
    1708  
    1709 !     we must renormalise pressures AGAIN, because exp(-phi) is nonlinear
    1710         ptot=0.0
     1697
     1698c periodicite de ps en longitude
    17111699        do j=1,jjp1
    1712            do i=1,iip1
    1713               ptot=ptot+ps(i,j)*aire(i,j)
    1714            enddo
    1715         enddo
    1716         do j=1,jjp1
    1717            do i=1,iip1
    1718               ps(i,j)=ps(i,j)*ptotn/ptot
    1719            enddo
    1720         enddo
    1721        
    1722 !     periodicity of surface ps in longitude
    1723         do j=1,jjp1
    1724            ps(1,j) = ps(iip1,j)
     1700          ps(1,j) = ps(iip1,j)
    17251701        end do
    1726        
    17271702      end if
    17281703
     1704         
    17291705c=======================================================================
    17301706c=======================================================================
    … …  
    17331709c    Initialisation de la physique / ecriture de newstartfi :
    17341710c=======================================================================
     1711
    17351712
    17361713      CALL inifilr
    … …  
    17641741      endif
    17651742
    1766 
    17671743C Calcul intermediaire
    17681744c
    … …  
    17981774      CALL caldyn0( itau,ucov,vcov,teta,ps,masse,pk,phis ,
    17991775     *                phi,w, pbaru,pbarv,day_ini+time )
    1800 c     CALL caldyn
    1801 c    $  ( itau,ucov,vcov,teta,ps,masse,pk,pkf,phis ,
    1802 c    $    phi,conser,du,dv,dteta,dp,w, pbaru,pbarv, day_ini )
    1803 
     1776
     1777         
    18041778      CALL dynredem0("restart.nc",day_ini,anneeref,phis,nqmx)
    18051779      CALL dynredem1("restart.nc",0.0,vcov,ucov,teta,q,nqmx,masse,ps)
Note: See TracChangeset for help on using the changeset viewer.