Index: trunk/mesoscale/LMDZ.MARS/README
===================================================================
--- trunk/mesoscale/LMDZ.MARS/README	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/README	(revision 54)
@@ -0,0 +1,1 @@
+link in_lmdz_mars_newphys/README
Index: trunk/mesoscale/LMDZ.MARS/compile
===================================================================
--- trunk/mesoscale/LMDZ.MARS/compile	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/compile	(revision 54)
@@ -0,0 +1,1 @@
+link in_lmdz_mars_newphys/compile
Index: trunk/mesoscale/LMDZ.MARS/create_make_gcm
===================================================================
--- trunk/mesoscale/LMDZ.MARS/create_make_gcm	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/create_make_gcm	(revision 54)
@@ -0,0 +1,253 @@
+#!/bin/sh
+#set -xv
+#
+# Ehouarn: added the following feature:
+#          create_make_gcm must be called with list of directories
+#          which contain header (or modules) which should be scanned
+#          by create_make_gcm when checking dependencies.
+#          Example:
+#          create_make_gcm dyn3d grid
+#          means that header/modules will be looked for in directories
+#          ., dyn3d and grid
+#
+machine=`hostname`
+os=`uname`
+gcm=`pwd`
+libf=$gcm/libf
+libo=$gcm/libo
+CRAY=0
+if [ "$machine" = "atlas" -o "$machine" = "etoile" -o "$machine" = "axis" ] ; then
+  CRAY=1
+fi
+XNEC=0
+if [ "$machine" = "rhodes" ] ; then
+  XNEC=1
+fi
+VPP=0
+if [ "$machine" = "nymphea0" ] ; then
+  VPP=1
+fi
+#
+echo "#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%"
+echo "# Definitions de Macros pour Make"
+echo "#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%"
+echo 
+echo "# create_make_gcm arguments: $*"
+echo
+echo "# Repertoires :"
+echo
+echo "GCM     = "$gcm
+if [ "$CRAY" = '0' ] ; then
+   echo "MACHINE = sun4"
+fi
+echo 'LIBF    = $(GCM)/libf'
+if [ "$CRAY" = '0' ] ; then
+#   echo 'LIBO    = $(GCM)/libo/$(MACHINE)'
+   echo 'LIBO    = $(LIBOGCM)/$(MACHINE)'
+else
+   echo 'LIBO    = $(GCM)/libo'
+fi
+#echo 'LOCAL_DIR=$(GCM)'
+#echo $localdir
+echo "LOCAL_DIR=`echo $localdir`"
+echo 'BIBIO    = $(LIBF)/bibio'
+echo 'AERONOMARS    = $(LIBF)/aeronomars'
+echo "FILTRE   = filtre"
+echo "PHYS  = "
+echo 'DYN  = $(DYNTYPE)'
+echo 'LIBPHY = $(LIBO)/libphy$(PHYS).a'
+echo 'DIRMAIN= $(DYN)$(DIM)d'
+echo 'RM=rm'
+echo
+echo "OPLINK = "
+echo
+echo '# Les differentes librairies pour l"edition des liens:'
+echo
+if [ "$XNEC" = '1' ] ; then
+  echo 'dyn3d      = $(LIBO)/libsxdyn3d.a $(LIBO)/libsx$(FILTRE).a'
+  echo 'dyn2d      = $(LIBO)/libsxdyn2d.a'
+  echo 'dyn1d      = $(LIBO)/libsxdyn1d.a'
+  echo 'L_DYN      = -lsxdyn$(DIM)d'
+  echo 'L_FILTRE   = -lsx$(FILTRE)'
+  echo 'L_PHY = -lsxphy$(PHYS) '
+  echo 'L_BIBIO    = -lsxbibio'
+  echo 'L_AERONOMARS    = -lsxaeronomars'
+  echo 'L_ADJNT    ='
+else
+  echo 'dyn3d      = $(LIBO)/libdyn3d.a $(LIBO)/lib$(FILTRE).a'
+  echo 'dyn2d      = $(LIBO)/libdyn2d.a'
+  echo 'dyn1d      = $(LIBO)/libdyn1d.a'
+  echo 'L_DYN      = -l$(DYN)$(DIM)d'
+  echo 'L_FILTRE   = -l$(FILTRE)'
+  echo 'L_PHY = -lphy$(PHYS) '
+  echo 'L_BIBIO    = -lbibio'
+  echo 'L_AERONOMARS    = -laeronomars'
+  echo 'L_ADJNT    ='
+fi
+echo
+echo "#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%"
+echo "# Option de compilation FORTRAN"
+echo "#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%"
+echo
+   echo 'COMPILE = $(F77) $(OPTIM) $(INCLUDE) -c'
+   echo 'COMPILE90 = $(F90) $(OPTIM90) $(INCLUDE) -c'
+   echo 'COMPTRU90 = $(F90) $(OPTIMTRU90) $(INCLUDE) -c'
+   echo "LINK    = $LINK"
+   echo "AR      = $AR"
+echo
+echo
+echo "#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%"
+echo "# Creation des differents executables"
+echo "#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%"
+echo
+echo "# Executables:"
+echo "# ------------"
+echo
+echo "PROG = code"
+echo
+echo 'main : $(DYN) bibio phys aeronomars $(OPTION_DEP) '
+echo '	cd $(LIBO) ; $(RANLIB) lib*.a ; cd $(GCM) ;\'
+echo '	cd $(LOCAL_DIR); \'
+echo '	$(COMPILE90) $(LIBF)/$(DIRMAIN)/$(PROG).F -o $(PROG).o ; \'
+if [ "$CRAY" = '0' ] ; then
+echo '	$(LINK) $(PROG).o -L$(LIBO) $(L_DYN) $(L_ADJNT) $(L_FILTRE) $(L_AERONOMARS) $(L_PHY) $(L_BIBIO) $(L_DYN) $(L_AERONOMARS) $(OPLINK) $(OPTION_LINK) -o $(LOCAL_DIR)/$(PROG).e ; $(RM) $(PROG).o '
+else
+echo '	$(LINK) $(PROG).o -L$(LIBO) $(L_DYN) $(L_ADJNT) $(L_FILTRE) $(L_AERONOMARS) $(L_PHY) $(L_BIBIO) $(L_DYN) $(L_AERONOMARS) $(OPLINK) $(OPTION_LINK) -o $(LOCAL_DIR)/$(PROG).e ; $(RM) $(PROG).o '
+fi
+echo
+echo 'dyn : $(LIBO)/libdyn$(DIM)d.a $(FILTRE)$(DIM)d'
+echo
+echo 'olddyn : $(LIBO)/libolddyn$(DIM)d.a $(FILTRE)$(DIM)d'
+echo
+echo 'phys : $(LIBPHY)'
+echo
+echo 'bibio : $(LIBO)/libbibio.a'
+echo
+echo 'aeronomars : $(LIBO)/libaeronomars.a'
+echo
+echo 'adjnt : $(LIBO)/libadjnt.a'
+echo
+echo '$(FILTRE)3d : $(LIBO)/lib$(FILTRE).a'
+echo
+echo '$(FILTRE)2d :'
+echo
+echo '$(FILTRE)1d :'
+echo
+echo "#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%"
+echo "# Contenu des differentes bibliotheques"
+echo "#%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%"
+echo
+echo
+cd $libf >/dev/null 2>&1
+for diri in ` ls `
+do
+   if [ -d $diri ] ; then
+   if [ "`ls $diri/*.F`" != "" ] || [ "`ls $diri/*.F90`" != "" ]  ; then 
+      cd $diri >/dev/null 2>&1
+      echo
+      listlib=""
+      for i in `ls *.F`
+      do
+         fili=`basename $i .F`
+         test=` (  head $i | grep '      PROGRAM' ) `
+         if [ "$test" = "" ] ; then 
+            listlib=$listlib" "$fili
+         fi
+      done
+      for i in `ls *.F90`
+      do
+         fili=`basename $i .F90`
+         test=` (  head $i | grep '      PROGRAM' ) `
+         if [ "$test" = "" ] ; then
+            listlib=$listlib" "$fili
+         fi 
+      done
+#
+      echo
+      echo
+      echo "#======================================================================="
+      echo "# Contenu de la bibliotheque correspondant au Directory "$diri
+      echo "#======================================================================="
+      echo
+      for fili in $listlib
+      do
+         echo '$(LIBO)/lib'$diri".a : " '$(LIBO)/lib'$diri".a("$fili".o)"
+         echo
+      done
+      echo '.PRECIOUS	: $(LIBO)/lib'$diri'.a'
+      echo
+      echo
+      echo "# Compilation des membres de la bibliotheque lib"$diri".a"
+      echo
+      for fili in $listlib
+      do
+         if [ -f $fili.F90 ] ; then
+           trufile=$fili.F90
+         else
+           trufile=$fili.F
+         fi
+         F90=0 ; egrep -i '^ *use ' $trufile > /dev/null 2>&1 && F90=1
+                 egrep -i '^ *module ' $trufile > /dev/null 2>&1 && F90=1
+                 egrep -i '#include*.inc ' $trufile > /dev/null 2>&1 && F90=1
+         str1='$(LIBO)/lib'$diri'.a('$fili'.o) : $(LIBF)/'$diri/$trufile
+         [ "$fili" = "chem.subs" ] && str1=$str1' $(LIBF)/'$diri/chem.mods.F
+         for stri in ` ( sed -n "/\#include/s/\#include//p" $trufile | sed 's/\"//g' ; egrep -i '^ *use ' $trufile | sed -e 's/,/ /' | awk ' { print $2 } ' ) `
+         do
+            stri=`echo $stri | tr [A-Z] [a-z]`
+            if [ -f $stri ] ; then
+               echo $str1 \\
+               str1='$(LIBF)/'$diri'/'$stri
+            else
+               if [ -f $stri.F ] || [ -f $stri.F90 ] ; then
+                  echo $str1 \\
+                  str1='$(LIBO)/lib'$diri'.a('$stri'.o)'
+               else
+#                  for dirinc in dyn3d grid bibio aeronomars phymars
+                  for dirinc in $*
+                  do
+                     if [ -f ../$dirinc/$stri ] ; then
+                        echo $str1 \\
+                        str1='$(LIBF)/'`cd .. ; ls */$stri | head -1`
+                     fi
+                     if [ -f ../$dirinc/$stri.F90 ] ; then
+                        echo $str1 \\
+                        str1='$(LIBO)/lib'$dirinc'.a('$stri'.o)'
+                     fi
+                  done
+               fi
+            fi
+         done
+         echo $str1
+	 if [ "$F90" -eq '0' ] ; then
+           echo '	cd $(LOCAL_DIR); \'
+	   echo '	$(COMPILE) $(LIBF)/'$diri'/'$trufile' ; \'
+	 else
+           echo '	cd $(LOCAL_DIR); \'
+           if [ -f $fili.F90 ] ; then
+	      echo '	$(COMPTRU90) $(LIBF)/'$diri'/'$trufile' ; \'
+           else
+	      echo '	$(COMPILE90) $(LIBF)/'$diri'/'$trufile' ; \'
+           fi
+           MODU=0; egrep -i '^ *module ' $trufile> /dev/null 2>&1 && MODU=1
+            if [ "$MODU" -eq '1' -a "$CRAY" != '1' ] ; then
+              if [ "$os" = 'UNIX_System_V' ] ; then
+                echo '	cp $(MOD_LOC_DIR)/*.$(MOD_SUFFIX) $(LIBO)/ ; \'
+              else
+                echo '	mv $(MOD_LOC_DIR)/'$fili'.$(MOD_SUFFIX) $(LIBO)/'$fili'.$(MOD_SUFFIX) ; \'
+#                echo '	if [ "$(MOD_LOC_DIR)" ne "$(LIBO)" ] ; then mv $(MOD_LOC_DIR)/'*'.$(MOD_SUFFIX) $(LIBO) ; fi ; \'
+              fi
+            fi
+	 fi
+	 if ( [ "$XNEC" -eq '1' ] || [ "$X6NEC" = '1' ] || [ "$X8BRODIE" = '1' ] ) ; then
+           echo '	sxar r $(LIBO)/libsx'$diri'.a '$fili'.o ; \'
+         fi
+         echo '	$(AR) r $(LIBO)/lib'$diri'.a '$fili'.o ; $(RM) '$fili'.o ; \'
+         echo '	cd $(GCM)'
+	 echo
+      done
+#	 
+      echo
+   cd $libf
+   fi
+   fi
+done
Index: trunk/mesoscale/LMDZ.MARS/deftank_gcm
===================================================================
--- trunk/mesoscale/LMDZ.MARS/deftank_gcm	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/deftank_gcm	(revision 54)
@@ -0,0 +1,1 @@
+link ../../mars/deftank
Index: trunk/mesoscale/LMDZ.MARS/in_lmdz_mars_newphys
===================================================================
--- trunk/mesoscale/LMDZ.MARS/in_lmdz_mars_newphys	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/in_lmdz_mars_newphys	(revision 54)
@@ -0,0 +1,1 @@
+link ../LMD_MM_MARS/SIMU/in_lmdz_mars_newphys
Index: trunk/mesoscale/LMDZ.MARS/libf/aeronomars
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/aeronomars	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/aeronomars	(revision 54)
@@ -0,0 +1,1 @@
+link ../libf_gcm/aeronomars
Index: trunk/mesoscale/LMDZ.MARS/libf/bibio
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/bibio	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/bibio	(revision 54)
@@ -0,0 +1,1 @@
+link ../libf_gcm/bibio
Index: trunk/mesoscale/LMDZ.MARS/libf/dyn3d
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/dyn3d	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/dyn3d	(revision 54)
@@ -0,0 +1,1 @@
+link ../libf_gcm/dyn3d
Index: trunk/mesoscale/LMDZ.MARS/libf/filtrez
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/filtrez	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/filtrez	(revision 54)
@@ -0,0 +1,1 @@
+link ../libf_gcm/filtrez
Index: trunk/mesoscale/LMDZ.MARS/libf/grid
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/grid	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/grid	(revision 54)
@@ -0,0 +1,1 @@
+link ../libf_gcm/grid
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/aerave.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/aerave.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/aerave.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/aerave.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/aerkind.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/aerkind.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/aerkind.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/aerkind.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/aeropacity.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/aeropacity.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/aeropacity.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/aeropacity.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/aeroptproperties.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/aeroptproperties.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/aeroptproperties.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/aeroptproperties.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/albedocaps.F90
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/albedocaps.F90	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/albedocaps.F90	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/albedocaps.F90
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/blackl.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/blackl.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/blackl.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/blackl.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/blendrad.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/blendrad.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/blendrad.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/blendrad.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/calldrag_noro.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/calldrag_noro.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/calldrag_noro.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/calldrag_noro.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/callkeys.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/callkeys.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/callkeys.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/callkeys.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/callradite.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/callradite.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/callradite.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/callradite.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/callsedim.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/callsedim.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/callsedim.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/callsedim.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/co2snow.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/co2snow.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/co2snow.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/co2snow.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/comcstfi.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/comcstfi.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/comcstfi.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/comcstfi.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/comdiurn.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/comdiurn.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/comdiurn.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/comdiurn.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/comg1d.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/comg1d.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/comg1d.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/comg1d.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/comgeomfi.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/comgeomfi.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/comgeomfi.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/comgeomfi.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/comorbit.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/comorbit.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/comorbit.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/comorbit.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/comsaison.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/comsaison.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/comsaison.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/comsaison.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/comsoil.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/comsoil.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/comsoil.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/comsoil.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/convadj.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/convadj.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/convadj.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/convadj.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/cvmgp.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/cvmgp.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/cvmgp.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/cvmgp.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/cvmgt.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/cvmgt.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/cvmgt.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/cvmgt.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/datafile.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/datafile.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/datafile.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/datafile.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/datareadnc.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/datareadnc.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/datareadnc.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/datareadnc.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/def_var.F90
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/def_var.F90	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/def_var.F90	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/def_var.F90
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/dimphys.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/dimphys.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/dimphys.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/dimphys.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/dimradmars.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/dimradmars.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/dimradmars.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/dimradmars.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/drag_noro.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/drag_noro.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/drag_noro.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/drag_noro.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/dustdevil.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/dustdevil.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/dustdevil.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/dustdevil.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/dustlift.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/dustlift.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/dustlift.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/dustlift.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/eofdump.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/eofdump.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/eofdump.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/eofdump.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/eofdump.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/eofdump.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/eofdump.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/eofdump.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/flusv.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/flusv.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/flusv.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/flusv.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/growthrate.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/growthrate.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/growthrate.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/growthrate.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/gwprofil.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/gwprofil.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/gwprofil.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/gwprofil.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/gwstress.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/gwstress.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/gwstress.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/gwstress.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/ineofdump.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/ineofdump.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/ineofdump.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/ineofdump.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/inifis.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/inifis.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/inifis.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/inifis.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/iniorbit.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/iniorbit.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/iniorbit.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/iniorbit.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/inistats.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/inistats.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/inistats.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/inistats.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/initracer.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/initracer.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/initracer.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/initracer.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/iniwrite.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/iniwrite.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/iniwrite.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/iniwrite.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/iniwritesoil.F90
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/iniwritesoil.F90	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/iniwritesoil.F90	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/iniwritesoil.F90
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/interp_line.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/interp_line.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/interp_line.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/interp_line.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/ismax.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/ismax.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/ismax.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/ismax.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/ismin.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/ismin.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/ismin.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/ismin.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/lwb.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/lwb.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/lwb.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/lwb.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/lwdiff.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/lwdiff.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/lwdiff.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/lwdiff.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/lwflux.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/lwflux.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/lwflux.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/lwflux.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/lwi.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/lwi.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/lwi.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/lwi.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/lwmain.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/lwmain.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/lwmain.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/lwmain.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/lwtt.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/lwtt.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/lwtt.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/lwtt.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/lwu.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/lwu.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/lwu.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/lwu.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/lwxb.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/lwxb.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/lwxb.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/lwxb.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/lwxd.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/lwxd.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/lwxd.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/lwxd.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/lwxn.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/lwxn.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/lwxn.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/lwxn.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/mkstat.F90
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/mkstat.F90	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/mkstat.F90	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/mkstat.F90
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/mucorr.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/mucorr.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/mucorr.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/mucorr.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/mufract.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/mufract.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/mufract.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/mufract.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/multipl.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/multipl.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/multipl.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/multipl.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/newcondens.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/newcondens.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/newcondens.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/newcondens.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/newsedim.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/newsedim.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/newsedim.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/newsedim.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/nirco2abs.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/nirco2abs.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/nirco2abs.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/nirco2abs.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/nltecool.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/nltecool.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/nltecool.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/nltecool.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/nltedata.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/nltedata.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/nltedata.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/nltedata.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/nlteparams.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/nlteparams.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/nlteparams.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/nlteparams.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/nlthermeq.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/nlthermeq.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/nlthermeq.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/nlthermeq.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/orbite.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/orbite.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/orbite.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/orbite.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/orodrag.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/orodrag.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/orodrag.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/orodrag.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/orosetup.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/orosetup.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/orosetup.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/orosetup.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/phyetat0.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/phyetat0.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/phyetat0.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/phyetat0.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/physdem1.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/physdem1.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/physdem1.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/physdem1.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/physiq.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/physiq.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/physiq.F	(revision 54)
@@ -0,0 +1,1485 @@
+      SUBROUTINE physiq(ngrid,nlayer,nq,
+     $            firstcall,lastcall,
+     $            pday,ptime,ptimestep,
+     $            pplev,pplay,pphi,
+     $            pu,pv,pt,pq,
+     $            pw,
+     $            pdu,pdv,pdt,pdq,pdpsrf,tracerdyn)
+
+
+      IMPLICIT NONE
+c=======================================================================
+c
+c   subject:
+c   --------
+c
+c   Organisation of the physical parametrisations of the LMD 
+c   martian atmospheric general circulation model.
+c
+c   The GCM can be run without or with tracer transport
+c   depending on the value of Logical "tracer" in file  "callphys.def"
+c   Tracers may be water vapor, ice OR chemical species OR dust particles
+c
+c   SEE comments in initracer.F about numbering of tracer species...
+c
+c   It includes:
+c
+c      1. Initialization:
+c      1.1 First call initializations
+c      1.2 Initialization for every call to physiq
+c      1.2.5 Compute mean mass and cp, R and thermal conduction coeff.
+c      2. Compute radiative transfer tendencies
+c         (longwave and shortwave) for CO2 and aerosols.
+c      3. Gravity wave and subgrid scale topography drag :
+c      4. Vertical diffusion (turbulent mixing):
+c      5. Convective adjustment
+c      6. Condensation and sublimation of carbon dioxide.
+c      7.  TRACERS :
+c       7a. water and water ice
+c       7b. call for photochemistry when tracers are chemical species
+c       7c. other scheme for tracer (dust) transport (lifting, sedimentation)
+c       7d. updates (CO2 pressure variations, surface budget)
+c      8. Contribution to tendencies due to thermosphere
+c      9. Surface and sub-surface temperature calculations
+c     10. Write outputs :
+c           - "startfi", "histfi" (if it's time)
+c           - Saving statistics (if "callstats = .true.")
+c           - Dumping eof (if "calleofdump = .true.")
+c           - Output any needed variables in "diagfi" 
+c     11. Diagnostic: mass conservation of tracers
+c 
+c   author: 
+c   ------- 
+c           Frederic Hourdin	15/10/93
+c           Francois Forget		1994
+c           Christophe Hourdin	02/1997 
+c           Subroutine completly rewritten by F.Forget (01/2000)
+c           Introduction of the photochemical module: S. Lebonnois (11/2002)
+c           Introduction of the thermosphere module: M. Angelats i Coll (2002)
+c           Water ice clouds: Franck Montmessin (update 06/2003)
+c           Radiatively active tracers: J.-B. Madeleine (10/2008-06/2009)
+c             Nb: See callradite.F for more information.
+c           
+c   arguments:
+c   ----------
+c
+c   input:
+c   ------
+c    ecri                  period (in dynamical timestep) to write output
+c    ngrid                 Size of the horizontal grid.
+c                          All internal loops are performed on that grid.
+c    nlayer                Number of vertical layers.
+c    nq                    Number of advected fields
+c    firstcall             True at the first call
+c    lastcall              True at the last call
+c    pday                  Number of days counted from the North. Spring
+c                          equinoxe.
+c    ptime                 Universal time (0<ptime<1): ptime=0.5 at 12:00 UT
+c    ptimestep             timestep (s)
+c    pplay(ngrid,nlayer)   Pressure at the middle of the layers (Pa)
+c    pplev(ngrid,nlayer+1) intermediate pressure levels (pa)
+c    pphi(ngrid,nlayer)    Geopotential at the middle of the layers (m2s-2)
+c    pu(ngrid,nlayer)      u component of the wind (ms-1)
+c    pv(ngrid,nlayer)      v component of the wind (ms-1)
+c    pt(ngrid,nlayer)      Temperature (K)
+c    pq(ngrid,nlayer,nq)   Advected fields
+c    pudyn(ngrid,nlayer)    \ 
+c    pvdyn(ngrid,nlayer)     \ Dynamical temporal derivative for the
+c    ptdyn(ngrid,nlayer)     / corresponding variables
+c    pqdyn(ngrid,nlayer,nq) /
+c    pw(ngrid,?)           vertical velocity
+c
+c   output:
+c   -------
+c
+c    pdu(ngrid,nlayermx)        \
+c    pdv(ngrid,nlayermx)         \  Temporal derivative of the corresponding
+c    pdt(ngrid,nlayermx)         /  variables due to physical processes.
+c    pdq(ngrid,nlayermx,nqmx)   /
+c    pdpsrf(ngrid)             /
+c    tracerdyn                 call tracer in dynamical part of GCM ?
+
+c
+c=======================================================================
+c
+c    0.  Declarations :
+c    ------------------
+
+#include "dimensions.h"
+#include "dimphys.h"
+#include "comgeomfi.h"
+#include "surfdat.h"
+#include "comsoil.h"
+#include "comdiurn.h"
+#include "callkeys.h"
+#include "comcstfi.h"
+#include "planete.h"
+#include "comsaison.h"
+#include "control.h"
+#include "dimradmars.h"
+#include "comg1d.h"
+#include "tracer.h"
+#include "nlteparams.h"
+
+#include "chimiedata.h"
+#include "watercap.h"
+#include "param.h"
+#include "param_v3.h"
+#include "conc.h"
+
+#include "netcdf.inc"
+
+
+
+c Arguments :
+c -----------
+
+c   inputs:
+c   -------
+      INTEGER ngrid,nlayer,nq
+      REAL ptimestep
+      REAL pplev(ngridmx,nlayer+1),pplay(ngridmx,nlayer)
+      REAL pphi(ngridmx,nlayer)
+      REAL pu(ngridmx,nlayer),pv(ngridmx,nlayer)
+      REAL pt(ngridmx,nlayer),pq(ngridmx,nlayer,nq)
+      REAL pw(ngridmx,nlayer) !Mars pvervel transmit par dyn3d
+      REAL zh(ngridmx,nlayermx)      ! potential temperature (K)
+      LOGICAL firstcall,lastcall
+
+      REAL pday
+      REAL ptime 
+      logical tracerdyn
+
+c   outputs:
+c   --------
+c     physical tendencies
+      REAL pdu(ngridmx,nlayer),pdv(ngridmx,nlayer)
+      REAL pdt(ngridmx,nlayer),pdq(ngridmx,nlayer,nq)
+      REAL pdpsrf(ngridmx) ! surface pressure tendency
+
+
+c Local saved variables:
+c ----------------------
+c     aerosol (dust or ice) extinction optical depth  at reference wavelength 
+c     "longrefvis" set in dimradmars.h , for one of the "naerkind"  kind of
+c      aerosol optical properties  :
+      REAL aerosol(ngridmx,nlayermx,naerkind)
+
+      INTEGER day_ini  ! Initial date of the run (sol since Ls=0) 
+      INTEGER icount     ! counter of calls to physiq during the run.
+      REAL tsurf(ngridmx)            ! Surface temperature (K)
+      REAL tsoil(ngridmx,nsoilmx)    ! sub-surface temperatures (K)
+      REAL co2ice(ngridmx)           ! co2 ice surface layer (kg.m-2)  
+      REAL albedo(ngridmx,2)         ! Surface albedo in each solar band
+      REAL emis(ngridmx)             ! Thermal IR surface emissivity
+      REAL dtrad(ngridmx,nlayermx)   ! Net atm. radiative heating rate (K.s-1)
+      REAL fluxrad_sky(ngridmx)      ! rad. flux from sky absorbed by surface (W.m-2)
+      REAL fluxrad(ngridmx)          ! Net radiative surface flux (W.m-2)
+      REAL capcal(ngridmx)           ! surface heat capacity (J m-2 K-1)
+      REAL fluxgrd(ngridmx)          ! surface conduction flux (W.m-2)
+      REAL qsurf(ngridmx,nqmx)       ! tracer on surface (e.g. kg.m-2)
+      REAL q2(ngridmx,nlayermx+1)    ! Turbulent Kinetic Energy 
+      INTEGER ig_vl1                 ! Grid Point near VL1   (for diagnostic) 
+
+c     Variables used by the water ice microphysical scheme:
+      REAL rice(ngridmx,nlayermx)    ! Water ice geometric mean radius (m)
+      REAL nuice(ngridmx,nlayermx)   ! Estimated effective variance
+                                     !   of the size distribution
+c     Albedo of deposited surface ice
+      REAL, PARAMETER :: alb_surfice = 0.4 ! 0.45
+
+      SAVE day_ini, icount
+      SAVE aerosol, tsurf,tsoil
+      SAVE co2ice,albedo,emis, q2
+      SAVE capcal,fluxgrd,dtrad,fluxrad,fluxrad_sky,qsurf
+      SAVE ig_vl1
+
+      REAL stephan   
+      DATA stephan/5.67e-08/  ! Stephan Boltzman constant
+      SAVE stephan
+
+c Local variables :
+c -----------------
+
+      REAL CBRT
+      EXTERNAL CBRT
+
+      CHARACTER*80 fichier 
+      INTEGER l,ig,ierr,igout,iq,i, tapphys
+
+      REAL fluxsurf_lw(ngridmx)      !incident LW (IR) surface flux (W.m-2)
+      REAL fluxsurf_sw(ngridmx,2)    !incident SW (solar) surface flux (W.m-2)
+      REAL fluxtop_lw(ngridmx)       !Outgoing LW (IR) flux to space (W.m-2)
+      REAL fluxtop_sw(ngridmx,2)     !Outgoing SW (solar) flux to space (W.m-2)
+      REAL tauref(ngridmx)           ! Reference column optical depth at 700 Pa
+                                     ! (used if active=F) 
+      REAL tau(ngridmx,naerkind)     ! Column dust optical depth at each point
+      REAL zls                       !  solar longitude (rad)
+      REAL zday                      ! date (time since Ls=0, in martian days)
+      REAL zzlay(ngridmx,nlayermx)   ! altitude at the middle of the layers
+      REAL zzlev(ngridmx,nlayermx+1) ! altitude at layer boundaries
+      REAL latvl1,lonvl1             ! Viking Lander 1 point (for diagnostic)
+
+c     Tendancies due to various processes:
+      REAL dqsurf(ngridmx,nqmx)
+      REAL zdtlw(ngridmx,nlayermx)     ! (K/s)
+      REAL zdtsw(ngridmx,nlayermx)     ! (K/s)
+      REAL cldtlw(ngridmx,nlayermx)     ! (K/s) LW heating rate for clear area
+      REAL cldtsw(ngridmx,nlayermx)     ! (K/s) SW heating rate for clear area
+      REAL zdtnirco2(ngridmx,nlayermx) ! (K/s)
+      REAL zdtnlte(ngridmx,nlayermx)   ! (K/s)
+      REAL zdtsurf(ngridmx)            ! (K/s)
+      REAL zdtcloud(ngridmx,nlayermx)
+      REAL zdvdif(ngridmx,nlayermx),zdudif(ngridmx,nlayermx)  ! (m.s-2)
+      REAL zdhdif(ngridmx,nlayermx), zdtsdif(ngridmx)         ! (K/s)
+      REAL zdvadj(ngridmx,nlayermx),zduadj(ngridmx,nlayermx)  ! (m.s-2)
+      REAL zdhadj(ngridmx,nlayermx)                           ! (K/s)
+      REAL zdtgw(ngridmx,nlayermx)                            ! (K/s)
+      REAL zdugw(ngridmx,nlayermx),zdvgw(ngridmx,nlayermx)    ! (m.s-2)
+      REAL zdtc(ngridmx,nlayermx),zdtsurfc(ngridmx)
+      REAL zdvc(ngridmx,nlayermx),zduc(ngridmx,nlayermx)
+
+      REAL zdqdif(ngridmx,nlayermx,nqmx), zdqsdif(ngridmx,nqmx)
+      REAL zdqsed(ngridmx,nlayermx,nqmx), zdqssed(ngridmx,nqmx)
+      REAL zdqdev(ngridmx,nlayermx,nqmx), zdqsdev(ngridmx,nqmx)
+      REAL zdqadj(ngridmx,nlayermx,nqmx)
+      REAL zdqc(ngridmx,nlayermx,nqmx)
+      REAL zdqcloud(ngridmx,nlayermx,nqmx)
+      REAL zdqscloud(ngridmx,nqmx)
+      REAL zdqchim(ngridmx,nlayermx,nqmx)
+      REAL zdqschim(ngridmx,nqmx)
+
+      REAL zdteuv(ngridmx,nlayermx)    ! (K/s)
+      REAL zdtconduc(ngridmx,nlayermx) ! (K/s)
+      REAL zdumolvis(ngridmx,nlayermx)
+      REAL zdvmolvis(ngridmx,nlayermx)
+      real zdqmoldiff(ngridmx,nlayermx,nqmx)
+
+c     Local variable for local intermediate calcul:
+      REAL zflubid(ngridmx)
+      REAL zplanck(ngridmx),zpopsk(ngridmx,nlayermx)
+      REAL zdum1(ngridmx,nlayermx)
+      REAL zdum2(ngridmx,nlayermx)
+      REAL ztim1,ztim2,ztim3, z1,z2
+      REAL ztime_fin
+      REAL zdh(ngridmx,nlayermx)
+      INTEGER length
+      PARAMETER (length=100)
+
+c local variables only used for diagnostic (output in file "diagfi" or "stats")
+c -----------------------------------------------------------------------------
+      REAL ps(ngridmx), zt(ngridmx,nlayermx)
+      REAL zu(ngridmx,nlayermx),zv(ngridmx,nlayermx)
+      REAL zq(ngridmx,nlayermx,nqmx)
+      REAL fluxtop_sw_tot(ngridmx), fluxsurf_sw_tot(ngridmx)
+      character*2 str2
+      character*5 str5
+      real zdtdif(ngridmx,nlayermx), zdtadj(ngridmx,nlayermx)
+      REAL ccn(ngridmx,nlayermx)   ! Cloud condensation nuclei
+                                   !   (particules kg-1)
+      real rdust(ngridmx,nlayermx) ! dust geometric mean radius (m)
+      real qtot1,qtot2 ! total aerosol mass
+      integer igmin, lmin
+      logical tdiag
+
+      real co2col(ngridmx)        ! CO2 column
+      REAL zplev(ngrid,nlayermx+1),zplay(ngrid,nlayermx)
+      REAL zstress(ngrid), cd
+      real hco2(nqmx),tmean, zlocal(nlayermx)
+      real rho(ngridmx,nlayermx)  ! density
+      real vmr(ngridmx,nlayermx)  ! volume mixing ratio
+      REAL mtot(ngridmx)          ! Total mass of water vapor (kg/m2)
+      REAL icetot(ngridmx)        ! Total mass of water ice (kg/m2)
+      REAL rave(ngridmx)          ! Mean water ice effective radius (m)
+      REAL opTES(ngridmx,nlayermx)! abs optical depth at 825 cm-1
+      REAL tauTES(ngridmx)        ! column optical depth at 825 cm-1
+      REAL Qabsice                ! Water ice absorption coefficient
+
+
+      REAL time_phys
+
+c=======================================================================
+
+c 1. Initialisation:
+c -----------------
+
+c  1.1   Initialisation only at first call
+c  ---------------------------------------
+      IF (firstcall) THEN
+
+c        variables set to 0
+c        ~~~~~~~~~~~~~~~~~~
+         call zerophys(ngrid*nlayer*naerkind,aerosol)
+         call zerophys(ngrid*nlayer,dtrad)
+         call zerophys(ngrid,fluxrad)
+
+c        read startfi 
+c        ~~~~~~~~~~~~
+
+! Read netcdf initial physical parameters.
+         CALL phyetat0 ("startfi.nc",0,0,
+     &         nsoilmx,nq,
+     &         day_ini,time_phys,
+     &         tsurf,tsoil,emis,q2,qsurf,co2ice)
+
+         if (pday.ne.day_ini) then
+           write(*,*) "PHYSIQ: ERROR: bad synchronization between ",
+     &                "physics and dynamics"
+           write(*,*) "dynamics day: ",pday
+           write(*,*) "physics day:  ",day_ini
+           stop
+         endif
+
+         write (*,*) 'In physiq day_ini =', day_ini
+
+c        Initialize albedo and orbital calculation
+c        ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+         CALL surfini(ngrid,co2ice,qsurf,albedo)
+         CALL iniorbit(aphelie,periheli,year_day,peri_day,obliquit)
+
+c        initialize soil 
+c        ~~~~~~~~~~~~~~~
+         IF (callsoil) THEN
+            CALL soil(ngrid,nsoilmx,firstcall,inertiedat,
+     s          ptimestep,tsurf,tsoil,capcal,fluxgrd)
+         ELSE
+            PRINT*,
+     &     'PHYSIQ WARNING! Thermal conduction in the soil turned off'
+            DO ig=1,ngrid
+               capcal(ig)=1.e5
+               fluxgrd(ig)=0.
+            ENDDO
+         ENDIF
+         icount=1
+
+
+c        initialize tracers
+c        ~~~~~~~~~~~~~~~~~~
+         tracerdyn=tracer
+         IF (tracer) THEN
+            CALL initracer(qsurf,co2ice)
+         ENDIF  ! end tracer
+
+c        Determining gridpoint near Viking Lander 1 (used for diagnostic only)
+c        ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+         if(ngrid.ne.1) then
+           latvl1= 22.27 
+           lonvl1= -47.94 
+           ig_vl1= 1+ int( (1.5-(latvl1-90.)*jjm/180.)  -2 )*iim +
+     &              int(1.5+(lonvl1+180)*iim/360.)
+           write(*,*) 'Viking Lander 1 GCM point: lat,lon',
+     &              lati(ig_vl1)*180/pi, long(ig_vl1)*180/pi
+         end if 
+
+c        Initialize thermospheric parameters
+c        ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+         if (callthermos) call param_read
+
+c        Initialize R and Cp as constant
+
+         if (.not.callthermos .and. .not.photochem) then
+                 do l=1,nlayermx
+                  do ig=1,ngridmx
+                   rnew(ig,l)=r
+                   cpnew(ig,l)=cpp
+                   mmean(ig,l)=mugaz
+                   enddo
+                  enddo  
+         endif         
+
+        IF (tracer.AND.water.AND.(ngridmx.NE.1)) THEN
+          write(*,*)"physiq: water_param Surface ice alb:",alb_surfice
+        ENDIF
+                   
+      ENDIF        !  (end of "if firstcall")
+
+c ---------------------------------------------------
+c 1.2   Initializations done at every physical timestep:
+c ---------------------------------------------------
+c
+      IF (ngrid.NE.ngridmx) THEN
+         PRINT*,'STOP in PHYSIQ'
+         PRINT*,'Probleme de dimensions :'
+         PRINT*,'ngrid     = ',ngrid
+         PRINT*,'ngridmx   = ',ngridmx
+         STOP
+      ENDIF
+
+c     Initialize various variables
+c     ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+      call zerophys(ngrid*nlayer, pdv)
+      call zerophys(ngrid*nlayer, pdu)
+      call zerophys(ngrid*nlayer, pdt)
+      call zerophys(ngrid*nlayer*nq, pdq)
+      call zerophys(ngrid, pdpsrf)
+      call zerophys(ngrid, zflubid)
+      call zerophys(ngrid, zdtsurf)
+      call zerophys(ngrid*nq, dqsurf)
+      igout=ngrid/2+1 
+
+
+      zday=pday+ptime ! compute time, in sols (and fraction thereof)
+
+c     Compute Solar Longitude (Ls) :
+c     ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+      if (season) then
+         call solarlong(zday,zls)
+      else
+         call solarlong(float(day_ini),zls)
+      end if
+
+c     Compute geopotential at interlayers
+c     ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+c     ponderation des altitudes au niveau des couches en dp/p
+
+      DO l=1,nlayer
+         DO ig=1,ngrid
+            zzlay(ig,l)=pphi(ig,l)/g
+         ENDDO
+      ENDDO
+      DO ig=1,ngrid
+         zzlev(ig,1)=0.
+         zzlev(ig,nlayer+1)=1.e7    ! dummy top of last layer above 10000 km...
+      ENDDO
+      DO l=2,nlayer
+         DO ig=1,ngrid
+            z1=(pplay(ig,l-1)+pplev(ig,l))/(pplay(ig,l-1)-pplev(ig,l))
+            z2=(pplev(ig,l)+pplay(ig,l))/(pplev(ig,l)-pplay(ig,l))
+            zzlev(ig,l)=(z1*zzlay(ig,l-1)+z2*zzlay(ig,l))/(z1+z2)
+         ENDDO
+      ENDDO
+
+
+!     Potential temperature calculation not the same in physiq and dynamic
+
+c     Compute potential temperature
+c     ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+      DO l=1,nlayer
+         DO ig=1,ngrid 
+            zpopsk(ig,l)=(pplay(ig,l)/pplev(ig,1))**rcp
+            zh(ig,l)=pt(ig,l)/zpopsk(ig,l)
+         ENDDO
+      ENDDO
+
+c-----------------------------------------------------------------------
+c    1.2.5 Compute mean mass, cp, and R
+c    --------------------------------
+
+      if(photochem.or.callthermos) then
+         call concentrations(pplay,pt,pdt,pq,pdq,ptimestep)
+      endif
+c-----------------------------------------------------------------------
+c    2. Compute radiative tendencies :
+c------------------------------------
+
+
+      IF (callrad) THEN
+         IF( MOD(icount-1,iradia).EQ.0) THEN
+
+c          Local Solar zenith angle
+c          ~~~~~~~~~~~~~~~~~~~~~~~~
+           CALL orbite(zls,dist_sol,declin)
+
+           IF(diurnal) THEN
+               ztim1=SIN(declin)
+               ztim2=COS(declin)*COS(2.*pi*(zday-.5))
+               ztim3=-COS(declin)*SIN(2.*pi*(zday-.5))
+
+               CALL solang(ngrid,sinlon,coslon,sinlat,coslat,
+     s         ztim1,ztim2,ztim3, mu0,fract)
+
+           ELSE
+               CALL mucorr(ngrid,declin, lati, mu0, fract,10000.,rad)
+           ENDIF
+
+c          NLTE cooling from CO2 emission
+c          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+           IF(callnlte) CALL nltecool(ngrid,nlayer,pplay,pt,zdtnlte)
+
+c          Find number of layers for LTE radiation calculations
+           IF(MOD(iphysiq*(icount-1),day_step).EQ.0)
+     &          CALL nlthermeq(ngrid,nlayer,pplev,pplay)
+
+c          Note: Dustopacity.F has been transferred to callradite.F
+         
+c          Call main radiative transfer scheme
+c          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+c          Transfer through CO2 (except NIR CO2 absorption)
+c            and aerosols (dust and water ice)
+
+c          Radiative transfer
+c          ------------------
+           CALL callradite(icount,ngrid,nlayer,nq,zday,zls,pq,albedo,
+     $     emis,mu0,pplev,pplay,pt,tsurf,fract,dist_sol,igout,
+     $     zdtlw,zdtsw,fluxsurf_lw,fluxsurf_sw,fluxtop_lw,fluxtop_sw,
+     &     tauref,tau,aerosol,ccn,rdust,rice,nuice)
+
+c          CO2 near infrared absorption
+c          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+           call zerophys(ngrid*nlayer,zdtnirco2)
+           if (callnirco2) then
+              call nirco2abs (ngrid,nlayer,pplay,dist_sol,
+     .                       mu0,fract,declin, zdtnirco2)
+           endif
+
+c          Radiative flux from the sky absorbed by the surface (W.m-2)
+c          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+           DO ig=1,ngrid
+               fluxrad_sky(ig)=emis(ig)*fluxsurf_lw(ig)
+     $         +fluxsurf_sw(ig,1)*(1.-albedo(ig,1))
+     $         +fluxsurf_sw(ig,2)*(1.-albedo(ig,2))
+           ENDDO
+
+
+c          Net atmospheric radiative heating rate (K.s-1)
+c          ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+           IF(callnlte) THEN
+              CALL blendrad(ngrid, nlayer, pplay,
+     &             zdtsw, zdtlw, zdtnirco2, zdtnlte, dtrad)
+           ELSE
+              DO l=1,nlayer
+                 DO ig=1,ngrid
+                    dtrad(ig,l)=zdtsw(ig,l)+zdtlw(ig,l)
+     &                          +zdtnirco2(ig,l)
+                  ENDDO
+              ENDDO
+           ENDIF
+
+
+
+        ENDIF ! of if(mod(icount-1,iradia).eq.0)
+
+c    Transformation of the radiative tendencies:
+c    -------------------------------------------
+
+c          Net radiative surface flux (W.m-2)
+c          ~~~~~~~~~~~~~~~~~~~~~~~~~~
+c
+           DO ig=1,ngrid
+               zplanck(ig)=tsurf(ig)*tsurf(ig)
+               zplanck(ig)=emis(ig)*
+     $         stephan*zplanck(ig)*zplanck(ig)
+               fluxrad(ig)=fluxrad_sky(ig)-zplanck(ig)
+           ENDDO
+
+
+         DO l=1,nlayer
+            DO ig=1,ngrid
+               pdt(ig,l)=pdt(ig,l)+dtrad(ig,l)
+            ENDDO
+         ENDDO
+
+      ENDIF ! of IF (callrad)
+
+c-----------------------------------------------------------------------
+c    3. Gravity wave and subgrid scale topography drag :
+c    -------------------------------------------------
+
+
+      IF(calllott)THEN
+
+        CALL calldrag_noro(ngrid,nlayer,ptimestep,
+     &                 pplay,pplev,pt,pu,pv,zdtgw,zdugw,zdvgw)
+ 
+        DO l=1,nlayer
+          DO ig=1,ngrid
+            pdv(ig,l)=pdv(ig,l)+zdvgw(ig,l)
+            pdu(ig,l)=pdu(ig,l)+zdugw(ig,l)
+            pdt(ig,l)=pdt(ig,l)+zdtgw(ig,l)
+          ENDDO
+        ENDDO
+      ENDIF
+
+c-----------------------------------------------------------------------
+c    4. Vertical diffusion (turbulent mixing):
+c    -----------------------------------------
+c
+      IF (calldifv) THEN
+
+
+         DO ig=1,ngrid
+            zflubid(ig)=fluxrad(ig)+fluxgrd(ig)
+         ENDDO
+
+         CALL zerophys(ngrid*nlayer,zdum1)
+         CALL zerophys(ngrid*nlayer,zdum2)
+         do l=1,nlayer
+            do ig=1,ngrid
+               zdh(ig,l)=pdt(ig,l)/zpopsk(ig,l)
+            enddo
+         enddo
+         
+c        Calling vdif (Martian version WITH CO2 condensation)
+         CALL vdifc(ngrid,nlayer,nq,co2ice,zpopsk,
+     $        ptimestep,capcal,lwrite,
+     $        pplay,pplev,zzlay,zzlev,z0,
+     $        pu,pv,zh,pq,tsurf,emis,qsurf,
+     $        zdum1,zdum2,zdh,pdq,zflubid,
+     $        zdudif,zdvdif,zdhdif,zdtsdif,q2,
+     &        zdqdif,zdqsdif)
+
+         DO l=1,nlayer
+            DO ig=1,ngrid
+               pdv(ig,l)=pdv(ig,l)+zdvdif(ig,l)
+               pdu(ig,l)=pdu(ig,l)+zdudif(ig,l)
+               pdt(ig,l)=pdt(ig,l)+zdhdif(ig,l)*zpopsk(ig,l)
+
+               zdtdif(ig,l)=zdhdif(ig,l)*zpopsk(ig,l) ! for diagnostic only
+
+            ENDDO
+         ENDDO
+
+         DO ig=1,ngrid
+            zdtsurf(ig)=zdtsurf(ig)+zdtsdif(ig)
+         ENDDO
+
+         if (tracer) then 
+           DO iq=1, nq
+            DO l=1,nlayer
+              DO ig=1,ngrid
+                 pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqdif(ig,l,iq) 
+              ENDDO
+            ENDDO
+           ENDDO
+           DO iq=1, nq
+              DO ig=1,ngrid
+                 dqsurf(ig,iq)=dqsurf(ig,iq) + zdqsdif(ig,iq)
+              ENDDO
+           ENDDO
+         end if ! of if (tracer)
+
+      ELSE    
+         DO ig=1,ngrid
+            zdtsurf(ig)=zdtsurf(ig)+
+     s      (fluxrad(ig)+fluxgrd(ig))/capcal(ig)
+         ENDDO
+      ENDIF ! of IF (calldifv)
+
+
+c-----------------------------------------------------------------------
+c   5. Dry convective adjustment:
+c   -----------------------------
+
+      IF(calladj) THEN
+
+         DO l=1,nlayer
+            DO ig=1,ngrid
+               zdh(ig,l)=pdt(ig,l)/zpopsk(ig,l)
+            ENDDO
+         ENDDO
+         CALL zerophys(ngrid*nlayer,zduadj)
+         CALL zerophys(ngrid*nlayer,zdvadj)
+         CALL zerophys(ngrid*nlayer,zdhadj)
+
+         CALL convadj(ngrid,nlayer,nq,ptimestep,
+     $                pplay,pplev,zpopsk,
+     $                pu,pv,zh,pq,
+     $                pdu,pdv,zdh,pdq,
+     $                zduadj,zdvadj,zdhadj,
+     $                zdqadj)
+
+         DO l=1,nlayer
+            DO ig=1,ngrid
+               pdu(ig,l)=pdu(ig,l)+zduadj(ig,l)
+               pdv(ig,l)=pdv(ig,l)+zdvadj(ig,l)
+               pdt(ig,l)=pdt(ig,l)+zdhadj(ig,l)*zpopsk(ig,l)
+
+               zdtadj(ig,l)=zdhadj(ig,l)*zpopsk(ig,l) ! for diagnostic only
+            ENDDO
+         ENDDO
+
+         if(tracer) then 
+           DO iq=1, nq
+            DO l=1,nlayer
+              DO ig=1,ngrid
+                 pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqadj(ig,l,iq) 
+              ENDDO
+            ENDDO
+           ENDDO
+         end if
+      ENDIF ! of IF(calladj)
+
+c-----------------------------------------------------------------------
+c   6. Carbon dioxide condensation-sublimation:
+c   -------------------------------------------
+
+      IF (callcond) THEN
+         CALL newcondens(ngrid,nlayer,nq,ptimestep,
+     $              capcal,pplay,pplev,tsurf,pt,
+     $              pphi,pdt,pdu,pdv,zdtsurf,pu,pv,pq,pdq,
+     $              co2ice,albedo,emis,
+     $              zdtc,zdtsurfc,pdpsrf,zduc,zdvc,zdqc,
+     $	            fluxsurf_sw,zls) 
+
+         DO l=1,nlayer
+           DO ig=1,ngrid
+             pdt(ig,l)=pdt(ig,l)+zdtc(ig,l)
+             pdv(ig,l)=pdv(ig,l)+zdvc(ig,l)
+             pdu(ig,l)=pdu(ig,l)+zduc(ig,l)
+           ENDDO
+         ENDDO
+         DO ig=1,ngrid
+            zdtsurf(ig) = zdtsurf(ig) + zdtsurfc(ig)
+         ENDDO
+
+         IF (tracer) THEN
+           DO iq=1, nq
+            DO l=1,nlayer
+              DO ig=1,ngrid
+                pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqc(ig,l,iq) 
+              ENDDO
+            ENDDO
+           ENDDO
+         ENDIF ! of IF (tracer)
+
+      ENDIF  ! of IF (callcond)
+
+c-----------------------------------------------------------------------
+c   7. Specific parameterizations for tracers 
+c:   -----------------------------------------
+
+      if (tracer) then 
+
+c   7a. Water and ice
+c     ---------------
+
+c        ---------------------------------------
+c        Water ice condensation in the atmosphere
+c        ----------------------------------------
+         IF (water) THEN
+
+           call watercloud(ngrid,nlayer,ptimestep,
+     &                pplev,pplay,pdpsrf,zzlev,zzlay, pt,pdt,
+     &                pq,pdq,zdqcloud,zdqscloud,zdtcloud,
+     &                nq,naerkind,tau,
+     &                ccn,rdust,rice,nuice)
+           if (activice) then
+c Temperature variation due to latent heat release
+           DO l=1,nlayer
+             DO ig=1,ngrid
+               pdt(ig,l)=pdt(ig,l)+zdtcloud(ig,l)
+             ENDDO
+           ENDDO
+           endif
+
+! increment water vapour and ice atmospheric tracers tendencies
+           IF (water) THEN
+             DO l=1,nlayer
+               DO ig=1,ngrid
+                 pdq(ig,l,igcm_h2o_vap)=pdq(ig,l,igcm_h2o_vap)+
+     &                                   zdqcloud(ig,l,igcm_h2o_vap)
+                 pdq(ig,l,igcm_h2o_ice)=pdq(ig,l,igcm_h2o_ice)+
+     &                                   zdqcloud(ig,l,igcm_h2o_ice)
+               ENDDO
+             ENDDO
+           ENDIF ! of IF (water) THEN
+! Increment water ice surface tracer tendency
+         DO ig=1,ngrid
+           dqsurf(ig,igcm_h2o_ice)=dqsurf(ig,igcm_h2o_ice)+
+     &                               zdqscloud(ig,igcm_h2o_ice)
+         ENDDO
+         
+         END IF  ! of IF (water)
+
+c   7b. Chemical species
+c     ------------------
+
+c        --------------
+c        photochemistry :
+c        --------------
+         IF (photochem .or. thermochem) then
+          call calchim(ptimestep,pplay,pplev,pt,pdt,dist_sol,mu0,
+     &      zzlay,zday,pq,pdq,rice,
+     &      zdqchim,zdqschim,zdqcloud,zdqscloud)
+!NB: Photochemistry includes condensation of H2O2
+
+           ! increment values of tracers:
+           DO iq=1,nq ! loop on all tracers; tendencies for non-chemistry
+                      ! tracers is zero anyways
+             DO l=1,nlayer
+               DO ig=1,ngrid
+                 pdq(ig,l,iq)=pdq(ig,l,iq)+zdqchim(ig,l,iq)
+               ENDDO
+             ENDDO
+           ENDDO ! of DO iq=1,nq
+           ! add condensation tendency for H2O2
+           if (igcm_h2o2.ne.0) then
+             DO l=1,nlayer
+               DO ig=1,ngrid
+                 pdq(ig,l,igcm_h2o2)=pdq(ig,l,igcm_h2o2)
+     &                                +zdqcloud(ig,l,igcm_h2o2)
+               ENDDO
+             ENDDO
+           endif
+
+           ! increment surface values of tracers:
+           DO iq=1,nq ! loop on all tracers; tendencies for non-chemistry
+                      ! tracers is zero anyways
+             DO ig=1,ngrid
+               dqsurf(ig,iq)=dqsurf(ig,iq)+zdqschim(ig,iq)
+             ENDDO
+           ENDDO ! of DO iq=1,nq
+           ! add condensation tendency for H2O2
+           if (igcm_h2o2.ne.0) then
+             DO ig=1,ngrid
+               dqsurf(ig,igcm_h2o2)=dqsurf(ig,igcm_h2o2)
+     &                                +zdqscloud(ig,igcm_h2o2)
+             ENDDO
+           endif
+
+         END IF  ! of IF (photochem.or.thermochem)
+
+c   7c. Aerosol particles
+c     -------------------
+
+c        ----------
+c        Dust devil :
+c        ----------
+         IF(callddevil) then 
+           call dustdevil(ngrid,nlayer,nq, pplev,pu,pv,pt, tsurf,q2,
+     &                zdqdev,zdqsdev)
+ 
+           if (dustbin.ge.1) then
+              do iq=1,nq
+                 DO l=1,nlayer
+                    DO ig=1,ngrid
+                       pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqdev(ig,l,iq)
+                    ENDDO
+                 ENDDO
+              enddo
+              do iq=1,nq
+                 DO ig=1,ngrid
+                    dqsurf(ig,iq)= dqsurf(ig,iq) + zdqsdev(ig,iq)
+                 ENDDO
+              enddo
+           endif  ! of if (dustbin.ge.1)
+
+         END IF ! of IF (callddevil)
+
+c        ------------- 
+c        Sedimentation :   acts also on water ice
+c        ------------- 
+         IF (sedimentation) THEN 
+           !call zerophys(ngrid*nlayer*nq, zdqsed)
+           zdqsed(1:ngrid,1:nlayer,1:nq)=0
+           !call zerophys(ngrid*nq, zdqssed)
+           zdqssed(1:ngrid,1:nq)=0
+
+           call callsedim(ngrid,nlayer, ptimestep,
+     &                pplev,zzlev, pt, rdust, rice,
+     &                pq, pdq, zdqsed, zdqssed,nq)
+
+           DO iq=1, nq
+             DO l=1,nlayer
+               DO ig=1,ngrid
+                    pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqsed(ig,l,iq)
+               ENDDO
+             ENDDO
+           ENDDO
+           DO iq=1, nq
+             DO ig=1,ngrid
+                dqsurf(ig,iq)= dqsurf(ig,iq) + zdqssed(ig,iq)
+             ENDDO
+           ENDDO
+         END IF   ! of IF (sedimentation)
+
+c   7d. Updates
+c     ---------
+
+        DO iq=1, nq
+          DO ig=1,ngrid
+
+c       ---------------------------------
+c       Updating tracer budget on surface
+c       ---------------------------------
+            qsurf(ig,iq)=qsurf(ig,iq)+ptimestep*dqsurf(ig,iq)
+
+          ENDDO  ! (ig)
+        ENDDO    ! (iq)
+
+      endif !  of if (tracer) 
+
+
+c-----------------------------------------------------------------------
+c   8. THERMOSPHERE CALCULATION
+c-----------------------------------------------------------------------
+
+      if (callthermos) then
+        call thermosphere(pplev,pplay,dist_sol,
+     $     mu0,ptimestep,ptime,zday,tsurf,zzlev,zzlay,
+     &     pt,pq,pu,pv,pdt,pdq,
+     $     zdteuv,zdtconduc,zdumolvis,zdvmolvis,zdqmoldiff)
+
+        DO l=1,nlayer
+          DO ig=1,ngrid
+            dtrad(ig,l)=dtrad(ig,l)+zdteuv(ig,l)
+            pdt(ig,l)=pdt(ig,l)+zdtconduc(ig,l)
+     &                         +zdteuv(ig,l)
+            pdv(ig,l)=pdv(ig,l)+zdvmolvis(ig,l)
+            pdu(ig,l)=pdu(ig,l)+zdumolvis(ig,l)
+            DO iq=1, nq
+              pdq(ig,l,iq)=pdq(ig,l,iq)+zdqmoldiff(ig,l,iq)
+            ENDDO
+          ENDDO
+        ENDDO
+
+      endif ! of if (callthermos)
+
+c-----------------------------------------------------------------------
+c   9. Surface  and sub-surface soil temperature
+c-----------------------------------------------------------------------
+c
+c
+c   9.1 Increment Surface temperature:
+c   ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+      DO ig=1,ngrid
+         tsurf(ig)=tsurf(ig)+ptimestep*zdtsurf(ig) 
+      ENDDO
+
+c  Prescribe a cold trap at south pole (except at high obliquity !!)
+c  Temperature at the surface is set there to be the temperature
+c  corresponding to equilibrium temperature between phases of CO2
+
+      IF (tracer.AND.water.AND.(ngridmx.NE.1)) THEN
+         if (caps.and.(obliquit.lt.27.)) then
+           ! NB: Updated surface pressure, at grid point 'ngrid', is
+           !     ps(ngrid)=pplev(ngrid,1)+pdpsrf(ngrid)*ptimestep
+           tsurf(ngrid)=1./(1./136.27-r/5.9e+5*alog(0.0095*
+     &                     (pplev(ngrid,1)+pdpsrf(ngrid)*ptimestep)))
+         endif
+c       -------------------------------------------------------------
+c       Change of surface albedo (set to 0.4) in case of ground frost
+c       everywhere except on the north permanent cap and in regions
+c       covered by dry ice. 
+c              ALWAYS PLACE these lines after newcondens !!!
+c       -------------------------------------------------------------
+         do ig=1,ngrid
+           if ((co2ice(ig).eq.0).and.
+     &        (qsurf(ig,igcm_h2o_ice).gt.0.005)) then
+              albedo(ig,1) = alb_surfice
+              albedo(ig,2) = alb_surfice
+           endif
+         enddo  ! of do ig=1,ngrid
+      ENDIF  ! of IF (tracer.AND.water.AND.(ngridmx.NE.1))
+
+c
+c   9.2 Compute soil temperatures and subsurface heat flux:
+c   ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+      IF (callsoil) THEN
+         CALL soil(ngrid,nsoilmx,.false.,inertiedat,
+     &          ptimestep,tsurf,tsoil,capcal,fluxgrd)
+      ENDIF
+
+c-----------------------------------------------------------------------
+c  10. Write output files
+c  ----------------------
+
+c    -------------------------------
+c    Dynamical fields incrementation
+c    -------------------------------
+c (FOR OUTPUT ONLY : the actual model integration is performed in the dynamics)
+      ! temperature, zonal and meridional wind
+      DO l=1,nlayer
+        DO ig=1,ngrid
+          zt(ig,l)=pt(ig,l)  + pdt(ig,l)*ptimestep
+          zu(ig,l)=pu(ig,l)  + pdu(ig,l)*ptimestep
+          zv(ig,l)=pv(ig,l)  + pdv(ig,l)*ptimestep
+        ENDDO
+      ENDDO
+
+      ! tracers
+      DO iq=1, nq
+        DO l=1,nlayer
+          DO ig=1,ngrid
+            zq(ig,l,iq)=pq(ig,l,iq) +pdq(ig,l,iq)*ptimestep
+          ENDDO
+        ENDDO
+      ENDDO
+
+      ! surface pressure
+      DO ig=1,ngrid
+          ps(ig)=pplev(ig,1) + pdpsrf(ig)*ptimestep
+      ENDDO
+
+      ! pressure
+      DO l=1,nlayer
+        DO ig=1,ngrid
+             zplev(ig,l)=pplev(ig,l)/pplev(ig,1)*ps(ig)
+             zplay(ig,l)=pplay(ig,l)/pplev(ig,1)*ps(ig)
+        ENDDO
+      ENDDO
+
+      ! Density 
+      DO l=1,nlayer
+         DO ig=1,ngrid
+            rho(ig,l) = zplay(ig,l)/(rnew(ig,l)*zt(ig,l))
+         ENDDO
+      ENDDO
+
+c    Compute surface stress : (NB: z0 is a common in planete.h)
+c     DO ig=1,ngrid
+c        cd = (0.4/log(zzlay(ig,1)/z0))**2
+c        zstress(ig) = rho(ig,1)*cd*(zu(ig,1)**2 + zv(ig,1)**2)
+c     ENDDO
+
+c     Sum of fluxes in solar spectral bands (for output only)
+      DO ig=1,ngrid
+	     fluxtop_sw_tot(ig)=fluxtop_sw(ig,1) + fluxtop_sw(ig,2)
+	     fluxsurf_sw_tot(ig)=fluxsurf_sw(ig,1) + fluxsurf_sw(ig,2)
+      ENDDO
+c ******* TEST ******************************************************
+      ztim1 = 999
+      DO l=1,nlayer
+        DO ig=1,ngrid
+           if (pt(ig,l).lt.ztim1) then
+               ztim1 = pt(ig,l)
+               igmin = ig
+               lmin = l 
+           end if
+        ENDDO
+      ENDDO
+      if(min(pt(igmin,lmin),zt(igmin,lmin)).lt.70.) then	
+        write(*,*) 'PHYSIQ: stability WARNING :'
+        write(*,*) 'pt, zt Tmin = ', pt(igmin,lmin), zt(igmin,lmin),
+     &              'ig l =', igmin, lmin
+      end if
+c *******************************************************************
+
+c     ---------------------
+c     Outputs to the screen 
+c     ---------------------
+
+      IF (lwrite) THEN
+         PRINT*,'Global diagnostics for the physics'
+         PRINT*,'Variables and their increments x and dx/dt * dt'
+         WRITE(*,'(a6,a10,2a15)') 'Ts','dTs','ps','dps'
+         WRITE(*,'(2f10.5,2f15.5)')
+     s   tsurf(igout),zdtsurf(igout)*ptimestep,
+     s   pplev(igout,1),pdpsrf(igout)*ptimestep
+         WRITE(*,'(a4,a6,5a10)') 'l','u','du','v','dv','T','dT'
+         WRITE(*,'(i4,6f10.5)') (l,
+     s   pu(igout,l),pdu(igout,l)*ptimestep,
+     s   pv(igout,l),pdv(igout,l)*ptimestep,
+     s   pt(igout,l),pdt(igout,l)*ptimestep,
+     s   l=1,nlayer)
+      ENDIF ! of IF (lwrite)
+
+      IF (ngrid.NE.1) THEN
+         print*,'Ls =',zls*180./pi,
+     &   ' tauref(700 Pa,lat=0) =',tauref(ngrid/2),
+     &   ' tau(Viking1) =',tau(ig_vl1,1)
+
+
+c        -------------------------------------------------------------------
+c        Writing NetCDF file  "RESTARTFI" at the end of the run
+c        -------------------------------------------------------------------
+c        Note: 'restartfi' is stored just before dynamics are stored
+c              in 'restart'. Between now and the writting of 'restart',
+c              there will have been the itau=itau+1 instruction and
+c              a reset of 'time' (lastacll = .true. when itau+1= itaufin)
+c              thus we store for time=time+dtvr
+
+         IF(lastcall) THEN
+            ztime_fin = ptime + ptimestep/(float(iphysiq)*daysec) 
+            write(*,*)'PHYSIQ: for physdem ztime_fin =',ztime_fin
+            call physdem1("restartfi.nc",long,lati,nsoilmx,nq,
+     .              ptimestep,pday,
+     .              ztime_fin,tsurf,tsoil,co2ice,emis,q2,qsurf,
+     .              area,albedodat,inertiedat,zmea,zstd,zsig,
+     .              zgam,zthe)
+         ENDIF
+
+c        -------------------------------------------------------------------
+c        Calculation of diagnostic variables written in both stats and
+c          diagfi files
+c        -------------------------------------------------------------------
+
+         if (tracer) then
+           if (water) then
+
+             call zerophys(ngrid,mtot)
+             call zerophys(ngrid,icetot)
+             call zerophys(ngrid,rave)
+             call zerophys(ngrid,tauTES)
+             do ig=1,ngrid 
+               do l=1,nlayermx
+                 mtot(ig) = mtot(ig) + 
+     &                      zq(ig,l,igcm_h2o_vap) * 
+     &                      (pplev(ig,l) - pplev(ig,l+1)) / g
+                 icetot(ig) = icetot(ig) + 
+     &                        zq(ig,l,igcm_h2o_ice) * 
+     &                        (pplev(ig,l) - pplev(ig,l+1)) / g
+                 rave(ig) = rave(ig) + 
+     &                      zq(ig,l,igcm_h2o_ice) *
+     &                      (pplev(ig,l) - pplev(ig,l+1)) / g * 
+     &                      rice(ig,l) * (1.+nuice_ref)
+c                Computing abs optical depth at 825 cm-1 in each
+c                  layer to simulate NEW TES retrieval
+                 Qabsice = min(
+     &             max(0.4e6*rice(ig,l)*(1.+nuice_ref)-0.05 ,0.),1.2
+     &                        )
+                 opTES(ig,l)= 0.75 * Qabsice * 
+     &             zq(ig,l,igcm_h2o_ice) *
+     &             (pplev(ig,l) - pplev(ig,l+1)) / g
+     &             / (rho_ice * rice(ig,l) * (1.+nuice_ref))
+                 tauTES(ig)=tauTES(ig)+ opTES(ig,l) 
+               enddo
+               rave(ig)=rave(ig)/max(icetot(ig),1.e-30)
+               if (icetot(ig)*1e3.lt.0.01) rave(ig)=0.
+             enddo
+
+           endif ! of if (water)
+         endif ! of if (tracer)
+
+c        -----------------------------------------------------------------
+c        WSTATS: Saving statistics
+c        -----------------------------------------------------------------
+c        ("stats" stores and accumulates 8 key variables in file "stats.nc"
+c        which can later be used to make the statistic files of the run:
+c        "stats")          only possible in 3D runs !
+
+         
+         IF (callstats) THEN
+
+           call wstats(ngrid,"ps","Surface pressure","Pa",2,ps)
+           call wstats(ngrid,"tsurf","Surface temperature","K",2,tsurf)
+           call wstats(ngrid,"co2ice","CO2 ice cover",
+     &                "kg.m-2",2,co2ice)
+           call wstats(ngrid,"fluxsurf_lw",
+     &                "Thermal IR radiative flux to surface","W.m-2",2,
+     &                fluxsurf_lw)
+           call wstats(ngrid,"fluxsurf_sw",
+     &                "Solar radiative flux to surface","W.m-2",2,
+     &                fluxsurf_sw_tot)
+           call wstats(ngrid,"fluxtop_lw",
+     &                "Thermal IR radiative flux to space","W.m-2",2,
+     &                fluxtop_lw)
+           call wstats(ngrid,"fluxtop_sw",
+     &                "Solar radiative flux to space","W.m-2",2,
+     &                fluxtop_sw_tot)
+           call wstats(ngrid,"temp","Atmospheric temperature","K",3,zt)
+           call wstats(ngrid,"u","Zonal (East-West) wind","m.s-1",3,zu)
+           call wstats(ngrid,"v","Meridional (North-South) wind",
+     &                "m.s-1",3,zv)
+           call wstats(ngrid,"w","Vertical (down-up) wind",
+     &                "m.s-1",3,pw)
+           call wstats(ngrid,"rho","Atmospheric density","none",3,rho)
+           call wstats(ngrid,"pressure","Pressure","Pa",3,pplay)
+c          call wstats(ngrid,"q2",
+c    &                "Boundary layer eddy kinetic energy",
+c    &                "m2.s-2",3,q2)
+c          call wstats(ngrid,"emis","Surface emissivity","w.m-1",2,
+c    &                emis)
+c          call wstats(ngrid,"ssurf","Surface stress","N.m-2",
+c    &                2,zstress)
+c          call wstats(ngrid,"sw_htrt","sw heat.rate",
+c    &                 "W.m-2",3,zdtsw)
+c          call wstats(ngrid,"lw_htrt","lw heat.rate",
+c    &                 "W.m-2",3,zdtlw)
+
+           if (tracer) then
+             if (water) then
+               vmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_vap)
+     &                  *mugaz/mmol(igcm_h2o_vap)
+               call wstats(ngrid,"vmr_h2ovapor",
+     &                    "H2O vapor volume mixing ratio","mol/mol",
+     &                    3,vmr)
+               vmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_ice)
+     &                  *mugaz/mmol(igcm_h2o_ice)
+               call wstats(ngrid,"vmr_h2oice",
+     &                    "H2O ice volume mixing ratio","mol/mol",
+     &                    3,vmr)
+               call wstats(ngrid,"h2o_ice_s",
+     &                    "surface h2o_ice","kg/m2",
+     &                    2,qsurf(1,igcm_h2o_ice))
+
+               call wstats(ngrid,"mtot",
+     &                    "total mass of water vapor","kg/m2",
+     &                    2,mtot)
+               call wstats(ngrid,"icetot",
+     &                    "total mass of water ice","kg/m2",
+     &                    2,icetot)
+               call wstats(ngrid,"reffice",
+     &                    "Mean reff","m",
+     &                    2,rave)
+c              call wstats(ngrid,"rice",
+c    &                    "Ice particle size","m",
+c    &                    3,rice)
+c              If activice is true, tauTES is computed in aeropacity.F;
+               if (.not.activice) then
+                 call wstats(ngrid,"tauTESap",
+     &                      "tau abs 825 cm-1","",
+     &                      2,tauTES)
+               endif
+
+             endif ! of if (water)
+
+             if (thermochem.or.photochem) then
+                do iq=1,nq
+                   if ((noms(iq).eq."o").or.(noms(iq).eq."co2").or.
+     .                (noms(iq).eq."co").or.(noms(iq).eq."n2").or.
+     .                (noms(iq).eq."h2").or.
+     .                (noms(iq).eq."o3")) then
+                        do l=1,nlayer
+                          do ig=1,ngrid
+                            vmr(ig,l)=zq(ig,l,iq)*mmean(ig,l)/mmol(iq)
+                          end do
+                        end do
+                        call wstats(ngrid,"vmr_"//trim(noms(iq)),
+     .                     "Volume mixing ratio","mol/mol",3,vmr)
+                   endif
+                enddo
+             endif ! of if (thermochem.or.photochem)
+
+           endif ! of if (tracer)
+
+           IF(lastcall) THEN
+             write (*,*) "Writing stats..."
+             call mkstats(ierr)
+           ENDIF
+
+         ENDIF !if callstats
+
+c        (Store EOF for Mars Climate database software)
+         IF (calleofdump) THEN
+            CALL eofdump(ngrid, nlayer, zu, zv, zt, rho, ps)
+         ENDIF
+
+c        ==========================================================
+c        WRITEDIAGFI: Outputs in netcdf file "DIAGFI", containing
+c          any variable for diagnostic (output with period
+c          "ecritphy", set in "run.def")
+c        ==========================================================
+c        WRITEDIAGFI can ALSO be called from any other subroutines
+c        for any variables !!
+        call WRITEDIAGFI(ngrid,"emis","Surface emissivity","w.m-1",2,
+     &                  emis)
+         call WRITEDIAGFI(ngrid,"tsurf","Surface temperature","K",2,
+     &                  tsurf)
+         call WRITEDIAGFI(ngrid,"ps","surface pressure","Pa",2,ps)
+        call WRITEDIAGFI(ngrid,"co2ice","co2 ice thickness","kg.m-2",2,
+     &                  co2ice)
+c         call WRITEDIAGFI(ngrid,"fluxsurf_lw","fluxsurf_lw","W.m-2",2,
+c     &                  fluxsurf_lw)
+c         call WRITEDIAGFI(ngrid,"fluxsurf_sw","fluxsurf_sw","W.m-2",2,
+c     &                  fluxsurf_sw_tot)
+c         call WRITEDIAGFI(ngrid,"fluxtop_lw","fluxtop_lw","W.m-2",2,
+c     &                  fluxtop_lw)
+c         call WRITEDIAGFI(ngrid,"fluxtop_sw","fluxtop_sw","W.m-2",2,
+c     &                  fluxtop_sw_tot)
+         call WRITEDIAGFI(ngrid,"temp","temperature","K",3,zt)
+c        call WRITEDIAGFI(ngrid,"tau","tau"," ",2,tau)
+        call WRITEDIAGFI(ngrid,"u","Zonal wind","m.s-1",3,zu)
+        call WRITEDIAGFI(ngrid,"v","Meridional wind","m.s-1",3,zv)
+        call WRITEDIAGFI(ngrid,"w","Vertical wind","m.s-1",3,pw)
+c         call WRITEDIAGFI(ngrid,"rho","density","none",3,rho)
+c        call WRITEDIAGFI(ngrid,"q2","q2","kg.m-3",3,q2)
+c        call WRITEDIAGFI(ngrid,'Teta','T potentielle','K',3,zh)
+c        call WRITEDIAGFI(ngrid,"pressure","Pressure","Pa",3,pplay)
+c        call WRITEDIAGFI(ngrid,"ssurf","Surface stress","N.m-2",2,
+c    &                  zstress)
+c        call WRITEDIAGFI(ngridmx,'sw_htrt','sw heat. rate',
+c    &                   'w.m-2',3,zdtsw)
+c        call WRITEDIAGFI(ngridmx,'lw_htrt','lw heat. rate',
+c    &                   'w.m-2',3,zdtlw)
+
+!!!!!!!!!!!!!!!!!!!!!!!!SOIL SOIL SOIL
+        call WRITEDIAGFI(ngrid,"tsoil","Soil temperature",
+     &                       "K",3,tsoil)
+        call WRITEDIAGFI(ngrid,"inertiedat","Soil inertia",
+     &                       "K",3,inertiedat)
+!!!!!!!!!!!!!!!!!!!!!!!!SOIL SOIL SOIL
+
+c        ----------------------------------------------------------
+c        Outputs of the CO2 cycle
+c        ----------------------------------------------------------
+
+         if (tracer.and.(igcm_co2.ne.0)) then
+!          call WRITEDIAGFI(ngrid,"co2_l1","co2 mix. ratio in 1st layer",
+!    &                     "kg/kg",2,zq(1,1,igcm_co2))
+!          call WRITEDIAGFI(ngrid,"co2","co2 mass mixing ratio",
+!    &                     "kg/kg",3,zq(1,1,igcm_co2))
+        
+         ! Compute co2 column
+         call zerophys(ngrid,co2col)
+         do l=1,nlayermx
+           do ig=1,ngrid
+             co2col(ig)=co2col(ig)+
+     &                  zq(ig,l,igcm_co2)*(pplev(ig,l)-pplev(ig,l+1))/g
+           enddo
+         enddo
+c         call WRITEDIAGFI(ngrid,"co2col","CO2 column","kg.m-2",2,
+c     &                  co2col)
+         endif ! of if (tracer.and.(igcm_co2.ne.0))
+
+c        ----------------------------------------------------------
+c        Outputs of the water cycle
+c        ----------------------------------------------------------
+         if (tracer) then
+           if (water) then
+
+            !!!! waterice = q01, voir readmeteo.F90
+            call WRITEDIAGFI(ngridmx,'q01',noms(iq),
+     &                      'kg/kg',3,
+     &                       zq(1:ngridmx,1:nlayermx,igcm_h2o_ice))
+            !!!! watervapor = q02, voir readmeteo.F90
+            call WRITEDIAGFI(ngridmx,'q02',noms(iq),
+     &                      'kg/kg',3,
+     &                       zq(1:ngridmx,1:nlayermx,igcm_h2o_vap))
+
+c            call WRITEDIAGFI(ngridmx,'qsurf'//str2,noms(iq),
+c     &                     'kg.m-2',2,qsurf(1,iq))
+
+
+c             CALL WRITEDIAGFI(ngridmx,'mtot',
+c     &                       'total mass of water vapor',
+c     &                       'kg/m2',2,mtot)
+c             CALL WRITEDIAGFI(ngridmx,'icetot',
+c     &                       'total mass of water ice',
+c     &                       'kg/m2',2,icetot)
+cc            vmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_ice)
+cc    &                *mugaz/mmol(igcm_h2o_ice)
+cc            call WRITEDIAGFI(ngridmx,'vmr_h2oice','h2o ice vmr',
+cc    &                       'mol/mol',3,vmr)
+c             CALL WRITEDIAGFI(ngridmx,'reffice',
+c     &                       'Mean reff',
+c     &                       'm',2,rave)
+cc            call WRITEDIAGFI(ngridmx,'rice','Ice particle size',
+cc    &                       'm',3,rice)
+cc            If activice is true, tauTES is computed in aeropacity.F;
+c             if (.not.activice) then
+c               CALL WRITEDIAGFI(ngridmx,'tauTESap',
+c     &                         'tau abs 825 cm-1',
+c     &                         '',2,tauTES)
+c             endif
+c             call WRITEDIAGFI(ngridmx,'h2o_ice_s',
+c     &                       'surface h2o_ice',
+c     &                       'kg.m-2',2,qsurf(1,igcm_h2o_ice))
+           endif !(water)
+
+
+           if (water.and..not.photochem) then
+             iq=nq
+c            write(str2(1:2),'(i2.2)') iq
+c            call WRITEDIAGFI(ngridmx,'dqs'//str2,'dqscloud',
+c    &                       'kg.m-2',2,zdqscloud(1,iq))
+c            call WRITEDIAGFI(ngridmx,'dqch'//str2,'var chim',
+c    &                       'kg/kg',3,zdqchim(1,1,iq))
+c            call WRITEDIAGFI(ngridmx,'dqd'//str2,'var dif',
+c    &                       'kg/kg',3,zdqdif(1,1,iq))
+c            call WRITEDIAGFI(ngridmx,'dqa'//str2,'var adj',
+c    &                       'kg/kg',3,zdqadj(1,1,iq))
+c            call WRITEDIAGFI(ngridmx,'dqc'//str2,'var c',
+c    &                       'kg/kg',3,zdqc(1,1,iq))
+           endif  !(water.and..not.photochem)
+         endif
+
+c        ----------------------------------------------------------
+c        Outputs of the dust cycle
+c        ----------------------------------------------------------
+
+c        call WRITEDIAGFI(ngridmx,'tauref',
+c    &                    'Dust ref opt depth','NU',2,tauref)
+
+         if (tracer.and.(dustbin.ne.0)) then
+c          call WRITEDIAGFI(ngridmx,'tau','taudust','SI',2,tau(1,1))
+           if (doubleq) then
+cc            call WRITEDIAGFI(ngridmx,'qsurf','qsurf',
+cc    &                       'kg.m-2',2,qsurf(1,1))
+cc            call WRITEDIAGFI(ngridmx,'Nsurf','N particles',
+cc    &                       'N.m-2',2,qsurf(1,2))
+cc            call WRITEDIAGFI(ngridmx,'dqsdev','ddevil lift',
+cc    &                       'kg.m-2.s-1',2,zdqsdev(1,1))
+cc            call WRITEDIAGFI(ngridmx,'dqssed','sedimentation',
+cc    &                       'kg.m-2.s-1',2,zdqssed(1,1))
+c             call WRITEDIAGFI(ngridmx,'reffdust','reffdust',
+c     &                        'm',3,rdust*ref_r0)
+             call WRITEDIAGFI(ngridmx,'dustq','Dust mass mr',
+     &                        'kg/kg',3,pq(1,1,igcm_dust_mass))
+            call WRITEDIAGFI(ngridmx,'dustN','Dust number',
+     &                        'part/kg',3,pq(1,1,igcm_dust_number))
+           else
+             do iq=1,dustbin
+               write(str2(1:2),'(i2.2)') iq
+               call WRITEDIAGFI(ngridmx,'q'//str2,'mix. ratio',
+     &                         'kg/kg',3,zq(1,1,iq))
+               call WRITEDIAGFI(ngridmx,'qsurf'//str2,'qsurf',
+     &                         'kg.m-2',2,qsurf(1,iq))
+             end do
+           endif ! (doubleq)
+c          if (submicron) then
+c            call WRITEDIAGFI(ngridmx,'dustsubm','subm mass mr',
+c    &                        'kg/kg',3,pq(1,1,igcm_dust_submicron))
+c          endif ! (submicron)
+         end if  ! (tracer.and.(dustbin.ne.0))
+
+c        ----------------------------------------------------------
+c        Output in netcdf file "diagsoil.nc" for subterranean
+c          variables (output every "ecritphy", as for writediagfi)
+c        ----------------------------------------------------------
+
+         ! Write soil temperature
+!        call writediagsoil(ngrid,"soiltemp","Soil temperature","K",
+!    &                     3,tsoil)
+         ! Write surface temperature
+!        call writediagsoil(ngrid,"tsurf","Surface temperature","K",
+!    &                     2,tsurf)
+
+c        ==========================================================
+c        END OF WRITEDIAGFI
+c        ==========================================================
+
+      ELSE     ! if(ngrid.eq.1)
+
+         print*,'Ls =',zls*180./pi,
+     &  '  tauref(700 Pa) =',tauref
+c      ----------------------------------------------------------------------
+c      Output in grads file "g1d" (ONLY when using testphys1d)
+c      (output at every X physical timestep)
+c      ----------------------------------------------------------------------
+c
+c        CALL writeg1d(ngrid,1,fluxsurf_lw,'Fs_ir','W.m-2')
+c         CALL writeg1d(ngrid,1,tsurf,'tsurf','K')
+c         CALL writeg1d(ngrid,1,ps,'ps','Pa')
+         
+c         CALL writeg1d(ngrid,nlayer,zt,'T','K')
+c        CALL writeg1d(ngrid,nlayer,pu,'u','m.s-1')
+c        CALL writeg1d(ngrid,nlayer,pv,'v','m.s-1')
+c        CALL writeg1d(ngrid,nlayer,pw,'w','m.s-1')
+
+! or output in diagfi.nc (for testphys1d)
+         call WRITEDIAGFI(ngridmx,'ps','Surface pressure','Pa',0,ps)
+         call WRITEDIAGFI(ngridmx,'temp','Temperature',
+     &                       'K',1,zt)
+
+         if(tracer) then
+c           CALL writeg1d(ngrid,1,tau,'tau','SI')
+            do iq=1,nq
+c              CALL writeg1d(ngrid,nlayer,zq(1,1,iq),noms(iq),'kg/kg') 
+               call WRITEDIAGFI(ngridmx,trim(noms(iq)),
+     &              trim(noms(iq)),'kg/kg',1,zq(1,1,iq))
+            end do
+         end if
+
+         zlocal(1)=-log(pplay(1,1)/pplev(1,1))* Rnew(1,1)*zt(1,1)/g
+
+         do l=2,nlayer-1
+            tmean=zt(1,l)
+            if(zt(1,l).ne.zt(1,l-1))
+     &        tmean=(zt(1,l)-zt(1,l-1))/log(zt(1,l)/zt(1,l-1))
+              zlocal(l)= zlocal(l-1)
+     &        -log(pplay(1,l)/pplay(1,l-1))*rnew(1,l)*tmean/g
+         enddo
+         zlocal(nlayer)= zlocal(nlayer-1)-
+     &                   log(pplay(1,nlayer)/pplay(1,nlayer-1))*
+     &                   rnew(1,nlayer)*tmean/g
+
+      END IF       ! if(ngrid.ne.1)
+
+      icount=icount+1
+      RETURN
+      END
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/planete.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/planete.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/planete.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/planete.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/profile.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/profile.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/profile.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/profile.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/readtesassim.F90
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/readtesassim.F90	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/readtesassim.F90	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/readtesassim.F90
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/scatter.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/scatter.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/scatter.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/scatter.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/scopyi.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/scopyi.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/scopyi.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/scopyi.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/soil.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/soil.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/soil.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/soil.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/soil_settings.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/soil_settings.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/soil_settings.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/soil_settings.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/solang.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/solang.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/solang.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/solang.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/solarlong.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/solarlong.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/solarlong.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/solarlong.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/statto.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/statto.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/statto.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/statto.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/suaer.F90
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/suaer.F90	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/suaer.F90	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/suaer.F90
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/sugwd.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/sugwd.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/sugwd.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/sugwd.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/sulw.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/sulw.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/sulw.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/sulw.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/surfdat.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/surfdat.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/surfdat.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/surfdat.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/surfini.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/surfini.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/surfini.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/surfini.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/swmain.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/swmain.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/swmain.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/swmain.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/swr_fouquart.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/swr_fouquart.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/swr_fouquart.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/swr_fouquart.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/swr_toon.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/swr_toon.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/swr_toon.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/swr_toon.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/swrayleigh.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/swrayleigh.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/swrayleigh.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/swrayleigh.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/tabfi.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/tabfi.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/tabfi.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/tabfi.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/testphys1d.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/testphys1d.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/testphys1d.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/testphys1d.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/tracer.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/tracer.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/tracer.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/tracer.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/updatereffrad.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/updatereffrad.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/updatereffrad.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/updatereffrad.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/vdif_cd.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/vdif_cd.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/vdif_cd.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/vdif_cd.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/vdif_kc.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/vdif_kc.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/vdif_kc.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/vdif_kc.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/vdifc.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/vdifc.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/vdifc.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/vdifc.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/vlz_fi.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/vlz_fi.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/vlz_fi.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/vlz_fi.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/watercap.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/watercap.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/watercap.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/watercap.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/watercloud.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/watercloud.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/watercloud.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/watercloud.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/watersat.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/watersat.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/watersat.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/watersat.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/writediagfi.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/writediagfi.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/writediagfi.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/writediagfi.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/writediagsoil.F90
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/writediagsoil.F90	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/writediagsoil.F90	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/writediagsoil.F90
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/writeg1d.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/writeg1d.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/writeg1d.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/writeg1d.F
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/wstats.F90
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/wstats.F90	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/wstats.F90	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/wstats.F90
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/yoegwd.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/yoegwd.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/yoegwd.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/yoegwd.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/yomaer.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/yomaer.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/yomaer.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/yomaer.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/yomlw.h
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/yomlw.h	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/yomlw.h	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/yomlw.h
Index: trunk/mesoscale/LMDZ.MARS/libf/phymars/zerophys.F
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf/phymars/zerophys.F	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf/phymars/zerophys.F	(revision 54)
@@ -0,0 +1,1 @@
+link ../../libf_gcm/phymars/zerophys.F
Index: trunk/mesoscale/LMDZ.MARS/libf_gcm
===================================================================
--- trunk/mesoscale/LMDZ.MARS/libf_gcm	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/libf_gcm	(revision 54)
@@ -0,0 +1,1 @@
+link ../../mars/libf
Index: trunk/mesoscale/LMDZ.MARS/makegcm
===================================================================
--- trunk/mesoscale/LMDZ.MARS/makegcm	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/makegcm	(revision 54)
@@ -0,0 +1,810 @@
+#!/bin/csh
+# $Header: /users/lmdz/cvsroot/LMDZ.3.3/makegcm,v 1.21 2001/07/04 08:41:44 lmdz Exp $
+#set verbose echo
+########################################################################
+# options par defaut pour la commande make
+########################################################################
+set dim="64x48x32"
+set physique=mars
+set phys="PHYS=$physique"
+set include='-I$(LIBF)/grid -I$(LIBF)/bibio -I. -I$(LIBF)/aeronomars'
+set ntrac = 1
+set filtre=filtrez
+set grille=reg
+set dyntype="dyn"
+########################################################################
+# path a changer contenant les sources et les objets du modele
+########################################################################
+
+###### VERSION LMDZ.3.3
+# set LMDGCM=$HOME/LMDZ.3.3
+#
+
+#### If you want you can set environment variables here
+#setenv LMDGCM /u/emlmd/LMDZ.MARS.BETA.work
+#setenv LIBOGCM /tmp15/emlmd/test_LMDZ.BETA/libo
+#if ( `uname -m` == "x86_64" ) then
+## 64 bit machines
+#  setenv NCDFLIB /donnees/emlmd/netcdf64-4.0.1_pgi/lib
+#  setenv NCDFINC /donnees/emlmd/netcdf64-4.0.1_pgi/include
+#else
+#  setenv NCDFLIB /donnees/emlmd/netcdf-4.0.1_pgi/lib
+#  setenv NCDFINC /donnees/emlmd/netcdf-4.0.1_pgi/include
+#endif
+####
+
+setenv localdir "`pwd`"
+set MODIPSL=0
+echo $localdir | grep modipsl >& /dev/null
+if ( ! $status ) then
+  set MODIPSL=1
+  setenv LMDGCM $localdir
+  cd ../..
+  setenv LIBOGCM "`pwd`/lib"
+  cd $localdir
+  if ( `hostname` == rhodes ) then
+    set NCDFINC=`grep sxnec ../../util/AA_make.gdef| grep NCDF_INC|sed -e "s/^.* =//"`
+    set NCDFLIB=`grep sxnec ../../util/AA_make.gdef| grep NCDF_LIB|sed -e 's/^.* =//'`
+  else
+    if ( `hostname` == nymphea0 ) then
+      set NCDFINC=`grep fjvpp ../../util/AA_make.gdef| grep NCDF_INC|sed -e "s/^.* =//"`
+      set NCDFLIB=`grep fjvpp ../../util/AA_make.gdef| grep NCDF_LIB|sed -e 's/^.* =//'`
+    else
+      echo 'Probleme de definition des variables NCDFINC et NCDFLIB'
+    endif 
+  endif 
+else
+  if ( ! $?LMDGCM ) then
+    echo You must initialize the variable LMDGCM in your environnement
+    echo for instance: "setenv LMDGCM /usr/myself/supergcm" in .cshrc
+    exit
+  endif
+  if ( ! $?LIBOGCM ) then
+    set LIBOGCM=$LMDGCM/libo
+  endif
+  if ( ! $?NCDFLIB ) then
+    echo You must initialize the variable NCDFLIB in your environnement
+    echo for instance: "setenv NCDFLIB /usr/myself/netcdf" in .cshrc
+    exit
+  endif
+  if ( ! $?NCDFINC ) then
+    echo You must initialize the variable NCDFINC in your environnement
+    echo for instance: "setenv NCDFINC /usr/myself/netcdf" in .cshrc
+    exit
+  endif
+endif
+set model=$LMDGCM
+set libo=$LIBOGCM
+
+########################################################################
+#  Les differentes platformes reconnues
+########################################################################
+
+set HP=0
+set IBM=0
+set SUN=0
+set VPP=0
+set CRAY=0
+set DEC=0
+set LINUX=0
+set NEC=0
+set XNEC=0
+if ( `uname` == HP-UX ) then
+   set machine=HP
+   set HP=1
+else if (`uname` == UNIX_System_V ) then
+   set machine=VPP
+   set VPP=1
+else if (`uname` == SunOS ) then
+   set machine=SUN
+   set SUN=1
+else if ( `uname` == AIX ) then
+   set machine=IBM
+   set IBM=1
+else if ( `uname` == OSF1 ) then
+   set machine=ALPHA
+   set DEC=1
+else if ( `uname` == Linux ) then
+   set machine=LINUX
+   set LINUX=1
+else if ( `hostname` == atlas || `hostname` == axis  || `hostname` == etoile ) then
+   set machine=CRAY
+   set CRAY=1
+else if ( `uname` == SUPER-UX ) then
+   set machine=NEC
+   set NEC=1
+else if ( `hostname` == rhodes) then
+   set machine=XNEC
+   set XNEC=1
+else
+   echo Vous travaillez sur une machine non prevue par le reglement
+   exit
+endif
+
+# create $libo directory if it doesn't exist
+if ( ! -d $libo )  then
+   mkdir $libo
+endif
+
+if $VPP then
+set netcdf=netcdf_v
+else 
+set netcdf=netcdf
+endif
+########################################################################
+#  Quelques initialisations de variables du shell.
+########################################################################
+
+set dyn=
+set opt_link=""
+set adjnt=""
+set opt_dep=""
+
+set optim90=""
+set oplink=""
+
+########################################################################
+#  Optimisations par defaut suivant les machines
+########################################################################
+
+echo "Optimisations par defaut suivant les machines"
+set libf=$model/libf
+#setenv localdir "LOCAL_DIR=`pwd`"
+#setenv localdir "`pwd`"
+cd $model
+if $CRAY then
+   set optim90="-Wp'-P' -DCRAY "'-p$(LIBO) -eiv '
+   set oplink="-Wl'-DSTACK=128 -f indef' -L$NCDFLIB -lnetcdf "
+   set mod_loc_dir=" "
+   set mod_suffix=" "
+else if $SUN then
+   set optim90=" -fast"
+   set optimtru90=" -fast -free"
+   set opt_link="-L$NCDFLIB -lnetcdf"
+   set mod_loc_dir=$localdir
+   set mod_suffix=mod
+else if $HP then
+else if $IBM then
+else if $VPP then
+   set optim90="$optim -X9 -w"
+   if $COUPLE then
+     set opt_link="-Wg,-c $IOIPSLDIR/liboasis2.4_mpi2.a /usr/lang/mpi2/lib64/libmpi.a /usr/lang/mpi2/lib64/libmp.a /usr/local/lib/lib64/libnetcdf_cc.a -L$IOIPSLDIR -lioipsl"
+     set oplink="-Wl,-t,-P,-dy "
+   else
+     set opt_link="-Wg,-c /usr/local/lib/lib64/libnetcdf_cc.a -L$IOIPSLDIR -lioipsl"
+     set oplink="-Wl,-t,-dy "
+   endif
+   set mod_loc_dir=$IOIPSLDIR
+   set mod_suffix=mod
+else if $DEC then
+else if $LINUX then
+#   set optim90=" -fast"
+#   set optimtru90=" -fast -c -Mfree "
+   set optim=" -O2 -Munroll -Mcache_align"
+   set optim90=" -O2 -Munroll -Mcache_align"
+   set optimtru90=" -O2 -Munroll -Mcache_align"
+#   set opt_link=" -Mfree -L/usr/local/pgi/linux86/lib -lpgf90 -lpgftnrtl -lpghpf -lpghpf2 -L$NCDFLIB -lnetcdf -Wl,-Bstatic "
+# Ehouarn: there is no /usr/local/pgi/linux86/lib, but it doesn't matter
+#          also changed -Wl,-Bstatic to -Bstatic
+##   set opt_link=" -Mfree -lpgf90 -lpgftnrtl -lpghpf -lpghpf2 -L$NCDFLIB -lnetcdf -Bstatic "
+   set opt_link=" -L$NCDFLIB -lnetcdf -Bstatic"
+   set mod_loc_dir=$localdir
+   set mod_suffix=mod
+else if $NEC then
+   set optim90=' -clear -C hopt -float0 -ew -P static -Wf,"-pvctl fullmsg noassume "'
+   set optimtru90=' -clear -f4 -C hopt -float0 -ew -P static -Wf,"-pvctl fullmsg noassume "'
+###   set opt_link=" -C hopt -float0 -ew -P static -L$IOIPSLDIR -lioipsl  -L/u/rech/psl/rpsl003/IOIPSL -lnetcdf_i8r8_v "
+   set opt_link=" -C hopt -float0 -ew -P static -L/SX/usr/local/lib  -lnetcdf_i8r8 "
+   set mod_loc_dir="."
+   set mod_suffix="mod"
+else if $XNEC then
+   set optim90=' -clear -R5 -C vopt -float0 -ew -P static -Wf,"-pvctl loopcnt=40000 fullmsg noassume "'
+   set optimtru90=' -clear -R5 -f4 -C vopt -float0 -ew -P static -Wf,"-pvctl loopcnt=40000 fullmsg noassume "'
+   if $MODIPSL then
+     if $COUPLE then
+       set opt_link="-L$IOIPSLDIR -lsxioipsl -loasis2.4_mpi2 -float0 -ew -P static -I$NCDFINC $NCDFLIB "
+     else
+       set opt_link="-L$IOIPSLDIR -lsxioipsl -float0 -ew -P static -I$NCDFINC $NCDFLIB "
+     endif
+     set mod_loc_dir="./"
+   else
+####     set opt_link=" -C hopt -float0 -ew -P static -L$IOIPSLDIR -lsxioipsl -L /u/rech/psl/rpsl003/IOIPSL -lnetcdf_i8r8_v -I$NCDFINC -L/SX/usr/local/lib -lnetcdf_i8r8"
+     set opt_link=" -C hopt -float0 -ew -P static -I$NCDFINC -L/SX/usr/local/lib -lnetcdf_i8r8"
+     set mod_loc_dir="."
+   endif
+   set mod_suffix="mod"
+
+endif
+
+set nomlib=$machine
+
+# Impose distinct name for 64 bit Linux machines (to avoid mixing 32 and 64)
+if (`uname` == Linux && `uname -m` == "x86_64") then
+  set nomlib=${machine}64
+endif
+
+########################################################################
+# lecture des options de mymake
+########################################################################
+
+top:
+if ($#argv > 0) then
+    switch ($1:q)
+
+    case -h:
+
+########################################################################
+# Manuel en ligne
+########################################################################
+more <<eod
+
+
+makegcm [Options] prog
+
+
+The makegcm script:
+-------------------
+
+1. compiles a series of subroutines located in the $LMDGCM/libf
+ sub-directories.
+ The objects are then stored in the libraries in $LIBOGCM.
+
+2. then, makegcm compiles program prog.f located by default in
+$LMDGCM/libf/dyn3d and makes the link with the libraries.
+
+Environment Variables '$LMDGCM' and '$LIBOGCM'
+ must be set as environment variables or directly
+ in the makegcm file.
+
+The makegcm command is used to control the different versions of the model
+ in parallel, compiled using the compilation options 
+ and the various dimensions, without having to recompile the whole model.
+
+The FORTRAN libraries are stored in directory $LIBOGCM.
+
+
+OPTIONS:
+--------
+
+The following options can either be defined by default by editing the
+makegcm "script", or in interactive mode:
+
+-d imxjmxlm  where im, jm, and lm are the number of longitudes,
+             latitudes and vertical layers respectively.
+
+-t ntrac   Selects the number of tracers present in the model
+
+             Options -d and -t overwrite file 
+             $LMDGCM/libf/grid/dimensions.h
+             which contains the 3 dimensions of the
+             horizontal grid 
+             im, jm, lm plus the number of tracers passively advected
+             by the dynamics ntrac,
+             in 4 PARAMETER FORTRAN format 
+             with a new file:
+             $LMDGCM/libf/grid/dimension/dimensions.im.jm.lm.tntrac
+             If the file does not exist already
+             it is created by the script
+             $LMDGCM/libf/grid/dimension/makdim
+
+-p PHYS    Selects the set of physical parameterizations
+           you want to compile the model with.
+           The model is then compiled using the physical
+           parameterization sources in directory:
+            $LMDGCM/libf/phyPHYS
+
+-g grille  Selects the grid type.
+           This option overwrites file
+           $LMDGCM/libf/grid/fxyprim.h
+           with file
+           $LMDGCM/libf/grid/fxy_grille.h
+           the grid can take the following values:
+           1. reg - the regular grid
+           2. sin - to obtain equidistant points in terms of sin(latitude)
+           3. new - to zoom into a part of the globe
+
+-O "compilation options" set of fortran compilation options to use
+
+-include path
+           Used if the subroutines contain #include files (ccp) that 
+           are located in directories that are not referenced by default.
+
+-adjnt     Compiles the adjoint model to the dynamical code.
+
+-olddyn    To compile GCM with "old dynamics"
+
+-filtre  filter
+           To select the longitudinal filter in the polar regions.
+           "filter" corresponds to the name of a directory located in
+           $LMDGCM/libf. The standard filter for the model is "filtrez"
+           which can be used for a regular grid and for a  
+           grid with longitudinal zoom.
+
+-link "-Ldir1 -lfile1 -Ldir2 -lfile2 ..."
+           Adds a link to FORTRAN libraries
+           libfile1.a, libfile2.a ... 
+           located in directories dir1, dir2 ...respectively
+           If dirn is a directory with an automatic path 
+           (/usr/lib ... for example) 
+           there is no need to specify  -Ldirn.
+
+eod
+exit
+
+########################################################################
+# Lecture des differentes options
+########################################################################
+
+    case -d:
+        set dim=$2 ; shift ; shift ; goto top
+                        
+    case -O:
+        set optim="$2" ; set optim90="$2" ; set optimtru90="$2" ; shift ; shift ; goto top
+
+     case -p
+        set physique="$2" ; set phys="PHYS=$physique" ; shift ; shift ; goto top
+
+     case -g
+        set grille="$2" ; shift ; shift ; goto top
+
+     case -t
+        set ntrac=$2 ; shift ; shift ; goto top
+
+     case -include
+        set include="$include -I$2" ; shift ; shift ; goto top
+
+     case -adjnt
+        set opt_dep="$opt_dep adjnt" ; set adjnt="-ladjnt -ldyn3d "
+
+     case -olddyn
+        set dyntype="olddyn" ; shift; goto top
+
+     case -filtre
+        set filtre=$2 ; shift ; shift ; goto top
+
+     case -link
+        set opt_link="$opt_link $2" ; shift ; shift ; goto top
+
+     case -debug
+        if $HP then
+           set optim90=" -g "
+        else if $SUN then
+           setenv PARALLEL 4
+           set optim90=" -g -C "
+           set optimtru90=" -g -C "
+        else if $CRAY then
+           set optim90="$optim90"" -G1 "
+        else if $LINUX then
+           set optim90="-g -Mbounds "
+        else 
+           echo "pas d option debug predefinie pour cette machine"
+           exit
+        endif
+        shift ; goto top
+
+     default
+        set code="$1" ; shift ; goto top
+
+   endsw
+endif
+
+echo "apres les opts dim $dim"
+
+########################################################################
+# cas special sans physique
+########################################################################
+if ( "$physique" == 'nophys' ) then
+   set phys="L_PHY= LIBPHY="
+endif
+
+########################################################################
+# choix du nombre de traceur par defaut si il n'a pas ete choisi,
+# suivant la physique
+########################################################################
+
+if ( $ntrac == 0  ) then
+    if ( "$physique" == 'nophys' ) then
+        set ntrac=1
+    else if ( "$physique" == 'lmd' ) then
+        set ntrac=2
+    else if ( "$physique" == 'lmd_test_li' ) then
+        set ntrac=2
+    else if ( "$physique" == 'ec' ) then
+        set ntrac=1
+    else
+        set ntrac = 1
+    endif
+endif
+
+########################################################################
+#subtilites sur le nom de la librairie
+########################################################################
+
+\rm tmp ; touch tmp
+\rm tmp90 ; touch tmp90
+foreach i ( $optim90 )
+   echo $i | sed -e 's/\"//g' -e "s/\'//g" -e 's/-//g' -e 's/://g' -e 's/=//g' -e 's/%//g' >> tmp
+end
+set suf=
+foreach i ( `sort tmp | uniq ` )
+   set suf=$suf$i
+end
+if ( ! $IBM ) then
+   set nomlib="$nomlib$suf"
+endif
+if ( $DEC ) then
+   set nomlib=DEC
+endif
+
+# dimension
+
+echo "dimension avant sed $dim"
+if ( $IBM ) then
+   set dim=`echo $dim | sed -en 's/[^0-9]/ /g'`
+   set dim_=`echo $dim | sed -en 's/[^0-9]/_/g'`
+else if ( $SUN || $XNEC ) then
+   set dim=`echo $dim | sed -e 's/[^0-9]/ /g'` 
+   set dim_=`echo $dim | sed -e 's/[^0-9]/_/g'`
+else
+#   set dim=`echo $dim | sed -n -e 's/[^0-9]/ /gp'`
+#   set dim_=`echo $dim | sed -n -e 's/[^0-9]/_/gp'`
+     set dim=`echo $dim | sed -e 's/[^0-9]/ /g'`
+     set dim_=`echo $dim | sed -e 's/[^0-9]/_/g'`
+endif
+
+# build final name of libraries directory:
+if ( "$dyntype" == "olddyn" ) then
+  set nomlib=${nomlib}_${dim_}_t${ntrac}_${grille}_olddyn
+else
+  set nomlib=${nomlib}_${dim_}_t${ntrac}_${grille}
+endif
+
+# Append 'physique' type, if it is not mars, to nomlib
+if ( "$physique" != "mars" ) then
+  set nomlib=${nomlib}_${physique}
+endif
+
+## M-A-F nomlib trop long sur CRAY pour ar
+if ( $CRAY ) then
+    set nomlib=F90_${dim_}_t${ntrac}
+endif
+if ( $NEC || $XNEC ) then
+    set nomlib=F90_${dim_}_t${ntrac}
+endif
+
+echo "calcul de la dimension"
+set dimc=`echo $dim | wc -w`
+
+if ( "$dimc" == "2" ) then
+  set include="$include "'-I$(LIBF)/dyn2d '
+  set dimh=$dim
+else
+  if ( "$dyntype" == "olddyn" ) then
+    set include="$include "'-I$(LIBF)/olddyn3d '
+  else
+    set include="$include "'-I$(LIBF)/dyn3d '
+  endif
+  set dimh=`echo $dim | awk ' { print $1 "." $2 } '`
+endif
+echo "dimc is $dimc"
+
+########################################################################
+# path pour les #include
+########################################################################
+
+set include="$include -I$NCDFINC "
+echo $include
+
+########################################################################
+# Gestion des dimensions du modele.
+# on cree ou remplace le fichier des dimensions/nombre de traceur
+########################################################################
+
+cd $libf/grid
+if ( -f dimensions.h ) then
+  echo "ATTENTION: vous etes sans doute en train de compiler le modele par ailleurs"
+  echo "Attendez que la premiere compilation soit terminee pour relancer la suivante."
+  echo "Si vous etes sur que vous ne compilez pas le modele par ailleurs,"
+  echo "vous pouvez continuer en repondant oui."
+  echo "Voulez-vous vraiment continuer?"
+  if ( $< == "oui" ) then
+    #remove old dimensions.h file
+    \rm -f $libf/grid/dimensions.h
+  else
+    exit
+  endif
+endif
+
+# Build the appropriate 'dimensions.h' file
+cd dimension
+makdim $ntrac $dim
+# echo contents of dimensions.h to standard output
+cat $libf/grid/dimensions.h
+
+cd $LMDGCM
+# set path to objects directory
+set libo=$libo/$nomlib
+# create objects directory, if it doesn't exist
+if ( ! -d $libo )  then
+   mkdir $libo
+   cd $model
+endif
+
+########################################################################
+# Differentes dynamiques (3d, 2d, 1d)
+########################################################################
+
+set dimension=`echo $dim | wc -w`
+echo dimension $dimension dim $dim
+if ( $dimension == 1 ) then
+  echo pas de dynamique
+  set dyn="L_DYN= DYN= L_FILTRE= DIRMAIN=phy$physique "
+endif
+endif
+
+cd $model
+if ( $dimension == 3 ) then
+  cd libf/grid
+  \rm fxyprim.h
+  cp -p fxy_${grille}.h fxyprim.h
+endif
+
+######################################################################
+#   Traitement special pour le nouveau rayonnement de Laurent Li.
+######################################################################
+
+if ( -f $libf/phy$physique/raddim.h ) then
+  if ( -f $libf/phy$physique/raddim.$dimh.h ) then
+    \rm $libf/phy$physique/raddim.h
+    cp -p $libf/phy$physique/raddim.$dimh.h $libf/phy$physique/raddim.h
+    echo $libf/phy$physique/raddim.$dimh.h 
+    cat $libf/phy$physique/raddim.$dimh.h 
+    cat $libf/phy$physique/raddim.h
+  else
+    echo On peut diminuer la taille de l executable en creant
+    echo le fichier $libf/phy$physique/raddim.$dimh.h
+    \cp -p $libf/phy$physique/raddim.defaut.h $libf/phy$physique/raddim.h
+  endif
+endif
+
+######################################################################
+# Gestion du filtre qui n'existe qu'en 3d.
+######################################################################
+
+# set filtre to 'oldfiltrez' if using -olddyn option
+if ( "$dyntype" == "olddyn" ) then
+  set filtre="oldfiltrez"
+endif
+
+if ( `expr $dimc \> 2` == 1 ) then
+   set filtre="FILTRE=$filtre"
+else
+   set filtre="FILTRE= L_FILTRE= "
+endif
+echo "MACRO FILTRE $filtre"
+
+echo "dimc $dimc"
+
+########################################################################
+#  Avant de lancer le make, on recree le makefile si necessaire
+########################################################################
+# c'est a dire dans 3 cas:
+# 1. si la liste des fichiers .F et .h a ete modifiee depuis la
+#    derniere creation du makefile
+# 2. si le fichier contenant cette liste "liste_des_sources"
+#    n'existe pas.
+# 3. Si le makefile n'existe pas.
+########################################################################
+
+cd $model
+find libf -name '*.[Fh]' -print >! tmp77
+#find libf -name '*.[Fh]' -exec egrep -i " *use *ioipsl" {} \; -print >! tmp90
+find libf -name '*.[Fh]90' -print >> tmp90
+
+if (    `diff tmp77 liste_des_sources_f77 | wc -w` \
+     || `diff tmp90 liste_des_sources_f90 | wc -w` \
+     || ! -f makefile \
+     || ! -f liste_des_sources_f90 \
+     || ! -f liste_des_sources_f77 ) then
+        echo "les fichiers suivants ont ete crees ou detruits"
+        echo "ou les fichiers suivants sont passes ou ne sont plus en Fortran 90"
+        diff liste_des_sources_f77 tmp77
+        diff liste_des_sources_f90 tmp90
+        \cp tmp77 liste_des_sources_f77
+        \cp tmp90 liste_des_sources_f90
+        echo "On recree le makefile"
+        if ("$dyntype" == "olddyn") then
+          ./create_make_gcm olddyn3d grid bibio aeronomars phy$physique >! tmp
+        else
+          ./create_make_gcm dyn3d grid bibio aeronomars phy$physique >! tmp
+        endif
+        \mv tmp makefile
+        echo "Nouveau makefile cree."
+endif
+
+########################################################################
+#  Execution de la comande make
+########################################################################
+
+echo PHYSIQUE $phys
+echo dynamique $dyn $dimension
+echo OPTIM90="$optim90" $filtre LIBO=$libo $dyn PHYS=$phys DIM=$dimc PROG=$code
+echo PATH pour les fichiers INCLUDE $include
+echo OPLINK="$oplink"
+
+if $HP then
+   set f77='fort77 +OP'
+   set f90='jensaisrien'
+   set opt_link="$opt_link -lm"
+else  if $VPP then
+   set f77=frt
+   set f90=$f77
+else if $CRAY then
+   set f77=f90
+   set f90=f90
+else if $LINUX then
+   set f77=pgf90
+   set f90=pgf90
+else if $SUN then
+   set f77=f90
+   set f90=f90
+else if $NEC then
+   set f77=f90
+   set f90=f90
+else if $XNEC then
+   set f77=sxmpif90
+   set f90=sxmpif90
+else
+   set f77=f77
+   set f90=f90
+endif
+
+cd $model
+
+if $VPP then
+  set make="gmake RANLIB=ls"
+else if $CRAY then
+  set make="make RANLIB=ls"
+else if $NEC then
+  set make="make RANLIB=ls"
+else if $LINUX then
+  set make="make -k RANLIB=ranlib"
+else if $XNEC then
+  set make="/usr/local/bin/gmake RANLIB=ls"
+  set make="/usr/freeware/bin/gmake RANLIB=ls"
+else
+  set make="make RANLIB=ranlib"
+endif
+
+
+
+#
+# les deux test suivants sont "temporaires" pour pallier des "faiblesses" du 
+# compilateur fortran Sun: f90: SC4.0 11 Sep 1995 FORTRAN 90 1.1
+#
+#if ($code == 'create_limit' && $SUN) then
+#   set link=f77
+#   set opt_link="-L$NCDFLIB -lnetcdf"
+#endif
+
+#if ($code == 'create_etat0' && $SUN) then
+#   if ( ! -f $libo/libdyn3d.a ) then
+#     echo "Priere de compiler gcm en premier pour des raisons d'optimisation"
+#     \rm $libf/grid/dimensions.h
+#     exit
+#   endif
+#   set optim90=" -dalign -fixed "
+#   set opt_link="-L$IOIPSLDIR -lioipsl  -L$NCDFLIB -lnetcdf"
+#   set link="$f90 $optim90"
+#   touch $LMDGCM/libf/dyn3d/startvar.F
+#   touch $LMDGCM/libf/dyn3d/etat0_netcdf.F
+#endif
+
+#
+# Encore un test temporaire: probleme de compilation sur VPP
+# l'optimisation de startvar se passe mal
+#
+#if ($code == 'create_etat0' && $VPP) then
+#   if ( ! -f $libo/libdyn3d.a ) then
+#     echo "Priere de compiler gcm en premier pour des raisons d'optimisation"
+#     \rm $libf/grid/dimensions.h                                             
+#     exit                       
+#   endif
+#   set optim90="$optim90"" -X9 -w"
+#   set opt_link="-L$NCDFLIB -lnetcdf"
+#   touch $LMDGCM/libf/dyn3d/startvar.F
+#   touch $LMDGCM/libf/dyn3d/etat0_netcdf.F
+#endif
+
+#
+# etat0_netcdf a besoin d'info de la physique
+# A revoir
+set include="$include"" -I$libf/phy$physique"
+
+
+#################################################################
+# Execution de la comande make... ENFIN!
+#################################################################
+
+if $VPP then
+  set optim90=" $optim90 -Am -M$libo"
+  set optimtru90="$optim90"
+else if $SUN then
+ set optim90=" $optim90 -M$libo"
+ set optimtru90=" $optimtru90 "
+else if $NEC then
+ set optim90=" $optim90 -I$libo "
+else if $XNEC then
+ set optim90=" $optim90 -I$libo "
+ set optimtru90=" $optimtru90 -I$libo "
+else if $LINUX then
+ set optim="$optim -module $libo"
+ set optim90="$optim90 -module $libo"
+ set optimtru90="$optimtru90 -module $libo"
+# Ehouarn: remove set mod_loc_dir def below; mod_loc_dir=$localdir (set above)
+# set mod_loc_dir=$libo
+endif
+
+set link="$f90 $optim90"
+
+set ar=ar
+
+if $XNEC then
+  set link="sxld $opt_link"
+  set link="$f90 "
+#  set ar=sxar
+endif
+
+
+cd $localdir
+
+echo $make -f $LMDGCM/makefile \
+OPTION_DEP="$opt_dep" OPTION_LINK="$opt_link" \
+OPTIM="$optim" \
+OPTIM90="$optim90" \
+OPTIMTRU90="$optimtru90" \
+INCLUDE="$include" \
+$filtre \
+LIBO=$libo \
+$dyn \
+$phys \
+DIM=$dimc \
+DYNTYPE="$dyntype" \
+L_ADJNT="$adjnt" \
+LOCAL_DIR="$localdir"  \
+F77="$f77" \
+F90="$f90" \
+OPLINK="$oplink" \
+LINK="$link" \
+GCM="$LMDGCM" \
+MOD_LOC_DIR=$mod_loc_dir \
+MOD_SUFFIX=$mod_suffix \
+AR=$ar \
+PROG=$code
+
+
+$make -f $LMDGCM/makefile \
+OPTION_DEP="$opt_dep" OPTION_LINK="$opt_link" \
+OPTIM="$optim" \
+OPTIM90="$optim90" \
+OPTIMTRU90="$optimtru90" \
+INCLUDE="$include" \
+$filtre \
+LIBO=$libo \
+$dyn \
+$phys \
+DIM=$dimc \
+DYNTYPE="$dyntype" \
+L_ADJNT="$adjnt" \
+LOCAL_DIR="$localdir"  \
+F77="$f77" \
+F90="$f90" \
+OPLINK="$oplink" \
+LINK="$link" \
+GCM="$LMDGCM" \
+MOD_LOC_DIR=$mod_loc_dir \
+MOD_SUFFIX=$mod_suffix \
+AR=$ar \
+PROG=$code
+
+\rm $libf/grid/dimensions.h
Index: trunk/mesoscale/LMDZ.MARS/makegcm_g95
===================================================================
--- trunk/mesoscale/LMDZ.MARS/makegcm_g95	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/makegcm_g95	(revision 54)
@@ -0,0 +1,812 @@
+#!/bin/csh
+# $Header: /users/lmdz/cvsroot/LMDZ.3.3/makegcm,v 1.21 2001/07/04 08:41:44 lmdz Exp $
+#set verbose echo
+########################################################################
+# options par defaut pour la commande make
+########################################################################
+set dim="64x48x32"
+set physique=mars
+set phys="PHYS=$physique"
+set include='-I$(LIBF)/grid -I$(LIBF)/bibio -I. -I$(LIBF)/aeronomars'
+set ntrac = 1
+set filtre=filtrez
+set grille=reg
+set dyntype="dyn"
+########################################################################
+# path a changer contenant les sources et les objets du modele
+########################################################################
+
+###### VERSION LMDZ.3.3
+# set LMDGCM=$HOME/LMDZ.3.3
+#
+
+#### If you want you can set environment variables here
+#setenv LMDGCM /u/emlmd/LMDZ.MARS.BETA.work
+#setenv LIBOGCM /tmp15/emlmd/libo
+#if ( `uname -m` == "x86_64" ) then
+# 64 bit machines
+#  setenv NCDFLIB /donnees/emlmd/netcdf64-4.0.1_g95/lib
+#  setenv NCDFINC /donnees/emlmd/netcdf64-4.0.1_g95/include
+#else
+#  setenv NCDFLIB /donnees/emlmd/netcdf-4.0.1_g95/lib
+#  setenv NCDFINC /donnees/emlmd/netcdf-4.0.1_g95/include
+#endif
+####
+
+setenv localdir "`pwd`"
+set MODIPSL=0
+echo $localdir | grep modipsl >& /dev/null
+if ( ! $status ) then
+  set MODIPSL=1
+  setenv LMDGCM $localdir
+  cd ../..
+  setenv LIBOGCM "`pwd`/lib"
+  cd $localdir
+  if ( `hostname` == rhodes ) then
+    set NCDFINC=`grep sxnec ../../util/AA_make.gdef| grep NCDF_INC|sed -e "s/^.* =//"`
+    set NCDFLIB=`grep sxnec ../../util/AA_make.gdef| grep NCDF_LIB|sed -e 's/^.* =//'`
+  else
+    if ( `hostname` == nymphea0 ) then
+      set NCDFINC=`grep fjvpp ../../util/AA_make.gdef| grep NCDF_INC|sed -e "s/^.* =//"`
+      set NCDFLIB=`grep fjvpp ../../util/AA_make.gdef| grep NCDF_LIB|sed -e 's/^.* =//'`
+    else
+      echo 'Probleme de definition des variables NCDFINC et NCDFLIB'
+    endif 
+  endif 
+else
+  if ( ! $?LMDGCM ) then
+    echo You must initialize the variable LMDGCM in your environnement
+    echo for instance: "setenv LMDGCM /usr/myself/supergcm" in .cshrc
+    exit
+  endif
+  if ( ! $?LIBOGCM ) then
+    set LIBOGCM=$LMDGCM/libo
+  endif
+  if ( ! $?NCDFLIB ) then
+    echo You must initialize the variable NCDFLIB in your environnement
+    echo for instance: "setenv NCDFLIB /usr/myself/netcdf" in .cshrc
+    exit
+  endif
+  if ( ! $?NCDFINC ) then
+    echo You must initialize the variable NCDFINC in your environnement
+    echo for instance: "setenv NCDFINC /usr/myself/netcdf" in .cshrc
+    exit
+  endif
+endif
+set model=$LMDGCM
+set libo=$LIBOGCM
+
+########################################################################
+#  Les differentes platformes reconnues
+########################################################################
+
+set HP=0
+set IBM=0
+set SUN=0
+set VPP=0
+set CRAY=0
+set DEC=0
+set LINUX=0
+set NEC=0
+set XNEC=0
+if ( `uname` == HP-UX ) then
+   set machine=HP
+   set HP=1
+else if (`uname` == UNIX_System_V ) then
+   set machine=VPP
+   set VPP=1
+else if (`uname` == SunOS ) then
+   set machine=SUN
+   set SUN=1
+else if ( `uname` == AIX ) then
+   set machine=IBM
+   set IBM=1
+else if ( `uname` == OSF1 ) then
+   set machine=ALPHA
+   set DEC=1
+else if ( `uname` == Linux ) then
+   set machine=LINUX
+   set LINUX=1
+else if ( `hostname` == atlas || `hostname` == axis  || `hostname` == etoile ) then
+   set machine=CRAY
+   set CRAY=1
+else if ( `uname` == SUPER-UX ) then
+   set machine=NEC
+   set NEC=1
+else if ( `hostname` == rhodes) then
+   set machine=XNEC
+   set XNEC=1
+else
+   echo Vous travaillez sur une machine non prevue par le reglement
+   exit
+endif
+
+# create $libo directory if it doesn't exist
+if ( ! -d $libo )  then
+   mkdir $libo
+endif
+
+if $VPP then
+set netcdf=netcdf_v
+else 
+set netcdf=netcdf
+endif
+########################################################################
+#  Quelques initialisations de variables du shell.
+########################################################################
+
+set dyn=
+set opt_link=""
+set adjnt=""
+set opt_dep=""
+
+set optim90=""
+set oplink=""
+
+########################################################################
+#  Optimisations par defaut suivant les machines
+########################################################################
+
+echo "Optimisations par defaut suivant les machines"
+set libf=$model/libf
+#setenv localdir "LOCAL_DIR=`pwd`"
+#setenv localdir "`pwd`"
+cd $model
+if $CRAY then
+   set optim90="-Wp'-P' -DCRAY "'-p$(LIBO) -eiv '
+   set oplink="-Wl'-DSTACK=128 -f indef' -L$NCDFLIB -lnetcdf "
+   set mod_loc_dir=" "
+   set mod_suffix=" "
+else if $SUN then
+   set optim90=" -fast"
+   set optimtru90=" -fast -free"
+   set opt_link="-L$NCDFLIB -lnetcdf"
+   set mod_loc_dir=$localdir
+   set mod_suffix=mod
+else if $HP then
+else if $IBM then
+else if $VPP then
+   set optim90="$optim -X9 -w"
+   if $COUPLE then
+     set opt_link="-Wg,-c $IOIPSLDIR/liboasis2.4_mpi2.a /usr/lang/mpi2/lib64/libmpi.a /usr/lang/mpi2/lib64/libmp.a /usr/local/lib/lib64/libnetcdf_cc.a -L$IOIPSLDIR -lioipsl"
+     set oplink="-Wl,-t,-P,-dy "
+   else
+     set opt_link="-Wg,-c /usr/local/lib/lib64/libnetcdf_cc.a -L$IOIPSLDIR -lioipsl"
+     set oplink="-Wl,-t,-dy "
+   endif
+   set mod_loc_dir=$IOIPSLDIR
+   set mod_suffix=mod
+else if $DEC then
+else if $LINUX then
+# Ehouarn 'g95' compiler
+   set optim="-O3 -fstatic -funroll-loops "
+   set optim90=" -O3 -fstatic -funroll-loops "
+   set optimtru90=" -O3 -fstatic -funroll-loops "
+#   set opt_link=" -Mfree -L/usr/local/pgi/linux86/lib -lpgf90 -lpgftnrtl -lpghpf -lpghpf2 -L$NCDFLIB -lnetcdf -Wl,-Bstatic "
+# Ehouarn: there is no /usr/local/pgi/linux86/lib, but it doesn't matter
+#          also changed -Wl,-Bstatic to -Bstatic
+   set opt_link=" -L$NCDFLIB -lnetcdf "
+   set mod_loc_dir=$localdir
+   set mod_suffix=mod
+else if $NEC then
+   set optim90=' -clear -C hopt -float0 -ew -P static -Wf,"-pvctl fullmsg noassume "'
+   set optimtru90=' -clear -f4 -C hopt -float0 -ew -P static -Wf,"-pvctl fullmsg noassume "'
+###   set opt_link=" -C hopt -float0 -ew -P static -L$IOIPSLDIR -lioipsl  -L/u/rech/psl/rpsl003/IOIPSL -lnetcdf_i8r8_v "
+   set opt_link=" -C hopt -float0 -ew -P static -L/SX/usr/local/lib  -lnetcdf_i8r8 "
+   set mod_loc_dir="."
+   set mod_suffix="mod"
+else if $XNEC then
+   set optim90=' -clear -R5 -C vopt -float0 -ew -P static -Wf,"-pvctl loopcnt=40000 fullmsg noassume "'
+   set optimtru90=' -clear -R5 -f4 -C vopt -float0 -ew -P static -Wf,"-pvctl loopcnt=40000 fullmsg noassume "'
+   if $MODIPSL then
+     if $COUPLE then
+       set opt_link="-L$IOIPSLDIR -lsxioipsl -loasis2.4_mpi2 -float0 -ew -P static -I$NCDFINC $NCDFLIB "
+     else
+       set opt_link="-L$IOIPSLDIR -lsxioipsl -float0 -ew -P static -I$NCDFINC $NCDFLIB "
+     endif
+     set mod_loc_dir="./"
+   else
+####     set opt_link=" -C hopt -float0 -ew -P static -L$IOIPSLDIR -lsxioipsl -L /u/rech/psl/rpsl003/IOIPSL -lnetcdf_i8r8_v -I$NCDFINC -L/SX/usr/local/lib -lnetcdf_i8r8"
+     set opt_link=" -C hopt -float0 -ew -P static -I$NCDFINC -L/SX/usr/local/lib -lnetcdf_i8r8"
+     set mod_loc_dir="."
+   endif
+   set mod_suffix="mod"
+
+endif
+
+# Ehouarn: add 'g95' string to pathname
+set nomlib=${machine}_g95
+
+# Impose distinct name for 64 bit Linux machines (to avoid mixing 32 and 64)
+if (`uname` == Linux && `uname -m` == "x86_64") then
+  set nomlib=${machine}64_g95
+endif
+
+########################################################################
+# lecture des options de mymake
+########################################################################
+
+top:
+if ($#argv > 0) then
+    switch ($1:q)
+
+    case -h:
+
+########################################################################
+# Manuel en ligne
+########################################################################
+more <<eod
+
+
+makegcm [Options] prog
+
+
+The makegcm script:
+-------------------
+
+1. compiles a series of subroutines located in the $LMDGCM/libf
+ sub-directories.
+ The objects are then stored in the libraries in $LIBOGCM.
+
+2. then, makegcm compiles program prog.f located by default in
+$LMDGCM/libf/dyn3d and makes the link with the libraries.
+
+Environment Variables '$LMDGCM' and '$LIBOGCM'
+ must be set as environment variables or directly
+ in the makegcm file.
+
+The makegcm command is used to control the different versions of the model
+ in parallel, compiled using the compilation options 
+ and the various dimensions, without having to recompile the whole model.
+
+The FORTRAN libraries are stored in directory $LIBOGCM.
+
+
+OPTIONS:
+--------
+
+The following options can either be defined by default by editing the
+makegcm "script", or in interactive mode:
+
+-d imxjmxlm  where im, jm, and lm are the number of longitudes,
+             latitudes and vertical layers respectively.
+
+-t ntrac   Selects the number of tracers present in the model
+
+             Options -d and -t overwrite file 
+             $LMDGCM/libf/grid/dimensions.h
+             which contains the 3 dimensions of the
+             horizontal grid 
+             im, jm, lm plus the number of tracers passively advected
+             by the dynamics ntrac,
+             in 4 PARAMETER FORTRAN format 
+             with a new file:
+             $LMDGCM/libf/grid/dimension/dimensions.im.jm.lm.tntrac
+             If the file does not exist already
+             it is created by the script
+             $LMDGCM/libf/grid/dimension/makdim
+
+-p PHYS    Selects the set of physical parameterizations
+           you want to compile the model with.
+           The model is then compiled using the physical
+           parameterization sources in directory:
+            $LMDGCM/libf/phyPHYS
+
+-g grille  Selects the grid type.
+           This option overwrites file
+           $LMDGCM/libf/grid/fxyprim.h
+           with file
+           $LMDGCM/libf/grid/fxy_grille.h
+           the grid can take the following values:
+           1. reg - the regular grid
+           2. sin - to obtain equidistant points in terms of sin(latitude)
+           3. new - to zoom into a part of the globe
+
+-O "compilation options" set of fortran compilation options to use
+
+-include path
+           Used if the subroutines contain #include files (ccp) that 
+           are located in directories that are not referenced by default.
+
+-adjnt     Compiles the adjoint model to the dynamical code.
+
+-olddyn    To compile GCM with "old dynamics"
+
+-filtre  filter
+           To select the longitudinal filter in the polar regions.
+           "filter" corresponds to the name of a directory located in
+           $LMDGCM/libf. The standard filter for the model is "filtrez"
+           which can be used for a regular grid and for a  
+           grid with longitudinal zoom.
+
+-link "-Ldir1 -lfile1 -Ldir2 -lfile2 ..."
+           Adds a link to FORTRAN libraries
+           libfile1.a, libfile2.a ... 
+           located in directories dir1, dir2 ...respectively
+           If dirn is a directory with an automatic path 
+           (/usr/lib ... for example) 
+           there is no need to specify  -Ldirn.
+
+eod
+exit
+
+########################################################################
+# Lecture des differentes options
+########################################################################
+
+    case -d:
+        set dim=$2 ; shift ; shift ; goto top
+                        
+    case -O:
+        set optim="$2" ; set optim90="$2" ; set optimtru90="$2" ; shift ; shift ; goto top
+
+     case -p
+        set physique="$2" ; set phys="PHYS=$physique" ; shift ; shift ; goto top
+
+     case -g
+        set grille="$2" ; shift ; shift ; goto top
+
+     case -t
+        set ntrac=$2 ; shift ; shift ; goto top
+
+     case -include
+        set include="$include -I$2" ; shift ; shift ; goto top
+
+     case -adjnt
+        set opt_dep="$opt_dep adjnt" ; set adjnt="-ladjnt -ldyn3d "
+
+     case -olddyn
+        set dyntype="olddyn" ; shift; goto top
+
+     case -filtre
+        set filtre=$2 ; shift ; shift ; goto top
+
+     case -link
+        set opt_link="$opt_link $2" ; shift ; shift ; goto top
+
+     case -debug
+        if $HP then
+           set optim90=" -g "
+        else if $SUN then
+           setenv PARALLEL 4
+           set optim90=" -g -C "
+           set optimtru90=" -g -C "
+        else if $CRAY then
+           set optim90="$optim90"" -G1 "
+        else if $LINUX then
+           set optim90="-g -Mbounds "
+        else 
+           echo "pas d option debug predefinie pour cette machine"
+           exit
+        endif
+        shift ; goto top
+
+     default
+        set code="$1" ; shift ; goto top
+
+   endsw
+endif
+
+echo "apres les opts dim $dim"
+
+########################################################################
+# cas special sans physique
+########################################################################
+if ( "$physique" == 'nophys' ) then
+   set phys="L_PHY= LIBPHY="
+endif
+
+########################################################################
+# choix du nombre de traceur par defaut si il n'a pas ete choisi,
+# suivant la physique
+########################################################################
+
+if ( $ntrac == 0  ) then
+    if ( "$physique" == 'nophys' ) then
+        set ntrac=1
+    else if ( "$physique" == 'lmd' ) then
+        set ntrac=2
+    else if ( "$physique" == 'lmd_test_li' ) then
+        set ntrac=2
+    else if ( "$physique" == 'ec' ) then
+        set ntrac=1
+    else
+        set ntrac = 1
+    endif
+endif
+
+########################################################################
+#subtilites sur le nom de la librairie
+########################################################################
+
+\rm tmp ; touch tmp
+\rm tmp90 ; touch tmp90
+foreach i ( $optim90 )
+   echo $i | sed -e 's/\"//g' -e "s/\'//g" -e 's/-//g' -e 's/://g' -e 's/=//g' -e 's/%//g' >> tmp
+end
+set suf=
+foreach i ( `sort tmp | uniq ` )
+   set suf=$suf$i
+end
+if ( ! $IBM ) then
+   set nomlib="$nomlib$suf"
+endif
+if ( $DEC ) then
+   set nomlib=DEC
+endif
+
+# dimension
+
+echo "dimension avant sed $dim"
+if ( $IBM ) then
+   set dim=`echo $dim | sed -en 's/[^0-9]/ /g'`
+   set dim_=`echo $dim | sed -en 's/[^0-9]/_/g'`
+else if ( $SUN || $XNEC ) then
+   set dim=`echo $dim | sed -e 's/[^0-9]/ /g'` 
+   set dim_=`echo $dim | sed -e 's/[^0-9]/_/g'`
+else
+#   set dim=`echo $dim | sed -n -e 's/[^0-9]/ /gp'`
+#   set dim_=`echo $dim | sed -n -e 's/[^0-9]/_/gp'`
+     set dim=`echo $dim | sed -e 's/[^0-9]/ /g'`
+     set dim_=`echo $dim | sed -e 's/[^0-9]/_/g'`
+endif
+
+# build final name of libraries directory:
+if ( "$dyntype" == "olddyn" ) then
+  set nomlib=${nomlib}_${dim_}_t${ntrac}_${grille}_olddyn
+else
+  set nomlib=${nomlib}_${dim_}_t${ntrac}_${grille}
+endif
+
+# Append 'physique' type, if it is not mars, to nomlib
+if ( "$physique" != "mars" ) then
+  set nomlib=${nomlib}_${physique}
+endif
+
+## M-A-F nomlib trop long sur CRAY pour ar
+if ( $CRAY ) then
+    set nomlib=F90_${dim_}_t${ntrac}
+endif
+if ( $NEC || $XNEC ) then
+    set nomlib=F90_${dim_}_t${ntrac}
+endif
+
+echo "calcul de la dimension"
+set dimc=`echo $dim | wc -w`
+
+if ( "$dimc" == "2" ) then
+  set include="$include "'-I$(LIBF)/dyn2d '
+  set dimh=$dim
+else
+  if ( "$dyntype" == "olddyn" ) then
+    set include="$include "'-I$(LIBF)/olddyn3d '
+  else
+    set include="$include "'-I$(LIBF)/dyn3d '
+  endif
+  set dimh=`echo $dim | awk ' { print $1 "." $2 } '`
+endif
+echo "dimc is $dimc"
+
+########################################################################
+# path pour les #include
+########################################################################
+
+set include="$include -I$NCDFINC "
+echo $include
+
+########################################################################
+# Gestion des dimensions du modele.
+# on cree ou remplace le fichier des dimensions/nombre de traceur
+########################################################################
+
+cd $libf/grid
+if ( -f dimensions.h ) then
+  echo "ATTENTION: vous etes sans doute en train de compiler le modele par ailleurs"
+  echo "Attendez que la premiere compilation soit terminee pour relancer la suivante."
+  echo "Si vous etes sur que vous ne compilez pas le modele par ailleurs,"
+  echo "vous pouvez continuer en repondant oui."
+  echo "Voulez-vous vraiment continuer?"
+  if ( $< == "oui" ) then
+    #remove old dimensions.h file
+    \rm -f $libf/grid/dimensions.h
+  else
+    exit
+  endif
+endif
+
+# Build the appropriate 'dimensions.h' file
+cd dimension
+makdim $ntrac $dim
+# echo contents of dimensions.h to standard output
+cat $libf/grid/dimensions.h
+
+cd $LMDGCM
+# set path to objects directory
+set libo=$libo/$nomlib
+# create objects directory, if it doesn't exist
+if ( ! -d $libo )  then
+   mkdir $libo
+   cd $model
+endif
+
+########################################################################
+# Differentes dynamiques (3d, 2d, 1d)
+########################################################################
+
+set dimension=`echo $dim | wc -w`
+echo dimension $dimension dim $dim
+if ( $dimension == 1 ) then
+  echo pas de dynamique
+  set dyn="L_DYN= DYN= L_FILTRE= DIRMAIN=phy$physique "
+endif
+endif
+
+cd $model
+if ( $dimension == 3 ) then
+  cd libf/grid
+  \rm fxyprim.h
+  cp -p fxy_${grille}.h fxyprim.h
+endif
+
+######################################################################
+#   Traitement special pour le nouveau rayonnement de Laurent Li.
+######################################################################
+
+if ( -f $libf/phy$physique/raddim.h ) then
+  if ( -f $libf/phy$physique/raddim.$dimh.h ) then
+    \rm $libf/phy$physique/raddim.h
+    cp -p $libf/phy$physique/raddim.$dimh.h $libf/phy$physique/raddim.h
+    echo $libf/phy$physique/raddim.$dimh.h 
+    cat $libf/phy$physique/raddim.$dimh.h 
+    cat $libf/phy$physique/raddim.h
+  else
+    echo On peut diminuer la taille de l executable en creant
+    echo le fichier $libf/phy$physique/raddim.$dimh.h
+    \cp -p $libf/phy$physique/raddim.defaut.h $libf/phy$physique/raddim.h
+  endif
+endif
+
+######################################################################
+# Gestion du filtre qui n'existe qu'en 3d.
+######################################################################
+
+# set filtre to 'oldfiltrez' if using -olddyn option
+if ( "$dyntype" == "olddyn" ) then
+  set filtre="oldfiltrez"
+endif
+
+if ( `expr $dimc \> 2` == 1 ) then
+   set filtre="FILTRE=$filtre"
+else
+   set filtre="FILTRE= L_FILTRE= "
+endif
+echo "MACRO FILTRE $filtre"
+
+echo "dimc $dimc"
+
+########################################################################
+#  Avant de lancer le make, on recree le makefile si necessaire
+########################################################################
+# c'est a dire dans 3 cas:
+# 1. si la liste des fichiers .F et .h a ete modifiee depuis la
+#    derniere creation du makefile
+# 2. si le fichier contenant cette liste "liste_des_sources"
+#    n'existe pas.
+# 3. Si le makefile n'existe pas.
+########################################################################
+
+cd $model
+find libf -name '*.[Fh]' -print >! tmp77
+#find libf -name '*.[Fh]' -exec egrep -i " *use *ioipsl" {} \; -print >! tmp90
+find libf -name '*.[Fh]90' -print >> tmp90
+
+if (    `diff tmp77 liste_des_sources_f77 | wc -w` \
+     || `diff tmp90 liste_des_sources_f90 | wc -w` \
+     || ! -f makefile \
+     || ! -f liste_des_sources_f90 \
+     || ! -f liste_des_sources_f77 ) then
+        echo "les fichiers suivants ont ete crees ou detruits"
+        echo "ou les fichiers suivants sont passes ou ne sont plus en Fortran 90"
+        diff liste_des_sources_f77 tmp77
+        diff liste_des_sources_f90 tmp90
+        \cp tmp77 liste_des_sources_f77
+        \cp tmp90 liste_des_sources_f90
+        echo "On recree le makefile"
+        if ("$dyntype" == "olddyn") then
+          ./create_make_gcm olddyn3d grid bibio aeronomars phy$physique >! tmp
+        else
+          ./create_make_gcm dyn3d grid bibio aeronomars phy$physique >! tmp
+        endif
+        \mv tmp makefile
+        echo "Nouveau makefile cree."
+endif
+
+########################################################################
+#  Execution de la comande make
+########################################################################
+
+echo PHYSIQUE $phys
+echo dynamique $dyn $dimension
+echo OPTIM90="$optim90" $filtre LIBO=$libo $dyn PHYS=$phys DIM=$dimc PROG=$code
+echo PATH pour les fichiers INCLUDE $include
+echo OPLINK="$oplink"
+
+if $HP then
+   set f77='fort77 +OP'
+   set f90='jensaisrien'
+   set opt_link="$opt_link -lm"
+else  if $VPP then
+   set f77=frt
+   set f90=$f77
+else if $CRAY then
+   set f77=f90
+   set f90=f90
+else if $LINUX then
+#   set f77=pgf90
+#   set f90=pgf90
+   set f77=g95
+   set f90=g95
+else if $SUN then
+   set f77=f90
+   set f90=f90
+else if $NEC then
+   set f77=f90
+   set f90=f90
+else if $XNEC then
+   set f77=sxmpif90
+   set f90=sxmpif90
+else
+   set f77=f77
+   set f90=f90
+endif
+
+cd $model
+
+if $VPP then
+  set make="gmake RANLIB=ls"
+else if $CRAY then
+  set make="make RANLIB=ls"
+else if $NEC then
+  set make="make RANLIB=ls"
+else if $LINUX then
+  set make="make -k RANLIB=ranlib"
+else if $XNEC then
+  set make="/usr/local/bin/gmake RANLIB=ls"
+  set make="/usr/freeware/bin/gmake RANLIB=ls"
+else
+  set make="make RANLIB=ranlib"
+endif
+
+
+
+#
+# les deux test suivants sont "temporaires" pour pallier des "faiblesses" du 
+# compilateur fortran Sun: f90: SC4.0 11 Sep 1995 FORTRAN 90 1.1
+#
+#if ($code == 'create_limit' && $SUN) then
+#   set link=f77
+#   set opt_link="-L$NCDFLIB -lnetcdf"
+#endif
+
+#if ($code == 'create_etat0' && $SUN) then
+#   if ( ! -f $libo/libdyn3d.a ) then
+#     echo "Priere de compiler gcm en premier pour des raisons d'optimisation"
+#     \rm $libf/grid/dimensions.h
+#     exit
+#   endif
+#   set optim90=" -dalign -fixed "
+#   set opt_link="-L$IOIPSLDIR -lioipsl  -L$NCDFLIB -lnetcdf"
+#   set link="$f90 $optim90"
+#   touch $LMDGCM/libf/dyn3d/startvar.F
+#   touch $LMDGCM/libf/dyn3d/etat0_netcdf.F
+#endif
+
+#
+# Encore un test temporaire: probleme de compilation sur VPP
+# l'optimisation de startvar se passe mal
+#
+#if ($code == 'create_etat0' && $VPP) then
+#   if ( ! -f $libo/libdyn3d.a ) then
+#     echo "Priere de compiler gcm en premier pour des raisons d'optimisation"
+#     \rm $libf/grid/dimensions.h                                             
+#     exit                       
+#   endif
+#   set optim90="$optim90"" -X9 -w"
+#   set opt_link="-L$NCDFLIB -lnetcdf"
+#   touch $LMDGCM/libf/dyn3d/startvar.F
+#   touch $LMDGCM/libf/dyn3d/etat0_netcdf.F
+#endif
+
+#
+# etat0_netcdf a besoin d'info de la physique
+# A revoir
+set include="$include"" -I$libf/phy$physique"
+
+
+#################################################################
+# Execution de la comande make... ENFIN!
+#################################################################
+
+if $VPP then
+  set optim90=" $optim90 -Am -M$libo"
+  set optimtru90="$optim90"
+else if $SUN then
+ set optim90=" $optim90 -M$libo"
+ set optimtru90=" $optimtru90 "
+else if $NEC then
+ set optim90=" $optim90 -I$libo "
+else if $XNEC then
+ set optim90=" $optim90 -I$libo "
+ set optimtru90=" $optimtru90 -I$libo "
+else if $LINUX then
+# Ehouarn : adapt to g95
+ set optim="$optim -fmod=$libo"
+ set optim90="$optim90 -fmod=$libo"
+ set optimtru90="$optimtru90 -fmod=$libo"
+# Ehouarn: remove set mod_loc_dir def below; mod_loc_dir=$localdir (set above)
+# set mod_loc_dir=$libo
+endif
+
+set link="$f90 $optim90"
+
+set ar=ar
+
+if $XNEC then
+  set link="sxld $opt_link"
+  set link="$f90 "
+#  set ar=sxar
+endif
+
+
+cd $localdir
+
+echo $make -f $LMDGCM/makefile \
+OPTION_DEP="$opt_dep" OPTION_LINK="$opt_link" \
+OPTIM="$optim" \
+OPTIM90="$optim90" \
+OPTIMTRU90="$optimtru90" \
+INCLUDE="$include" \
+$filtre \
+LIBO=$libo \
+$dyn \
+$phys \
+DIM=$dimc \
+DYNTYPE="$dyntype" \
+L_ADJNT="$adjnt" \
+LOCAL_DIR="$localdir"  \
+F77="$f77" \
+F90="$f90" \
+OPLINK="$oplink" \
+LINK="$link" \
+GCM="$LMDGCM" \
+MOD_LOC_DIR=$mod_loc_dir \
+MOD_SUFFIX=$mod_suffix \
+AR=$ar \
+PROG=$code
+
+
+$make -f $LMDGCM/makefile \
+OPTION_DEP="$opt_dep" OPTION_LINK="$opt_link" \
+OPTIM="$optim" \
+OPTIM90="$optim90" \
+OPTIMTRU90="$optimtru90" \
+INCLUDE="$include" \
+$filtre \
+LIBO=$libo \
+$dyn \
+$phys \
+DIM=$dimc \
+DYNTYPE="$dyntype" \
+L_ADJNT="$adjnt" \
+LOCAL_DIR="$localdir"  \
+F77="$f77" \
+F90="$f90" \
+OPLINK="$oplink" \
+LINK="$link" \
+GCM="$LMDGCM" \
+MOD_LOC_DIR=$mod_loc_dir \
+MOD_SUFFIX=$mod_suffix \
+AR=$ar \
+PROG=$code
+
+\rm $libf/grid/dimensions.h
Index: trunk/mesoscale/LMDZ.MARS/myGCM/DEFS_JB/callphys.def
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/DEFS_JB/callphys.def	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/DEFS_JB/callphys.def	(revision 54)
@@ -0,0 +1,117 @@
+##General options
+##~~~~~~~~~~~~~~~
+#Run with or without tracer transport ?
+tracer=.true.
+
+#Diurnal cycle ?  if diurnal=False, diurnal averaged solar heating
+diurnal=.true.
+
+#Seasonal cycle ? if season=False, Ls stays constant, to value set in "start"
+season = .true.   
+
+#write some more output on the screen ?
+lwrite = .false.    
+
+#Save statistics in file "stats.nc" ?
+stats =.false.
+
+#Save EOF profiles in file "profiles" for Climate Database?
+calleofdump = .false.
+
+## Dust scenario. Used if the dust is prescribed (i.e. if tracer=F or active=F)
+## ~~~~~~~~~~~~~
+#  =1 Dust opt.deph read in startfi; =2 Viking scenario; =3 MGS scenario,
+#  =4 Mars Year 24 from TES assimilation (same as =24 for now)
+#  =24 Mars Year 24 from TES assimilation (ie: MCD reference case)
+#  =25 Mars Year 25 from TES assimilation (ie: a year with a global dust storm)
+#  =26 Mars Year 26 from TES assimilation
+iaervar = 26 
+# Dust opacity at 700 Pa (when constant, i.e. for the iaervar=1 case)
+tauvis=0.2
+# Dust vertical distribution: 
+# (=0: old distrib. (Pollack90), =1: top set by "topdustref",
+#  =2: Viking scenario; =3 MGS scenario)
+iddist  = 3
+# Dust top altitude (km). (Matters only if iddist=1)
+topdustref = 55.
+
+## Physical Parameterizations :
+## ~~~~~~~~~~~~~~~~~~~~~~~~~~
+# call radiative transfer ?
+callrad   = .true.
+# call NLTE radiative schemes ?   matters only if callrad=T
+callnlte = .true.
+# call CO2 NIR absorption ?   matters only if callrad=T
+callnirco2 = .true.
+# call turbulent vertical diffusion ?
+calldifv  = .true.
+# call convective adjustment ?
+calladj  = .true.
+# call CO2 condensation ?
+callcond =.true.
+# call thermal conduction in the soil ?
+callsoil = .true.
+# call Lott's gravity wave/subgrid topography scheme ?
+calllott  = .true.
+
+## Radiative transfer options :
+## ~~~~~~~~~~~~~~~~~~~~~~~~~~
+# the rad.transfer is computed every "iradia" physical timestep
+iradia    = 1
+# Output of the exchange coefficient mattrix ? for diagnostic only
+callg2d  = .false.
+# Rayleigh scattering : (should be .false. for now)
+rayleigh  = .false.
+
+## Tracer (dust water, ice and/or chemical species) options (used if tracer=T):
+## ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+# DUST: Transported dust ? (if >0, use 'dustbin' dust bins)
+dustbin    = 2
+# DUST: Radiatively active dust ? (matters if dustbin>0)
+active  = .true.
+# DUST: use mass and number mixing ratios to predict dust size ?
+doubleq   = .true.
+# DUST: use a small population of dust particules (submicron dust)?
+submicron = .false.
+# DUST: lifted by GCM surface winds ?
+lifting = .true.
+# DUST: lifted by dust devils ?
+callddevil = .false.
+# DUST: Scavenging by CO2 snowfall ?
+scavenging = .false.
+# DUST/WATERICE: Gravitationnal sedimentation ?
+sedimentation = .true.
+# WATERICE: Radiatively active transported atmospheric water ice ?
+activice   = .true.
+# WATER: Compute water cycle
+water = .true.
+# WATER: current permanent caps at both poles. True IS RECOMMENDED
+#        (with .true., North cap is a source of water and South pole
+#         is a cold trap)
+caps  = .true.
+# PHOTOCHEMISTRY: include chemical species 
+photochem  = .false.
+
+## Thermospheric options (relevant if tracer=T) :
+##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+# call thermosphere ?
+callthermos = .false.
+# WATER: included without cycle (only if water=.false.)
+thermoswater = .false.
+# call thermal conduction ? (only if callthermos=.true.)
+callconduct = .false.
+# call EUV heating ? (only if callthermos=.true.)
+calleuv=.false.
+# call molecular viscosity ? (only if callthermos=.true.)
+callmolvis = .false.
+# call molecular diffusion ? (only if callthermos=.true.)
+callmoldiff = .false.
+# call thermospheric photochemistry ? (only if callthermos=.true.)
+thermochem = .false.
+# date for solar flux calculation: (1985 < date < 2002)
+## (Solar min=1996.4 ave=1993.4 max=1990.6)
+solarcondate = 1993.4
+
+## albedoes of CO2 ice caps from TES ?
+TESicealbedo=.true.
+
Index: trunk/mesoscale/LMDZ.MARS/myGCM/DEFS_JB/run.def
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/DEFS_JB/run.def	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/DEFS_JB/run.def	(revision 54)
@@ -0,0 +1,112 @@
+#
+#----------------------------------------------------------------------- 
+#Parametres de controle du run:                                          
+#------------------------------                                          
+
+# Nombre de jours d'integration                                         
+     nday=3
+
+# nombre de pas par jour (multiple de iperiod) ( ici pour  dt = 1 min )      
+ day_step = 960
+
+# periode pour le pas Matsuno (en pas)                                  
+  iperiod=5
+
+# periode de sortie des variables de controle (en pas)                  
+  iconser=120
+
+# periode d'ecriture du fichier histoire (en jour)                      
+    iecri=100
+
+# periode de stockage fichier histmoy (en jour)                         
+ periodav=60.
+
+# periode de la dissipation (en pas)                                    
+  idissip=1
+
+# choix de l'operateur de dissipation (star ou  non star )              
+ lstardis=.true.
+
+# avec ou sans coordonnee hybrides
+ hybrid=.true.
+
+# nombre d'iterations de l'operateur de dissipation   gradiv            
+nitergdiv=1
+
+# nombre d'iterations de l'operateur de dissipation  nxgradrot          
+nitergrot=2
+
+# nombre d'iterations de l'operateur de dissipation  divgrad            
+   niterh=2
+
+# temps de dissipation des plus petites long.d ondes pour u,v (gradiv)  
+ tetagdiv= 1800.
+
+# temps de dissipation des plus petites long.d ondes pour u,v(nxgradrot)
+ tetagrot=2500.
+
+# temps de dissipation des plus petites long.d ondes pour  h ( divgrad) 
+ tetatemp=2500.
+
+# coefficient pour gamdissip                                            
+  coefdis=0.
+
+# choix du shema d'integration temporelle (Matsuno ou Matsuno-leapfrog) 
+  purmats=.false.
+
+# avec ou sans physique                                                 
+   physic=.true.
+
+# periode de la physique (en pas)                                       
+  iphysiq=20
+
+# choix d'une grille reguliere                                          
+  grireg=.true.
+
+# frequence (en pas) de l'ecriture du fichier diagfi               
+ ecritphy=40
+
+# longitude en degres du centre du zoom                                 
+   clon=63.
+
+# latitude en degres du centre du zoom                                  
+   clat=0.
+
+# facteur de grossissement du zoom,selon longitude                      
+  grossismx=1.
+
+# facteur de grossissement du zoom ,selon latitude                      
+ grossismy=1.
+
+#  Fonction  f(y)  hyperbolique  si = .true.  , sinon  sinusoidale         
+  fxyhypb=.false.
+
+# extension en longitude  de la zone du zoom  ( fraction de la zone totale)
+   dzoomx= 0.
+
+# extension en latitude de la zone  du zoom  ( fraction de la zone totale)
+   dzoomy=0.
+
+#  raideur du zoom en  X
+    taux=2.
+
+#  raideur du zoom en  Y
+    tauy=2.
+
+#  Fonction  f(y) avec y = Sin(latit.) si = .TRUE. ,  Sinon  y = latit.
+  ysinus= .false.
+
+# Avec sponge layer
+  callsponge  = .true.
+ 
+# Sponge:  mode0(u=v=0), mode1(u=umoy,v=0), mode2(u=umoy,v=vmoy)
+  mode_sponge= 2
+
+# Sponge:  hauteur de sponge (km)
+  hsponge= 90
+
+# Sponge:  tetasponge (secondes)
+  tetasponge = 50000
+
+# some definitions for the physics, in file 'callphys.def'
+INCLUDEDEF=callphys.def
Index: trunk/mesoscale/LMDZ.MARS/myGCM/DEFS_JB/traceur.def
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/DEFS_JB/traceur.def	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/DEFS_JB/traceur.def	(revision 54)
@@ -0,0 +1,6 @@
+5
+co2
+dust_mass
+dust_number
+h2o_vap
+h2o_ice
Index: trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def	(revision 54)
@@ -0,0 +1,1 @@
+link DEFS_JB/callphys.def
Index: trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def.csttau
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def.csttau	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def.csttau	(revision 54)
@@ -0,0 +1,1 @@
+link ../in_lmdz_mars_newphys/myGCM/callphys.def.csttau
Index: trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def.water
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def.water	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def.water	(revision 54)
@@ -0,0 +1,1 @@
+link ../in_lmdz_mars_newphys/myGCM/callphys.def.water
Index: trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def_sav
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def_sav	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/callphys.def_sav	(revision 54)
@@ -0,0 +1,1 @@
+link ../in_lmdz_mars_newphys/myGCM/callphys.def
Index: trunk/mesoscale/LMDZ.MARS/myGCM/diagfi.nc
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/diagfi.nc	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/diagfi.nc	(revision 54)
@@ -0,0 +1,1 @@
+link ../../TMPDIR/GCMINI/diagfi.nc
Index: trunk/mesoscale/LMDZ.MARS/myGCM/gcm.e
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/gcm.e	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/gcm.e	(revision 54)
@@ -0,0 +1,1 @@
+link ../gcm.e
Index: trunk/mesoscale/LMDZ.MARS/myGCM/launch_gcm
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/launch_gcm	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/launch_gcm	(revision 54)
@@ -0,0 +1,1 @@
+link ../in_lmdz_mars_newphys/myGCM/launch_gcm
Index: trunk/mesoscale/LMDZ.MARS/myGCM/run.def
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/run.def	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/run.def	(revision 54)
@@ -0,0 +1,1 @@
+link DEFS_JB/run.def
Index: trunk/mesoscale/LMDZ.MARS/myGCM/run.def.64x48x25
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/run.def.64x48x25	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/run.def.64x48x25	(revision 54)
@@ -0,0 +1,1 @@
+link ../in_lmdz_mars_newphys/myGCM/run.def.64x48x25
Index: trunk/mesoscale/LMDZ.MARS/myGCM/run.def_sav
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/run.def_sav	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/run.def_sav	(revision 54)
@@ -0,0 +1,1 @@
+link ../in_lmdz_mars_newphys/myGCM/run.def
Index: trunk/mesoscale/LMDZ.MARS/myGCM/run_mcd_3days
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/run_mcd_3days	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/run_mcd_3days	(revision 54)
@@ -0,0 +1,1 @@
+link ../in_lmdz_mars_newphys/myGCM/run_mcd_3days
Index: trunk/mesoscale/LMDZ.MARS/myGCM/start.nc
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/start.nc	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/start.nc	(revision 54)
@@ -0,0 +1,1 @@
+link DEFS_JB/start40.nc
Index: trunk/mesoscale/LMDZ.MARS/myGCM/startfi.nc
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/startfi.nc	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/startfi.nc	(revision 54)
@@ -0,0 +1,1 @@
+link DEFS_JB/startfi40.nc
Index: trunk/mesoscale/LMDZ.MARS/myGCM/traceur.def
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/traceur.def	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/traceur.def	(revision 54)
@@ -0,0 +1,1 @@
+link DEFS_JB/traceur.def
Index: trunk/mesoscale/LMDZ.MARS/myGCM/traceur.def_sav
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/traceur.def_sav	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/traceur.def_sav	(revision 54)
@@ -0,0 +1,1 @@
+link ../in_lmdz_mars_newphys/myGCM/traceur.def
Index: trunk/mesoscale/LMDZ.MARS/myGCM/what_I_need.info
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/what_I_need.info	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/what_I_need.info	(revision 54)
@@ -0,0 +1,1 @@
+link ../in_lmdz_mars_newphys/myGCM/what_I_need.info
Index: trunk/mesoscale/LMDZ.MARS/myGCM/z2sig.def
===================================================================
--- trunk/mesoscale/LMDZ.MARS/myGCM/z2sig.def	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/myGCM/z2sig.def	(revision 54)
@@ -0,0 +1,1 @@
+link ../deftank_gcm/z2sig.def
Index: trunk/mesoscale/LMDZ.MARS/writediagfi.info
===================================================================
--- trunk/mesoscale/LMDZ.MARS/writediagfi.info	(revision 54)
+++ trunk/mesoscale/LMDZ.MARS/writediagfi.info	(revision 54)
@@ -0,0 +1,1 @@
+link in_lmdz_mars_newphys/writediagfi.info
Index: trunk/mesoscale/LMD_LES_MARS/modif_mars/module_first_rk_step_part1.F
===================================================================
--- trunk/mesoscale/LMD_LES_MARS/modif_mars/module_first_rk_step_part1.F	(revision 51)
+++ trunk/mesoscale/LMD_LES_MARS/modif_mars/module_first_rk_step_part1.F	(revision 54)
@@ -616,4 +616,8 @@
      &        ,MARS_CICE=grid%mars_cice,MARS_EMISS=grid%mars_emiss        &
      &        ,MARS_TSOIL=grid%mars_tsoil                                 &
+#ifdef NEWPHYS
+     &        ,MARS_ISOIL=grid%mars_isoil       &
+     &        ,MARS_DSOIL=grid%mars_dsoil       &
+#endif
      &        ,MARS_GW=grid%mars_gw                                       &
      &        ,NUM_SOIL_LAYERS=config_flags%num_soil_layers               &
Index: trunk/mesoscale/LMD_LES_MARS/modif_mars/module_model_constants.F
===================================================================
--- trunk/mesoscale/LMD_LES_MARS/modif_mars/module_model_constants.F	(revision 51)
+++ trunk/mesoscale/LMD_LES_MARS/modif_mars/module_model_constants.F	(revision 54)
@@ -19,14 +19,14 @@
 
 #if ( NMM_CORE == 1 )
-   REAL    , PARAMETER :: r_d          = 192.  ! gas constant m2 s-2 K-1
-   REAL    , PARAMETER :: cp           = 844.6   ! r= 8.314511E+0 *1000.E+0/mugaz
+   REAL    , PARAMETER :: r_d          = 191.  ! gas constant m2 s-2 K-1
+   REAL    , PARAMETER :: cp           = 744.5   ! r= 8.314511E+0 *1000.E+0/mugaz
 #else
-   REAL    , PARAMETER :: r_d          = 192.
-   REAL    , PARAMETER :: cp           = 844.6
+   REAL    , PARAMETER :: r_d          = 191.
+   REAL    , PARAMETER :: cp           = 744.5
 #endif
 
    REAL    , PARAMETER :: r_v          = 461.6  ! gas constant for water vapor
    REAL    , PARAMETER :: cv           = cp-r_d
-   REAL    , PARAMETER :: cpv          = 4.4*r_v  ! cp / R = 4.4 ?? ! rcp = r/cp  ~0.256793
+   REAL    , PARAMETER :: cpv          = 3.9*r_v  ! cp / R = 4.4 ?? ! rcp = r/cp  ~0.256793
    REAL    , PARAMETER :: cvv          = cpv-r_v
    REAL    , PARAMETER :: cvpm         = -cv/cp
@@ -73,4 +73,7 @@
    REAL    , PARAMETER :: wdtemisiceN   = 0.4  ! time scale for snow metamorphism (north) !2  
    REAL    , PARAMETER :: wdtemisiceS   = 0.4  ! time scale for snow metamorphism (south) !2 
+#ifdef NEWPHYS
+   REAL    , PARAMETER :: wvolcapa      = 1.e6       ! volumetric capacity of soil (new soil model)
+#endif
 !
 !-------------------------------
Index: trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/compile
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/compile	(revision 51)
+++ trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/compile	(revision 54)
@@ -24,4 +24,11 @@
 #nz=32
 ##################
+tracers=5
+nx=64
+ny=48
+nz=25
+##################
+
+
 
 
Index: trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/callphys.def
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/callphys.def	(revision 51)
+++ trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/callphys.def	(revision 54)
@@ -3,5 +3,4 @@
 #Run with or without tracer transport ?
 tracer=.false.
-##true???
 
 #Diurnal cycle ?  if diurnal=False, diurnal averaged solar heating
Index: trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/callphys.def.csttau
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/callphys.def.csttau	(revision 54)
+++ trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/callphys.def.csttau	(revision 54)
@@ -0,0 +1,122 @@
+##General options
+##~~~~~~~~~~~~~~~
+#Run with or without tracer transport ?
+tracer=.true.
+
+#Diurnal cycle ?  if diurnal=False, diurnal averaged solar heating
+diurnal=.true.
+
+#Seasonal cycle ? if season=False, Ls stays constant, to value set in "start"
+season = .true.   
+
+#write some more output on the screen ?
+lwrite = .false.    
+
+#Save statistics in file "stats.nc" ?
+stats =.false.
+
+#Save EOF profiles in file "profiles" for Climate Database?
+calleofdump = .false.
+
+## Dust scenario. Used if the dust is prescribed (i.e. if tracer=F or active=F)
+## ~~~~~~~~~~~~~
+#  =1 Dust opt.deph read in startfi; =2 Viking scenario; =3 MGS scenario,
+#  =4 Mars Year 24 from TES assimilation (same as =24 for now)
+#  =24 Mars Year 24 from TES assimilation (ie: MCD reference case)
+#  =25 Mars Year 25 from TES assimilation (ie: a year with a global dust storm)
+#  =26 Mars Year 26 from TES assimilation
+iaervar = 1 
+# Dust opacity at 700 Pa (when constant, i.e. for the iaervar=1 case)
+tauvis=0.5
+# Dust vertical distribution: 
+# (=0: old distrib. (Pollack90), =1: top set by "topdustref",
+#  =2: Viking scenario; =3 MGS scenario)
+iddist  = 3
+# Dust top altitude (km). (Matters only if iddist=1)
+topdustref = 55.
+
+## Physical Parameterizations :
+## ~~~~~~~~~~~~~~~~~~~~~~~~~~
+# call radiative transfer ?
+callrad   = .true.
+# call NLTE radiative schemes ?   matters only if callrad=T
+callnlte = .true.
+# call CO2 NIR absorption ?   matters only if callrad=T
+callnirco2 = .true.
+# call turbulent vertical diffusion ?
+calldifv  = .true.
+# call convective adjustment ?
+calladj  = .true.
+# call CO2 condensation ?
+callcond =.true.
+# call thermal conduction in the soil ?
+callsoil = .true.
+# call Lott's gravity wave/subgrid topography scheme ?
+calllott  = .true.
+
+## Radiative transfer options :
+## ~~~~~~~~~~~~~~~~~~~~~~~~~~
+# the rad.transfer is computed every "iradia" physical timestep
+iradia    = 1
+# Output of the exchange coefficient mattrix ? for diagnostic only
+callg2d  = .false.
+# Rayleigh scattering : (should be .false. for now)
+rayleigh  = .false.
+
+## Tracer (dust water, ice and/or chemical species) options (used if tracer=T):
+## ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+# DUST: Transported dust ? (if >0, use 'dustbin' dust bins)
+dustbin    = 0
+# DUST: Radiatively active dust ? (matters if dustbin>0)
+active  = .false.
+# DUST: use mass and number mixing ratios to predict dust size ?
+# (must also have dustbin=1)
+doubleq   = .false.
+# DUST: lifted by GCM surface winds ?
+lifting = .false.
+# DUST: lifted by dust devils ?
+callddevil = .false.
+# DUST: Scavenging by CO2 snowfall ?
+scavenging = .false.
+# DUST/WATERICE: Gravitationnal sedimentation ?
+sedimentation = .false.
+# WATERICE: Radiatively active transported atmospheric water ice ?
+activice   = .false.
+# WATER: Compute water cycle
+water = .false.
+# WATER: current permanent caps at both poles. True IS RECOMMENDED
+#        (with .true., North cap is a source of water and South pole
+#         is a cold trap)
+caps  = .true.
+# PHOTOCHEMISTRY: include chemical species 
+photochem  = .false.
+
+## Thermospheric options (relevant if tracer=T) :
+##~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+# call thermosphere ?
+callthermos = .false.
+# WATER: included without cycle (only if water=.false.)
+thermoswater = .false.
+# call thermal conduction ? (only if callthermos=.true.)
+callconduct = .false.
+# call EUV heating ? (only if callthermos=.true.)
+calleuv=.false.
+# call molecular viscosity ? (only if callthermos=.true.)
+callmolvis = .false.
+# call molecular diffusion ? (only if callthermos=.true.)
+callmoldiff = .false.
+# call thermospheric photochemistry ? (only if callthermos=.true.)
+thermochem = .false.
+# date for solar flux calculation: (1985 < date < 2002)
+## (Solar min=1996.4 ave=1993.4 max=1990.6)
+solarcondate = 1993.4
+       
+
+## albedoes of CO2 ice caps from TES ?
+TESicealbedo=.true.
+## Coefficient for Northern cap albedoes
+TESice_Ncoef=1.4
+## Coefficient for Southern cap albedoes
+TESice_Scoef=1.4
+
+
Index: trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/run.def
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/run.def	(revision 54)
+++ trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/run.def	(revision 54)
@@ -0,0 +1,112 @@
+#
+#----------------------------------------------------------------------- 
+#Parametres de controle du run:                                          
+#------------------------------                                          
+
+# Nombre de jours d'integration                                         
+     nday=20
+
+# nombre de pas par jour (multiple de iperiod) ( ici pour  dt = 1 min )      
+ day_step = 960
+
+# periode pour le pas Matsuno (en pas)                                  
+  iperiod=5
+
+# periode de sortie des variables de controle (en pas)                  
+  iconser=120
+
+# periode d'ecriture du fichier histoire (en jour)                      
+    iecri=100
+
+# periode de stockage fichier histmoy (en jour)                         
+ periodav=60.
+
+# periode de la dissipation (en pas)                                    
+  idissip=5
+
+# choix de l'operateur de dissipation (star ou  non star )              
+ lstardis=.true.
+
+# avec ou sans coordonnee hybrides
+ hybrid=.true.
+
+# nombre d'iterations de l'operateur de dissipation   gradiv            
+nitergdiv=1
+
+# nombre d'iterations de l'operateur de dissipation  nxgradrot          
+nitergrot=2
+
+# nombre d'iterations de l'operateur de dissipation  divgrad            
+   niterh=2
+
+# temps de dissipation des plus petites long.d ondes pour u,v (gradiv)  
+ tetagdiv= 3500.
+
+# temps de dissipation des plus petites long.d ondes pour u,v(nxgradrot)
+ tetagrot=4500.
+
+# temps de dissipation des plus petites long.d ondes pour  h ( divgrad) 
+ tetatemp=4500.
+
+# coefficient pour gamdissip                                            
+  coefdis=0.
+
+# choix du shema d'integration temporelle (Matsuno ou Matsuno-leapfrog) 
+  purmats=.false.
+
+# avec ou sans physique                                                 
+   physic=.true.
+
+# periode de la physique (en pas)                                       
+  iphysiq=20
+
+# choix d'une grille reguliere                                          
+  grireg=.true.
+
+# frequence (en pas) de l'ecriture du fichier diagfi               
+ ecritphy=40
+
+# longitude en degres du centre du zoom                                 
+   clon=63.
+
+# latitude en degres du centre du zoom                                  
+   clat=0.
+
+# facteur de grossissement du zoom,selon longitude                      
+  grossismx=1.
+
+# facteur de grossissement du zoom ,selon latitude                      
+ grossismy=1.
+
+#  Fonction  f(y)  hyperbolique  si = .true.  , sinon  sinusoidale         
+  fxyhypb=.false.
+
+# extension en longitude  de la zone du zoom  ( fraction de la zone totale)
+   dzoomx= 0.
+
+# extension en latitude de la zone  du zoom  ( fraction de la zone totale)
+   dzoomy=0.
+
+#  raideur du zoom en  X
+    taux=2.
+
+#  raideur du zoom en  Y
+    tauy=2.
+
+#  Fonction  f(y) avec y = Sin(latit.) si = .TRUE. ,  Sinon  y = latit.
+  ysinus= .false.
+
+# Avec sponge layer
+  callsponge  = .true.
+ 
+# Sponge:  mode0(u=v=0), mode1(u=umoy,v=0), mode2(u=umoy,v=vmoy)
+  mode_sponge= 2
+
+# Sponge:  hauteur de sponge (km)
+  hsponge= 130
+
+# Sponge:  tetasponge (secondes)
+  tetasponge = 50000
+
+# some definitions for the physics, in file 'callphys.def'
+INCLUDEDEF=callphys.def
Index: trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/traceur.def
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/traceur.def	(revision 54)
+++ trunk/mesoscale/LMD_MM_MARS/SIMU/in_lmdz_mars_newphys/myGCM/traceur.def	(revision 54)
@@ -0,0 +1,2 @@
+1
+co2
Index: trunk/mesoscale/LMD_MM_MARS/SIMU/runmeso
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SIMU/runmeso	(revision 51)
+++ trunk/mesoscale/LMD_MM_MARS/SIMU/runmeso	(revision 54)
@@ -165,4 +165,5 @@
       0) tracers='1' ;;
       1) tracers='2' ;;
