Changeset 1056 for trunk/LMDZ.TITAN/libf/phytitan/calchim.F

Timestamp:

Oct 7, 2013, 6:42:03 PM (11 years ago)

Author:

slebonnois

Message:

SL: Titan runs ! see DOC/chantiers/commit_importants.log

File:

• trunk/LMDZ.TITAN/libf/phytitan/calchim.F (modified) (11 diffs)

trunk/LMDZ.TITAN/libf/phytitan/calchim.F

@@ -1 @@
-      SUBROUTINE calchim(ny,qy_c,nomqy_c,declin_rad,ls_rad,dtchim,
+      SUBROUTINE calchim(nlon,ny,qy_c,nomqy_c,declin_rad,ls_rad,dtchim,
      .                   ctemp,cplay,cplev,
      .                   dqyc)
@@ -9 +9 @@
 c     Auteur: S. Lebonnois,  01/2000 | 09/2003 
 c            adaptation pour Titan 3D: 02/2009
+c            adaptation pour // : 04/2013
 c
 c-------------------------------------------------
 c
       use dimphy
+      use common_mod, only:utilaer,maer,prodaer,csn,csh,psurfhaze,
+     .                     NLEV,NC,ND,NR
+      USE comgeomphy,  only: rlatd
+      use moyzon_mod,  only: klat
       implicit none
 #include "dimensions.h"
@@ -19 +24 @@
 #include "YOMCST.h"
 
-#include "titan_for.h"  
-!!!  doit etre en accord avec titan.h
-#include "aerprod.h"
-
 c    Arguments
 c    ---------
 
+      INTEGER      nlon                   ! nb of horiz points
       INTEGER      ny                     ! nb de composes (nqmax-nmicro) 
-      REAL         qy_c(jjm+1,klev,NC)    ! Especes chimiques apres adv.+diss.
+      REAL         qy_c(nlon,klev,NC)     ! Especes chimiques apres adv.+diss.
       character*10 nomqy_c(NC+1)          ! Noms des especes chimiques
       REAL         declin_rad,ls_rad      ! declinaison et long solaire en radians
       REAL         dtchim                 ! pas de temps chimie
-      REAL         ctemp(jjm+1,klev)      ! Temperature
-      REAL         cplay(jjm+1,klev)      ! pression (Pa)
-      REAL         cplev(jjm+1,klev)      ! pression intercouches (Pa)
-      
-      REAL         dqyc(jjm+1,klev,NC)    ! Tendances especes chimiques 
+      REAL         ctemp(nlon,klev)      ! Temperature
+      REAL         cplay(nlon,klev)      ! pression (Pa)
+      REAL         cplev(nlon,klev)      ! pression intercouches (Pa)
+      
+      REAL         dqyc(nlon,klev,NC)    ! Tendances especes chimiques 
       
 c    Local variables :
 c    -----------------
+
+      integer i,j,l,ic,jm1
+
 c variables envoyees dans la chimie: double precision
 
-      integer i,j,l,ic
       REAL  temp_c(klev),press_c(klev)     ! T,p(mbar) a 1 lat donnee
       REAL  declin_c                       ! declinaison en degres
@@ -64 +68 @@
       
       REAL  mass(NC),duree
-      REAL  tablefluxtop(NC,jjm+1,5)
+      REAL  tablefluxtop(NC,jjp1,5)
       REAL  botCH4
       DATA  botCH4/0.05/
@@ -90 +94 @@
 c ************************************
 
-        allocate(krpd(15,ND+1,klev,jjm+1),krate(klev,NR))
+        allocate(krpd(15,ND+1,klev,jjp1),krate(klev,NR))
 
 c Verification dimension verticale: coherence titan_for.h et klev
@@ -110 +114 @@
         endif
 
-c calcul de temp_c, densites et press_c a l'equateur:
-c --------------------------------------------------
-
-        print*,'pression, densites et temp a l equateur (chimie):'
+c calcul de temp_c, densites et press_c au milieu de l'ensemble des points:
+c ----------------------------------------------------------------------
+
+        print*,'pression, densites et temp (chimie):'
         print*,'level, press_c, nb, temp_c'
         DO l=1,klev
 c     temp_c (K):
-         temp_c(l)  = ctemp(jjm/2+1,l)
+         temp_c(l)  = ctemp(nlon/2+1,l)
 c     press_c (mbar):
-         press_c(l) = cplay(jjm/2+1,l)/100.
+         press_c(l) = cplay(nlon/2+1,l)/100.
 c     nb (cm-3):
          nb(l) = 1.e-4*press_c(l) / (RKBOL*temp_c(l))
@@ -393 +397 @@
 c                BOUCLE SUR LES LATITUDES
 c
-      DO j=1,jjp1 
-      
+      DO j=1,nlon 
+      
+      if (j.eq.1) then 
+         jm1=1
+      else
+         jm1=j-1
+      endif
+
+      if((j.eq.1).or.(klat(j).ne.klat(jm1))) then
+
 c***********************************************************************
 c***********************************************************************
@@ -459 +471 @@
 c!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
 
-       if (firstcal.and.(j.eq.1)) then
-         print*,'Alt, densites et temp au pole (chimie):'
-         print*,'level, z_bas, nb, temp_c'
-         do l=1,klev
-          print*,l,rinter(l)-RA/1000.,nb(l),temp_c(l)
-         enddo
-       endif
-
-       if (firstcal.and.(j.eq.jjm/2)) then
-c         print*,'g,mugaz'
-c         print*,g,mugaz
-         print*,'Alt, densites et temp a l equateur (chimie):'
-         print*,'level, z_bas, nb, temp_c'
-         do l=1,klev
-          print*,l,rinter(l)-RA/1000.,nb(l),temp_c(l)
-         enddo
-       endif
-       
 c-----------------------------------------------------------------------
 c
@@ -507 +501 @@
 c   --------------------
        
-       call gptitan(jjp1,rinter,temp_c,nb,
+       call gptitan(rinter,temp_c,nb,
      $              nomqy_c,cqy,fluxtop,
-     $              declin_c,duree,(j-1),mass,
+     $              declin_c,duree,(klat(j)-1),mass,
      $              botCH4,krpd,krate,reactif,
      $              nom_prod,nom_perte,prod,perte,
@@ -515 +509 @@
      $              htoh2,surfhaze)
        
-c       if ( j.eq.jjm/2 ) 
-c    $        print*,cqy(1,1),cqy(klev,1),cqy(1,2),cqy(klev,2)
-c       if ( j.eq.jjm/2 ) 
-c    $  print*,qy_c(j,1,1),qy_c(j,klev,1),qy_c(j,1,2),qy_c(j,klev,2)
-
-c       stop
-
 c   Tendances composition 
 c   ---------------------
@@ -551 +538 @@
 c***********************************************************************
 c***********************************************************************
-c
+
 c              FIN: BOUCLE SUR LES LATITUDES 
-c
+
+      else      ! same latitude, we don't do calculations again
+        dqyc(j,:,:) = dqyc(jm1,:,:)
+        if (aerprod.eq.1) then
+          prodaer(j,:,:) = prodaer(jm1,:,:) 
+          maer(j,:,:)    = maer(jm1,:,:)
+          csn(j,:,:)     = csn(jm1,:,:)
+          csh(j,:,:)     = csh(jm1,:,:)
+        endif
+      endif
+
       ENDDO