Index: trunk/LMDZ.MARS/libf/phymars/dyn1d/init_testphys1d_mod.F90
===================================================================
--- trunk/LMDZ.MARS/libf/phymars/dyn1d/init_testphys1d_mod.F90	(revision 3186)
+++ trunk/LMDZ.MARS/libf/phymars/dyn1d/init_testphys1d_mod.F90	(revision 3188)
@@ -7,5 +7,5 @@
 SUBROUTINE init_testphys1d(start1Dname,startfiname,startfiles_1D,therestart1D,therestartfi,ngrid,nlayer,odpref, &
                            nq,q,time,psurf,u,v,temp,ndt,ptif,pks,dttestphys,zqsat,dq,dqdyn,day0,day,gru,grv,w,  &
-                           play,plev,latitude,longitude,cell_area,atm_wat_profile,atm_wat_tau)
+                           play,plev,latitude,longitude,cell_area,atm_wat_profile,atm_wat_tau,CO2cond_ps)
 
 use ioipsl_getincom,          only: getin ! To use 'getin'
@@ -78,4 +78,5 @@
 real, dimension(1),                  intent(out) :: latitude, longitude, cell_area
 real,                                intent(out) :: atm_wat_profile, atm_wat_tau ! Force atmospheric water profiles
+real,                                intent(out) :: CO2cond_ps                 ! Relaxation coefficient for psurf
 
 !=======================================================================
@@ -308,4 +309,14 @@
 write(*,*) " psurf = ",psurf
 
+! Coefficient to control the surface pressure change
+! To be defined in "callphys.def" so that the PEM can read it asa well
+CO2cond_ps = 1. ! Default value
+write(*,*) 'Coefficient to control the surface pressure change?'
+call getin("CO2cond_ps",CO2cond_ps)
+if (CO2cond_ps < 0. .or. CO2cond_ps > 1.) then
+    write(*,*) 'Value for ''CO2cond_ps'' is not between 0 and 1.'
+    error stop 'Please, specify a correct value!'
+endif
+
 ! Reference pressures
 pa = 20.     ! transition pressure (for hybrid coord.)
Index: trunk/LMDZ.MARS/libf/phymars/dyn1d/testphys1d.F90
===================================================================
--- trunk/LMDZ.MARS/libf/phymars/dyn1d/testphys1d.F90	(revision 3186)
+++ trunk/LMDZ.MARS/libf/phymars/dyn1d/testphys1d.F90	(revision 3188)
@@ -71,6 +71,7 @@
 real                                :: dttestphys    ! testphys1d timestep
 real, dimension(nlayer)             :: play          ! Pressure at the middle of the layers (Pa)
-real, dimension(nlayer + 1)         :: plev          ! intermediate pressure levels (pa)
+real, dimension(nlayer + 1)         :: plev          ! intermediate pressure levels (Pa)
 real                                :: psurf         ! Surface pressure
+real                                :: CO2cond_ps    ! Coefficient to control the surface pressure change
 real, dimension(nlayer)             :: u, v          ! zonal, meridional wind
 real                                :: gru, grv      ! prescribed "geostrophic" background wind
@@ -149,5 +150,5 @@
 call init_testphys1d('start1D.txt','startfi.nc',startfiles_1D,therestart1D,therestartfi,ngrid,nlayer,odpref, &
                      nq,q,time,psurf,u,v,temp,ndt,ptif,pks,dttestphys,zqsat,dq,dqdyn,day0,day,gru,grv,w,     &
-                     play,plev,latitude,longitude,cell_area,atm_wat_profile,atm_wat_tau)
+                     play,plev,latitude,longitude,cell_area,atm_wat_profile,atm_wat_tau,CO2cond_ps)
 
 ! Write a "startfi" file
@@ -261,5 +262,8 @@
     ! Compute pressure for next time step
     ! -----------------------------------
-    psurf = psurf + dttestphys*dpsurf(1) ! surface pressure change
+    ! CO2cond_ps is a coefficient to control the surface pressure change
+    ! If CO2cond_ps = 0, psurf is kept constant; If CO2cond_ps = 1, psurf varies normally
+    ! CO2cond_ps = 300*400/144.37e6 = 8.31e-4 (ratio of polar cap surface over planetary surface) is a typical value for tests
+    psurf = psurf + CO2cond_ps*dttestphys*dpsurf(1) ! Surface pressure change
     plev = ap + psurf*bp
     play = aps + psurf*bps