+     11) tracers='5' ;;
       *) tracers='1' ;;
   esac
Index: trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/Registry/Registry.EM
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/Registry/Registry.EM	(revision 51)
+++ trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/Registry/Registry.EM	(revision 54)
@@ -113,4 +113,10 @@
 state  real  QH2O_ICE  ikjftb  scalar  1  -  i01rhusdf=(bdy_interp:dt) "QH2O_ICE"   "Water ice mixing ratio"     "kg kg-1"
 state  real  QDUST     ikjftb  scalar  1  -  i01rhusdf=(bdy_interp:dt) "QDUST"      "Dust mixing ratio"          "kg kg-1"
+####
+state  real  Qco2        ikjftb  scalar  1  -  i01rhusdf=(bdy_interp:dt)  "Qco2"         "CO2 mixing ratio"          "kg kg-1"
+#state  real  Qdust_mass  ikjftb  scalar  1  -  i01rhusdf=(bdy_interp:dt)  "Qdust_mass"   "dust_mass mixing ratio"    "kg kg-1"
+state  real  Qdust_number  ikjftb  scalar  1  -  i01rhusdf=(bdy_interp:dt)  "Qdust_number"   "dust_number mixing ratio"    "kg kg-1"
+#state  real  Qh2o_vap  ikjftb  scalar  1  -  i01rhusdf=(bdy_interp:dt)  "Qh2o_vap"   "h2o_vap mixing ratio"    "kg kg-1"
+#state  real  Qh2o_ice  ikjftb  scalar  1  -  i01rhusdf=(bdy_interp:dt)  "Qh2o_ice"   "h2o_ice mixing ratio"    "kg kg-1"
 ####
 ####
@@ -1324,5 +1330,7 @@
 package   nowater      mars==0                      -              moist:qv
 package   water        mars==1                      -              moist:qv;scalar:qh2o,qh2o_ice
-package   dust         mars==2                      -              moist:qv;scalar:qdust 
+package   dust         mars==2                      -              moist:qv;scalar:qdust
+#package   newwater     mars==11                     -              scalar:qco2,qh2o,qh2o_ice
+package   newwater     mars==11                     -              scalar:qco2,qh2o,qh2o_ice,qdust,qdust_number
 ##### MARS OPTIONS
 ##### MARS OPTIONS
Index: trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/dyn_em/module_initialize_real.F
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/dyn_em/module_initialize_real.F	(revision 51)
+++ trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/dyn_em/module_initialize_real.F	(revision 54)
@@ -1253,5 +1253,5 @@
                             its , ite , jts , jte , kts , kte )
 
-if (config_flags%mars == 1) then
+if ( (config_flags%mars == 1) .OR. (config_flags%mars == 11) ) then
 !if (size(scalar(0,0,0,:)) > 2) then
 !! autre possibilite: activer les flags dans metgrid.tbl
Index: trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/dyn_em/solve_em.F
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/dyn_em/solve_em.F	(revision 51)
+++ trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/dyn_em/solve_em.F	(revision 54)
@@ -2692,5 +2692,7 @@
 
 !!!!****MARS: pass water vapor at boundaries (and all other tracers if nested domains)
-IF ( ( (config_flags%mars .eq. 1) .AND. (is .EQ. 2) ) .OR. config_flags%nested ) THEN
+IF (  ( (config_flags%mars .eq. 1) .AND. (is .EQ. 2) ) &
+ .OR. ( (config_flags%mars .eq. 11) .AND. (is .EQ. 2) ) &
+ .OR. config_flags%nested ) THEN
          CALL relax_bdy_scalar ( scalar_tend(ims,kms,jms,is),            &
                                  scalar(ims,kms,jms,is),  grid%em_mut,         &
@@ -2750,5 +2752,7 @@
 !!!!****MARS: if parent domain and any other tracer than water vapor, no bdy and simple zero flux bdy condition
 IF (  config_flags%specified  ) THEN
-IF ( (config_flags%mars .ne. 1) .OR. ((config_flags%mars .eq. 1) .and. (is .ne. 2)) ) THEN
+IF (    (config_flags%mars .ne. 1) &
+   .OR. ((config_flags%mars .eq. 11) .and. (is .ne. 2)) &
+   .OR. ((config_flags%mars .eq. 1 ) .and. (is .ne. 2)) ) THEN
            CALL flow_dep_bdy  ( scalar(ims,kms,jms,is),     &
                                 grid%em_ru_m, grid%em_rv_m, config_flags,   &
Index: trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/main/real_em.F
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/main/real_em.F	(revision 51)
+++ trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/main/real_em.F	(revision 54)
@@ -679,5 +679,5 @@
 !                    ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
 !!!!MARS : pb a regler avec P_QH2O qui n est pas reconnu si QH2O n est pas dans active dans le registry 
-IF ( config_flags%mars .eq. 1 ) THEN
+IF ( (config_flags%mars .eq. 1) .OR. (config_flags%mars == 11) ) THEN
       CALL couple ( grid%em_mu_2 , grid%em_mub , qbdy3dtemp1 , grid%scalar(:,:,:,2) , 't' , grid%msft , &
                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
@@ -742,5 +742,5 @@
 !                                                                 ips , ipe , jps , jpe , kps , kpe )
 !!!!MARS
-IF ( config_flags%mars .eq. 1 ) THEN
+IF ( (config_flags%mars .eq. 1) .OR. (config_flags%mars == 11) ) THEN
       CALL stuff_bdy     ( qbdy3dtemp1 , grid%scalar_b(:,:,:,:,2)   , 'T' , ijds , ijde , spec_bdy_width      , &
                                                                  ids , ide , jds , jde , kds , kde , &
@@ -808,5 +808,5 @@
 !                    ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
 !!!!MARS
-IF ( config_flags%mars .eq. 1 ) THEN
+IF ( (config_flags%mars .eq. 1) .OR. (config_flags%mars == 11) ) THEN
       CALL couple ( grid%em_mu_2 , grid%em_mub , qbdy3dtemp2 , grid%scalar(:,:,:,2) , 't' , grid%msft , &
                     ids, ide, jds, jde, kds, kde, ims, ime, jms, jme, kms, kme, ips, ipe, jps, jpe, kps, kpe )
@@ -874,5 +874,5 @@
 !                                                            ips , ipe , jps , jpe , kps , kpe )
 !!!!MARS
-IF ( config_flags%mars .eq. 1 ) THEN
+IF ( (config_flags%mars .eq. 1) .OR. (config_flags%mars == 11) ) THEN
       CALL stuff_bdytend ( qbdy3dtemp2 , qbdy3dtemp1 , REAL(interval_seconds) , grid%scalar_bt(:,:,:,:,2) , 'T' , &
                                                             ijds , ijde , spec_bdy_width      , &
@@ -1028,5 +1028,5 @@
 !                                                                    ips , ipe , jps , jpe , kps , kpe )
 !!!!MARS
-IF ( config_flags%mars .eq. 1 ) THEN
+IF ( (config_flags%mars .eq. 1) .OR. (config_flags%mars == 11) ) THEN
         CALL stuff_bdy     ( qbdy3dtemp1 , grid%scalar_b(:,:,:,:,2)   , 'T', ijds , ijde , spec_bdy_width      , &
                                                                     ids , ide , jds , jde , kds , kde , &
Index: trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd/libf/phymars/meso_physiq.F
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd/libf/phymars/meso_physiq.F	(revision 51)
+++ trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/mars_lmd/libf/phymars/meso_physiq.F	(revision 54)
@@ -873,22 +873,22 @@
 
 
-c-----------------------------------------------------------------------
-c    3. Gravity wave and subgrid scale topography drag :
-c    -------------------------------------------------
-
-
-      IF(calllott)THEN
-
-        CALL calldrag_noro(ngrid,nlayer,ptimestep,
-     &                 pplay,pplev,pt,pu,pv,zdtgw,zdugw,zdvgw)
- 
-        DO l=1,nlayer
-          DO ig=1,ngrid
-            pdv(ig,l)=pdv(ig,l)+zdvgw(ig,l)
-            pdu(ig,l)=pdu(ig,l)+zdugw(ig,l)
-            pdt(ig,l)=pdt(ig,l)+zdtgw(ig,l)
-          ENDDO
-        ENDDO
-      ENDIF
+!c-----------------------------------------------------------------------
+!c    3. Gravity wave and subgrid scale topography drag :
+!c    -------------------------------------------------
+!
+!
+!      IF(calllott)THEN
+!
+!        CALL calldrag_noro(ngrid,nlayer,ptimestep,
+!     &                 pplay,pplev,pt,pu,pv,zdtgw,zdugw,zdvgw)
+! 
+!        DO l=1,nlayer
+!          DO ig=1,ngrid
+!            pdv(ig,l)=pdv(ig,l)+zdvgw(ig,l)
+!            pdu(ig,l)=pdu(ig,l)+zdugw(ig,l)
+!            pdt(ig,l)=pdt(ig,l)+zdtgw(ig,l)
+!          ENDDO
+!        ENDDO
+!      ENDIF
 
 c-----------------------------------------------------------------------
@@ -1153,44 +1153,44 @@
 c     ------------------
 
-c        --------------
-c        photochemistry :
-c        --------------
-         IF(photochem .or. thermochem) then
-          call calchim(ptimestep,pplay,pplev,pt,pdt,dist_sol,mu0,
-     .      zzlay,zday,pq,pdq,zdqchim,zdqschim,zdqcloud,zdqscloud)
-           
-c Photochemistry includes condensation of H2O2
-
-           do iq=nqchem_min,nq
-            if (noms(iq).eq."h2o2") then
-             DO l=1,nlayer
-               DO ig=1,ngrid
-                    pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqchim(ig,l,iq)
-                    pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqcloud(ig,l,iq)
-               ENDDO
-             ENDDO
-            else
-             DO l=1,nlayer
-               DO ig=1,ngrid
-                    pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqchim(ig,l,iq)
-               ENDDO
-             ENDDO
-            endif
-           ENDDO
-           do iq=nqchem_min,nq
-            if (noms(iq).eq."h2o2") then
-               DO ig=1,ngrid
-                    dqsurf(ig,iq)= dqsurf(ig,iq) + zdqschim(ig,iq)
-                    dqsurf(ig,iq)= dqsurf(ig,iq) + zdqscloud(ig,iq)
-               ENDDO
-            else
-               DO ig=1,ngrid
-                    dqsurf(ig,iq)= dqsurf(ig,iq) + zdqschim(ig,iq)
-               ENDDO
-            endif
-           ENDDO
-
-         END IF  ! (photochem.or.thermochem)
-c        print*,'photochem ok'
+!c        --------------
+!c        photochemistry :
+!c        --------------
+!         IF(photochem .or. thermochem) then
+!          call calchim(ptimestep,pplay,pplev,pt,pdt,dist_sol,mu0,
+!     .      zzlay,zday,pq,pdq,zdqchim,zdqschim,zdqcloud,zdqscloud)
+!           
+!c Photochemistry includes condensation of H2O2
+!
+!           do iq=nqchem_min,nq
+!            if (noms(iq).eq."h2o2") then
+!             DO l=1,nlayer
+!               DO ig=1,ngrid
+!                    pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqchim(ig,l,iq)
+!                    pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqcloud(ig,l,iq)
+!               ENDDO
+!             ENDDO
+!            else
+!             DO l=1,nlayer
+!               DO ig=1,ngrid
+!                    pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqchim(ig,l,iq)
+!               ENDDO
+!             ENDDO
+!            endif
+!           ENDDO
+!           do iq=nqchem_min,nq
+!            if (noms(iq).eq."h2o2") then
+!               DO ig=1,ngrid
+!                    dqsurf(ig,iq)= dqsurf(ig,iq) + zdqschim(ig,iq)
+!                    dqsurf(ig,iq)= dqsurf(ig,iq) + zdqscloud(ig,iq)
+!               ENDDO
+!            else
+!               DO ig=1,ngrid
+!                    dqsurf(ig,iq)= dqsurf(ig,iq) + zdqschim(ig,iq)
+!               ENDDO
+!            endif
+!           ENDDO
+!
+!         END IF  ! (photochem.or.thermochem)
+!c        print*,'photochem ok'
 
 c   7c. Aerosol particles
@@ -1274,32 +1274,33 @@
 
 
-c-----------------------------------------------------------------------
-c   8.5 THERMOSPHERE CALCULATION
-c-----------------------------------------------------------------------
-
-      if (callthermos) then
-        call thermosphere(pplev,pplay,dist_sol,
-     $     mu0,ptimestep,ptime,zday,tsurf,zzlev,zzlay,
-     &     pt,pq,pu,pv,pdt,pdq,
-     $     zdteuv,zdtconduc,zdumolvis,zdvmolvis,zdqmoldiff)
-c           do iq=nqchem_min,nq
-c           write(*,*) 'thermo iq,pq',iq,pq(690,1,iq)
-c           enddo
-
-        DO l=1,nlayer
-          DO ig=1,ngrid
-            dtrad(ig,l)=dtrad(ig,l)+zdteuv(ig,l)
-            pdt(ig,l)=pdt(ig,l)+zdtconduc(ig,l)
-     &                         +zdteuv(ig,l)
-            pdv(ig,l)=pdv(ig,l)+zdvmolvis(ig,l)
-            pdu(ig,l)=pdu(ig,l)+zdumolvis(ig,l)
-            DO iq=1, nq
-              pdq(ig,l,iq)=pdq(ig,l,iq)+zdqmoldiff(ig,l,iq)
-            ENDDO
-          ENDDO
-        ENDDO
-
-
-      endif
+!c-----------------------------------------------------------------------
+!c   8.5 THERMOSPHERE CALCULATION
+!c-----------------------------------------------------------------------
+!
+!      if (callthermos) then
+!        call thermosphere(pplev,pplay,dist_sol,
+!     $     mu0,ptimestep,ptime,zday,tsurf,zzlev,zzlay,
+!     &     pt,pq,pu,pv,pdt,pdq,
+!     $     zdteuv,zdtconduc,zdumolvis,zdvmolvis,zdqmoldiff)
+!c           do iq=nqchem_min,nq
+!c           write(*,*) 'thermo iq,pq',iq,pq(690,1,iq)
+!c           enddo
+!
+!        DO l=1,nlayer
+!          DO ig=1,ngrid
+!            dtrad(ig,l)=dtrad(ig,l)+zdteuv(ig,l)
+!            pdt(ig,l)=pdt(ig,l)+zdtconduc(ig,l)
+!     &                         +zdteuv(ig,l)
+!            pdv(ig,l)=pdv(ig,l)+zdvmolvis(ig,l)
+!            pdu(ig,l)=pdu(ig,l)+zdumolvis(ig,l)
+!            DO iq=1, nq
+!              pdq(ig,l,iq)=pdq(ig,l,iq)+zdqmoldiff(ig,l,iq)
+!            ENDDO
+!          ENDDO
+!        ENDDO
+!
+!
+!      endif
+
 c-----------------------------------------------------------------------
 c   8. Surface  and sub-surface soil temperature
Index: trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/phys/module_lmd_driver.F
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/phys/module_lmd_driver.F	(revision 51)
+++ trunk/mesoscale/LMD_MM_MARS/SRC/WRFV2/phys/module_lmd_driver.F	(revision 54)
@@ -697,15 +697,20 @@
    tsoil_val = tsoil_val*0. + tsurf_val
   ENDIF
+#ifdef NEWPHYS
+  isoil_val(:)=MARS_ISOIL(i,:,j)
+  dsoil_val(:)=MARS_DSOIL(i,:,j)
+#endif
 ELSE
   IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** IDEALIZED SIMULATION tsoil is set to tsurf'
   do k=1,nsoil
    tsoil_val(k) = tsurf_val
+#ifdef NEWPHYS
+   IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** IDEALIZED SIMULATION isoil and dsoil standard'
+   isoil_val(k) = inertiedat_val
+   dsoil_val(k) = sqrt(887.75/3.14)*((2.**(k-0.5))-1.) * inertiedat_val / wvolcapa
+   IF ( (i == ips) .AND. (j == jps) ) PRINT *,'** Mars ** ISOIL DSOIL are ',isoil_val(k), dsoil_val(k)
+#endif
   enddo
 ENDIF
-
-#ifdef NEWPHYS
-isoil_val(:)=MARS_ISOIL(i,:,j)
-dsoil_val(:)=MARS_DSOIL(i,:,j)
-#endif
 
 !-------------------!
Index: trunk/mesoscale/LMD_MM_MARS/makemeso
===================================================================
--- trunk/mesoscale/LMD_MM_MARS/makemeso	(revision 51)
+++ trunk/mesoscale/LMD_MM_MARS/makemeso	(revision 54)
@@ -52,6 +52,7 @@
 justphys=0
 debug=0
+fresh_start=0
 phys=""
-while getopts "drc:njhgp" options; do
+while getopts "drc:njhgpf" options; do
   case $options in
    d ) donotcompile=1;;     ## just to check the compile folder
@@ -61,5 +62,6 @@
    j ) justphys=1;;         ## just compile LMD physics
    g ) debug=1;;            ## debug mode
-   p ) phys="newphys_";;
+   p ) phys="newphys_";;    ## with new physics
+   f ) fresh_start=1;;      ## a fresh start
    h ) echo "
 # Use:    
@@ -83,4 +85,6 @@
 #
 # makemeso -p                ## with new LMD physics
+#
+# makemeso -f                ## fresh start [clean -a]
 " ; exit ;;
   esac
@@ -185,26 +189,34 @@
     if [[ "${phys}" == "newphys_" ]]
     then
-    sed s+"mars_lmd"+"mars_lmd_new"+g copy_model > copy_model_tmp
-    chmod 755 copy_model_tmp
-    ./copy_model_tmp
-    \rm copy_model_tmp
-    mv WRFV2 ${conf_wrf}/
-    cd ${conf_wrf}/WRFV2
-    ln -sf mars_lmd_new mars_lmd
-    ln -sf meso_callkeys_newphys.h meso_callkeys.h
-    cd Registry ; Registry.bash ; cd ..
+     sed s+"mars_lmd"+"mars_lmd_new"+g copy_model > copy_model_tmp
+     chmod 755 copy_model_tmp
+     ./copy_model_tmp
+     \rm copy_model_tmp
+     mv WRFV2 ${conf_wrf}/
+     cd ${conf_wrf}/WRFV2
+     ln -sf mars_lmd_new mars_lmd
+     ln -sf meso_callkeys_newphys.h meso_callkeys.h
+     cd Registry ; Registry.bash ; cd ..
     else
-    ./copy_model
-    mv WRFV2 ${conf_wrf}/
-    cd ${conf_wrf}/WRFV2
-    cd Registry ; Registry.bash ; cd .. 
-    fi
+     ./copy_model
+     mv WRFV2 ${conf_wrf}/
+     cd ${conf_wrf}/WRFV2
+     cd Registry ; Registry.bash ; cd .. 
+     fi
   else
     cd ${conf_wrf}/WRFV2 
-    echo Did you modify anything in the Registry or clean ? y for yes, any key for no ; read answer
-    case ${answer} in
-      y) cd Registry ; mv Registry Registry.bak ; Registry.bash ; cd .. ;;
-      *) answer='no' ;;
-    esac 
+    if [ ${fresh_start} -eq 1 ] 
+    then
+      echo '*** FRESH START, I clean everything'
+      clean -a > /dev/null 2> /dev/null
+      cd Registry ; mv Registry Registry.bak ; Registry.bash ; cd ..
+      answer='y' ## a voir... pas si sur
+    else
+      echo Did you modify anything in the Registry or clean ? y for yes, any key for no ; read answer
+      case ${answer} in
+        y) cd Registry ; mv Registry Registry.bak ; Registry.bash ; cd .. ;;
+        *) answer='no' ;;
+      esac
+    fi 
   fi
   # summary
