Changeset 227 for codes/icosagcm/branches/SATURN_DYNAMICO
- Timestamp:
- 07/16/14 18:05:01 (10 years ago)
- Location:
- codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf
- Files:
-
- 1 added
- 2 deleted
- 85 edited
Legend:
- Unmodified
- Added
- Removed
-
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/dyn3d_common/cpdet_mod.F90
r222 r227 142 142 ! Parallel version of t2tpot, for an arbitrary number of columns 143 143 USE control_mod, only : planet_type 144 USE parallel_lmdz, only : OMP_CHUNK 144 145 IMPLICIT none 145 146 … … 179 180 ! Parallel version of t2tpot, over the full dynamics (scalar) grid 180 181 ! (more efficient than multiple calls to t2tpot_p() with slices of data) 181 USE parallel_lmdz, only : jj_begin,jj_end 182 USE parallel_lmdz, only : jj_begin,jj_end,OMP_CHUNK 182 183 USE control_mod, only : planet_type 183 184 IMPLICIT none … … 223 224 ! Parallel version of tpot2t, for an arbitrary number of columns 224 225 USE control_mod, only : planet_type 226 USE parallel_lmdz, only : OMP_CHUNK 225 227 IMPLICIT none 226 228 ! for cpp, nu_venus and t0_venus: … … 259 261 ! Parallel version of tpot2t, over the full dynamics (scalar) grid 260 262 ! (more efficient than multiple calls to tpot2t_p() with slices of data) 261 USE parallel_lmdz, only : jj_begin,jj_end 263 USE parallel_lmdz, only : jj_begin,jj_end,OMP_CHUNK 262 264 USE control_mod, only : planet_type 263 265 IMPLICIT none -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/dyn3dpar/gcm.F
r222 r227 599 599 c write(78,*) 'q',q 600 600 601 c$OMP PARALLEL DEFAULT(SHARED) COPYIN(/temps/,/logici/,/logicl/) 601 !c$OMP PARALLEL DEFAULT(SHARED) COPYIN(/temps/,/logici/,/logicl/) 602 !variable temps no longer exists 603 c$OMP PARALLEL DEFAULT(SHARED) 604 c$OMP1 COPYIN(/temps_r/,/temps_i/,/temps_c/,/logici/,/logicl/) 602 605 CALL leapfrog_p(ucov,vcov,teta,ps,masse,phis,q, 603 606 . time_0) -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/dyn3dpar/sponge_mod_p.F90
r222 r227 31 31 32 32 USE Write_Field_p 33 use parallel_lmdz, only: pole_sud,pole_nord,jj_begin,jj_end 33 use parallel_lmdz, only: pole_sud,pole_nord,jj_begin,jj_end,OMP_CHUNK 34 34 implicit none 35 35 #include "dimensions.h" -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/aeropacity.F90
r222 r227 30 30 ! pq Aerosol mixing ratio 31 31 ! reffrad(ngrid,nlayer,naerkind) Aerosol effective radius 32 ! QREFvis3d(ngrid,nlayer mx,naerkind) \ 3d extinction coefficients33 ! QREFir3d(ngrid,nlayer mx,naerkind) / at reference wavelengths32 ! QREFvis3d(ngrid,nlayer,naerkind) \ 3d extinction coefficients 33 ! QREFir3d(ngrid,nlayer,naerkind) / at reference wavelengths 34 34 ! 35 35 ! Output … … 40 40 !======================================================================= 41 41 42 #include "dimensions.h"43 #include "dimphys.h"42 !#include "dimensions.h" 43 !#include "dimphys.h" 44 44 #include "callkeys.h" 45 45 #include "comcstfi.h" 46 #include "comvert.h"46 !#include "comvert.h" 47 47 48 48 INTEGER,INTENT(IN) :: ngrid ! number of atmospheric columns … … 54 54 REAL,INTENT(OUT) :: aerosol(ngrid,nlayer,naerkind) ! aerosol optical depth 55 55 REAL,INTENT(IN) :: reffrad(ngrid,nlayer,naerkind) ! aerosol effective radius 56 REAL,INTENT(IN) :: QREFvis3d(ngrid,nlayer mx,naerkind) ! extinction coefficient in the visible57 REAL,INTENT(IN) :: QREFir3d(ngrid,nlayer mx,naerkind)56 REAL,INTENT(IN) :: QREFvis3d(ngrid,nlayer,naerkind) ! extinction coefficient in the visible 57 REAL,INTENT(IN) :: QREFir3d(ngrid,nlayer,naerkind) 58 58 REAL,INTENT(OUT):: tau_col(ngrid) !column integrated visible optical depth 59 59 ! BENJAMIN MODIFS 60 real,intent(in) :: cloudfrac(ngrid,nlayer mx) ! cloud fraction60 real,intent(in) :: cloudfrac(ngrid,nlayer) ! cloud fraction 61 61 real,intent(out) :: totcloudfrac(ngrid) ! total cloud fraction 62 62 logical,intent(in) :: clearsky … … 67 67 68 68 LOGICAL,SAVE :: firstcall=.true. 69 !$OMP THREADPRIVATE(firstcall) 69 70 REAL CBRT 70 71 EXTERNAL CBRT … … 72 73 INTEGER,SAVE :: i_co2ice=0 ! co2 ice 73 74 INTEGER,SAVE :: i_h2oice=0 ! water ice 75 !$OMP THREADPRIVATE(i_co2ice,i_h2oice) 74 76 CHARACTER(LEN=20) :: tracername ! to temporarily store text 75 77 … … 135 137 iaer=iaero_co2 136 138 ! 1. Initialization 137 aerosol(1:ngrid,1:nlayer mx,iaer)=0.0139 aerosol(1:ngrid,1:nlayer,iaer)=0.0 138 140 ! 2. Opacity calculation 139 141 if (noaero) then ! aerosol set to zero 140 aerosol(1:ngrid,1:nlayer mx,iaer)=0.0142 aerosol(1:ngrid,1:nlayer,iaer)=0.0 141 143 elseif (aerofixco2.or.(i_co2ice.eq.0)) then ! CO2 ice cloud prescribed 142 aerosol(1:ngrid,1:nlayer mx,iaer)=1.e-9144 aerosol(1:ngrid,1:nlayer,iaer)=1.e-9 143 145 !aerosol(1:ngrid,12,iaer)=4.0 ! single cloud layer option 144 146 else … … 172 174 iaer=iaero_h2o 173 175 ! 1. Initialization 174 aerosol(1:ngrid,1:nlayer mx,iaer)=0.0176 aerosol(1:ngrid,1:nlayer,iaer)=0.0 175 177 ! 2. Opacity calculation 176 178 if (aerofixh2o.or.(i_h2oice.eq.0).or.clearsky) then 177 aerosol(1:ngrid,1:nlayer mx,iaer) =1.e-9179 aerosol(1:ngrid,1:nlayer,iaer) =1.e-9 178 180 179 181 ! put cloud at cloudlvl … … 226 228 227 229 if(CLFvarying)then 228 call totalcloudfrac(ngrid,n q,cloudfrac,totcloudfrac,pplev,pq,aerosol(1:ngrid,1:nlayermx,iaer))230 call totalcloudfrac(ngrid,nlayer,nq,cloudfrac,totcloudfrac,pplev,pq,aerosol(1,1,iaer)) 229 231 do ig=1, ngrid 230 232 do l=1,nlayer-1 ! to stop the rad tran bug … … 253 255 iaer=iaero_dust 254 256 ! 1. Initialization 255 aerosol(1:ngrid,1:nlayer mx,iaer)=0.0257 aerosol(1:ngrid,1:nlayer,iaer)=0.0 256 258 257 259 topdust=30.0 ! km (used to be 10.0 km) LK … … 264 266 ! Typical mixing ratio profile 265 267 266 zp=(p reff/pplay(ig,l))**(70./topdust)268 zp=(pplev(ig,1)/pplay(ig,l))**(70./topdust) 267 269 expfactor=max(exp(0.007*(1.-max(zp,1.))),1.e-3) 268 270 … … 277 279 ! Rescaling each layer to reproduce the choosen (or assimilated) 278 280 ! dust extinction opacity at visible reference wavelength, which 279 ! is scaled to the "preff" reference surface pressure available in comvert.h 280 ! and stored in startfi.nc 281 ! is scaled to the surface pressure pplev(ig,1) 281 282 282 283 taudusttmp(1:ngrid)=0. … … 291 292 aerosol(ig,l,iaer) = max(1E-20, & 292 293 dusttau & 293 * pplev(ig,1) / p reff&294 * pplev(ig,1) / pplev(ig,1) & 294 295 * aerosol(ig,l,iaer) & 295 296 / taudusttmp(ig)) … … 307 308 308 309 ! 1. Initialization 309 aerosol(1:ngrid,1:nlayer mx,iaer)=0.0310 aerosol(1:ngrid,1:nlayer,iaer)=0.0 310 311 311 312 … … 317 318 ! Typical mixing ratio profile 318 319 319 zp=(p reff/pplay(ig,l))**(70./30) !emulating topdust320 zp=(pplev(ig,1)/pplay(ig,l))**(70./30) !emulating topdust 320 321 expfactor=max(exp(0.007*(1.-max(zp,1.))),1.e-3) 321 322 … … 335 336 aerosol(ig,l,iaer) = max(1E-20, & 336 337 1 & 337 * pplev(ig,1) / p reff&338 * pplev(ig,1) / pplev(ig,1) & 338 339 * aerosol(ig,l,iaer) & 339 340 / tauh2so4tmp(ig)) … … 367 368 iaer=iaero_back2lay 368 369 ! 1. Initialization 369 aerosol(1:ngrid,1:nlayer mx,iaer)=0.0370 aerosol(1:ngrid,1:nlayer,iaer)=0.0 370 371 ! 2. Opacity calculation 371 372 DO ig=1,ngrid -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/aeroptproperties.F90
r222 r227 32 32 ! ============================================================== 33 33 34 #include "dimensions.h"35 #include "dimphys.h"34 !#include "dimensions.h" 35 !#include "dimphys.h" 36 36 #include "callkeys.h" 37 37 … … 60 60 INTEGER :: grid_i,grid_j 61 61 ! Intermediate variable 62 REAL :: var_tmp,var3d_tmp(ngrid,nlayer mx)62 REAL :: var_tmp,var3d_tmp(ngrid,nlayer) 63 63 ! Bilinear interpolation factors 64 64 REAL :: kx,ky,k1,k2,k3,k4 … … 67 67 ! Pi! 68 68 REAL,SAVE :: pi 69 !$OMP THREADPRIVATE(pi) 69 70 ! Variables used by the Gauss-Legendre integration: 70 71 INTEGER radius_id,gausind … … 87 88 0.08327674160932,0.06267204829828, & 88 89 0.04060142982019,0.01761400714091/ 90 !$OMP THREADPRIVATE(radgaus,weightgaus) 89 91 ! Indices 90 92 INTEGER :: i,j,k,l,m,iaer,idomain … … 102 104 ! Grid used to remember which calculation is done 103 105 LOGICAL,SAVE :: checkgrid(refftabsize,nuefftabsize,naerkind,2) = .false. 106 !$OMP THREADPRIVATE(refftab,nuefftab,logvratgrid,vratgrid,checkgrid) 104 107 ! Optical properties of the grid (VISIBLE) 105 108 REAL,SAVE :: qsqrefVISgrid(refftabsize,nuefftabsize,L_NSPECTV,naerkind) … … 108 111 REAL,SAVE :: omegVISgrid(refftabsize,nuefftabsize,L_NSPECTV,naerkind) 109 112 REAL,SAVE :: gVISgrid(refftabsize,nuefftabsize,L_NSPECTV,naerkind) 113 !$OMP THREADPRIVATE(qsqrefVISgrid,qextVISgrid,qscatVISgrid,omegVISgrid,gVISgrid) 110 114 ! Optical properties of the grid (INFRARED) 111 115 REAL,SAVE :: qsqrefIRgrid(refftabsize,nuefftabsize,L_NSPECTI,naerkind) … … 114 118 REAL,SAVE :: omegIRgrid(refftabsize,nuefftabsize,L_NSPECTI,naerkind) 115 119 REAL,SAVE :: gIRgrid(refftabsize,nuefftabsize,L_NSPECTI,naerkind) 120 !$OMP THREADPRIVATE(qsqrefIRgrid,qextIRgrid,qscatIRgrid,omegIRgrid,gIRgrid) 116 121 ! Optical properties of the grid (REFERENCE WAVELENGTHS) 117 122 REAL,SAVE :: qrefVISgrid(refftabsize,nuefftabsize,naerkind) … … 121 126 REAL,SAVE :: omegrefVISgrid(refftabsize,nuefftabsize,naerkind) 122 127 REAL,SAVE :: omegrefIRgrid(refftabsize,nuefftabsize,naerkind) 128 !$OMP THREADPRIVATE(qrefVISgrid,qscatrefVISgrid,qrefIRgrid,qscatrefIRgrid,omegrefVISgrid,& 129 !$OMP omegrefIRgrid) 123 130 ! Firstcall 124 131 LOGICAL,SAVE :: firstcall = .true. 132 !$OMP THREADPRIVATE(firstcall) 125 133 ! Variables used by the Gauss-Legendre integration: 126 134 REAL,SAVE :: normd(refftabsize,nuefftabsize,naerkind,2) 127 135 REAL,SAVE :: dista(refftabsize,nuefftabsize,naerkind,2,ngau) 128 136 REAL,SAVE :: distb(refftabsize,nuefftabsize,naerkind,2,ngau) 137 !$OMP THREADPRIVATE(normd,dista,distb) 129 138 130 139 REAL,SAVE :: radGAUSa(ngau,naerkind,2) 131 140 REAL,SAVE :: radGAUSb(ngau,naerkind,2) 141 !$OMP THREADPRIVATE(radGAUSa,radGAUSb) 132 142 133 143 REAL,SAVE :: qsqrefVISa(L_NSPECTV,ngau,naerkind) … … 141 151 REAL,SAVE :: gVISa(L_NSPECTV,ngau,naerkind) 142 152 REAL,SAVE :: gVISb(L_NSPECTV,ngau,naerkind) 153 !$OMP THREADPRIVATE(qsqrefVISa,qrefVISa,qsqrefVISb,qrefVISb,omegVISa, & 154 !$OMP omegrefVISa,omegVISb,omegrefVISb,gVISa,gVISb) 143 155 144 156 REAL,SAVE :: qsqrefIRa(L_NSPECTI,ngau,naerkind) … … 152 164 REAL,SAVE :: gIRa(L_NSPECTI,ngau,naerkind) 153 165 REAL,SAVE :: gIRb(L_NSPECTI,ngau,naerkind) 166 !$OMP THREADPRIVATE(qsqrefIRa,qrefIRa,qsqrefIRb,qrefIRb,omegIRa,omegrefIRa,& 167 !$OMP omegIRb,omegrefIRb,gIRa,gIRb) 154 168 155 169 REAL :: radiusm … … 161 175 INTEGER :: ngrid,nlayer 162 176 ! Aerosol effective radius used for radiative transfer (meter) 163 REAL :: reffrad(ngrid,nlayermx,naerkind)177 REAL,INTENT(IN) :: reffrad(ngrid,nlayer,naerkind) 164 178 ! Aerosol effective variance used for radiative transfer (n.u.) 165 REAL :: nueffrad(ngrid,nlayermx,naerkind)179 REAL,INTENT(IN) :: nueffrad(ngrid,nlayer,naerkind) 166 180 167 181 ! Outputs 168 182 ! ------- 169 183 170 REAL :: QVISsQREF3d(ngrid,nlayermx,L_NSPECTV,naerkind)171 REAL :: omegaVIS3d(ngrid,nlayermx,L_NSPECTV,naerkind)172 REAL :: gVIS3d(ngrid,nlayermx,L_NSPECTV,naerkind)173 174 REAL :: QIRsQREF3d(ngrid,nlayermx,L_NSPECTI,naerkind)175 REAL :: omegaIR3d(ngrid,nlayermx,L_NSPECTI,naerkind)176 REAL :: gIR3d(ngrid,nlayermx,L_NSPECTI,naerkind)177 178 REAL :: QREFvis3d(ngrid,nlayermx,naerkind)179 REAL :: QREFir3d(ngrid,nlayermx,naerkind)180 181 REAL :: omegaREFvis3d(ngrid,nlayer mx,naerkind)182 REAL :: omegaREFir3d(ngrid,nlayer mx,naerkind)184 REAL,INTENT(OUT) :: QVISsQREF3d(ngrid,nlayer,L_NSPECTV,naerkind) 185 REAL,INTENT(OUT) :: omegaVIS3d(ngrid,nlayer,L_NSPECTV,naerkind) 186 REAL,INTENT(OUT) :: gVIS3d(ngrid,nlayer,L_NSPECTV,naerkind) 187 188 REAL,INTENT(OUT) :: QIRsQREF3d(ngrid,nlayer,L_NSPECTI,naerkind) 189 REAL,INTENT(OUT) :: omegaIR3d(ngrid,nlayer,L_NSPECTI,naerkind) 190 REAL,INTENT(OUT) :: gIR3d(ngrid,nlayer,L_NSPECTI,naerkind) 191 192 REAL,INTENT(OUT) :: QREFvis3d(ngrid,nlayer,naerkind) 193 REAL,INTENT(OUT) :: QREFir3d(ngrid,nlayer,naerkind) 194 195 REAL :: omegaREFvis3d(ngrid,nlayer,naerkind) 196 REAL :: omegaREFir3d(ngrid,nlayer,naerkind) 183 197 184 198 DO iaer = 1, naerkind ! Loop on aerosol kind … … 724 738 k2*omegrefIRgrid(grid_i+1,1,iaer) 725 739 ENDIF ! -------------------------------- 726 ENDDO !nlayer mx740 ENDDO !nlayer 727 741 ENDDO !ngrid 728 742 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/aerosol_mod.F90
r222 r227 16 16 ! two-layer simple aerosol model 17 17 integer :: iaero_back2lay = 0 18 !$OMP THREADPRIVATE(iaero_co2,iaero_h2o,iaero_dust,iaero_h2so4,noaero,iaero_back2lay) 18 19 19 20 !================================================================== -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/ave_stelspec.F90
r222 r227 38 38 integer ifine 39 39 40 real,allocatable :: lam(:),stel_f(:)40 real,allocatable,save :: lam(:),stel_f(:) !read by master 41 41 real band,lamm,lamp 42 42 real dl … … 116 116 End Select 117 117 118 !$OMP MASTER 118 119 allocate(lam(Nfine),stel_f(Nfine)) 119 120 … … 153 154 close(111) 154 155 endif 156 !$OMP END MASTER 157 !$OMP BARRIER 155 158 156 159 ! sum data by band … … 173 176 174 177 STELLAR(1:L_NSPECTV)=STELLAR(1:L_NSPECTV)/sum(STELLAR(1:L_NSPECTV)) 175 deallocate(lam,stel_f) 176 178 !$OMP BARRIER 179 !$OMP MASTER 180 if (allocated(lam)) deallocate(lam) 181 if (allocated(stel_f)) deallocate(stel_f) 182 !$OMP END MASTER 183 !$OMP BARRIER 177 184 endif 178 185 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/bilinearbig.F90
r222 r227 14 14 real*8 f2d_arr(nX,nY) 15 15 real*8,save :: x,y 16 !$OMP THREADPRIVATE(x,y) 16 17 17 18 integer strlen -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/calc_cpp_mugaz.F90
r222 r227 20 20 implicit none 21 21 22 #include "dimensions.h"23 #include "dimphys.h"22 !#include "dimensions.h" 23 !#include "dimphys.h" 24 24 #include "comcstfi.h" 25 25 #include "callkeys.h" -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/calcenergy_kcm.F90
r222 r227 1 subroutine calcenergy_kcm( Tsurf,T,Play,Plev,Qsurf,Q,muvar,Eatmtot)1 subroutine calcenergy_kcm(nlayer,Tsurf,T,Play,Plev,Qsurf,Q,muvar,Eatmtot) 2 2 3 3 … … 11 11 ! ---------------------------------------------------------------- 12 12 13 #include "dimensions.h"14 #include "dimphys.h"13 !#include "dimensions.h" 14 !#include "dimphys.h" 15 15 #include "comcstfi.h" 16 16 !#include "callkeys.h" … … 19 19 20 20 ! inputs 21 real Tsurf,T(1:nlayermx) 22 real Play(1:nlayermx),Plev(1:nlayermx+1) 23 real Qsurf,Q(1:nlayermx) 24 real muvar(1,nlayermx+1) 21 integer,intent(in) :: nlayer 22 real Tsurf,T(1:nlayer) 23 real Play(1:nlayer),Plev(1:nlayer+1) 24 real Qsurf,Q(1:nlayer) 25 real muvar(1,nlayer+1) 25 26 26 27 ! internal … … 29 30 real s_c,rho_v,L 30 31 double precision p_v, s_v, nul 31 real VMR(1:nlayer mx)32 real VMR(1:nlayer) 32 33 33 34 ! output … … 41 42 42 43 43 do il=1,nlayer mx44 do il=1,nlayer 44 45 VMR(il)=Q(il)*(muvar(1,il+1)/mH2O) 45 46 end do … … 47 48 48 49 Eatmtot = 0.0 49 do il=1,nlayer mx50 do il=1,nlayer 50 51 51 52 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/callcorrk.F90
r222 r227 12 12 use watercommon_h 13 13 use datafile_mod, only: datadir 14 use ioipsl_getincom 14 ! use ioipsl_getincom 15 use ioipsl_getincom_p 15 16 use gases_h 16 17 use radii_mod, only : su_aer_radii,co2_reffrad,h2o_reffrad,dust_reffrad,h2so4_reffrad,back2lay_reffrad … … 35 36 !================================================================== 36 37 37 #include "dimphys.h"38 !#include "dimphys.h" 38 39 #include "comcstfi.h" 39 40 #include "callkeys.h" … … 42 43 ! Declaration of the arguments (INPUT - OUTPUT) on the LMD GCM grid 43 44 ! Layer #1 is the layer near the ground. 44 ! Layer #nlayer mxis the layer at the top.45 ! Layer #nlayer is the layer at the top. 45 46 46 47 INTEGER,INTENT(IN) :: ngrid ! number of atmospheric columns … … 52 53 REAL,INTENT(IN) :: emis(ngrid) ! LW emissivity 53 54 real,intent(in) :: mu0(ngrid) ! cosine of sun incident angle 54 REAL,INTENT(IN) :: pplev(ngrid,nlayer mx+1) ! inter-layer pressure (Pa)55 REAL,INTENT(IN) :: pplay(ngrid,nlayer mx) ! mid-layer pressure (Pa)56 REAL,INTENT(IN) :: pt(ngrid,nlayer mx) ! air temperature (K)55 REAL,INTENT(IN) :: pplev(ngrid,nlayer+1) ! inter-layer pressure (Pa) 56 REAL,INTENT(IN) :: pplay(ngrid,nlayer) ! mid-layer pressure (Pa) 57 REAL,INTENT(IN) :: pt(ngrid,nlayer) ! air temperature (K) 57 58 REAL,INTENT(IN) :: tsurf(ngrid) ! surface temperature (K) 58 59 REAL,INTENT(IN) :: fract(ngrid) ! fraction of day 59 60 REAL,INTENT(IN) :: dist_star ! distance star-planet (AU) 60 REAL,INTENT(OUT) :: aerosol(ngrid,nlayer mx,naerkind) ! aerosol tau (kg/kg)61 real,intent(in) :: muvar(ngrid,nlayer mx+1)62 REAL,INTENT(OUT) :: dtlw(ngrid,nlayer mx) ! heating rate (K/s) due to LW63 REAL,INTENT(OUT) :: dtsw(ngrid,nlayer mx) ! heating rate (K/s) due to SW61 REAL,INTENT(OUT) :: aerosol(ngrid,nlayer,naerkind) ! aerosol tau (kg/kg) 62 real,intent(in) :: muvar(ngrid,nlayer+1) 63 REAL,INTENT(OUT) :: dtlw(ngrid,nlayer) ! heating rate (K/s) due to LW 64 REAL,INTENT(OUT) :: dtsw(ngrid,nlayer) ! heating rate (K/s) due to SW 64 65 REAL,INTENT(OUT) :: fluxsurf_lw(ngrid) ! incident LW flux to surf (W/m2) 65 66 REAL,INTENT(OUT) :: fluxsurf_sw(ngrid) ! incident SW flux to surf (W/m2) … … 71 72 REAL,INTENT(OUT) :: tau_col(ngrid) ! diagnostic from aeropacity 72 73 ! for H2O cloud fraction in aeropacity 73 real,intent(in) :: cloudfrac(ngrid,nlayer mx)74 real,intent(in) :: cloudfrac(ngrid,nlayer) 74 75 real,intent(out) :: totcloudfrac(ngrid) 75 76 logical,intent(in) :: clearsky … … 79 80 ! Globally varying aerosol optical properties on GCM grid 80 81 ! Not needed everywhere so not in radcommon_h 81 REAL :: QVISsQREF3d(ngrid,nlayer mx,L_NSPECTV,naerkind)82 REAL :: omegaVIS3d(ngrid,nlayer mx,L_NSPECTV,naerkind)83 REAL :: gVIS3d(ngrid,nlayer mx,L_NSPECTV,naerkind)84 85 REAL :: QIRsQREF3d(ngrid,nlayer mx,L_NSPECTI,naerkind)86 REAL :: omegaIR3d(ngrid,nlayer mx,L_NSPECTI,naerkind)87 REAL :: gIR3d(ngrid,nlayer mx,L_NSPECTI,naerkind)88 89 ! REAL :: omegaREFvis3d(ngrid,nlayer mx,naerkind)90 ! REAL :: omegaREFir3d(ngrid,nlayer mx,naerkind) ! not sure of the point of these...82 REAL :: QVISsQREF3d(ngrid,nlayer,L_NSPECTV,naerkind) 83 REAL :: omegaVIS3d(ngrid,nlayer,L_NSPECTV,naerkind) 84 REAL :: gVIS3d(ngrid,nlayer,L_NSPECTV,naerkind) 85 86 REAL :: QIRsQREF3d(ngrid,nlayer,L_NSPECTI,naerkind) 87 REAL :: omegaIR3d(ngrid,nlayer,L_NSPECTI,naerkind) 88 REAL :: gIR3d(ngrid,nlayer,L_NSPECTI,naerkind) 89 90 ! REAL :: omegaREFvis3d(ngrid,nlayer,naerkind) 91 ! REAL :: omegaREFir3d(ngrid,nlayer,naerkind) ! not sure of the point of these... 91 92 92 93 REAL,ALLOCATABLE,SAVE :: reffrad(:,:,:) ! aerosol effective radius (m) 93 94 REAL,ALLOCATABLE,SAVE :: nueffrad(:,:,:) ! aerosol effective variance 95 !$OMP THREADPRIVATE(reffrad,nueffrad) 94 96 95 97 !----------------------------------------------------------------------- … … 128 130 logical global1d 129 131 save szangle,global1d 132 !$OMP THREADPRIVATE(szangle,global1d) 130 133 real*8 taugsurf(L_NSPECTV,L_NGAUSS-1) 131 134 real*8 taugsurfi(L_NSPECTI,L_NGAUSS-1) … … 141 144 real*8,save :: GIAER(L_LEVELS+1,L_NSPECTI,naerkind) 142 145 143 !REAL :: QREFvis3d(ngrid,nlayer mx,naerkind)144 !REAL :: QREFir3d(ngrid,nlayer mx,naerkind)146 !REAL :: QREFvis3d(ngrid,nlayer,naerkind) 147 !REAL :: QREFir3d(ngrid,nlayer,naerkind) 145 148 !save QREFvis3d, QREFir3d 146 149 real, dimension(:,:,:), save, allocatable :: QREFvis3d 147 150 real, dimension(:,:,:), save, allocatable :: QREFir3d 151 !$OMP THREADPRIVATE(QXVAER,QSVAER,GVAER,QXIAER,QSIAER,GIAER,QREFvis3d,QREFir3d) 148 152 149 153 … … 177 181 178 182 ! included by RW for runaway greenhouse 1D study 179 real vtmp(nlayer mx)183 real vtmp(nlayer) 180 184 REAL*8 muvarrad(L_LEVELS) 181 185 … … 195 199 !!! ALLOCATED instances are necessary because of CLFvarying 196 200 !!! strategy to call callcorrk twice in physiq... 197 IF(.not.ALLOCATED(QREFvis3d)) ALLOCATE(QREFvis3d(ngrid,nlayer mx,naerkind))198 IF(.not.ALLOCATED(QREFir3d)) ALLOCATE(QREFir3d(ngrid,nlayer mx,naerkind))201 IF(.not.ALLOCATED(QREFvis3d)) ALLOCATE(QREFvis3d(ngrid,nlayer,naerkind)) 202 IF(.not.ALLOCATED(QREFir3d)) ALLOCATE(QREFir3d(ngrid,nlayer,naerkind)) 199 203 ! Effective radius and variance of the aerosols 200 204 IF(.not.ALLOCATED(reffrad)) allocate(reffrad(ngrid,nlayer,naerkind)) … … 207 211 call abort 208 212 endif 209 call su_aer_radii(ngrid, reffrad,nueffrad)213 call su_aer_radii(ngrid,nlayer,reffrad,nueffrad) 210 214 211 215 … … 214 218 ! set up correlated k 215 219 print*, "callcorrk: Correlated-k data base folder:",trim(datadir) 216 call getin ("corrkdir",corrkdir)220 call getin_p("corrkdir",corrkdir) 217 221 print*, "corrkdir = ",corrkdir 218 222 write( tmp1, '(i3)' ) L_NSPECTI … … 238 242 PRINT*, 'Simulate global averaged conditions ?' 239 243 global1d = .false. ! default value 240 call getin ("global1d",global1d)244 call getin_p("global1d",global1d) 241 245 write(*,*) "global1d = ",global1d 242 246 ! Test of incompatibility: … … 251 255 PRINT *,'(assumed for averaged solar flux S/4)' 252 256 szangle=60.0 ! default value 253 call getin ("szangle",szangle)257 call getin_p("szangle",szangle) 254 258 write(*,*) "szangle = ",szangle 255 259 endif … … 266 270 267 271 if ((iaer.eq.iaero_co2).and.tracer.and.(igcm_co2_ice.gt.0)) then ! treat condensed co2 particles. 268 call co2_reffrad(ngrid,n q,pq,reffrad(1,1,iaero_co2))269 print*,'Max. CO2 ice particle size = ',maxval(reffrad(1:ngrid,1:nlayer mx,iaer))/1.e-6,' um'270 print*,'Min. CO2 ice particle size = ',minval(reffrad(1:ngrid,1:nlayer mx,iaer))/1.e-6,' um'272 call co2_reffrad(ngrid,nlayer,nq,pq,reffrad(1,1,iaero_co2)) 273 print*,'Max. CO2 ice particle size = ',maxval(reffrad(1:ngrid,1:nlayer,iaer))/1.e-6,' um' 274 print*,'Min. CO2 ice particle size = ',minval(reffrad(1:ngrid,1:nlayer,iaer))/1.e-6,' um' 271 275 end if 272 276 if ((iaer.eq.iaero_h2o).and.water) then ! treat condensed water particles. to be generalized for other aerosols 273 call h2o_reffrad(ngrid, pq(1,1,igcm_h2o_ice),pt, &277 call h2o_reffrad(ngrid,nlayer,pq(1,1,igcm_h2o_ice),pt, & 274 278 reffrad(1,1,iaero_h2o),nueffrad(1,1,iaero_h2o)) 275 print*,'Max. H2O cloud particle size = ',maxval(reffrad(1:ngrid,1:nlayer mx,iaer))/1.e-6,' um'276 print*,'Min. H2O cloud particle size = ',minval(reffrad(1:ngrid,1:nlayer mx,iaer))/1.e-6,' um'279 print*,'Max. H2O cloud particle size = ',maxval(reffrad(1:ngrid,1:nlayer,iaer))/1.e-6,' um' 280 print*,'Min. H2O cloud particle size = ',minval(reffrad(1:ngrid,1:nlayer,iaer))/1.e-6,' um' 277 281 endif 278 282 if(iaer.eq.iaero_dust)then 279 call dust_reffrad(ngrid, reffrad(1,1,iaero_dust))283 call dust_reffrad(ngrid,nlayer,reffrad(1,1,iaero_dust)) 280 284 print*,'Dust particle size = ',reffrad(1,1,iaer)/1.e-6,' um' 281 285 endif 282 286 if(iaer.eq.iaero_h2so4)then 283 call h2so4_reffrad(ngrid, reffrad(1,1,iaero_h2so4))287 call h2so4_reffrad(ngrid,nlayer,reffrad(1,1,iaero_h2so4)) 284 288 print*,'H2SO4 particle size =',reffrad(1,1,iaer)/1.e-6,' um' 285 289 endif … … 318 322 do iaer=1,naerkind 319 323 DO nw=1,L_NSPECTV 320 do l=1,nlayer mx321 322 temp1=QVISsQREF3d(ig,nlayer mx+1-l,nw,iaer) &323 *QREFvis3d(ig,nlayer mx+1-l,iaer)324 325 temp2=QVISsQREF3d(ig,max(nlayer mx-l,1),nw,iaer) &326 *QREFvis3d(ig,max(nlayer mx-l,1),iaer)324 do l=1,nlayer 325 326 temp1=QVISsQREF3d(ig,nlayer+1-l,nw,iaer) & 327 *QREFvis3d(ig,nlayer+1-l,iaer) 328 329 temp2=QVISsQREF3d(ig,max(nlayer-l,1),nw,iaer) & 330 *QREFvis3d(ig,max(nlayer-l,1),iaer) 327 331 328 332 qxvaer(2*l,nw,iaer) = temp1 329 333 qxvaer(2*l+1,nw,iaer)=(temp1+temp2)/2 330 334 331 temp1=temp1*omegavis3d(ig,nlayer mx+1-l,nw,iaer)332 temp2=temp2*omegavis3d(ig,max(nlayer mx-l,1),nw,iaer)335 temp1=temp1*omegavis3d(ig,nlayer+1-l,nw,iaer) 336 temp2=temp2*omegavis3d(ig,max(nlayer-l,1),nw,iaer) 333 337 334 338 qsvaer(2*l,nw,iaer) = temp1 335 339 qsvaer(2*l+1,nw,iaer)=(temp1+temp2)/2 336 340 337 temp1=gvis3d(ig,nlayer mx+1-l,nw,iaer)338 temp2=gvis3d(ig,max(nlayer mx-l,1),nw,iaer)341 temp1=gvis3d(ig,nlayer+1-l,nw,iaer) 342 temp2=gvis3d(ig,max(nlayer-l,1),nw,iaer) 339 343 340 344 gvaer(2*l,nw,iaer) = temp1 … … 344 348 345 349 qxvaer(1,nw,iaer)=qxvaer(2,nw,iaer) 346 qxvaer(2*nlayer mx+1,nw,iaer)=0.350 qxvaer(2*nlayer+1,nw,iaer)=0. 347 351 348 352 qsvaer(1,nw,iaer)=qsvaer(2,nw,iaer) 349 qsvaer(2*nlayer mx+1,nw,iaer)=0.353 qsvaer(2*nlayer+1,nw,iaer)=0. 350 354 351 355 gvaer(1,nw,iaer)=gvaer(2,nw,iaer) 352 gvaer(2*nlayer mx+1,nw,iaer)=0.356 gvaer(2*nlayer+1,nw,iaer)=0. 353 357 354 358 end do … … 356 360 ! longwave 357 361 DO nw=1,L_NSPECTI 358 do l=1,nlayer mx359 360 temp1=QIRsQREF3d(ig,nlayer mx+1-l,nw,iaer) &361 *QREFir3d(ig,nlayer mx+1-l,iaer)362 363 temp2=QIRsQREF3d(ig,max(nlayer mx-l,1),nw,iaer) &364 *QREFir3d(ig,max(nlayer mx-l,1),iaer)362 do l=1,nlayer 363 364 temp1=QIRsQREF3d(ig,nlayer+1-l,nw,iaer) & 365 *QREFir3d(ig,nlayer+1-l,iaer) 366 367 temp2=QIRsQREF3d(ig,max(nlayer-l,1),nw,iaer) & 368 *QREFir3d(ig,max(nlayer-l,1),iaer) 365 369 366 370 qxiaer(2*l,nw,iaer) = temp1 367 371 qxiaer(2*l+1,nw,iaer)=(temp1+temp2)/2 368 372 369 temp1=temp1*omegair3d(ig,nlayer mx+1-l,nw,iaer)370 temp2=temp2*omegair3d(ig,max(nlayer mx-l,1),nw,iaer)373 temp1=temp1*omegair3d(ig,nlayer+1-l,nw,iaer) 374 temp2=temp2*omegair3d(ig,max(nlayer-l,1),nw,iaer) 371 375 372 376 qsiaer(2*l,nw,iaer) = temp1 373 377 qsiaer(2*l+1,nw,iaer)=(temp1+temp2)/2 374 378 375 temp1=gir3d(ig,nlayer mx+1-l,nw,iaer)376 temp2=gir3d(ig,max(nlayer mx-l,1),nw,iaer)379 temp1=gir3d(ig,nlayer+1-l,nw,iaer) 380 temp2=gir3d(ig,max(nlayer-l,1),nw,iaer) 377 381 378 382 giaer(2*l,nw,iaer) = temp1 … … 382 386 383 387 qxiaer(1,nw,iaer)=qxiaer(2,nw,iaer) 384 qxiaer(2*nlayer mx+1,nw,iaer)=0.388 qxiaer(2*nlayer+1,nw,iaer)=0. 385 389 386 390 qsiaer(1,nw,iaer)=qsiaer(2,nw,iaer) 387 qsiaer(2*nlayer mx+1,nw,iaer)=0.391 qsiaer(2*nlayer+1,nw,iaer)=0. 388 392 389 393 giaer(1,nw,iaer)=giaer(2,nw,iaer) 390 giaer(2*nlayer mx+1,nw,iaer)=0.394 giaer(2*nlayer+1,nw,iaer)=0. 391 395 392 396 end do … … 481 485 elseif(varfixed)then 482 486 483 do l=1,nlayer mx! here we will assign fixed water vapour profiles globally487 do l=1,nlayer ! here we will assign fixed water vapour profiles globally 484 488 RH = satval * ((pplay(ig,l)/pplev(ig,1) - 0.02) / 0.98) 485 489 if(RH.lt.0.0) RH=0.0 … … 507 511 ! call watersat(Ttemp,ptemp,qsat) 508 512 509 qvar(2*nlayer mx+1)= RH * qsat ! ~realistic profile (e.g. 80% saturation at ground)513 qvar(2*nlayer+1)= RH * qsat ! ~realistic profile (e.g. 80% saturation at ground) 510 514 511 515 else … … 541 545 542 546 muvarrad(1) = muvarrad(2) 543 muvarrad(2*nlayer mx+1)=muvar(ig,1)547 muvarrad(2*nlayer+1)=muvar(ig,1) 544 548 545 549 print*,'Recalculating qvar with VARIABLE epsi for kastprof' … … 579 583 580 584 muvarrad(1) = muvarrad(2) 581 muvarrad(2*nlayer mx+1)=muvar(ig,1)585 muvarrad(2*nlayer+1)=muvar(ig,1) 582 586 endif 583 587 … … 607 611 608 612 tlevrad(1) = tlevrad(2) 609 tlevrad(2*nlayer mx+1)=tsurf(ig)613 tlevrad(2*nlayer+1)=tsurf(ig) 610 614 611 615 tmid(1) = tlevrad(2) … … 866 870 IF( ALLOCATED( gasi ) ) DEALLOCATE( gasi ) 867 871 IF( ALLOCATED( gasv ) ) DEALLOCATE( gasv ) 872 !$OMP BARRIER 873 !$OMP MASTER 868 874 IF( ALLOCATED( pgasref ) ) DEALLOCATE( pgasref ) 869 875 IF( ALLOCATED( tgasref ) ) DEALLOCATE( tgasref ) 870 876 IF( ALLOCATED( wrefvar ) ) DEALLOCATE( wrefvar ) 871 877 IF( ALLOCATED( pfgasref ) ) DEALLOCATE( pfgasref ) 878 !$OMP END MASTER 879 !$OMP BARRIER 872 880 IF ( ALLOCATED(reffrad)) DEALLOCATE(reffrad) 873 881 IF ( ALLOCATED(nueffrad)) DEALLOCATE(nueffrad) -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/callkeys.h
r222 r227 117 117 real J2 118 118 real MassPlanet 119 120 logical :: iscallphys=.false.!existence of callphys.def -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/callsedim.F
r222 r227 28 28 c ------------- 29 29 30 #include "dimensions.h"31 #include "dimphys.h"30 !#include "dimensions.h" 31 !#include "dimphys.h" 32 32 #include "comcstfi.h" 33 33 #include "callkeys.h" … … 65 65 real epaisseur (ngrid,nlay) ! Layer thickness (m) 66 66 real wq(ngrid,nlay+1) ! displaced tracer mass (kg.m-2) 67 c real dens(ngrid,nlay ermx) ! Mean density of the ice part. accounting for dust core67 c real dens(ngrid,nlay) ! Mean density of the ice part. accounting for dust core 68 68 69 69 70 70 LOGICAL,SAVE :: firstcall=.true. 71 !$OMP THREADPRIVATE(firstcall) 71 72 72 73 c ** un petit test de coherence … … 125 126 if (water.and.(iq.eq.igcm_h2o_ice)) then 126 127 ! compute radii for h2o_ice 127 call h2o_reffrad(ngrid, zqi(1,1,igcm_h2o_ice),zt,128 call h2o_reffrad(ngrid,nlay,zqi(1,1,igcm_h2o_ice),zt, 128 129 & reffrad(1,1,iaero_h2o),nueffrad(1,1,iaero_h2o)) 129 130 ! call sedimentation for h2o_ice -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/comdiurn_h.F90
r222 r227 6 6 real, allocatable, dimension(:) :: sinlon, coslon, sinlat, coslat 7 7 logical :: ldiurn 8 !$OMP THREADPRIVATE(sinlon,coslon,sinlat,coslat,ldiurn) !ldiurn is unused 8 9 9 10 end module comdiurn_h -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/comgeomfi_h.F90
r222 r227 6 6 REAL,ALLOCATABLE,DIMENSION(:) :: long,lati,area 7 7 REAL :: totarea, totarea_planet 8 !$OMP THREADPRIVATE(long,lati,area,totarea) 8 9 9 10 end module comgeomfi_h -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/comsaison_h.F90
r222 r227 9 9 10 10 real, allocatable, dimension(:) :: mu0,fract 11 !$OMP THREADPRIVATE(isaison,callsais,dist_star,declin,mu0,fract) 11 12 12 13 end module comsaison_h -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/comsoil_h.F90
r222 r227 14 14 ! in physdem (or set via tabfi, or initialized in 15 15 ! soil_settings.F) 16 !$OMP THREADPRIVATE(layer,mlayer,inertiedat,volcapa) 16 17 17 18 contains -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/condense_cloud.F90
r222 r227 32 32 ! ptsrf(ngrid) Surface temperature 33 33 ! 34 ! pdt(ngrid,nlayer mx) Time derivative before condensation/sublimation of pt34 ! pdt(ngrid,nlayer) Time derivative before condensation/sublimation of pt 35 35 ! pdtsrf(ngrid) Time derivative before condensation/sublimation of ptsrf 36 36 ! pqsurf(ngrid,nq) Sedimentation flux at the surface (kg.m-2.s-1) … … 38 38 ! Outputs 39 39 ! ------- 40 ! pdpsrf(ngrid) 41 ! pdtc(ngrid,nlayer mx) / to the time derivatives of Ps, pt, and ptsrf42 ! pdtsrfc(ngrid) 40 ! pdpsrf(ngrid) \ Contribution of condensation/sublimation 41 ! pdtc(ngrid,nlayer) / to the time derivatives of Ps, pt, and ptsrf 42 ! pdtsrfc(ngrid) / 43 43 ! 44 44 ! Both … … 56 56 !================================================================== 57 57 58 #include "dimensions.h"59 #include "dimphys.h"58 !#include "dimensions.h" 59 !#include "dimphys.h" 60 60 #include "comcstfi.h" 61 #include "comvert.h"61 !#include "comvert.h" 62 62 #include "callkeys.h" 63 63 … … 99 99 100 100 REAL reffrad(ngrid,nlayer) ! radius (m) of the co2 ice particles 101 REAL*8 zt(ngrid,nlayer mx)102 REAL zq(ngrid,nlayer mx,nq)101 REAL*8 zt(ngrid,nlayer) 102 REAL zq(ngrid,nlayer,nq) 103 103 REAL zcpi 104 REAL ztcond (ngrid,nlayer mx)105 REAL ztnuc (ngrid,nlayer mx)104 REAL ztcond (ngrid,nlayer) 105 REAL ztnuc (ngrid,nlayer) 106 106 REAL ztcondsol(ngrid) 107 107 REAL zdiceco2(ngrid) 108 REAL zcondicea(ngrid,nlayer mx), zcondices(ngrid)108 REAL zcondicea(ngrid,nlayer), zcondices(ngrid) 109 109 REAL zfallice(ngrid), Mfallice(ngrid) 110 REAL zmflux(nlayer mx+1)111 REAL zu(nlayer mx),zv(nlayermx)110 REAL zmflux(nlayer+1) 111 REAL zu(nlayer),zv(nlayer) 112 112 REAL ztsrf(ngrid) 113 REAL ztc(nlayer mx), ztm(nlayermx+1)114 REAL zum(nlayer mx+1) , zvm(nlayermx+1)113 REAL ztc(nlayer), ztm(nlayer+1) 114 REAL zum(nlayer+1) , zvm(nlayer+1) 115 115 LOGICAL condsub(ngrid) 116 116 REAL subptimestep 117 117 Integer Ntime 118 real masse (ngrid,nlayer mx), w(ngrid,nlayermx,nq)119 real wq(ngrid,nlayer mx+1)118 real masse (ngrid,nlayer), w(ngrid,nlayer,nq) 119 real wq(ngrid,nlayer+1) 120 120 real vstokes,reff 121 121 122 122 ! Special diagnostic variables 123 real tconda1(ngrid,nlayer mx)124 real tconda2(ngrid,nlayer mx)125 real zdtsig (ngrid,nlayer mx)126 real zdt (ngrid,nlayer mx)123 real tconda1(ngrid,nlayer) 124 real tconda2(ngrid,nlayer) 125 real zdtsig (ngrid,nlayer) 126 real zdt (ngrid,nlayer) 127 127 128 128 !----------------------------------------------------------------------- … … 133 133 REAL,SAVE :: cpice=1000. 134 134 REAL,SAVE,ALLOCATABLE,DIMENSION(:) :: emisref 135 !$OMP THREADPRIVATE(latcond,ccond,cpice,emisref) 135 136 136 137 LOGICAL,SAVE :: firstcall=.true. 138 !$OMP THREADPRIVATE(firstcall) 137 139 REAL,EXTERNAL :: SSUM 138 140 … … 140 142 141 143 INTEGER,SAVE :: i_co2ice=0 ! co2 ice 144 !$OMP THREADPRIVATE(i_co2ice) 142 145 CHARACTER(LEN=20) :: tracername ! to temporarily store text 143 146 … … 205 208 ! zcondices(ngrid) condensation rate on the ground (kg/m2/s) 206 209 ! zfallice(ngrid) flux of ice falling on surface (kg/m2/s) 207 ! pdtc(ngrid,nlayer mx) dT/dt due to phase changes (K/s)210 ! pdtc(ngrid,nlayer) dT/dt due to phase changes (K/s) 208 211 209 212 … … 301 304 302 305 ! sedimentation computed from radius computed from q in module radii_mod 303 call co2_reffrad(ngrid,n q,zq,reffrad)306 call co2_reffrad(ngrid,nlayer,nq,zq,reffrad) 304 307 305 308 do ilay=1,nlayer … … 321 324 ! Computing q after sedimentation 322 325 323 call vlz_fi(ngrid, zq(1,1,i_co2ice),2.,masse,w(1,1,i_co2ice),wq)326 call vlz_fi(ngrid,nlayer,zq(1,1,i_co2ice),2.,masse,w(1,1,i_co2ice),wq) 324 327 325 328 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/convadj.F
r222 r227 29 29 ! ------------ 30 30 31 #include "dimensions.h"32 #include "dimphys.h"31 !#include "dimensions.h" 32 !#include "dimphys.h" 33 33 #include "comcstfi.h" 34 34 #include "callkeys.h" … … 57 57 INTEGER jcnt, jadrs(ngrid) 58 58 59 REAL sig(nlay ermx+1),sdsig(nlayermx),dsig(nlayermx)60 REAL zu(ngrid,nlay ermx),zv(ngrid,nlayermx)61 REAL zh(ngrid,nlay ermx)62 REAL zu2(ngrid,nlay ermx),zv2(ngrid,nlayermx)63 REAL zh2(ngrid,nlay ermx), zhc(ngrid,nlayermx)59 REAL sig(nlay+1),sdsig(nlay),dsig(nlay) 60 REAL zu(ngrid,nlay),zv(ngrid,nlay) 61 REAL zh(ngrid,nlay) 62 REAL zu2(ngrid,nlay),zv2(ngrid,nlay) 63 REAL zh2(ngrid,nlay), zhc(ngrid,nlay) 64 64 REAL zhm,zsm,zdsm,zum,zvm,zalpha,zhmc 65 65 … … 67 67 INTEGER iq,ico2 68 68 save ico2 69 REAL zq(ngrid,nlayermx,nq), zq2(ngrid,nlayermx,nq) 69 !$OMP THREADPRIVATE(ico2) 70 REAL zq(ngrid,nlay,nq), zq2(ngrid,nlay,nq) 70 71 REAL zqm(nq),zqco2m 71 72 real m_co2, m_noco2, A , B 72 73 save A, B 74 !$OMP THREADPRIVATE(A,B) 73 75 74 76 real mtot1, mtot2 , mm1, mm2 … … 77 79 save firstcall 78 80 data firstcall/.true./ 81 !$OMP THREADPRIVATE(firstcall) 79 82 80 83 ! for conservation test -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/dimphy.F90
r222 r227 2 2 3 3 INTEGER,SAVE :: klon ! number of atmospheric columns (for this OpenMP subgrid) 4 INTEGER,SAVE :: klev ! number of atmospheric layers 5 INTEGER,SAVE :: klevp1 ! number of atmospheric layers+1 6 INTEGER,SAVE :: klevm1 ! number of atmospheric layers-1 4 INTEGER,SAVE :: klev ! number of atmospheric layers, read by master 5 INTEGER,SAVE :: klevp1 ! number of atmospheric layers+1, read by master 6 INTEGER,SAVE :: klevm1 ! number of atmospheric layers-1, read by master 7 7 ! INTEGER,SAVE :: kflev 8 8 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/evap.F
r222 r227 1 subroutine evap(ngrid,n q,dtime,pt, pq, pdq, pdt,1 subroutine evap(ngrid,nlayer,nq,dtime,pt, pq, pdq, pdt, 2 2 $ dqevap,dtevap, qevap, tevap) 3 3 … … 7 7 implicit none 8 8 9 #include "dimensions.h"10 #include "dimphys.h"9 !#include "dimensions.h" 10 !#include "dimphys.h" 11 11 #include "comcstfi.h" 12 12 … … 24 24 !================================================================== 25 25 26 INTEGER ngrid,n q26 INTEGER ngrid,nlayer,nq 27 27 28 28 ! Arguments: 29 REAL pt(ngrid,nlayer mx)30 REAL pq(ngrid,nlayer mx,nq)31 REAL pdt(ngrid,nlayer mx)32 REAL pdq(ngrid,nlayer mx,nq)33 REAL dqevap(ngrid,nlayer mx)34 REAL dtevap(ngrid,nlayer mx)35 REAL qevap(ngrid,nlayer mx,nq)29 REAL pt(ngrid,nlayer) 30 REAL pq(ngrid,nlayer,nq) 31 REAL pdt(ngrid,nlayer) 32 REAL pdq(ngrid,nlayer,nq) 33 REAL dqevap(ngrid,nlayer) 34 REAL dtevap(ngrid,nlayer) 35 REAL qevap(ngrid,nlayer,nq) 36 36 REAL dtime 37 37 38 38 ! Local: 39 REAL tevap(ngrid,nlayer mx)39 REAL tevap(ngrid,nlayer) 40 40 REAL zlvdcp 41 41 REAL zlsdcp … … 47 47 ! 48 48 49 DO l=1,nlayer mx49 DO l=1,nlayer 50 50 DO ig=1,ngrid 51 51 qevap(ig,l,igcm_h2o_vap)=pq(ig,l,igcm_h2o_vap) … … 57 57 ENDDO 58 58 59 DO l = 1, nlayer mx59 DO l = 1, nlayer 60 60 DO ig = 1, ngrid 61 61 zlvdcp=RLVTT/RCPD!/(1.0+RVTMP2*qevap(ig,l,igcm_h2o_vap)) … … 77 77 ENDDO 78 78 79 RETURN80 79 END -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/forceWCfn.F
r222 r227 1 subroutine forceWCfn(ngrid,n q,pplev,pt,dq,dqs)1 subroutine forceWCfn(ngrid,nlayer,nq,pplev,pt,dq,dqs) 2 2 3 3 USE tracer_h … … 18 18 !================================================================== 19 19 20 #include "dimensions.h"21 #include "dimphys.h"20 !#include "dimensions.h" 21 !#include "dimphys.h" 22 22 #include "comcstfi.h" 23 23 24 INTEGER ngrid,n q24 INTEGER ngrid,nlayer,nq 25 25 26 26 real masse, Wtot, Wdiff 27 27 28 real pplev(ngrid,nlayer mx+1)28 real pplev(ngrid,nlayer+1) 29 29 real pt(ngrid) 30 30 31 31 real dqs(ngrid,nq) 32 real dq(ngrid,nlayer mx,nq)32 real dq(ngrid,nlayer,nq) 33 33 34 34 integer iq, ig, ilay … … 37 37 do ig=1,ngrid 38 38 Wtot = 0.0 39 do ilay=1,nlayer mx39 do ilay=1,nlayer 40 40 masse = (pplev(ig,ilay) - pplev(ig,ilay+1))/g 41 41 Wtot = Wtot + masse*dq(ig,ilay,iq) -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/gases_h.F90
r222 r227 14 14 integer :: ngasmx 15 15 integer :: vgas 16 character*20,allocatable,DIMENSION(:) :: gnom ! name of the gas 16 character*20,allocatable,DIMENSION(:) :: gnom ! name of the gas, read by master 17 17 real,allocatable,DIMENSION(:) :: gfrac 18 18 … … 31 31 integer :: igas_C2H2 32 32 integer :: igas_C2H6 33 !!$OMP THREADPRIVATE(ngasmx,vgas,gnom,gfrac,& 34 ! !$OMP igas_H2,igas_He,igas_H2O,igas_CO2,igas_CO,igas_N2,& 35 ! !$OMP igas_O2,igas_SO2,igas_H2S,igas_CH4,igas_NH3,igas_C2H2,igas_C2H6) 33 36 34 37 end module gases_h -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/hydrol.F90
r222 r227 4 4 pctsrf_sic,sea_ice) 5 5 6 use ioipsl_getincom 6 ! use ioipsl_getincom 7 use ioipsl_getincom_p 7 8 use watercommon_h, only: T_h2O_ice_liq, RLFTT, rhowater, mx_eau_sol 8 9 USE surfdat_h … … 39 40 !================================================================== 40 41 41 #include "dimensions.h"42 #include "dimphys.h"42 !#include "dimensions.h" 43 !#include "dimphys.h" 43 44 #include "comcstfi.h" 44 45 #include "callkeys.h" … … 50 51 real albedoice 51 52 save albedoice 53 !$OMP THREADPRIVATE(albedoice) 52 54 53 55 real snowlayer … … 59 61 logical,save :: activerunoff ! enable simple runoff scheme? 60 62 logical,save :: oceanalbvary ! ocean albedo varies with the diurnal cycle? 63 !$OMP THREADPRIVATE(oceanbulkavg,activerunoff,oceanalbvary) 61 64 62 65 ! Arguments … … 66 69 real totalrunoff, tsea, oceanarea 67 70 save oceanarea 71 !$OMP THREADPRIVATE(runoff,oceanarea) 68 72 69 73 real ptimestep … … 96 100 97 101 integer, save :: ivap, iliq, iice 102 !$OMP THREADPRIVATE(ivap,iliq,iice) 98 103 99 104 logical, save :: firstcall 105 !$OMP THREADPRIVATE(firstcall) 100 106 101 107 data firstcall /.true./ … … 108 114 write(*,*)"Activate runnoff into oceans?" 109 115 activerunoff=.false. 110 call getin ("activerunoff",activerunoff)116 call getin_p("activerunoff",activerunoff) 111 117 write(*,*)" activerunoff = ",activerunoff 112 118 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/ini_archive.F
r222 r227 40 40 41 41 #include "dimensions.h" 42 #include "dimphys.h"42 !#include "dimphys.h" 43 43 #include "paramet.h" 44 44 #include "comconst.h" -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/inifis.F
r222 r227 10 10 use comsoil_h, only: ini_comsoil_h 11 11 use control_mod, only: ecritphy 12 use planete_mod, only: nres 12 13 use planetwide_mod, only: planetwide_sumval 13 14 … … 49 50 use datafile_mod, only: datadir 50 51 ! to use 'getin' 51 USE ioipsl_getincom 52 ! USE ioipsl_getincom 53 USE ioipsl_getincom_p 52 54 IMPLICIT NONE 53 #include "dimensions.h"54 #include "dimphys.h"55 #include "planete.h"55 !#include "dimensions.h" 56 !#include "dimphys.h" 57 !#include "planete.h" 56 58 #include "comcstfi.h" 57 59 #include "callkeys.h" … … 78 80 real psurf,pN2 ! added by RW for Gliese 581d N2+CO2 79 81 82 !$OMP MASTER 80 83 rad=prad 81 84 daysec=pdaysec … … 88 91 avocado = 6.02214179e23 ! added by RW 89 92 90 ! -------------------------------------------------------- 91 ! The usual Tests 92 ! -------------- 93 IF (nlayer.NE.nlayermx) THEN 94 PRINT*,'STOP in inifis' 95 PRINT*,'Probleme de dimensions :' 96 PRINT*,'nlayer = ',nlayer 97 PRINT*,'nlayermx = ',nlayermx 98 STOP 99 ENDIF 93 !$OMP END MASTER 94 !$OMP BARRIER 100 95 101 96 ! read in 'ecritphy' (frequency of calls to physics, in dynamical steps) 102 97 ! (also done in dyn3d/defrun_new but not in LMDZ.COMMON) 103 call getin ("ecritphy",ecritphy)98 call getin_p("ecritphy",ecritphy) 104 99 105 100 ! -------------------------------------------------------------- … … 107 102 ! -------------------------------------------------------------- 108 103 104 !$OMP MASTER 109 105 ! check that 'callphys.def' file is around 110 106 OPEN(99,file='callphys.def',status='old',form='formatted' 111 107 & ,iostat=ierr) 112 108 CLOSE(99) 109 IF(ierr.EQ.0) iscallphys=.true. !iscallphys initialised as false in callkeys.h 110 !$OMP END MASTER 111 !$OMP BARRIER 113 112 114 IF(ierr.EQ.0) THEN 113 !!! IF(ierr.EQ.0) THEN 114 IF(iscallphys) THEN 115 115 PRINT* 116 116 PRINT* … … 121 121 write(*,*) "Directory where external input files are:" 122 122 ! default 'datadir' is set in "datadir_mod" 123 call getin ("datadir",datadir) ! default path123 call getin_p("datadir",datadir) ! default path 124 124 write(*,*) " datadir = ",trim(datadir) 125 125 126 126 write(*,*) "Run with or without tracer transport ?" 127 127 tracer=.false. ! default value 128 call getin ("tracer",tracer)128 call getin_p("tracer",tracer) 129 129 write(*,*) " tracer = ",tracer 130 130 … … 132 132 & " due to tracer evaporation/condensation?" 133 133 mass_redistrib=.false. ! default value 134 call getin ("mass_redistrib",mass_redistrib)134 call getin_p("mass_redistrib",mass_redistrib) 135 135 write(*,*) " mass_redistrib = ",mass_redistrib 136 136 … … 138 138 write(*,*) "(if diurnal=false, diurnal averaged solar heating)" 139 139 diurnal=.true. ! default value 140 call getin ("diurnal",diurnal)140 call getin_p("diurnal",diurnal) 141 141 write(*,*) " diurnal = ",diurnal 142 142 … … 145 145 & "set in 'start'" 146 146 season=.true. ! default value 147 call getin ("season",season)147 call getin_p("season",season) 148 148 write(*,*) " season = ",season 149 149 150 150 write(*,*) "Tidally resonant rotation ?" 151 151 tlocked=.false. ! default value 152 call getin ("tlocked",tlocked)152 call getin_p("tlocked",tlocked) 153 153 write(*,*) "tlocked = ",tlocked 154 154 155 155 write(*,*) "Saturn ring shadowing ?" 156 156 rings_shadow = .false. 157 call getin ("rings_shadow", rings_shadow)157 call getin_p("rings_shadow", rings_shadow) 158 158 write(*,*) "rings_shadow = ", rings_shadow 159 159 160 160 write(*,*) "Compute latitude-dependent gravity field?" 161 161 oblate = .false. 162 call getin ("oblate", oblate)162 call getin_p("oblate", oblate) 163 163 write(*,*) "oblate = ", oblate 164 164 165 165 write(*,*) "Flattening of the planet (a-b)/a " 166 166 flatten = 0.0 167 call getin ("flatten", flatten)167 call getin_p("flatten", flatten) 168 168 write(*,*) "flatten = ", flatten 169 169 … … 171 171 write(*,*) "Needed if oblate=.true.: J2" 172 172 J2 = 0.0 173 call getin ("J2", J2)173 call getin_p("J2", J2) 174 174 write(*,*) "J2 = ", J2 175 175 176 176 write(*,*) "Needed if oblate=.true.: Planet mass (*1e24 kg)" 177 177 MassPlanet = 0.0 178 call getin ("MassPlanet", MassPlanet)178 call getin_p("MassPlanet", MassPlanet) 179 179 write(*,*) "MassPlanet = ", MassPlanet 180 180 181 181 write(*,*) "Needed if oblate=.true.: Planet mean radius (m)" 182 182 Rmean = 0.0 183 call getin ("Rmean", Rmean)183 call getin_p("Rmean", Rmean) 184 184 write(*,*) "Rmean = ", Rmean 185 185 … … 193 193 write(*,*) "Tidal resonance ratio ?" 194 194 nres=0 ! default value 195 call getin ("nres",nres)195 call getin_p("nres",nres) 196 196 write(*,*) "nres = ",nres 197 197 198 198 write(*,*) "Write some extra output to the screen ?" 199 199 lwrite=.false. ! default value 200 call getin ("lwrite",lwrite)200 call getin_p("lwrite",lwrite) 201 201 write(*,*) " lwrite = ",lwrite 202 202 203 203 write(*,*) "Save statistics in file stats.nc ?" 204 204 callstats=.true. ! default value 205 call getin ("callstats",callstats)205 call getin_p("callstats",callstats) 206 206 write(*,*) " callstats = ",callstats 207 207 208 208 write(*,*) "Test energy conservation of model physics ?" 209 209 enertest=.false. ! default value 210 call getin ("enertest",enertest)210 call getin_p("enertest",enertest) 211 211 write(*,*) " enertest = ",enertest 212 212 213 213 write(*,*) "Check to see if cpp values used match gases.def ?" 214 214 check_cpp_match=.true. ! default value 215 call getin ("check_cpp_match",check_cpp_match)215 call getin_p("check_cpp_match",check_cpp_match) 216 216 write(*,*) " check_cpp_match = ",check_cpp_match 217 217 218 218 write(*,*) "call radiative transfer ?" 219 219 callrad=.true. ! default value 220 call getin ("callrad",callrad)220 call getin_p("callrad",callrad) 221 221 write(*,*) " callrad = ",callrad 222 222 223 223 write(*,*) "call correlated-k radiative transfer ?" 224 224 corrk=.true. ! default value 225 call getin ("corrk",corrk)225 call getin_p("corrk",corrk) 226 226 write(*,*) " corrk = ",corrk 227 227 228 228 write(*,*) "prohibit calculations outside corrk T grid?" 229 229 strictboundcorrk=.true. ! default value 230 call getin ("strictboundcorrk",strictboundcorrk)230 call getin_p("strictboundcorrk",strictboundcorrk) 231 231 write(*,*) "strictboundcorrk = ",strictboundcorrk 232 232 … … 234 234 & "(matters only if callrad=T)" 235 235 callgasvis=.false. ! default value 236 call getin ("callgasvis",callgasvis)236 call getin_p("callgasvis",callgasvis) 237 237 write(*,*) " callgasvis = ",callgasvis 238 238 … … 240 240 & "(matters only if callrad=T)" 241 241 continuum=.true. ! default value 242 call getin ("continuum",continuum)242 call getin_p("continuum",continuum) 243 243 write(*,*) " continuum = ",continuum 244 244 245 245 write(*,*) "use analytic function for H2O continuum ?" 246 246 H2Ocont_simple=.false. ! default value 247 call getin ("H2Ocont_simple",H2Ocont_simple)247 call getin_p("H2Ocont_simple",H2Ocont_simple) 248 248 write(*,*) " H2Ocont_simple = ",H2Ocont_simple 249 249 250 250 write(*,*) "call turbulent vertical diffusion ?" 251 251 calldifv=.true. ! default value 252 call getin ("calldifv",calldifv)252 call getin_p("calldifv",calldifv) 253 253 write(*,*) " calldifv = ",calldifv 254 254 255 255 write(*,*) "use turbdiff instead of vdifc ?" 256 256 UseTurbDiff=.true. ! default value 257 call getin ("UseTurbDiff",UseTurbDiff)257 call getin_p("UseTurbDiff",UseTurbDiff) 258 258 write(*,*) " UseTurbDiff = ",UseTurbDiff 259 259 260 260 write(*,*) "call convective adjustment ?" 261 261 calladj=.true. ! default value 262 call getin ("calladj",calladj)262 call getin_p("calladj",calladj) 263 263 write(*,*) " calladj = ",calladj 264 264 265 265 write(*,*) "call CO2 condensation ?" 266 266 co2cond=.false. ! default value 267 call getin ("co2cond",co2cond)267 call getin_p("co2cond",co2cond) 268 268 write(*,*) " co2cond = ",co2cond 269 269 ! Test of incompatibility … … 275 275 write(*,*) "CO2 supersaturation level ?" 276 276 co2supsat=1.0 ! default value 277 call getin ("co2supsat",co2supsat)277 call getin_p("co2supsat",co2supsat) 278 278 write(*,*) " co2supsat = ",co2supsat 279 279 280 280 write(*,*) "Radiative timescale for Newtonian cooling ?" 281 281 tau_relax=30. ! default value 282 call getin ("tau_relax",tau_relax)282 call getin_p("tau_relax",tau_relax) 283 283 write(*,*) " tau_relax = ",tau_relax 284 284 tau_relax=tau_relax*24*3600 ! convert Earth days --> seconds … … 286 286 write(*,*)"call thermal conduction in the soil ?" 287 287 callsoil=.true. ! default value 288 call getin ("callsoil",callsoil)288 call getin_p("callsoil",callsoil) 289 289 write(*,*) " callsoil = ",callsoil 290 290 … … 292 292 & " physical timestep" 293 293 iradia=1 ! default value 294 call getin ("iradia",iradia)294 call getin_p("iradia",iradia) 295 295 write(*,*)" iradia = ",iradia 296 296 297 297 write(*,*)"Rayleigh scattering ?" 298 298 rayleigh=.false. 299 call getin ("rayleigh",rayleigh)299 call getin_p("rayleigh",rayleigh) 300 300 write(*,*)" rayleigh = ",rayleigh 301 301 302 302 write(*,*) "Use blackbody for stellar spectrum ?" 303 303 stelbbody=.false. ! default value 304 call getin ("stelbbody",stelbbody)304 call getin_p("stelbbody",stelbbody) 305 305 write(*,*) " stelbbody = ",stelbbody 306 306 307 307 write(*,*) "Stellar blackbody temperature ?" 308 308 stelTbb=5800.0 ! default value 309 call getin ("stelTbb",stelTbb)309 call getin_p("stelTbb",stelTbb) 310 310 write(*,*) " stelTbb = ",stelTbb 311 311 312 312 write(*,*)"Output mean OLR in 1D?" 313 313 meanOLR=.false. 314 call getin ("meanOLR",meanOLR)314 call getin_p("meanOLR",meanOLR) 315 315 write(*,*)" meanOLR = ",meanOLR 316 316 317 317 write(*,*)"Output spectral OLR in 3D?" 318 318 specOLR=.false. 319 call getin ("specOLR",specOLR)319 call getin_p("specOLR",specOLR) 320 320 write(*,*)" specOLR = ",specOLR 321 321 322 322 write(*,*)"Operate in kastprof mode?" 323 323 kastprof=.false. 324 call getin ("kastprof",kastprof)324 call getin_p("kastprof",kastprof) 325 325 write(*,*)" kastprof = ",kastprof 326 326 327 327 write(*,*)"Uniform absorption in radiative transfer?" 328 328 graybody=.false. 329 call getin ("graybody",graybody)329 call getin_p("graybody",graybody) 330 330 write(*,*)" graybody = ",graybody 331 331 … … 333 333 write(*,*) "Use slab-ocean ?" 334 334 ok_slab_ocean=.false. ! default value 335 call getin ("ok_slab_ocean",ok_slab_ocean)335 call getin_p("ok_slab_ocean",ok_slab_ocean) 336 336 write(*,*) "ok_slab_ocean = ",ok_slab_ocean 337 337 338 338 write(*,*) "Use slab-sea-ice ?" 339 339 ok_slab_sic=.true. ! default value 340 call getin ("ok_slab_sic",ok_slab_sic)340 call getin_p("ok_slab_sic",ok_slab_sic) 341 341 write(*,*) "ok_slab_sic = ",ok_slab_sic 342 342 343 343 write(*,*) "Use heat transport for the ocean ?" 344 344 ok_slab_heat_transp=.true. ! default value 345 call getin ("ok_slab_heat_transp",ok_slab_heat_transp)345 call getin_p("ok_slab_heat_transp",ok_slab_heat_transp) 346 346 write(*,*) "ok_slab_heat_transp = ",ok_slab_heat_transp 347 347 … … 357 357 write(*,*)"Stratospheric temperature for kastprof mode?" 358 358 Tstrat=167.0 359 call getin ("Tstrat",Tstrat)359 call getin_p("Tstrat",Tstrat) 360 360 write(*,*)" Tstrat = ",Tstrat 361 361 362 362 write(*,*)"Remove lower boundary?" 363 363 nosurf=.false. 364 call getin ("nosurf",nosurf)364 call getin_p("nosurf",nosurf) 365 365 write(*,*)" nosurf = ",nosurf 366 366 … … 375 375 . "... matters only if callsoil=F" 376 376 intheat=0. 377 call getin ("intheat",intheat)377 call getin_p("intheat",intheat) 378 378 write(*,*)" intheat = ",intheat 379 379 380 380 write(*,*)"Use Newtonian cooling for radiative transfer?" 381 381 newtonian=.false. 382 call getin ("newtonian",newtonian)382 call getin_p("newtonian",newtonian) 383 383 write(*,*)" newtonian = ",newtonian 384 384 … … 399 399 write(*,*)"Test physics timescale in 1D?" 400 400 testradtimes=.false. 401 call getin ("testradtimes",testradtimes)401 call getin_p("testradtimes",testradtimes) 402 402 write(*,*)" testradtimes = ",testradtimes 403 403 … … 411 411 write(*,*)"Default planetary temperature?" 412 412 tplanet=215.0 413 call getin ("tplanet",tplanet)413 call getin_p("tplanet",tplanet) 414 414 write(*,*)" tplanet = ",tplanet 415 415 416 416 write(*,*)"Which star?" 417 417 startype=1 ! default value = Sol 418 call getin ("startype",startype)418 call getin_p("startype",startype) 419 419 write(*,*)" startype = ",startype 420 420 421 421 write(*,*)"Value of stellar flux at 1 AU?" 422 422 Fat1AU=1356.0 ! default value = Sol today 423 call getin ("Fat1AU",Fat1AU)423 call getin_p("Fat1AU",Fat1AU) 424 424 write(*,*)" Fat1AU = ",Fat1AU 425 425 … … 429 429 write(*,*)"Varying H2O cloud fraction?" 430 430 CLFvarying=.false. ! default value 431 call getin ("CLFvarying",CLFvarying)431 call getin_p("CLFvarying",CLFvarying) 432 432 write(*,*)" CLFvarying = ",CLFvarying 433 433 434 434 write(*,*)"Value of fixed H2O cloud fraction?" 435 435 CLFfixval=1.0 ! default value 436 call getin ("CLFfixval",CLFfixval)436 call getin_p("CLFfixval",CLFfixval) 437 437 write(*,*)" CLFfixval = ",CLFfixval 438 438 439 439 write(*,*)"fixed radii for Cloud particles?" 440 440 radfixed=.false. ! default value 441 call getin ("radfixed",radfixed)441 call getin_p("radfixed",radfixed) 442 442 write(*,*)" radfixed = ",radfixed 443 443 … … 448 448 write(*,*)"Number mixing ratio of CO2 ice particles:" 449 449 Nmix_co2=1.e6 ! default value 450 call getin ("Nmix_co2",Nmix_co2)450 call getin_p("Nmix_co2",Nmix_co2) 451 451 write(*,*)" Nmix_co2 = ",Nmix_co2 452 452 453 453 ! write(*,*)"Number of radiatively active aerosols:" 454 454 ! naerkind=0. ! default value 455 ! call getin ("naerkind",naerkind)455 ! call getin_p("naerkind",naerkind) 456 456 ! write(*,*)" naerkind = ",naerkind 457 457 458 458 write(*,*)"Opacity of dust (if used):" 459 459 dusttau=0. ! default value 460 call getin ("dusttau",dusttau)460 call getin_p("dusttau",dusttau) 461 461 write(*,*)" dusttau = ",dusttau 462 462 463 463 write(*,*)"Radiatively active CO2 aerosols?" 464 464 aeroco2=.false. ! default value 465 call getin ("aeroco2",aeroco2)465 call getin_p("aeroco2",aeroco2) 466 466 write(*,*)" aeroco2 = ",aeroco2 467 467 468 468 write(*,*)"Fixed CO2 aerosol distribution?" 469 469 aerofixco2=.false. ! default value 470 call getin ("aerofixco2",aerofixco2)470 call getin_p("aerofixco2",aerofixco2) 471 471 write(*,*)" aerofixco2 = ",aerofixco2 472 472 473 473 write(*,*)"Radiatively active water ice?" 474 474 aeroh2o=.false. ! default value 475 call getin ("aeroh2o",aeroh2o)475 call getin_p("aeroh2o",aeroh2o) 476 476 write(*,*)" aeroh2o = ",aeroh2o 477 477 478 478 write(*,*)"Fixed H2O aerosol distribution?" 479 479 aerofixh2o=.false. ! default value 480 call getin ("aerofixh2o",aerofixh2o)480 call getin_p("aerofixh2o",aerofixh2o) 481 481 write(*,*)" aerofixh2o = ",aerofixh2o 482 482 483 483 write(*,*)"Radiatively active sulfuric acid aersols?" 484 484 aeroh2so4=.false. ! default value 485 call getin ("aeroh2so4",aeroh2so4)485 call getin_p("aeroh2so4",aeroh2so4) 486 486 write(*,*)" aeroh2so4 = ",aeroh2so4 487 487 … … 490 490 write(*,*)"Radiatively active two-layer aersols?" 491 491 aeroback2lay=.false. ! default value 492 call getin ("aeroback2lay",aeroback2lay)492 call getin_p("aeroback2lay",aeroback2lay) 493 493 write(*,*)" aeroback2lay = ",aeroback2lay 494 494 … … 496 496 & "in the tropospheric layer (visible)" 497 497 obs_tau_col_tropo=8.D0 498 call getin ("obs_tau_col_tropo",obs_tau_col_tropo)498 call getin_p("obs_tau_col_tropo",obs_tau_col_tropo) 499 499 write(*,*)" obs_tau_col_tropo = ",obs_tau_col_tropo 500 500 … … 502 502 & "in the stratospheric layer (visible)" 503 503 obs_tau_col_strato=0.08D0 504 call getin ("obs_tau_col_strato",obs_tau_col_strato)504 call getin_p("obs_tau_col_strato",obs_tau_col_strato) 505 505 write(*,*)" obs_tau_col_strato = ",obs_tau_col_strato 506 506 507 507 write(*,*)"TWOLAY AEROSOL: pres_bottom_tropo? in pa" 508 508 pres_bottom_tropo=66000.0 509 call getin ("pres_bottom_tropo",pres_bottom_tropo)509 call getin_p("pres_bottom_tropo",pres_bottom_tropo) 510 510 write(*,*)" pres_bottom_tropo = ",pres_bottom_tropo 511 511 512 512 write(*,*)"TWOLAY AEROSOL: pres_top_tropo? in pa" 513 513 pres_top_tropo=18000.0 514 call getin ("pres_top_tropo",pres_top_tropo)514 call getin_p("pres_top_tropo",pres_top_tropo) 515 515 write(*,*)" pres_top_tropo = ",pres_top_tropo 516 516 517 517 write(*,*)"TWOLAY AEROSOL: pres_bottom_strato? in pa" 518 518 pres_bottom_strato=2000.0 519 call getin ("pres_bottom_strato",pres_bottom_strato)519 call getin_p("pres_bottom_strato",pres_bottom_strato) 520 520 write(*,*)" pres_bottom_strato = ",pres_bottom_strato 521 521 522 522 write(*,*)"TWOLAY AEROSOL: pres_top_strato? in pa" 523 523 pres_top_strato=100.0 524 call getin ("pres_top_strato",pres_top_strato)524 call getin_p("pres_top_strato",pres_top_strato) 525 525 write(*,*)" pres_top_strato = ",pres_top_strato 526 526 … … 528 528 & "tropospheric layer, in meters" 529 529 size_tropo=2.e-6 530 call getin ("size_tropo",size_tropo)530 call getin_p("size_tropo",size_tropo) 531 531 write(*,*)" size_tropo = ",size_tropo 532 532 … … 534 534 & "stratospheric layer, in meters" 535 535 size_strato=1.e-7 536 call getin ("size_strato",size_strato)536 call getin_p("size_strato",size_strato) 537 537 write(*,*)" size_strato = ",size_strato 538 538 … … 541 541 write(*,*)"Cloud pressure level (with kastprof only):" 542 542 cloudlvl=0. ! default value 543 call getin ("cloudlvl",cloudlvl)543 call getin_p("cloudlvl",cloudlvl) 544 544 write(*,*)" cloudlvl = ",cloudlvl 545 545 … … 547 547 Tstrat=167.0 548 548 varactive=.false. 549 call getin ("varactive",varactive)549 call getin_p("varactive",varactive) 550 550 write(*,*)" varactive = ",varactive 551 551 552 552 write(*,*)"Is the variable gas species distribution set?" 553 553 varfixed=.false. 554 call getin ("varfixed",varfixed)554 call getin_p("varfixed",varfixed) 555 555 write(*,*)" varfixed = ",varfixed 556 556 557 557 write(*,*)"What is the saturation % of the variable species?" 558 558 satval=0.8 559 call getin ("satval",satval)559 call getin_p("satval",satval) 560 560 write(*,*)" satval = ",satval 561 561 … … 570 570 write(*,*) "Gravitationnal sedimentation ?" 571 571 sedimentation=.false. ! default value 572 call getin ("sedimentation",sedimentation)572 call getin_p("sedimentation",sedimentation) 573 573 write(*,*) " sedimentation = ",sedimentation 574 574 575 575 write(*,*) "Compute water cycle ?" 576 576 water=.false. ! default value 577 call getin ("water",water)577 call getin_p("water",water) 578 578 write(*,*) " water = ",water 579 579 … … 587 587 write(*,*) "Include water condensation ?" 588 588 watercond=.false. ! default value 589 call getin ("watercond",watercond)589 call getin_p("watercond",watercond) 590 590 write(*,*) " watercond = ",watercond 591 591 … … 599 599 write(*,*) "Include water precipitation ?" 600 600 waterrain=.false. ! default value 601 call getin ("waterrain",waterrain)601 call getin_p("waterrain",waterrain) 602 602 write(*,*) " waterrain = ",waterrain 603 603 604 604 write(*,*) "Include surface hydrology ?" 605 605 hydrology=.false. ! default value 606 call getin ("hydrology",hydrology)606 call getin_p("hydrology",hydrology) 607 607 write(*,*) " hydrology = ",hydrology 608 608 609 609 write(*,*) "Evolve surface water sources ?" 610 610 sourceevol=.false. ! default value 611 call getin ("sourceevol",sourceevol)611 call getin_p("sourceevol",sourceevol) 612 612 write(*,*) " sourceevol = ",sourceevol 613 613 614 614 write(*,*) "Ice evolution timestep ?" 615 615 icetstep=100.0 ! default value 616 call getin ("icetstep",icetstep)616 call getin_p("icetstep",icetstep) 617 617 write(*,*) " icetstep = ",icetstep 618 618 619 619 write(*,*) "Snow albedo ?" 620 620 albedosnow=0.5 ! default value 621 call getin ("albedosnow",albedosnow)621 call getin_p("albedosnow",albedosnow) 622 622 write(*,*) " albedosnow = ",albedosnow 623 623 624 624 write(*,*) "Maximum ice thickness ?" 625 625 maxicethick=2.0 ! default value 626 call getin ("maxicethick",maxicethick)626 call getin_p("maxicethick",maxicethick) 627 627 write(*,*) " maxicethick = ",maxicethick 628 628 629 629 write(*,*) "Freezing point of seawater ?" 630 630 Tsaldiff=-1.8 ! default value 631 call getin ("Tsaldiff",Tsaldiff)631 call getin_p("Tsaldiff",Tsaldiff) 632 632 write(*,*) " Tsaldiff = ",Tsaldiff 633 633 634 634 write(*,*) "Does user want to force cpp and mugaz?" 635 635 force_cpp=.false. ! default value 636 call getin ("force_cpp",force_cpp)636 call getin_p("force_cpp",force_cpp) 637 637 write(*,*) " force_cpp = ",force_cpp 638 638 … … 640 640 mugaz = -99999. 641 641 PRINT *,'MEAN MOLECULAR MASS in g mol-1 ?' 642 call getin ("mugaz",mugaz)642 call getin_p("mugaz",mugaz) 643 643 IF (mugaz.eq.-99999.) THEN 644 644 PRINT *, "mugaz must be set if force_cpp = T" … … 649 649 cpp = -99999. 650 650 PRINT *,'SPECIFIC HEAT CAPACITY in J K-1 kg-1 ?' 651 call getin ("cpp",cpp)651 call getin_p("cpp",cpp) 652 652 IF (cpp.eq.-99999.) THEN 653 653 PRINT *, "cpp must be set if force_cpp = T" … … 711 711 ENDDO 712 712 713 !$OMP MASTER 713 714 pi=2.*asin(1.) ! NB: pi is a common in comcstfi.h 715 !$OMP END MASTER 716 !$OMP BARRIER 714 717 715 718 ! allocate "comsoil_h" arrays -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/iniorbit.F
r222 r227 1 1 SUBROUTINE iniorbit 2 2 $ (papoastr,pperiastr,pyear_day,pperi_day,pobliq) 3 4 USE planete_mod, only: apoastr, periastr, year_day, obliquit, 5 & peri_day, e_elips, p_elips, timeperi 3 6 IMPLICIT NONE 4 7 … … 19 22 c ---------- 20 23 c - Doit etre appele avant d'utiliser orbite. 21 c - initialise une partie du common planete .h24 c - initialise une partie du common planete_mod 22 25 c 23 26 c Arguments: … … 35 38 c ------------- 36 39 37 #include "planete.h"40 !#include "planete.h" 38 41 #include "comcstfi.h" 39 42 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/iniphysiq.F90
r222 r227 17 17 rlatd ! latitudes 18 18 use infotrac, only : nqtot ! number of advected tracers 19 use planete_mod, only: ini_planete_mod 19 20 20 21 implicit none 22 include "dimensions.h" 23 include "comvert.h" 21 24 22 25 real,intent(in) :: prad ! radius of the planet (m) … … 58 61 ENDIF 59 62 60 !$OMP PARALLEL PRIVATE(ibegin,iend) 61 !$OMP+SHARED(parea,pcu,pcv,plon,plat)63 !$OMP PARALLEL PRIVATE(ibegin,iend) & 64 !$OMP SHARED(parea,pcu,pcv,plon,plat) 62 65 63 66 offset=0 … … 67 70 rlond(1:klon_omp)=plon(offset+klon_omp_begin:offset+klon_omp_end) 68 71 rlatd(1:klon_omp)=plat(offset+klon_omp_begin:offset+klon_omp_end) 72 73 ! copy over preff , ap() and bp() 74 call ini_planete_mod(nlayer,preff,ap,bp) 69 75 70 76 ! copy some fundamental parameters to physics -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/initracer.F
r222 r227 24 24 25 25 #include "dimensions.h" 26 #include "dimphys.h"26 !#include "dimphys.h" 27 27 #include "comcstfi.h" 28 28 #include "callkeys.h" -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/iniwrite_specIR.F
r222 r227 34 34 #include "netcdf.inc" 35 35 #include "serre.h" 36 #include"dimphys.h"36 !#include"dimphys.h" 37 37 38 38 c Arguments: -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/iniwrite_specVI.F
r222 r227 34 34 #include "netcdf.inc" 35 35 #include "serre.h" 36 #include"dimphys.h"36 !#include"dimphys.h" 37 37 38 38 c Arguments: -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/interpolateH2H2.F90
r222 r227 43 43 logical firstcall 44 44 45 save wn_arr, temp_arr, abs_arr 45 save wn_arr, temp_arr, abs_arr !read by master 46 46 47 47 character*100 dt_file … … 72 72 dt_file=TRIM(datadir)//'/continuum_data/H2-H2_norm_2011.cia' 73 73 74 !$OMP MASTER 74 75 open(33,file=dt_file,form='formatted',status='old',iostat=ios) 75 76 if (ios.ne.0) then ! file not found … … 102 103 endif 103 104 close(33) 105 !$OMP END MASTER 106 !$OMP BARRIER 104 107 105 108 print*,'interpolateH2H2: At wavenumber ',wn,' cm^-1' -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/interpolateH2He.F90
r222 r227 45 45 logical firstcall 46 46 47 save wn_arr, temp_arr, abs_arr 47 save wn_arr, temp_arr, abs_arr !read by master 48 48 49 49 character*100 dt_file … … 73 73 ! 1.1 Open the ASCII files 74 74 dt_file=TRIM(datadir)//'/continuum_data/H2-He_norm_2011.cia' 75 75 76 !$OMP MASTER 76 77 open(33,file=dt_file,form='formatted',status='old',iostat=ios) 77 78 if (ios.ne.0) then ! file not found … … 104 105 endif 105 106 close(33) 107 !$OMP END MASTER 108 !$OMP BARRIER 106 109 107 110 print*,'interpolateH2He: At wavenumber ',wn,' cm^-1' -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/interpolateH2Ocont_CKD.F90
r222 r227 43 43 logical firstcall 44 44 45 save wn_arr, temp_arr, abs_arrS, abs_arrF 45 save wn_arr, temp_arr, abs_arrS, abs_arrF !read by master 46 46 47 47 character*100 dt_file … … 57 57 ! 1.1 Open the ASCII files 58 58 59 !$OMP MASTER 59 60 ! nu array 60 61 dt_file=TRIM(datadir)//'/continuum_data/H2O_CONT_NU.dat' … … 129 130 print*,' H2O pressure ',presS,' Pa' 130 131 print*,' air pressure ',presF,' Pa' 132 !$OMP END MASTER 133 !$OMP BARRIER 131 134 132 135 endif -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/interpolateN2H2.F90
r222 r227 44 44 logical firstcall 45 45 46 save wn_arr, temp_arr, abs_arr 46 save wn_arr, temp_arr, abs_arr !read by master 47 47 48 48 character*100 dt_file … … 72 72 dt_file=TRIM(datadir)//'/continuum_data/N2-H2_2011.cia' 73 73 74 !$OMP MASTER 74 75 open(33,file=dt_file,form='formatted',status='old',iostat=ios) 75 76 if (ios.ne.0) then ! file not found … … 102 103 endif 103 104 close(33) 105 !$OMP END MASTER 106 !$OMP BARRIER 104 107 105 108 print*,'interpolateN2H2: At wavenumber ',wn,' cm^-1' -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/interpolateN2N2.F90
r222 r227 43 43 logical firstcall 44 44 45 save wn_arr, temp_arr, abs_arr 45 save wn_arr, temp_arr, abs_arr !read by master 46 46 47 47 character*100 dt_file … … 70 70 dt_file=TRIM(datadir)//'/continuum_data/N2-N2_2011.cia' 71 71 72 !$OMP MASTER 72 73 open(33,file=dt_file,form='formatted',status='old',iostat=ios) 73 74 if (ios.ne.0) then ! file not found … … 100 101 endif 101 102 close(33) 103 !$OMP END MASTER 104 !$OMP BARRIER 102 105 103 106 print*,'interpolateN2N2: At wavenumber ',wn,' cm^-1' -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/iostart.F90
r222 r227 5 5 INTEGER,SAVE :: nid_start ! NetCDF file identifier for startfi.nc file 6 6 INTEGER,SAVE :: nid_restart ! NetCDF file identifier for restartfi.nc file 7 !$OMP THREADPRIVATE(nid_start,nid_restart) 7 8 8 9 ! restartfi.nc file dimension identifiers: (see open_restartphy()) … … 16 17 INTEGER,SAVE :: idim8 ! "ocean_layers" dimension 17 18 INTEGER,SAVE :: timeindex ! current time index (for time-dependent fields) 19 !$OMP THREADPRIVATE(idim1,idim2,idim3,idim4,idim5,idim6,idim7,timeindex) 18 20 INTEGER,PARAMETER :: length=100 ! size of tab_cntrl array 19 21 … … 473 475 INTEGER :: ierr 474 476 LOGICAL,SAVE :: already_created=.false. 477 !$OMP THREADPRIVATE(already_created) 475 478 476 479 IF (is_master) THEN … … 956 959 INTEGER :: idim1d 957 960 logical,save :: firsttime=.true. 961 !$OMP THREADPRIVATE(firsttime) 958 962 959 963 IF (is_master) THEN -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/kcm1d.F90
r222 r227 7 7 use ioipsl_getincom, only: getin 8 8 use comsaison_h, only: mu0, fract, dist_star 9 use planete_mod 9 10 ! use control_mod 10 11 … … 30 31 31 32 #include "dimensions.h" 32 #include "dimphys.h"33 !#include "dimphys.h" 33 34 #include "callkeys.h" 34 35 #include "comcstfi.h" 35 #include "planete.h"36 !#include "planete.h" 36 37 !#include "control.h" 37 38 … … 42 43 integer nlayer,nlevel,nq 43 44 integer ilay,ilev,iq,iw,iter 44 real play( nlayermx) ! Pressure at the middle of the layers [Pa]45 real zlay( nlayermx) ! Altitude at middle of the layers [km]46 real plev( nlayermx+1) ! Intermediate pressure levels [Pa]47 real temp( nlayermx) ! temperature at the middle of the layers [K]45 real play(llm) ! Pressure at the middle of the layers [Pa] 46 real zlay(llm) ! Altitude at middle of the layers [km] 47 real plev(llm+1) ! Intermediate pressure levels [Pa] 48 real temp(llm) ! temperature at the middle of the layers [K] 48 49 real,allocatable :: q(:,:) ! tracer mixing ratio [kg/kg] 49 50 real,allocatable :: vmr(:,:) ! tracer mixing ratio [mol/mol] … … 53 54 real emis, albedo 54 55 55 real muvar( nlayermx+1)56 57 real dtsw( nlayermx) ! heating rate (K/s) due to SW58 real dtlw( nlayermx) ! heating rate (K/s) due to LW56 real muvar(llm+1) 57 58 real dtsw(llm) ! heating rate (K/s) due to SW 59 real dtlw(llm) ! heating rate (K/s) due to LW 59 60 real fluxsurf_lw ! incident LW flux to surf (W/m2) 60 61 real fluxtop_lw ! outgoing LW flux to space (W/m2) … … 64 65 65 66 ! not used 66 real reffrad( nlayermx,naerkind)67 real nueffrad( nlayermx,naerkind)68 real cloudfrac( nlayermx)67 real reffrad(llm,naerkind) 68 real nueffrad(llm,naerkind) 69 real cloudfrac(llm) 69 70 real totcloudfrac 70 71 real tau_col 71 72 72 73 real dTstrat 73 real aerosol( nlayermx,naerkind) ! aerosol tau (kg/kg)74 real aerosol(llm,naerkind) ! aerosol tau (kg/kg) 74 75 real OLR_nu(1,L_NSPECTI) 75 76 real OSR_nu(1,L_NSPECTV) … … 93 94 94 95 95 nlayer= nlayermx96 nlayer=llm 96 97 nlevel=nlayer+1 97 98 … … 225 226 ! allocate arrays which depend on number of tracers 226 227 allocate(nametrac(nq)) 227 allocate(q(nlayer mx,nq))228 allocate(vmr(nlayer mx,nq))228 allocate(q(nlayer,nq)) 229 allocate(vmr(nlayer,nq)) 229 230 allocate(qsurf(nq)) 230 231 … … 268 269 269 270 ! Create vertical profiles 270 call kcmprof_fn( psurf,qsurf(1),tsurf, &271 call kcmprof_fn(nlayer,psurf,qsurf(1),tsurf, & 271 272 tstrat,play,plev,zlay,temp,q(:,1),muvar(1)) 272 273 … … 316 317 ! Calculate total atmospheric energy 317 318 Eatmtot=0.0 318 ! call calcenergy_kcm( tsurf,temp,play,plev,qsurf,&319 ! call calcenergy_kcm(nlayer,tsurf,temp,play,plev,qsurf,& 319 320 ! q(:,1),muvar,Eatmtot) 320 321 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/kcmprof_fn.F90
r222 r227 1 subroutine kcmprof_fn( psurf_rcm,qsurf_rcm,Tsurf_rcm,Tstra_rcm,P_rcm,Pl_rcm,z_rcm,T_rcm,q_rcm,m_rcm)1 subroutine kcmprof_fn(nlayer,psurf_rcm,qsurf_rcm,Tsurf_rcm,Tstra_rcm,P_rcm,Pl_rcm,z_rcm,T_rcm,q_rcm,m_rcm) 2 2 3 3 use params_h … … 12 12 ! ---------------------------------------------------------------- 13 13 14 #include "dimensions.h"15 #include "dimphys.h"14 !#include "dimensions.h" 15 !#include "dimphys.h" 16 16 #include "comcstfi.h" 17 17 #include "callkeys.h" … … 21 21 22 22 ! rcm inputs 23 integer nlayer 23 24 real Tsurf_rcm,Tstra_rcm 24 25 25 26 ! rcm outputs 26 27 real psurf_rcm,qsurf_rcm 27 real P_rcm(1:nlayer mx)28 real Pl_rcm(1:nlayer mx+1)29 real z_rcm(1:nlayer mx)30 real T_rcm(1:nlayer mx),q_rcm(1:nlayermx)31 real m_rcm(1:nlayer mx+1)28 real P_rcm(1:nlayer) 29 real Pl_rcm(1:nlayer+1) 30 real z_rcm(1:nlayer) 31 real T_rcm(1:nlayer),q_rcm(1:nlayer) 32 real m_rcm(1:nlayer+1) 32 33 33 34 ! rcm for interpolation (should really use log coords?) … … 161 162 162 163 ! define fine pressure grid 163 dlogp_rcm = -(log(psurf_rcm)-log(ptop))/nlayer mx164 dlogp_rcm = -(log(psurf_rcm)-log(ptop))/nlayer 164 165 165 166 P_rcm(1) = psurf_rcm*exp(dlogp_rcm) 166 do ilay_rcm=1,nlayer mx-1167 do ilay_rcm=1,nlayer-1 167 168 P_rcm(ilay_rcm+1) = P_rcm(ilay_rcm)*exp(dlogp_rcm) 168 169 enddo 169 170 170 171 Pl_rcm(1) = psurf_rcm 171 do ilev_rcm=2,nlayer mx172 do ilev_rcm=2,nlayer 172 173 ! log-linear interpolation 173 174 Pl_rcm(ilev_rcm) = exp( log( P_rcm(ilev_rcm)*P_rcm(ilev_rcm-1) )/2 ) … … 318 319 do ilay=2,nlay 319 320 320 if(ilay_rcm.le.nlayer mx)then321 if(ilay_rcm.le.nlayer)then 321 322 ! interpolate rcm variables 322 323 … … 349 350 350 351 ifinal_rcm=ilay_rcm-1 351 if(ifinal_rcm.lt.nlayer mx)then352 if(ifinal_rcm.lt.nlayer)then 352 353 if(verbose)then 353 354 print*,'Interpolation in kcmprof stopped at layer',ilay_rcm,'!' 354 355 endif 355 356 356 do ilay_rcm=ifinal_rcm+1,nlayer mx357 do ilay_rcm=ifinal_rcm+1,nlayer 357 358 358 359 z_rcm(ilay_rcm) = z_rcm(ilay_rcm-1) … … 364 365 endif 365 366 366 do ilay=2,nlayer mx367 do ilay=2,nlayer 367 368 if(T_rcm(ilay).lt.Ttop)then 368 369 T_rcm(ilay)=Ttop … … 373 374 if(co2cond)then 374 375 print*,'CO2 condensation haircut - assumes CO2-dominated atmosphere!' 375 do ilay=2,nlayer mx376 do ilay=2,nlayer 376 377 if(P_rcm(ilay).lt.518000.)then 377 378 TCO2cond = (-3167.8)/(log(.01*P_rcm(ilay))-23.23) ! Fanale's formula -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/largescale.F90
r222 r227 1 subroutine largescale(ngrid,n q,ptimestep, pplev, pplay, pt, pq,&2 1 subroutine largescale(ngrid,nlayer,nq,ptimestep, pplev, pplay, & 2 pt, pq, pdt, pdq, pdtlsc, pdqvaplsc, pdqliqlsc, rneb) 3 3 4 4 5 5 ! to use 'getin' 6 use ioipsl_getincom 6 ! use ioipsl_getincom 7 use ioipsl_getincom_p 7 8 use watercommon_h, only : RLVTT, RCPD, RVTMP2, & 8 9 T_h2O_ice_clouds,T_h2O_ice_liq,Psat_waterDP,Lcpdqsat_waterDP … … 23 24 !================================================================== 24 25 25 #include "dimensions.h"26 #include "dimphys.h"26 !#include "dimensions.h" 27 !#include "dimphys.h" 27 28 #include "comcstfi.h" 28 29 29 30 #include "callkeys.h" 30 31 31 INTEGER ngrid,n q32 INTEGER ngrid,nlayer,nq 32 33 33 34 ! Arguments 34 35 REAL ptimestep ! intervalle du temps (s) 35 REAL pplev(ngrid,nlayer mx+1) ! pression a inter-couche36 REAL pplay(ngrid,nlayer mx) ! pression au milieu de couche37 REAL pt(ngrid,nlayer mx) ! temperature (K)38 REAL pq(ngrid,nlayer mx,nq) ! tracer mixing ratio (kg/kg)39 REAL pdt(ngrid,nlayer mx) ! physical temperature tenedency (K/s)40 REAL pdq(ngrid,nlayer mx,nq)! physical tracer tenedency (K/s)41 REAL pdtlsc(ngrid,nlayer mx) ! incrementation de la temperature (K)42 REAL pdqvaplsc(ngrid,nlayer mx) ! incrementation de la vapeur d'eau43 REAL pdqliqlsc(ngrid,nlayer mx) ! incrementation de l'eau liquide44 REAL rneb(ngrid,nlayer mx) ! fraction nuageuse36 REAL pplev(ngrid,nlayer+1) ! pression a inter-couche 37 REAL pplay(ngrid,nlayer) ! pression au milieu de couche 38 REAL pt(ngrid,nlayer) ! temperature (K) 39 REAL pq(ngrid,nlayer,nq) ! tracer mixing ratio (kg/kg) 40 REAL pdt(ngrid,nlayer) ! physical temperature tenedency (K/s) 41 REAL pdq(ngrid,nlayer,nq)! physical tracer tenedency (K/s) 42 REAL pdtlsc(ngrid,nlayer) ! incrementation de la temperature (K) 43 REAL pdqvaplsc(ngrid,nlayer) ! incrementation de la vapeur d'eau 44 REAL pdqliqlsc(ngrid,nlayer) ! incrementation de l'eau liquide 45 REAL rneb(ngrid,nlayer) ! fraction nuageuse 45 46 46 47 47 48 ! Options du programme 48 49 REAL, SAVE :: ratqs ! determine largeur de la distribution de vapeur 50 !$OMP THREADPRIVATE(ratqs) 49 51 50 52 ! Variables locales … … 63 65 64 66 ! evaporation calculations 65 REAL dqevap(ngrid,nlayer mx),dtevap(ngrid,nlayermx)66 REAL qevap(ngrid,nlayer mx,nq)67 REAL tevap(ngrid,nlayer mx)67 REAL dqevap(ngrid,nlayer),dtevap(ngrid,nlayer) 68 REAL qevap(ngrid,nlayer,nq) 69 REAL tevap(ngrid,nlayer) 68 70 69 71 DOUBLE PRECISION zx_q(ngrid) 70 72 LOGICAL,SAVE :: firstcall=.true. 73 !$OMP THREADPRIVATE(firstcall) 71 74 72 75 … … 75 78 write(*,*) "value for ratqs? " 76 79 ratqs=0.2 ! default value 77 call getin ("ratqs",ratqs)80 call getin_p("ratqs",ratqs) 78 81 write(*,*) " ratqs = ",ratqs 79 82 … … 83 86 ! GCM -----> subroutine variables, initialisation of outputs 84 87 85 pdtlsc(1:ngrid,1:nlayer mx) = 0.086 pdqvaplsc(1:ngrid,1:nlayer mx) = 0.087 pdqliqlsc(1:ngrid,1:nlayer mx) = 0.088 rneb(1:ngrid,1:nlayer mx) = 0.088 pdtlsc(1:ngrid,1:nlayer) = 0.0 89 pdqvaplsc(1:ngrid,1:nlayer) = 0.0 90 pdqliqlsc(1:ngrid,1:nlayer) = 0.0 91 rneb(1:ngrid,1:nlayer) = 0.0 89 92 Lcp=RLVTT/RCPD 90 93 91 94 92 95 ! Evaporate cloud water/ice 93 call evap(ngrid,n q,ptimestep,pt,pq,pdq,pdt,dqevap,dtevap,qevap,tevap)96 call evap(ngrid,nlayer,nq,ptimestep,pt,pq,pdq,pdt,dqevap,dtevap,qevap,tevap) 94 97 ! note: we use qevap but not tevap in largescale/moistadj 95 98 ! otherwise is a big mess … … 97 100 98 101 ! Boucle verticale (du haut vers le bas) 99 DO k = nlayer mx, 1, -1102 DO k = nlayer, 1, -1 100 103 101 104 zt(1:ngrid)=pt(1:ngrid,k)+(pdt(1:ngrid,k)+dtevap(1:ngrid,k))*ptimestep … … 171 174 pdtlsc(1:ngrid,k) = pdqliqlsc(1:ngrid,k)*real(Lcp) 172 175 173 Enddo ! k= nlayer mx, 1, -1176 Enddo ! k= nlayer, 1, -1 174 177 175 178 176 return177 179 end -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/lect_start_archive.F
r222 r227 1 SUBROUTINE lect_start_archive(date,tsurf,tsoil,emis,q2, 1 SUBROUTINE lect_start_archive(ngrid,nlayer, 2 & date,tsurf,tsoil,emis,q2, 2 3 & t,ucov,vcov,ps,h,phisold_newgrid, 3 4 & q,qsurf,surfith,nid, … … 30 31 31 32 #include "dimensions.h" 32 #include "dimphys.h"33 !#include "dimphys.h" 33 34 !#include "planete.h" 34 35 #include "paramet.h" … … 47 48 c======================================================================= 48 49 50 INTEGER,INTENT(IN) :: ngrid, nlayer 51 49 52 c Old variables dimensions (from file) 50 53 c------------------------------------ … … 99 102 c variable physique 100 103 c------------------ 101 REAL tsurf(ngrid mx) ! surface temperature102 REAL tsoil(ngrid mx,nsoilmx) ! soil temperature103 REAL co2ice(ngrid mx) ! CO2 ice layer104 REAL emis(ngrid mx)105 REAL q2(ngrid mx,nlayermx+1),qsurf(ngridmx,nqtot)106 REAL tslab(ngrid mx,noceanmx)107 REAL rnat(ngrid mx),pctsrf_sic(ngridmx)108 REAL tsea_ice(ngrid mx),sea_ice(ngridmx)109 c REAL phisfi(ngrid mx)104 REAL tsurf(ngrid) ! surface temperature 105 REAL tsoil(ngrid,nsoilmx) ! soil temperature 106 REAL co2ice(ngrid) ! CO2 ice layer 107 REAL emis(ngrid) 108 REAL q2(ngrid,llm+1),qsurf(ngrid,nqtot) 109 REAL tslab(ngrid,noceanmx) 110 REAL rnat(ngrid),pctsrf_sic(ngrid) 111 REAL tsea_ice(ngrid),sea_ice(ngrid) 112 c REAL phisfi(ngrid) 110 113 111 114 INTEGER i,j,l … … 322 325 allocate(varp1 (imold+1,jmold+1,llm+1)) 323 326 324 write(*,*) 'q2',ngrid mx,nlayermx+1327 write(*,*) 'q2',ngrid,llm+1 325 328 write(*,*) 'q2S',iip1,jjp1,llm+1 326 329 write(*,*) 'q2old',imold+1,jmold+1,lmold+1 … … 1052 1055 call interp_horiz (tsurfold,tsurfs,imold,jmold,iim,jjm,1, 1053 1056 & rlonuold,rlatvold,rlonu,rlatv) 1054 call gr_dyn_fi (1,iim+1,jjm+1,ngrid mx,tsurfs,tsurf)1057 call gr_dyn_fi (1,iim+1,jjm+1,ngrid,tsurfs,tsurf) 1055 1058 c write(44,*) 'tsurf', tsurf 1056 1059 … … 1059 1062 ! & imold,jmold,iim,jjm,nsoilmx, 1060 1063 ! & rlonuold,rlatvold,rlonu,rlatv) 1061 ! call gr_dyn_fi (nsoilmx,iim+1,jjm+1,ngrid mx,tsoils,tsoil)1064 ! call gr_dyn_fi (nsoilmx,iim+1,jjm+1,ngrid,tsoils,tsoil) 1062 1065 c write(45,*) 'tsoil',tsoil 1063 1066 … … 1065 1068 call interp_horiz (emisold,emiss,imold,jmold,iim,jjm,1, 1066 1069 & rlonuold,rlatvold,rlonu,rlatv) 1067 call gr_dyn_fi (1,iim+1,jjm+1,ngrid mx,emiss,emis)1070 call gr_dyn_fi (1,iim+1,jjm+1,ngrid,emiss,emis) 1068 1071 c write(46,*) 'emis',emis 1069 1072 … … 1201 1204 1202 1205 ! Reshape inertiedatS to scalar grid as inertiedat 1203 call gr_dyn_fi (nsoilmx,iim+1,jjm+1,ngrid mx,1206 call gr_dyn_fi (nsoilmx,iim+1,jjm+1,ngrid, 1204 1207 & inertiedatS,inertiedat) 1205 1208 … … 1287 1290 1288 1291 ! Reshape tsoilS to scalar grid as tsoil 1289 call gr_dyn_fi (nsoilmx,iim+1,jjm+1,ngrid mx,tsoilS,tsoil)1292 call gr_dyn_fi (nsoilmx,iim+1,jjm+1,ngrid,tsoilS,tsoil) 1290 1293 1291 1294 c----------------------------------------------------------------------- … … 1294 1297 call interp_horiz (tslabold,tslabs,imold,jmold,iim,jjm,noceanmx, 1295 1298 & rlonuold,rlatvold,rlonu,rlatv) 1296 call gr_dyn_fi (1,iim+1,jjm+1,ngrid mx,tslabs,tslab)1299 call gr_dyn_fi (1,iim+1,jjm+1,ngrid,tslabs,tslab) 1297 1300 1298 1301 call interp_horiz (rnatold,rnats,imold,jmold,iim,jjm,1, 1299 1302 & rlonuold,rlatvold,rlonu,rlatv) 1300 call gr_dyn_fi (1,iim+1,jjm+1,ngrid mx,rnats,rnat)1303 call gr_dyn_fi (1,iim+1,jjm+1,ngrid,rnats,rnat) 1301 1304 1302 1305 call interp_horiz (pctsrf_sicold,pctsrf_sics,imold,jmold,iim, 1303 1306 & jjm,1,rlonuold,rlatvold,rlonu,rlatv) 1304 call gr_dyn_fi (1,iim+1,jjm+1,ngrid mx,pctsrf_sics,pctsrf_sic)1307 call gr_dyn_fi (1,iim+1,jjm+1,ngrid,pctsrf_sics,pctsrf_sic) 1305 1308 1306 1309 call interp_horiz (tsea_iceold,tsea_ices,imold,jmold,iim,jjm,1, 1307 1310 & rlonuold,rlatvold,rlonu,rlatv) 1308 call gr_dyn_fi (1,iim+1,jjm+1,ngrid mx,tsea_ices,tsea_ice)1311 call gr_dyn_fi (1,iim+1,jjm+1,ngrid,tsea_ices,tsea_ice) 1309 1312 1310 1313 call interp_horiz (sea_iceold,sea_ices,imold,jmold,iim,jjm,1, 1311 1314 & rlonuold,rlatvold,rlonu,rlatv) 1312 call gr_dyn_fi (1,iim+1,jjm+1,ngrid mx,sea_ices,sea_ice)1315 call gr_dyn_fi (1,iim+1,jjm+1,ngrid,sea_ices,sea_ice) 1313 1316 1314 1317 c----------------------------------------------------------------------- … … 1334 1337 & rlonuold,rlatvold,rlonu,rlatv) 1335 1338 write (*,*) 'lect_start_archive: q2s ', q2s (1,2,1) ! INFO 1336 call gr_dyn_fi (llm+1,iim+1,jjm+1,ngrid mx,q2s,q2)1339 call gr_dyn_fi (llm+1,iim+1,jjm+1,ngrid,q2s,q2) 1337 1340 write (*,*) 'lect_start_archive: q2 ', q2 (1,2) ! INFO 1338 1341 c write(47,*) 'q2',q2 … … 1384 1387 enddo 1385 1388 1386 call gr_dyn_fi (nqtot,iim+1,jjm+1,ngrid mx,qsurfs,qsurf)1389 call gr_dyn_fi (nqtot,iim+1,jjm+1,ngrid,qsurfs,qsurf) 1387 1390 1388 1391 c traceurs 3D … … 1435 1438 enddo 1436 1439 1437 ! call gr_dyn_fi (1,iim+1,jjm+1,ngrid mx,co2ices,co2ice)1440 ! call gr_dyn_fi (1,iim+1,jjm+1,ngrid,co2ices,co2ice) 1438 1441 ! no need to transfer "co2ice" any more; it is in qsurf(igcm_co2_ice) 1439 1442 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/mass_redistribution.F90
r222 r227 9 9 USE comgeomfi_h 10 10 USE tracer_h 11 USE planete_mod, only: bp 11 12 12 13 IMPLICIT NONE … … 50 51 ! ------------------ 51 52 ! 52 #include "dimensions.h"53 #include "dimphys.h"53 !#include "dimensions.h" 54 !#include "dimphys.h" 54 55 #include "comcstfi.h" 55 #include "comvert.h"56 #include "paramet.h"56 !#include "comvert.h" 57 !#include "paramet.h" 57 58 #include "callkeys.h" 58 59 … … 84 85 85 86 ! vertical reorganization of sigma levels 86 REAL zzu(nlayer mx),zzv(nlayermx)87 REAL zzq(nlayer mx,nq),zzt(nlayermx)87 REAL zzu(nlayer),zzv(nlayer) 88 REAL zzq(nlayer,nq),zzt(nlayer) 88 89 ! Dummy variables 89 90 INTEGER n,l,ig,iq 90 REAL zdtsig(ngrid,nlayer mx)91 REAL zmass(ngrid,nlayer mx),zzmass(nlayermx),w(nlayermx+1)92 REAL zdmass_sum(ngrid,nlayer mx+1)93 REAL zmflux(nlayer mx+1)94 REAL zq1(nlayer mx)91 REAL zdtsig(ngrid,nlayer) 92 REAL zmass(ngrid,nlayer),zzmass(nlayer),w(nlayer+1) 93 REAL zdmass_sum(ngrid,nlayer+1) 94 REAL zmflux(nlayer+1) 95 REAL zq1(nlayer) 95 96 REAL ztsrf(ngrid) 96 REAL ztm(nlayermx+1) 97 REAL zum(nlayermx+1) , zvm(nlayermx+1) 98 REAL zqm(nlayermx+1,nq),zqm1(nlayermx+1) 97 REAL ztm(nlayer+1) 98 REAL zum(nlayer+1) , zvm(nlayer+1) 99 REAL zqm(nlayer+1,nq),zqm1(nlayer+1) 100 REAL sigma(nlayer+1) 99 101 100 102 ! local saved variables 101 103 LOGICAL, SAVE :: firstcall=.true. 104 !$OMP THREADPRIVATE(firstcall) 102 105 103 106 !---------------------------------------------------------------------- … … 127 130 128 131 DO ig=1,ngrid 129 zdmass_sum(ig,nlayer mx+1)=0.130 DO l = nlayer mx, 1, -1132 zdmass_sum(ig,nlayer+1)=0. 133 DO l = nlayer, 1, -1 131 134 zmass(ig,l) = (pplev(ig,l)-pplev(ig,l+1))/glat(ig) 132 135 zdmass_sum(ig,l)= zdmass_sum(ig,l+1)+pdmassmr(ig,l) … … 181 184 ! Correction to account for redistribution between sigma or hybrid 182 185 ! layers when changing surface pressure 183 ! zzx quantities have dimension (nlayer mx) to speed up calculation186 ! zzx quantities have dimension (nlayer) to speed up calculation 184 187 ! ************************************************************* 185 188 186 189 DO ig=1,ngrid 187 zzt(1:nlayer mx) = pt(ig,1:nlayermx) + pdt(ig,1:nlayermx) * ptimestep188 zzu(1:nlayer mx) = pu(ig,1:nlayermx) + pdu(ig,1:nlayermx) * ptimestep189 zzv(1:nlayer mx) = pv(ig,1:nlayermx) + pdv(ig,1:nlayermx) * ptimestep190 zzq(1:nlayer mx,1:nq)=pq(ig,1:nlayermx,1:nq)+pdq(ig,1:nlayermx,1:nq)*ptimestep ! must add the water that has fallen???190 zzt(1:nlayer) = pt(ig,1:nlayer) + pdt(ig,1:nlayer) * ptimestep 191 zzu(1:nlayer) = pu(ig,1:nlayer) + pdu(ig,1:nlayer) * ptimestep 192 zzv(1:nlayer) = pv(ig,1:nlayer) + pdv(ig,1:nlayer) * ptimestep 193 zzq(1:nlayer,1:nq)=pq(ig,1:nlayer,1:nq)+pdq(ig,1:nlayer,1:nq)*ptimestep ! must add the water that has fallen??? 191 194 192 195 ! Mass fluxes of air through the sigma levels (kg.m-2.s-1) (>0 when up) 193 196 ! """""""""""""""""""""""""""""""""""""""""""""""""""""""""""""" 194 195 197 zmflux(1) = zmassboil(ig) 198 sigma(1)=1 196 199 DO l=1,nlayer 200 ! Ehouarn: shouldn't we rather compute sigma levels than use bp()? 201 ! sigma(l+1)=pplev(ig,l+1)/pplev(ig,1) 202 ! zmflux(l+1) = zmflux(l) + pdmassmr(ig,l) - & 203 ! (sigma(l)-sigma(l+1))*(zdmass_sum(ig,1)+zmflux(1)) 204 ! if (abs(zmflux(l+1)).lt.1E-13.OR.sigma(l+1).eq.0.) zmflux(l+1)=0. 205 ! Ehouarn: but for now leave things as before 197 206 zmflux(l+1) = zmflux(l) + pdmassmr(ig,l) - (bp(l)-bp(l+1))*(zdmass_sum(ig,1)+zmflux(1)) 198 207 ! zmflux set to 0 if very low to avoid: top layer is disappearing in v1ld … … 202 211 ! Mass of each layer 203 212 ! ------------------ 204 zzmass(1:nlayer mx)=zmass(ig,1:nlayermx)*(1.+pdpsrfmr(ig)*ptimestep/pplev(ig,1))213 zzmass(1:nlayer)=zmass(ig,1:nlayer)*(1.+pdpsrfmr(ig)*ptimestep/pplev(ig,1)) 205 214 206 215 … … 219 228 220 229 ! Van Leer scheme: 221 w(1:nlayer mx+1)=-zmflux(1:nlayermx+1)*ptimestep230 w(1:nlayer+1)=-zmflux(1:nlayer+1)*ptimestep 222 231 call vl1d(zzt,2.,zzmass,w,ztm) 223 232 call vl1d(zzu ,2.,zzmass,w,zum) 224 233 call vl1d(zzv ,2.,zzmass,w,zvm) 225 234 do iq=1,nq 226 zq1(1:nlayer mx)=zzq(1:nlayermx,iq)235 zq1(1:nlayer)=zzq(1:nlayer,iq) 227 236 zqm1(1)=zqm(1,iq) 228 237 ! print*,iq -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/mod_grid_phy_lmdz.F90
r222 r227 12 12 INTEGER,SAVE :: nbp_lev ! == llm 13 13 INTEGER,SAVE :: klon_glo 14 !$OMP THREADPRIVATE(nbp_lon,nbp_lat,nbp_lev,klon_glo) 14 15 15 16 INTERFACE grid1dTo2d_glo -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/mod_phys_lmdz_mpi_data.F90
r222 r227 16 16 INTEGER,SAVE :: klon_mpi_end 17 17 INTEGER,SAVE :: klon_mpi 18 !!$OMP THREADPRIVATE(ii_begin,ii_end,jj_begin,jj_end,jj_nb,ij_begin,& 19 ! !$OMP ij_end,ij_nb,klon_mpi_begin,klon_mpi_end,klon_mpi) 18 20 19 21 INTEGER,SAVE,ALLOCATABLE,DIMENSION(:) :: jj_para_nb … … 31 33 INTEGER,SAVE,ALLOCATABLE,DIMENSION(:) :: klon_mpi_para_begin 32 34 INTEGER,SAVE,ALLOCATABLE,DIMENSION(:) :: klon_mpi_para_end 35 !!$OMP THREADPRIVATE(jj_para_nb,jj_para_begin,jj_para_end,ii_para_begin,ii_para_end,& 36 ! !$OMP ij_para_nb,ij_para_begin,ij_para_end,klon_mpi_para_nb,klon_mpi_para_begin,& 37 ! !$OMP klon_mpi_para_end) 33 38 34 39 … … 38 43 LOGICAL,SAVE :: is_mpi_root 39 44 LOGICAL,SAVE :: is_using_mpi 45 !!$OMP THREADPRIVATE(mpi_rank,mpi_size,mpi_root,is_mpi_root,is_using_mpi) 40 46 41 47 … … 43 49 LOGICAL,SAVE :: is_south_pole 44 50 INTEGER,SAVE :: COMM_LMDZ_PHY 51 !!$OMP THREADPRIVATE(is_north_pole,is_south_pole,COMM_LMDZ_PHY) 45 52 46 53 CONTAINS -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/mod_phys_lmdz_para.F90
r222 r227 13 13 14 14 !$OMP THREADPRIVATE(klon_loc,is_master) 15 !$OMP THREADPRIVATE(is_sequential,is_parallel) 15 16 16 17 CONTAINS -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/moistadj.F90
r222 r227 1 subroutine moistadj(ngrid, n q, pt, pq, pdq, pplev, pplay, pdtmana, pdqmana, ptimestep, rneb)1 subroutine moistadj(ngrid, nlayer, nq, pt, pq, pdq, pplev, pplay, pdtmana, pdqmana, ptimestep, rneb) 2 2 3 3 use watercommon_h, only: T_h2O_ice_liq, RLVTT, RCPD, RCPV, Psat_water, Lcpdqsat_water 4 USE tracer_h 4 USE tracer_h, only: igcm_h2o_vap, igcm_h2o_ice 5 5 6 6 implicit none … … 20 20 !===================================================================== 21 21 22 #include "dimensions.h"23 #include "dimphys.h"22 !#include "dimensions.h" 23 !#include "dimphys.h" 24 24 #include "comcstfi.h" 25 25 26 INTEGER ngrid, nq 27 28 REAL pt(ngrid,nlayermx) ! temperature (K) 29 REAL pq(ngrid,nlayermx,nq) ! tracer (kg/kg) 30 REAL pdq(ngrid,nlayermx,nq) 31 32 REAL pdqmana(ngrid,nlayermx,nq) ! tendency of tracers (kg/kg.s-1) 33 REAL pdtmana(ngrid,nlayermx) ! temperature increment 26 INTEGER,INTENT(IN) :: ngrid, nlayer, nq 27 28 REAL,INTENT(IN) :: pt(ngrid,nlayer) ! temperature (K) 29 REAL,INTENT(IN) :: pq(ngrid,nlayer,nq) ! tracer (kg/kg) 30 REAL,INTENT(IN) :: pdq(ngrid,nlayer,nq) 31 REAL,INTENT(IN) :: pplev(ngrid,nlayer+1) ! inter-layer pressure (Pa) 32 REAL,INTENT(IN) :: pplay(ngrid,nlayer) ! mid-layer pressure (Pa) 33 REAL,INTENT(IN) :: ptimestep ! physics timestep (s) 34 REAL,INTENT(OUT) :: pdqmana(ngrid,nlayer,nq) ! tracer tendencies (kg/kg.s-1) 35 REAL,INTENT(OUT) :: pdtmana(ngrid,nlayer) ! temperature increment(K/s) 36 REAL,INTENT(OUT) :: rneb(ngrid,nlayer) ! cloud fraction 34 37 35 38 ! local variables 36 REAL zt(ngrid,nlayermx) ! temperature (K) 37 REAL zq(ngrid,nlayermx) ! humidite specifique (kg/kg) 38 REAL pplev(ngrid,nlayermx+1) ! pression a inter-couche (Pa) 39 REAL pplay(ngrid,nlayermx) ! pression au milieu de couche (Pa) 40 41 REAL d_t(ngrid,nlayermx) ! temperature increment 42 REAL d_q(ngrid,nlayermx) ! incrementation pour vapeur d'eau 43 REAL d_ql(ngrid,nlayermx) ! incrementation pour l'eau liquide 44 REAL rneb(ngrid,nlayermx) ! cloud fraction 45 REAL ptimestep 39 REAL zt(ngrid,nlayer) ! temperature (K) 40 REAL zq(ngrid,nlayer) ! humidite specifique (kg/kg) 41 42 REAL d_t(ngrid,nlayer) ! temperature increment 43 REAL d_q(ngrid,nlayer) ! incrementation pour vapeur d'eau 44 REAL d_ql(ngrid,nlayer) ! incrementation pour l'eau liquide 46 45 47 46 ! REAL t_coup … … 55 54 INTEGER k1, k1p, k2, k2p 56 55 LOGICAL itest(ngrid) 57 REAL delta_q(ngrid, nlayer mx)58 DOUBLE PRECISION :: cp_new_t(nlayer mx), v_cptt(ngrid,nlayermx)59 REAL cp_delta_t(nlayer mx)60 DOUBLE PRECISION :: v_cptj(nlayer mx), v_cptjk1, v_ssig56 REAL delta_q(ngrid, nlayer) 57 DOUBLE PRECISION :: cp_new_t(nlayer), v_cptt(ngrid,nlayer) 58 REAL cp_delta_t(nlayer) 59 DOUBLE PRECISION :: v_cptj(nlayer), v_cptjk1, v_ssig 61 60 REAL v_p, v_t, v_zqs,v_cptt2,v_pratio,v_dlnpsat 62 REAL zqs(ngrid,nlayer mx), zdqs(ngrid,nlayermx),zpsat(ngrid,nlayermx),zdlnpsat(ngrid,nlayermx)61 REAL zqs(ngrid,nlayer), zdqs(ngrid,nlayer),zpsat(ngrid,nlayer),zdlnpsat(ngrid,nlayer) 63 62 REAL zq1(ngrid), zq2(ngrid) 64 DOUBLE PRECISION :: gamcpdz(ngrid,2:nlayer mx)63 DOUBLE PRECISION :: gamcpdz(ngrid,2:nlayer) 65 64 DOUBLE PRECISION :: zdp, zdpm 66 65 … … 68 67 REAL zflo ! flotabilite 69 68 70 DOUBLE PRECISION :: local_q(ngrid,nlayer mx),local_t(ngrid,nlayermx)69 DOUBLE PRECISION :: local_q(ngrid,nlayer),local_t(ngrid,nlayer) 71 70 72 71 REAL zdelta, zcor, zcvm5 … … 78 77 INTEGER,SAVE :: i_h2o=0 ! water vapour 79 78 INTEGER,SAVE :: i_ice=0 ! water ice 80 81 LOGICAL firstcall 82 SAVE firstcall 83 84 DATA firstcall /.TRUE./ 79 !$OMP THREADPRIVATE(i_h2o,i_ice) 80 81 LOGICAL,SAVE :: firstcall=.TRUE. 82 !$OMP THREADPRIVATE(firstcall) 85 83 86 84 IF (firstcall) THEN … … 96 94 97 95 ! GCM -----> subroutine variables 98 zq(1:ngrid,1:nlayer mx) = pq(1:ngrid,1:nlayermx,i_h2o)+ pdq(1:ngrid,1:nlayermx,i_h2o)*ptimestep99 zt(1:ngrid,1:nlayer mx) = pt(1:ngrid,1:nlayermx)100 pdqmana(1:ngrid,1:nlayer mx,1:nq)=0.0101 102 DO k = 1, nlayer mx96 zq(1:ngrid,1:nlayer) = pq(1:ngrid,1:nlayer,i_h2o)+ pdq(1:ngrid,1:nlayer,i_h2o)*ptimestep 97 zt(1:ngrid,1:nlayer) = pt(1:ngrid,1:nlayer) 98 pdqmana(1:ngrid,1:nlayer,1:nq)=0.0 99 100 DO k = 1, nlayer 103 101 DO i = 1, ngrid 104 102 if(zq(i,k).lt.0.)then … … 108 106 ENDDO 109 107 110 local_q(1:ngrid,1:nlayer mx) = zq(1:ngrid,1:nlayermx)111 local_t(1:ngrid,1:nlayer mx) = zt(1:ngrid,1:nlayermx)112 rneb(1:ngrid,1:nlayer mx) = 0.0113 d_ql(1:ngrid,1:nlayer mx) = 0.0114 d_t(1:ngrid,1:nlayer mx) = 0.0115 d_q(1:ngrid,1:nlayer mx) = 0.0108 local_q(1:ngrid,1:nlayer) = zq(1:ngrid,1:nlayer) 109 local_t(1:ngrid,1:nlayer) = zt(1:ngrid,1:nlayer) 110 rneb(1:ngrid,1:nlayer) = 0.0 111 d_ql(1:ngrid,1:nlayer) = 0.0 112 d_t(1:ngrid,1:nlayer) = 0.0 113 d_q(1:ngrid,1:nlayer) = 0.0 116 114 117 115 ! Calculate v_cptt 118 DO k = 1, nlayer mx116 DO k = 1, nlayer 119 117 DO i = 1, ngrid 120 118 v_cptt(i,k) = RCPD * local_t(i,k) … … 128 126 129 127 ! Calculate Gamma * Cp * dz: (gamma is the critical gradient) 130 DO k = 2, nlayer mx128 DO k = 2, nlayer 131 129 DO i = 1, ngrid 132 130 zdp = pplev(i,k)-pplev(i,k+1) … … 159 157 810 CONTINUE ! look for k1, the base of the column 160 158 k2 = k2 + 1 161 IF (k2 .GT. nlayer mx) GOTO 9999159 IF (k2 .GT. nlayer) GOTO 9999 162 160 zflo = v_cptt(i,k2-1) - v_cptt(i,k2) - gamcpdz(i,k2) 163 161 zsat=(local_q(i,k2-1)-zqs(i,k2-1))*(pplev(i,k2-1)-pplev(i,k2)) & … … 169 167 170 168 820 CONTINUE !! look for k2, the top of the column 171 IF (k2 .EQ. nlayer mx) GOTO 821169 IF (k2 .EQ. nlayer) GOTO 821 172 170 k2p = k2 + 1 173 171 zsat=zsat+(pplev(i,k2p)-pplev(i,k2p+1))*(local_q(i,k2p)-zqs(i,k2p)) … … 227 225 ! ENDDO 228 226 229 DO k = 2, nlayer mx227 DO k = 2, nlayer 230 228 zdpm = pplev(i,k-1) - pplev(i,k) 231 229 zdp = pplev(i,k) - pplev(i,k+1) … … 272 270 ! a l'endroit ou la vapeur d'eau est diminuee par l'ajustement): 273 271 274 DO k = 1, nlayer mx272 DO k = 1, nlayer 275 273 DO i = 1, ngrid 276 274 IF (itest(i)) THEN … … 291 289 ENDIF 292 290 ENDDO 293 DO k = 1, nlayer mx291 DO k = 1, nlayer 294 292 DO i = 1, ngrid 295 293 IF (itest(i)) THEN … … 300 298 ENDDO 301 299 ENDDO 302 DO k = 1, nlayer mx300 DO k = 1, nlayer 303 301 DO i = 1, ngrid 304 302 IF (itest(i)) THEN … … 308 306 ENDDO 309 307 310 DO k = 1, nlayer mx308 DO k = 1, nlayer 311 309 DO i = 1, ngrid 312 310 local_q(i, k) = MAX(local_q(i, k), seuil_vap) … … 314 312 ENDDO 315 313 316 DO k = 1, nlayer mx314 DO k = 1, nlayer 317 315 DO i = 1, ngrid 318 316 d_t(i,k) = local_t(i,k) - zt(i,k) … … 322 320 323 321 ! now subroutine -----> GCM variables 324 DO k = 1, nlayer mx322 DO k = 1, nlayer 325 323 DO i = 1, ngrid 326 324 … … 333 331 334 332 335 RETURN336 333 END -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/newsedim.F
r222 r227 16 16 ! ------------ 17 17 18 #include "dimensions.h"19 #include "dimphys.h"18 !#include "dimensions.h" 19 !#include "dimphys.h" 20 20 #include "comcstfi.h" 21 21 … … 44 44 45 45 LOGICAL,SAVE :: firstcall=.true. 46 !$OMP THREADPRIVATE(firstcall) 46 47 47 48 c Traceurs : 48 49 c ~~~~~~~~ 49 real traversee (ngrid,nlay ermx)50 real vstokes(ngrid,nlay ermx)51 real w(ngrid,nlay ermx)50 real traversee (ngrid,nlay) 51 real vstokes(ngrid,nlay) 52 real w(ngrid,nlay) 52 53 real ptop, dztop, Ep, Stra 53 54 … … 62 63 c local and saved variable 63 64 real,save :: a,b 65 !$OMP THREADPRIVATE(a,b) 64 66 65 67 c ** un petit test de coherence … … 191 193 end do 192 194 193 call vlz_fi(ngrid, pqi,2.,masse,w,wq)195 call vlz_fi(ngrid,nlay,pqi,2.,masse,w,wq) 194 196 c write(*,*) ' newsed: wq(6), wq(7), q(6)', 195 197 c & wq(1,6),wq(1,7),pqi(1,6) 196 198 197 RETURN198 199 END -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/newstart.F
r222 r227 23 23 use datafile_mod, only: datadir 24 24 ! to use 'getin' 25 USE ioipsl_getincom, only: getin 25 ! USE ioipsl_getincom, only: getin 26 USE ioipsl_getincom_p, only: getin_p 26 27 use control_mod, only: day_step, iphysiq, anneeref 27 28 use phyredem, only: physdem0, physdem1 … … 32 33 33 34 #include "dimensions.h" 34 #include "dimphys.h" 35 #include "planete.h" 35 !#include "dimphys.h" 36 integer, parameter :: ngridmx = (2+(jjm-1)*iim - 1/jjm) 37 !#include "planete.h" 36 38 #include "paramet.h" 37 39 #include "comconst.h" … … 98 100 real emisread ! added by RW 99 101 REAL,ALLOCATABLE :: qsurf(:,:) 100 REAL q2(ngridmx, nlayermx+1)102 REAL q2(ngridmx,llm+1) 101 103 ! REAL rnaturfi(ngridmx) 102 104 real alb(iip1,jjp1),albfi(ngridmx) ! albedos … … 171 173 172 174 ! added by BC for cloud fraction setup 173 REAL hice(ngridmx),cloudfrac(ngridmx, nlayermx)175 REAL hice(ngridmx),cloudfrac(ngridmx,llm) 174 176 REAL totalfrac(ngridmx) 175 177 … … 341 343 write(*,*) 'Reading file STARTFI' 342 344 fichnom = 'startfi.nc' 343 CALL phyetat0 (ngridmx, fichnom,tab0,Lmodif,nsoilmx,nqtot,344 . day_ini,time,345 CALL phyetat0 (ngridmx,llm,fichnom,tab0,Lmodif,nsoilmx, 346 . nqtot,day_ini,time, 345 347 . tsurf,tsoil,emis,q2,qsurf, !) ! temporary modif by RDW 346 348 . cloudfrac,totalfrac,hice,rnat,pctsrf_sic,tslab,tsea_ice, … … 480 482 ! First get the correct value of datadir, if not already done: 481 483 ! default 'datadir' is set in "datafile_mod" 482 call getin ("datadir",datadir)484 call getin_p("datadir",datadir) 483 485 write(*,*) 'Available surface data files are:' 484 486 filestring='ls '//trim(datadir)//' | grep .nc' … … 534 536 535 537 write(*,*) 'Reading file START_ARCHIVE' 536 CALL lect_start_archive(date,tsurf,tsoil,emis,q2, 538 CALL lect_start_archive(ngridmx,llm, 539 & date,tsurf,tsoil,emis,q2, 537 540 & t,ucov,vcov,ps,teta,phisold_newgrid,q,qsurf, 538 541 & surfith,nid, … … 956 959 ucov(1:iip1,1:jjp1,1:llm)=0. 957 960 vcov(1:iip1,1:jjm,1:llm)=0. 958 q2(1:ngridmx,1: nlayermx+1)=0.961 q2(1:ngridmx,1:llm+1)=0. 959 962 else 960 963 write(*,*)'problem reading file ',trim(txt),' !' … … 1296 1299 ucov(1:iip1,1:jjp1,1:llm)=0 1297 1300 vcov(1:iip1,1:jjm,1:llm)=0 1298 q2(1:ngridmx,1: nlayermx+1)=01301 q2(1:ngridmx,1:llm+1)=0 1299 1302 1300 1303 c radequi : Radiative equilibrium profile of temperatures and no winds … … 1326 1329 ucov(1:iip1,1:jjp1,1:llm)=0 1327 1330 vcov(1:iip1,1:jjm,1:llm)=0 1328 q2(1:ngridmx,1: nlayermx+1)=01331 q2(1:ngridmx,1:llm+1)=0 1329 1332 1330 1333 c coldstart : T set 1K above CO2 frost point and no winds … … 1362 1365 ucov(1:iip1,1:jjp1,1:llm)=0 1363 1366 vcov(1:iip1,1:jjm,1:llm)=0 1364 q2(1:ngridmx,1: nlayermx+1)=01367 q2(1:ngridmx,1:llm+1)=0 1365 1368 1366 1369 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/newtrelax.F90
r222 r227 1 subroutine newtrelax(ngrid, mu0,sinlat,popsk,temp,pplay,pplev,dtrad,firstcall)1 subroutine newtrelax(ngrid,nlayer,mu0,sinlat,popsk,temp,pplay,pplev,dtrad,firstcall) 2 2 3 3 implicit none 4 4 5 #include "dimensions.h"6 #include "dimphys.h"5 !#include "dimensions.h" 6 !#include "dimphys.h" 7 7 #include "comcstfi.h" 8 8 #include "callkeys.h" … … 21 21 !================================================================== 22 22 23 integer ngrid24 23 25 24 ! Input 26 real mu0(ngrid) ! cosine of sun incident angle 27 real sinlat(ngrid) ! sine of latitude 28 real temp(ngrid,nlayermx) ! temperature at each layer (K) 29 real pplay(ngrid,nlayermx) ! pressure at each layer (Pa) 30 real pplev(ngrid,nlayermx+1) ! pressure at each level (Pa) 31 real popsk(ngrid,nlayermx) ! pot. T to T converter 25 integer,intent(in) :: ngrid, nlayer 26 logical,intent(in) :: firstcall 27 real,intent(in) :: mu0(ngrid) ! cosine of sun incident angle 28 real,intent(in) :: sinlat(ngrid) ! sine of latitude 29 real,intent(in) :: temp(ngrid,nlayer) ! temperature at each layer (K) 30 real,intent(in) :: pplay(ngrid,nlayer) ! pressure at each layer (Pa) 31 real,intent(in) :: pplev(ngrid,nlayer+1) ! pressure at each level (Pa) 32 real,intent(in) :: popsk(ngrid,nlayer) ! pot. T to T converter 32 33 33 34 ! Output 34 real dtrad(ngrid,nlayermx)35 real,intent(out) :: dtrad(ngrid,nlayer) 35 36 36 37 ! Internal 37 38 real Trelax_V, Trelax_H 38 39 real,allocatable,dimension(:,:),save :: Trelax 40 !$OMP THREADPRIVATE(Trelax) 39 41 40 42 real T_trop ! relaxation temperature at tropopause (K) … … 44 46 real sig, f_sig, sig_trop 45 47 integer l,ig 46 logical firstcall47 48 48 49 … … 53 54 if(firstcall) then 54 55 55 ALLOCATE(Trelax(ngrid,nlayer mx))56 ALLOCATE(Trelax(ngrid,nlayer)) 56 57 57 58 print*,'-----------------------------------------------------' … … 66 67 T_surf = 126. + 239.*mu0(ig) 67 68 T_trop = 140. + 52.*mu0(ig) 68 do l=1,nlayer mx69 do l=1,nlayer 69 70 70 71 if(mu0(ig).le.0.0)then ! night side … … 86 87 sig_trop=(T_trop/T_surf)**(1./rcp) 87 88 88 do l=1,nlayer mx89 do l=1,nlayer 89 90 do ig=1,ngrid 90 91 … … 109 110 endif 110 111 111 firstcall=.false.112 113 112 endif 114 113 115 114 ! Calculate radiative forcing 116 do l=1,nlayer mx115 do l=1,nlayer 117 116 do ig=1,ngrid 118 117 dtrad(ig,l) = -(temp(ig,l) - Trelax(ig,l)) / tau_relax … … 129 128 !call writediagfi(ngrid,'ThetaZ','stellar zenith angle','deg',2,mu0) 130 129 131 return132 130 end subroutine newtrelax -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/ocean_slab_mod.F90
r222 r227 15 15 16 16 17 #include "dimensions.h"18 #include "dimphys.h"17 !#include "dimensions.h" 18 !#include "dimphys.h" 19 19 #include "comcstfi.h" 20 20 #include "callkeys.h" … … 26 26 27 27 !LOGICAL, PRIVATE, SAVE :: ok_slab_sic,ok_slab_heaT_h2O_ice_liqBS 28 ! $OMP THREADPRIVATE(ok_slab_sic,ok_slab_heaT_h2O_ice_liqBS)28 !!$OMP THREADPRIVATE(ok_slab_sic,ok_slab_heaT_h2O_ice_liqBS) 29 29 !INTEGER, PRIVATE, SAVE :: slab_ekman, slab_cadj 30 ! $OMP THREADPRIVATE(slab_ekman,slab_cadj)30 !!$OMP THREADPRIVATE(slab_ekman,slab_cadj) 31 31 INTEGER, PRIVATE, SAVE :: lmt_pas, julien, idayvrai 32 32 !$OMP THREADPRIVATE(lmt_pas,julien,idayvrai) … … 67 67 68 68 69 #include "dimensions.h"70 #include "dimphys.h"69 !#include "dimensions.h" 70 !#include "dimphys.h" 71 71 #include "comcstfi.h" 72 72 #include "callkeys.h" … … 210 210 use slab_ice_h 211 211 212 #include "dimensions.h"213 #include "dimphys.h"212 !#include "dimensions.h" 213 !#include "dimphys.h" 214 214 #include "comcstfi.h" 215 215 #include "callkeys.h" … … 451 451 use slab_ice_h 452 452 453 #include "dimensions.h"454 #include "comcstfi.h"453 !#include "dimensions.h" 454 !#include "comcstfi.h" 455 455 ! INCLUDE "iniprint.h" 456 456 #include "callkeys.h" … … 534 534 ! Transport diffusif 535 535 ! IF (ok_soil_hdiff) THEN 536 CALL divgrad_phy(n oceanmx,tmp_tslab_loc,dtdiff_loc)536 CALL divgrad_phy(ngrid,noceanmx,tmp_tslab_loc,dtdiff_loc) 537 537 dtdiff_loc=dtdiff_loc*soil_hdiff*50./SUM(slabh)!*100. 538 538 ! dtdiff_loc(:,1)=dtdiff_loc(:,1)*soil_hdiff*50./SUM(slabh)*0.8 … … 544 544 ! Calcul de transport par Ekman 545 545 546 CALL slab_ekman2( taux_slab,tauy_slab,tslab,dtekman_loc)546 CALL slab_ekman2(ngrid,taux_slab,tauy_slab,tslab,dtekman_loc) 547 547 548 548 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/orbite.F
r222 r227 1 1 subroutine orbite(pls,pdist_star,pdecli) 2 3 use planete_mod, only: p_elips, e_elips, timeperi, obliquit 2 4 implicit none 3 5 … … 23 25 c ------------- 24 26 25 #include "planete.h"27 !#include "planete.h" 26 28 #include "comcstfi.h" 27 29 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/phyetat0.F90
r222 r227 1 subroutine phyetat0 (ngrid, fichnom,tab0,Lmodif,nsoil,nq, &1 subroutine phyetat0 (ngrid,nlayer,fichnom,tab0,Lmodif,nsoil,nq, & 2 2 day_ini,time,tsurf,tsoil, & 3 3 emis,q2,qsurf,cloudfrac,totcloudfrac,hice, & … … 20 20 !====================================================================== 21 21 #include "netcdf.inc" 22 #include "dimensions.h"23 #include "dimphys.h"24 #include "planete.h"22 !#include "dimensions.h" 23 !#include "dimphys.h" 24 !#include "planete.h" 25 25 #include "comcstfi.h" 26 26 … … 33 33 ! inputs: 34 34 integer,intent(in) :: ngrid 35 integer,intent(in) :: nlayer 35 36 character*(*),intent(in) :: fichnom ! "startfi.nc" file 36 37 integer,intent(in) :: tab0 … … 45 46 real,intent(out) :: tsoil(ngrid,nsoil) ! soil temperature 46 47 real,intent(out) :: emis(ngrid) ! surface emissivity 47 real,intent(out) :: q2(ngrid, llm+1) !48 real,intent(out) :: q2(ngrid,nlayer+1) ! 48 49 real,intent(out) :: qsurf(ngrid,nq) ! tracers on surface 49 50 ! real co2ice(ngrid) ! co2 ice cover 50 real,intent(out) :: cloudfrac(ngrid,nlayer mx)51 real,intent(out) :: cloudfrac(ngrid,nlayer) 51 52 real,intent(out) :: hice(ngrid), totcloudfrac(ngrid) 52 53 real,intent(out) :: pctsrf_sic(ngrid),tslab(ngrid,noceanmx) -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/phyredem.F90
r222 r227 17 17 put_var, put_field, length 18 18 use mod_grid_phy_lmdz, only : klon_glo 19 use planete_mod, only: year_day, periastr, apoastr, peri_day, & 20 obliquit, z0, lmixmin, emin_turb 19 21 20 22 implicit none 21 #include "planete.h"23 !#include "planete.h" 22 24 #include "comcstfi.h" 23 25 character(len=*), intent(in) :: filename … … 147 149 !#include "temps.h" 148 150 #include "comcstfi.h" 149 #include "planete.h"151 !#include "planete.h" 150 152 #include "callkeys.h" 151 153 !====================================================================== -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/physiq.F90
r222 r227 5 5 pplev,pplay,pphi, & 6 6 pu,pv,pt,pq, & 7 pw,&7 flxw, & 8 8 pdu,pdv,pdt,pdq,pdpsrf,tracerdyn) 9 9 … … 12 12 use gases_h, only: gnom, gfrac 13 13 use radcommon_h, only: sigma, eclipse, glat, grav 14 use radii_mod, only: h2o_reffrad, co2_reffrad , h2o_cloudrad14 use radii_mod, only: h2o_reffrad, co2_reffrad 15 15 use aerosol_mod, only: iaero_co2, iaero_h2o 16 16 use surfdat_h, only: phisfi, albedodat, zmea, zstd, zsig, zgam, zthe, & … … 32 32 use planetwide_mod, only: planetwide_minval,planetwide_maxval,planetwide_sumval 33 33 use mod_phys_lmdz_para, only : is_master 34 34 use planete_mod, only: apoastr, periastr, year_day, peri_day, & 35 obliquit, nres, z0 35 36 36 37 implicit none … … 97 98 ! ptdyn(ngrid,nlayer) / corresponding variables 98 99 ! pqdyn(ngrid,nlayer,nq) / 99 ! pw(ngrid,?) vertical velocity100 ! flxw(ngrid,nlayer) vertical mass flux (kg/s) at layer lower boundary 100 101 ! 101 102 ! output 102 103 ! ------ 103 104 ! 104 ! pdu(ngrid,nlayer mx) \105 ! pdv(ngrid,nlayer mx) \ Temporal derivative of the corresponding106 ! pdt(ngrid,nlayer mx) / variables due to physical processes.107 ! pdq(ngrid,nlayer mx) /105 ! pdu(ngrid,nlayer) \ 106 ! pdv(ngrid,nlayer) \ Temporal derivative of the corresponding 107 ! pdt(ngrid,nlayer) / variables due to physical processes. 108 ! pdq(ngrid,nlayer) / 108 109 ! pdpsrf(ngrid) / 109 110 ! tracerdyn call tracer in dynamical part of GCM ? … … 131 132 ! ------------------ 132 133 133 #include "dimensions.h"134 #include "dimphys.h"134 !#include "dimensions.h" 135 !#include "dimphys.h" 135 136 #include "callkeys.h" 136 137 #include "comcstfi.h" 137 #include "planete.h"138 !#include "planete.h" 138 139 !#include "control.h" 139 140 #include "netcdf.inc" … … 160 161 real,intent(in) :: pt(ngrid,nlayer) ! temperature (K) 161 162 real,intent(in) :: pq(ngrid,nlayer,nq) ! tracers (.../kg_of_air) 162 real,intent(in) :: pw(ngrid,nlayer) ! vertical velocity (m/s)163 163 real,intent(in) :: flxw(ngrid,nlayer) ! vertical mass flux (ks/s) 164 ! at lower boundary of layer 164 165 165 166 … … 179 180 ! "longrefvis" set in dimradmars.h , for one of the "naerkind" kind of 180 181 ! aerosol optical properties: 181 ! real aerosol(ngrid,nlayer mx,naerkind)182 ! real aerosol(ngrid,nlayer,naerkind) 182 183 ! this is now internal to callcorrk and hence no longer needed here 183 184 … … 187 188 real, dimension(:,:),allocatable,save :: tsoil ! sub-surface temperatures (K) 188 189 real, dimension(:),allocatable,save :: albedo ! Surface albedo 190 !$OMP THREADPRIVATE(tsurf,tsoil,albedo) 189 191 190 192 real,dimension(:),allocatable,save :: albedo0 ! Surface albedo 191 193 real,dimension(:),allocatable,save :: rnat ! added by BC 194 !$OMP THREADPRIVATE(albedo0,rnat) 192 195 193 196 real,dimension(:),allocatable,save :: emis ! Thermal IR surface emissivity … … 199 202 real,dimension(:,:),allocatable,save :: qsurf ! tracer on surface (e.g. kg.m-2) 200 203 real,dimension(:,:),allocatable,save :: q2 ! Turbulent Kinetic Energy 204 !$OMP THREADPRIVATE(emis,dtrad,fluxrad_sky,fluxrad,capcal,fluxgrd,qsurf,q2) 201 205 202 206 save day_ini, icount 207 !$OMP THREADPRIVATE(day_ini,icount) 203 208 204 209 ! Local variables : … … 207 212 ! aerosol (dust or ice) extinction optical depth at reference wavelength 208 213 ! for the "naerkind" optically active aerosols: 209 real aerosol(ngrid,nlayermx,naerkind) 210 real zh(ngrid,nlayermx) ! potential temperature (K) 214 real aerosol(ngrid,nlayer,naerkind) 215 real zh(ngrid,nlayer) ! potential temperature (K) 216 real pw(ngrid,nlayer) ! vertical velocity (m/s) (>0 when downwards) 211 217 212 218 character*80 fichier … … 225 231 real,dimension(:,:),allocatable,save :: zdtsw ! (K/s) 226 232 real,dimension(:),allocatable,save :: sensibFlux ! turbulent flux given by the atm to the surface 233 !$OMP THREADPRIVATE(fluxsurf_lw,fluxsurf_sw,fluxtop_lw,fluxabs_sw,fluxtop_dn,fluxdyn,OLR_nu,OSR_nu,& 234 !$OMP zdtlw,zdtsw,sensibFlux) 227 235 228 236 real zls ! solar longitude (rad) 229 237 real zday ! date (time since Ls=0, in martian days) 230 real zzlay(ngrid,nlayer mx) ! altitude at the middle of the layers231 real zzlev(ngrid,nlayer mx+1) ! altitude at layer boundaries238 real zzlay(ngrid,nlayer) ! altitude at the middle of the layers 239 real zzlev(ngrid,nlayer+1) ! altitude at layer boundaries 232 240 real latvl1,lonvl1 ! Viking Lander 1 point (for diagnostic) 233 241 234 242 ! Tendencies due to various processes: 235 243 real dqsurf(ngrid,nq) 236 real cldtlw(ngrid,nlayer mx) ! (K/s) LW heating rate for clear areas237 real cldtsw(ngrid,nlayer mx) ! (K/s) SW heating rate for clear areas244 real cldtlw(ngrid,nlayer) ! (K/s) LW heating rate for clear areas 245 real cldtsw(ngrid,nlayer) ! (K/s) SW heating rate for clear areas 238 246 real zdtsurf(ngrid) ! (K/s) 239 real dtlscale(ngrid,nlayer mx)240 real zdvdif(ngrid,nlayer mx),zdudif(ngrid,nlayermx) ! (m.s-2)241 real zdhdif(ngrid,nlayer mx), zdtsdif(ngrid) ! (K/s)242 real zdtdif(ngrid,nlayer mx) ! (K/s)243 real zdvadj(ngrid,nlayer mx),zduadj(ngrid,nlayermx) ! (m.s-2)244 real zdhadj(ngrid,nlayer mx) ! (K/s)245 real zdtgw(ngrid,nlayer mx) ! (K/s)246 real zdtmr(ngrid,nlayer mx)247 real zdugw(ngrid,nlayer mx),zdvgw(ngrid,nlayermx) ! (m.s-2)248 real zdtc(ngrid,nlayer mx),zdtsurfc(ngrid)249 real zdvc(ngrid,nlayer mx),zduc(ngrid,nlayermx)250 real zdumr(ngrid,nlayer mx),zdvmr(ngrid,nlayermx)247 real dtlscale(ngrid,nlayer) 248 real zdvdif(ngrid,nlayer),zdudif(ngrid,nlayer) ! (m.s-2) 249 real zdhdif(ngrid,nlayer), zdtsdif(ngrid) ! (K/s) 250 real zdtdif(ngrid,nlayer) ! (K/s) 251 real zdvadj(ngrid,nlayer),zduadj(ngrid,nlayer) ! (m.s-2) 252 real zdhadj(ngrid,nlayer) ! (K/s) 253 real zdtgw(ngrid,nlayer) ! (K/s) 254 real zdtmr(ngrid,nlayer) 255 real zdugw(ngrid,nlayer),zdvgw(ngrid,nlayer) ! (m.s-2) 256 real zdtc(ngrid,nlayer),zdtsurfc(ngrid) 257 real zdvc(ngrid,nlayer),zduc(ngrid,nlayer) 258 real zdumr(ngrid,nlayer),zdvmr(ngrid,nlayer) 251 259 real zdtsurfmr(ngrid) 252 260 253 real zdmassmr(ngrid,nlayer mx),zdpsrfmr(ngrid)261 real zdmassmr(ngrid,nlayer),zdpsrfmr(ngrid) 254 262 real zdmassmr_col(ngrid) 255 263 256 real zdqdif(ngrid,nlayer mx,nq), zdqsdif(ngrid,nq)257 real zdqevap(ngrid,nlayer mx)258 real zdqsed(ngrid,nlayer mx,nq), zdqssed(ngrid,nq)259 real zdqdev(ngrid,nlayer mx,nq), zdqsdev(ngrid,nq)260 real zdqadj(ngrid,nlayer mx,nq)261 real zdqc(ngrid,nlayer mx,nq)262 real zdqmr(ngrid,nlayer mx,nq),zdqsurfmr(ngrid,nq)263 real zdqlscale(ngrid,nlayer mx,nq)264 real zdqdif(ngrid,nlayer,nq), zdqsdif(ngrid,nq) 265 real zdqevap(ngrid,nlayer) 266 real zdqsed(ngrid,nlayer,nq), zdqssed(ngrid,nq) 267 real zdqdev(ngrid,nlayer,nq), zdqsdev(ngrid,nq) 268 real zdqadj(ngrid,nlayer,nq) 269 real zdqc(ngrid,nlayer,nq) 270 real zdqmr(ngrid,nlayer,nq),zdqsurfmr(ngrid,nq) 271 real zdqlscale(ngrid,nlayer,nq) 264 272 real zdqslscale(ngrid,nq) 265 real zdqchim(ngrid,nlayer mx,nq)273 real zdqchim(ngrid,nlayer,nq) 266 274 real zdqschim(ngrid,nq) 267 275 268 real zdteuv(ngrid,nlayer mx) ! (K/s)269 real zdtconduc(ngrid,nlayer mx) ! (K/s)270 real zdumolvis(ngrid,nlayer mx)271 real zdvmolvis(ngrid,nlayer mx)272 real zdqmoldiff(ngrid,nlayer mx,nq)276 real zdteuv(ngrid,nlayer) ! (K/s) 277 real zdtconduc(ngrid,nlayer) ! (K/s) 278 real zdumolvis(ngrid,nlayer) 279 real zdvmolvis(ngrid,nlayer) 280 real zdqmoldiff(ngrid,nlayer,nq) 273 281 274 282 ! Local variables for local calculations: 275 283 real zflubid(ngrid) 276 real zplanck(ngrid),zpopsk(ngrid,nlayer mx)277 real zdum1(ngrid,nlayer mx)278 real zdum2(ngrid,nlayer mx)284 real zplanck(ngrid),zpopsk(ngrid,nlayer) 285 real zdum1(ngrid,nlayer) 286 real zdum2(ngrid,nlayer) 279 287 real ztim1,ztim2,ztim3, z1,z2 280 288 real ztime_fin 281 real zdh(ngrid,nlayer mx)289 real zdh(ngrid,nlayer) 282 290 real gmplanet 283 291 real taux(ngrid),tauy(ngrid) … … 288 296 ! local variables only used for diagnostics (output in file "diagfi" or "stats") 289 297 ! ------------------------------------------------------------------------------ 290 real ps(ngrid), zt(ngrid,nlayer mx)291 real zu(ngrid,nlayer mx),zv(ngrid,nlayermx)292 real zq(ngrid,nlayer mx,nq)298 real ps(ngrid), zt(ngrid,nlayer) 299 real zu(ngrid,nlayer),zv(ngrid,nlayer) 300 real zq(ngrid,nlayer,nq) 293 301 character*2 str2 294 302 character*5 str5 295 real zdtadj(ngrid,nlayer mx)296 real zdtdyn(ngrid,nlayer mx)303 real zdtadj(ngrid,nlayer) 304 real zdtdyn(ngrid,nlayer) 297 305 real,allocatable,dimension(:,:),save :: ztprevious 298 real reff(ngrid,nlayermx) ! effective dust radius (used if doubleq=T) 306 !$OMP THREADPRIVATE(ztprevious) 307 real reff(ngrid,nlayer) ! effective dust radius (used if doubleq=T) 299 308 real qtot1,qtot2 ! total aerosol mass 300 309 integer igmin, lmin 301 310 logical tdiag 302 311 303 real zplev(ngrid,nlayer mx+1),zplay(ngrid,nlayermx)312 real zplev(ngrid,nlayer+1),zplay(ngrid,nlayer) 304 313 real zstress(ngrid), cd 305 real hco2(nq), tmean, zlocal(nlayer mx)306 real vmr(ngrid,nlayer mx) ! volume mixing ratio314 real hco2(nq), tmean, zlocal(nlayer) 315 real vmr(ngrid,nlayer) ! volume mixing ratio 307 316 308 317 real time_phys … … 310 319 ! reinstated by RW for diagnostic 311 320 real,allocatable,dimension(:),save :: tau_col 321 !$OMP THREADPRIVATE(tau_col) 312 322 313 323 ! included by RW to reduce insanity of code … … 318 328 319 329 ! included by RW for H2O precipitation 320 real zdtrain(ngrid,nlayer mx)321 real zdqrain(ngrid,nlayer mx,nq)330 real zdtrain(ngrid,nlayer) 331 real zdqrain(ngrid,nlayer,nq) 322 332 real zdqsrain(ngrid) 323 333 real zdqssnow(ngrid) … … 325 335 326 336 ! included by RW for H2O Manabe scheme 327 real dtmoist(ngrid,nlayer mx)328 real dqmoist(ngrid,nlayer mx,nq)329 330 real qvap(ngrid,nlayer mx)331 real dqvaplscale(ngrid,nlayer mx)332 real dqcldlscale(ngrid,nlayer mx)333 real rneb_man(ngrid,nlayer mx)334 real rneb_lsc(ngrid,nlayer mx)337 real dtmoist(ngrid,nlayer) 338 real dqmoist(ngrid,nlayer,nq) 339 340 real qvap(ngrid,nlayer) 341 real dqvaplscale(ngrid,nlayer) 342 real dqcldlscale(ngrid,nlayer) 343 real rneb_man(ngrid,nlayer) 344 real rneb_lsc(ngrid,nlayer) 335 345 336 346 ! included by RW to account for surface cooling by evaporation … … 339 349 340 350 ! to test energy conservation (RW+JL) 341 real mass(ngrid,nlayer mx),massarea(ngrid,nlayermx)351 real mass(ngrid,nlayer),massarea(ngrid,nlayer) 342 352 real dEtot, dEtots, AtmToSurf_TurbFlux 343 353 real,save :: dEtotSW, dEtotsSW, dEtotLW, dEtotsLW 344 real dEzRadsw(ngrid,nlayermx),dEzRadlw(ngrid,nlayermx),dEzdiff(ngrid,nlayermx) 354 !$OMP THREADPRIVATE(dEtotSW, dEtotsSW, dEtotLW, dEtotsLW) 355 real dEzRadsw(ngrid,nlayer),dEzRadlw(ngrid,nlayer),dEzdiff(ngrid,nlayer) 345 356 real dEdiffs(ngrid),dEdiff(ngrid) 346 real madjdE(ngrid), lscaledE(ngrid),madjdEz(ngrid,nlayer mx), lscaledEz(ngrid,nlayermx)357 real madjdE(ngrid), lscaledE(ngrid),madjdEz(ngrid,nlayer), lscaledEz(ngrid,nlayer) 347 358 !JL12 conservation test for mean flow kinetic energy has been disabled temporarily 348 359 real dtmoist_max,dtmoist_min … … 351 362 352 363 ! included by BC for evaporation 353 real qevap(ngrid,nlayer mx,nq)354 real tevap(ngrid,nlayer mx)355 real dqevap1(ngrid,nlayer mx)356 real dtevap1(ngrid,nlayer mx)364 real qevap(ngrid,nlayer,nq) 365 real tevap(ngrid,nlayer) 366 real dqevap1(ngrid,nlayer) 367 real dtevap1(ngrid,nlayer) 357 368 358 369 ! included by BC for hydrology 359 370 real,allocatable,save :: hice(:) 371 !$OMP THREADPRIVATE(hice) 360 372 361 373 ! included by RW to test water conservation (by routine) … … 365 377 logical watertest 366 378 save watertest 379 !$OMP THREADPRIVATE(watertest) 367 380 368 381 ! included by RW for RH diagnostic 369 real qsat(ngrid,nlayer mx), RH(ngrid,nlayermx), H2Omaxcol(ngrid),psat_tmp382 real qsat(ngrid,nlayer), RH(ngrid,nlayer), H2Omaxcol(ngrid),psat_tmp 370 383 371 384 ! included by RW for hydrology … … 378 391 ! included by BC for double radiative transfer call 379 392 logical clearsky 380 real zdtsw1(ngrid,nlayer mx)381 real zdtlw1(ngrid,nlayer mx)393 real zdtsw1(ngrid,nlayer) 394 real zdtlw1(ngrid,nlayer) 382 395 real fluxsurf_lw1(ngrid) 383 396 real fluxsurf_sw1(ngrid) … … 392 405 real,allocatable,dimension(:,:),save :: cloudfrac 393 406 real,allocatable,dimension(:),save :: totcloudfrac 407 !$OMP THREADPRIVATE(cloudfrac,totcloudfrac) 394 408 395 409 ! included by RW for vdifc water conservation test … … 400 414 ! double precision qsurf_hist(ngrid,nq) 401 415 real,allocatable,dimension(:,:),save :: qsurf_hist 416 !$OMP THREADPRIVATE(qsurf_hist) 402 417 403 418 ! included by RW for temp convadj conservation test 404 real playtest(nlayer mx)405 real plevtest(nlayer mx)406 real ttest(nlayer mx)407 real qtest(nlayer mx)419 real playtest(nlayer) 420 real plevtest(nlayer) 421 real ttest(nlayer) 422 real qtest(nlayer) 408 423 integer igtest 409 424 410 425 ! included by RW for runway greenhouse 1D study 411 real muvar(ngrid,nlayer mx+1)426 real muvar(ngrid,nlayer+1) 412 427 413 428 ! included by RW for variable H2O particle sizes 414 429 real,allocatable,dimension(:,:,:),save :: reffrad ! aerosol effective radius (m) 415 430 real,allocatable,dimension(:,:,:),save :: nueffrad ! aerosol effective radius variance 416 ! real :: nueffrad_dummy(ngrid,nlayermx,naerkind) !! AS. This is temporary. Check below why. 417 real :: reffh2oliq(ngrid,nlayermx) ! liquid water particles effective radius (m) 418 real :: reffh2oice(ngrid,nlayermx) ! water ice particles effective radius (m) 419 ! real reffH2O(ngrid,nlayermx) 431 !$OMP THREADPRIVATE(reffrad,nueffrad) 432 ! real :: nueffrad_dummy(ngrid,nlayer,naerkind) !! AS. This is temporary. Check below why. 433 real :: reffh2oliq(ngrid,nlayer) ! liquid water particles effective radius (m) 434 real :: reffh2oice(ngrid,nlayer) ! water ice particles effective radius (m) 435 ! real reffH2O(ngrid,nlayer) 420 436 real reffcol(ngrid,naerkind) 421 437 … … 427 443 integer num_run 428 444 logical,save :: ice_update 445 !$OMP THREADPRIVATE(ice_initial,ice_min,ice_update) 429 446 430 447 ! included by MS to compute the daily average of rings shadowing … … 445 462 real, dimension(:),allocatable,save :: zmasq 446 463 integer, dimension(:),allocatable,save ::knindex 464 !$OMP THREADPRIVATE(pctsrf_sic,tslab,tsea_ice,sea_ice,zmasq,knindex) 447 465 448 466 real :: tsurf2(ngrid) … … 469 487 ALLOCATE(rnat(ngrid)) 470 488 ALLOCATE(emis(ngrid)) 471 ALLOCATE(dtrad(ngrid,nlayer mx))489 ALLOCATE(dtrad(ngrid,nlayer)) 472 490 ALLOCATE(fluxrad_sky(ngrid)) 473 491 ALLOCATE(fluxrad(ngrid)) … … 475 493 ALLOCATE(fluxgrd(ngrid)) 476 494 ALLOCATE(qsurf(ngrid,nq)) 477 ALLOCATE(q2(ngrid,nlayer mx+1))478 ALLOCATE(ztprevious(ngrid,nlayer mx))479 ALLOCATE(cloudfrac(ngrid,nlayer mx))495 ALLOCATE(q2(ngrid,nlayer+1)) 496 ALLOCATE(ztprevious(ngrid,nlayer)) 497 ALLOCATE(cloudfrac(ngrid,nlayer)) 480 498 ALLOCATE(totcloudfrac(ngrid)) 481 499 ALLOCATE(hice(ngrid)) 482 500 ALLOCATE(qsurf_hist(ngrid,nq)) 483 ALLOCATE(reffrad(ngrid,nlayer mx,naerkind))484 ALLOCATE(nueffrad(ngrid,nlayer mx,naerkind))501 ALLOCATE(reffrad(ngrid,nlayer,naerkind)) 502 ALLOCATE(nueffrad(ngrid,nlayer,naerkind)) 485 503 ALLOCATE(ice_initial(ngrid)) 486 504 ALLOCATE(ice_min(ngrid)) … … 494 512 ALLOCATE(OSR_nu(ngrid,L_NSPECTV)) 495 513 ALLOCATE(sensibFlux(ngrid)) 496 ALLOCATE(zdtlw(ngrid,nlayer mx))497 ALLOCATE(zdtsw(ngrid,nlayer mx))514 ALLOCATE(zdtlw(ngrid,nlayer)) 515 ALLOCATE(zdtsw(ngrid,nlayer)) 498 516 ALLOCATE(tau_col(ngrid)) 499 517 ALLOCATE(pctsrf_sic(ngrid)) … … 541 559 ! read startfi (initial parameters) 542 560 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 543 call phyetat0(ngrid, "startfi.nc",0,0,nsoilmx,nq, &561 call phyetat0(ngrid,nlayer,"startfi.nc",0,0,nsoilmx,nq, & 544 562 day_ini,time_phys,tsurf,tsoil,emis,q2,qsurf, & 545 563 cloudfrac,totcloudfrac,hice, & … … 587 605 588 606 589 CALL init_masquv( zmasq)607 CALL init_masquv(ngrid,zmasq) 590 608 591 609 … … 628 646 ice_update=.false. 629 647 if(sourceevol)then 648 !$OMP MASTER 630 649 open(128,file='num_run',form='formatted', & 631 650 status="old",iostat=ierr) … … 637 656 read(128,*) num_run 638 657 close(128) 658 !$OMP END MASTER 659 !$OMP BARRIER 639 660 640 661 if(num_run.ne.0.and.mod(num_run,2).eq.0)then … … 721 742 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 722 743 723 pdu(1:ngrid,1:nlayer mx) = 0.0724 pdv(1:ngrid,1:nlayer mx) = 0.0744 pdu(1:ngrid,1:nlayer) = 0.0 745 pdv(1:ngrid,1:nlayer) = 0.0 725 746 if ( .not.nearco2cond ) then 726 pdt(1:ngrid,1:nlayer mx) = 0.0747 pdt(1:ngrid,1:nlayer) = 0.0 727 748 endif 728 pdq(1:ngrid,1:nlayer mx,1:nq) = 0.0749 pdq(1:ngrid,1:nlayer,1:nq) = 0.0 729 750 pdpsrf(1:ngrid) = 0.0 730 751 zflubid(1:ngrid) = 0.0 … … 773 794 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 774 795 775 zzlay(1:ngrid,1:nlayer mx)=pphi(1:ngrid,1:nlayermx)776 do l=1,nlayer mx796 zzlay(1:ngrid,1:nlayer)=pphi(1:ngrid,1:nlayer) 797 do l=1,nlayer 777 798 zzlay(1:ngrid,l)= zzlay(1:ngrid,l)/glat(1:ngrid) 778 799 enddo … … 802 823 enddo 803 824 825 ! Compute vertical velocity (m/s) from vertical mass flux 826 ! w = F / (rho*area) and rho = r*T/P 827 ! but first linearly interpolate mass flux to mid-layers 828 do l=1,nlayer-1 829 pw(1:ngrid,l)=0.5*(flxw(1:ngrid,l)+flxw(1:ngrid,l+1)) 830 enddo 831 pw(1:ngrid,nlayer)=0.5*flxw(1:ngrid,nlayer) ! since flxw(nlayer+1)=0 832 do l=1,nlayer 833 pw(1:ngrid,l)=(pw(1:ngrid,l)*pplay(1:ngrid,l)) / & 834 (r*pt(1:ngrid,l)*area(1:ngrid)) 835 enddo 804 836 805 837 !----------------------------------------------------------------------- … … 892 924 endif 893 925 if(water) then 894 muvar(1:ngrid,1:nlayer mx)=mugaz/(1.e0+(1.e0/epsi-1.e0)*pq(1:ngrid,1:nlayermx,igcm_h2o_vap))895 muvar(1:ngrid,nlayer mx+1)=mugaz/(1.e0+(1.e0/epsi-1.e0)*pq(1:ngrid,nlayermx,igcm_h2o_vap))926 muvar(1:ngrid,1:nlayer)=mugaz/(1.e0+(1.e0/epsi-1.e0)*pq(1:ngrid,1:nlayer,igcm_h2o_vap)) 927 muvar(1:ngrid,nlayer+1)=mugaz/(1.e0+(1.e0/epsi-1.e0)*pq(1:ngrid,nlayer,igcm_h2o_vap)) 896 928 ! take into account water effect on mean molecular weight 897 929 else 898 muvar(1:ngrid,1:nlayer mx+1)=mugaz930 muvar(1:ngrid,1:nlayer+1)=mugaz 899 931 endif 900 932 … … 950 982 tau_col(ig) = ntf*tau_col1(ig) + tf*tau_col(ig) 951 983 952 zdtlw(ig,1:nlayer mx) = ntf*zdtlw1(ig,1:nlayermx) + tf*zdtlw(ig,1:nlayermx)953 zdtsw(ig,1:nlayer mx) = ntf*zdtsw1(ig,1:nlayermx) + tf*zdtsw(ig,1:nlayermx)984 zdtlw(ig,1:nlayer) = ntf*zdtlw1(ig,1:nlayer) + tf*zdtlw(ig,1:nlayer) 985 zdtsw(ig,1:nlayer) = ntf*zdtsw1(ig,1:nlayer) + tf*zdtsw(ig,1:nlayer) 954 986 955 987 OSR_nu(ig,1:L_NSPECTV) = ntf*OSR_nu1(ig,1:L_NSPECTV) + tf*OSR_nu(ig,1:L_NSPECTV) … … 976 1008 977 1009 ! Net atmospheric radiative heating rate (K.s-1) 978 dtrad(1:ngrid,1:nlayer mx)=zdtsw(1:ngrid,1:nlayermx)+zdtlw(1:ngrid,1:nlayermx)1010 dtrad(1:ngrid,1:nlayer)=zdtsw(1:ngrid,1:nlayer)+zdtlw(1:ngrid,1:nlayer) 979 1011 980 1012 elseif(newtonian)then … … 982 1014 ! b) Call Newtonian cooling scheme 983 1015 ! ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 984 call newtrelax(ngrid, mu0,sinlat,zpopsk,pt,pplay,pplev,dtrad,firstcall)1016 call newtrelax(ngrid,nlayer,mu0,sinlat,zpopsk,pt,pplay,pplev,dtrad,firstcall) 985 1017 986 1018 zdtsurf(1:ngrid) = +(pt(1:ngrid,1)-tsurf(1:ngrid))/ptimestep … … 1001 1033 1002 1034 1003 dtrad(1:ngrid,1:nlayer mx)=0.01035 dtrad(1:ngrid,1:nlayer)=0.0 1004 1036 ! hence no atmospheric radiative heating 1005 1037 … … 1015 1047 zplanck(1:ngrid)=emis(1:ngrid)*sigma*zplanck(1:ngrid)*zplanck(1:ngrid) 1016 1048 fluxrad(1:ngrid)=fluxrad_sky(1:ngrid)-zplanck(1:ngrid) 1017 pdt(1:ngrid,1:nlayer mx)=pdt(1:ngrid,1:nlayermx)+dtrad(1:ngrid,1:nlayermx)1049 pdt(1:ngrid,1:nlayer)=pdt(1:ngrid,1:nlayer)+dtrad(1:ngrid,1:nlayer) 1018 1050 1019 1051 !------------------------- … … 1048 1080 zflubid(1:ngrid)=fluxrad(1:ngrid)+fluxgrd(1:ngrid) 1049 1081 1050 zdum1(1:ngrid,1:nlayer mx)=0.01051 zdum2(1:ngrid,1:nlayer mx)=0.01082 zdum1(1:ngrid,1:nlayer)=0.0 1083 zdum2(1:ngrid,1:nlayer)=0.0 1052 1084 1053 1085 … … 1066 1098 else 1067 1099 1068 zdh(1:ngrid,1:nlayer mx)=pdt(1:ngrid,1:nlayermx)/zpopsk(1:ngrid,1:nlayermx)1100 zdh(1:ngrid,1:nlayer)=pdt(1:ngrid,1:nlayer)/zpopsk(1:ngrid,1:nlayer) 1069 1101 1070 1102 call vdifc(ngrid,nlayer,nq,rnat,zpopsk, & … … 1077 1109 taux,tauy,lastcall) 1078 1110 1079 zdtdif(1:ngrid,1:nlayer mx)=zdhdif(1:ngrid,1:nlayermx)*zpopsk(1:ngrid,1:nlayermx) ! for diagnostic only1080 zdqevap(1:ngrid,1:nlayer mx)=0.1111 zdtdif(1:ngrid,1:nlayer)=zdhdif(1:ngrid,1:nlayer)*zpopsk(1:ngrid,1:nlayer) ! for diagnostic only 1112 zdqevap(1:ngrid,1:nlayer)=0. 1081 1113 1082 1114 end if !turbdiff 1083 1115 1084 pdv(1:ngrid,1:nlayer mx)=pdv(1:ngrid,1:nlayermx)+zdvdif(1:ngrid,1:nlayermx)1085 pdu(1:ngrid,1:nlayer mx)=pdu(1:ngrid,1:nlayermx)+zdudif(1:ngrid,1:nlayermx)1086 pdt(1:ngrid,1:nlayer mx)=pdt(1:ngrid,1:nlayermx)+zdtdif(1:ngrid,1:nlayermx)1116 pdv(1:ngrid,1:nlayer)=pdv(1:ngrid,1:nlayer)+zdvdif(1:ngrid,1:nlayer) 1117 pdu(1:ngrid,1:nlayer)=pdu(1:ngrid,1:nlayer)+zdudif(1:ngrid,1:nlayer) 1118 pdt(1:ngrid,1:nlayer)=pdt(1:ngrid,1:nlayer)+zdtdif(1:ngrid,1:nlayer) 1087 1119 zdtsurf(1:ngrid)=zdtsurf(1:ngrid)+zdtsdif(1:ngrid) 1088 1120 … … 1094 1126 1095 1127 if (tracer) then 1096 pdq(1:ngrid,1:nlayer mx,1:nq)=pdq(1:ngrid,1:nlayermx,1:nq)+ zdqdif(1:ngrid,1:nlayermx,1:nq)1128 pdq(1:ngrid,1:nlayer,1:nq)=pdq(1:ngrid,1:nlayer,1:nq)+ zdqdif(1:ngrid,1:nlayer,1:nq) 1097 1129 dqsurf(1:ngrid,1:nq)=dqsurf(1:ngrid,1:nq) + zdqsdif(1:ngrid,1:nq) 1098 1130 end if ! of if (tracer) … … 1171 1203 if(calladj) then 1172 1204 1173 zdh(1:ngrid,1:nlayer mx) = pdt(1:ngrid,1:nlayermx)/zpopsk(1:ngrid,1:nlayermx)1174 zduadj(1:ngrid,1:nlayer mx)=0.01175 zdvadj(1:ngrid,1:nlayer mx)=0.01176 zdhadj(1:ngrid,1:nlayer mx)=0.01205 zdh(1:ngrid,1:nlayer) = pdt(1:ngrid,1:nlayer)/zpopsk(1:ngrid,1:nlayer) 1206 zduadj(1:ngrid,1:nlayer)=0.0 1207 zdvadj(1:ngrid,1:nlayer)=0.0 1208 zdhadj(1:ngrid,1:nlayer)=0.0 1177 1209 1178 1210 … … 1184 1216 zdqadj) 1185 1217 1186 pdu(1:ngrid,1:nlayer mx) = pdu(1:ngrid,1:nlayermx) + zduadj(1:ngrid,1:nlayermx)1187 pdv(1:ngrid,1:nlayer mx) = pdv(1:ngrid,1:nlayermx) + zdvadj(1:ngrid,1:nlayermx)1188 pdt(1:ngrid,1:nlayer mx) = pdt(1:ngrid,1:nlayermx) + zdhadj(1:ngrid,1:nlayermx)*zpopsk(1:ngrid,1:nlayermx)1189 zdtadj(1:ngrid,1:nlayer mx) = zdhadj(1:ngrid,1:nlayermx)*zpopsk(1:ngrid,1:nlayermx) ! for diagnostic only1218 pdu(1:ngrid,1:nlayer) = pdu(1:ngrid,1:nlayer) + zduadj(1:ngrid,1:nlayer) 1219 pdv(1:ngrid,1:nlayer) = pdv(1:ngrid,1:nlayer) + zdvadj(1:ngrid,1:nlayer) 1220 pdt(1:ngrid,1:nlayer) = pdt(1:ngrid,1:nlayer) + zdhadj(1:ngrid,1:nlayer)*zpopsk(1:ngrid,1:nlayer) 1221 zdtadj(1:ngrid,1:nlayer) = zdhadj(1:ngrid,1:nlayer)*zpopsk(1:ngrid,1:nlayer) ! for diagnostic only 1190 1222 1191 1223 if(tracer) then 1192 pdq(1:ngrid,1:nlayer mx,1:nq) = pdq(1:ngrid,1:nlayermx,1:nq) + zdqadj(1:ngrid,1:nlayermx,1:nq)1224 pdq(1:ngrid,1:nlayer,1:nq) = pdq(1:ngrid,1:nlayer,1:nq) + zdqadj(1:ngrid,1:nlayer,1:nq) 1193 1225 end if 1194 1226 … … 1240 1272 zdqc) 1241 1273 1242 pdt(1:ngrid,1:nlayer mx)=pdt(1:ngrid,1:nlayermx)+zdtc(1:ngrid,1:nlayermx)1243 pdv(1:ngrid,1:nlayer mx)=pdv(1:ngrid,1:nlayermx)+zdvc(1:ngrid,1:nlayermx)1244 pdu(1:ngrid,1:nlayer mx)=pdu(1:ngrid,1:nlayermx)+zduc(1:ngrid,1:nlayermx)1274 pdt(1:ngrid,1:nlayer)=pdt(1:ngrid,1:nlayer)+zdtc(1:ngrid,1:nlayer) 1275 pdv(1:ngrid,1:nlayer)=pdv(1:ngrid,1:nlayer)+zdvc(1:ngrid,1:nlayer) 1276 pdu(1:ngrid,1:nlayer)=pdu(1:ngrid,1:nlayer)+zduc(1:ngrid,1:nlayer) 1245 1277 zdtsurf(1:ngrid) = zdtsurf(1:ngrid) + zdtsurfc(1:ngrid) 1246 1278 1247 pdq(1:ngrid,1:nlayer mx,1:nq)=pdq(1:ngrid,1:nlayermx,1:nq)+ zdqc(1:ngrid,1:nlayermx,1:nq)1279 pdq(1:ngrid,1:nlayer,1:nq)=pdq(1:ngrid,1:nlayer,1:nq)+ zdqc(1:ngrid,1:nlayer,1:nq) 1248 1280 ! Note: we do not add surface co2ice tendency 1249 1281 ! because qsurf(:,igcm_co2_ice) is updated in condens_co2cloud … … 1283 1315 ! ---------------- 1284 1316 1285 dqmoist(1:ngrid,1:nlayer mx,1:nq)=0.1286 dtmoist(1:ngrid,1:nlayer mx)=0.1287 1288 call moistadj(ngrid,n q,pt,pq,pdq,pplev,pplay,dtmoist,dqmoist,ptimestep,rneb_man)1289 1290 pdq(1:ngrid,1:nlayer mx,igcm_h2o_vap) = pdq(1:ngrid,1:nlayermx,igcm_h2o_vap) &1291 +dqmoist(1:ngrid,1:nlayer mx,igcm_h2o_vap)1292 pdq(1:ngrid,1:nlayer mx,igcm_h2o_ice) =pdq(1:ngrid,1:nlayermx,igcm_h2o_ice) &1293 +dqmoist(1:ngrid,1:nlayer mx,igcm_h2o_ice)1294 pdt(1:ngrid,1:nlayer mx) = pdt(1:ngrid,1:nlayermx)+dtmoist(1:ngrid,1:nlayermx)1317 dqmoist(1:ngrid,1:nlayer,1:nq)=0. 1318 dtmoist(1:ngrid,1:nlayer)=0. 1319 1320 call moistadj(ngrid,nlayer,nq,pt,pq,pdq,pplev,pplay,dtmoist,dqmoist,ptimestep,rneb_man) 1321 1322 pdq(1:ngrid,1:nlayer,igcm_h2o_vap) = pdq(1:ngrid,1:nlayer,igcm_h2o_vap) & 1323 +dqmoist(1:ngrid,1:nlayer,igcm_h2o_vap) 1324 pdq(1:ngrid,1:nlayer,igcm_h2o_ice) =pdq(1:ngrid,1:nlayer,igcm_h2o_ice) & 1325 +dqmoist(1:ngrid,1:nlayer,igcm_h2o_ice) 1326 pdt(1:ngrid,1:nlayer) = pdt(1:ngrid,1:nlayer)+dtmoist(1:ngrid,1:nlayer) 1295 1327 1296 1328 !------------------------- … … 1321 1353 ! Large scale condensation/evaporation 1322 1354 ! -------------------------------- 1323 call largescale(ngrid,n q,ptimestep,pplev,pplay,pt,pq,pdt,pdq,dtlscale,dqvaplscale,dqcldlscale,rneb_lsc)1324 1325 pdt(1:ngrid,1:nlayer mx) = pdt(1:ngrid,1:nlayermx)+dtlscale(1:ngrid,1:nlayermx)1326 pdq(1:ngrid,1:nlayer mx,igcm_h2o_vap) = pdq(1:ngrid,1:nlayermx,igcm_h2o_vap)+dqvaplscale(1:ngrid,1:nlayermx)1327 pdq(1:ngrid,1:nlayer mx,igcm_h2o_ice) = pdq(1:ngrid,1:nlayermx,igcm_h2o_ice)+dqcldlscale(1:ngrid,1:nlayermx)1355 call largescale(ngrid,nlayer,nq,ptimestep,pplev,pplay,pt,pq,pdt,pdq,dtlscale,dqvaplscale,dqcldlscale,rneb_lsc) 1356 1357 pdt(1:ngrid,1:nlayer) = pdt(1:ngrid,1:nlayer)+dtlscale(1:ngrid,1:nlayer) 1358 pdq(1:ngrid,1:nlayer,igcm_h2o_vap) = pdq(1:ngrid,1:nlayer,igcm_h2o_vap)+dqvaplscale(1:ngrid,1:nlayer) 1359 pdq(1:ngrid,1:nlayer,igcm_h2o_ice) = pdq(1:ngrid,1:nlayer,igcm_h2o_ice)+dqcldlscale(1:ngrid,1:nlayer) 1328 1360 1329 1361 !------------------------- … … 1363 1395 if(waterrain)then 1364 1396 1365 zdqrain(1:ngrid,1:nlayer mx,1:nq) = 0.01397 zdqrain(1:ngrid,1:nlayer,1:nq) = 0.0 1366 1398 zdqsrain(1:ngrid) = 0.0 1367 1399 zdqssnow(1:ngrid) = 0.0 1368 1400 1369 call rain(ngrid,n q,ptimestep,pplev,pplay,pt,pdt,pq,pdq, &1401 call rain(ngrid,nlayer,nq,ptimestep,pplev,pplay,pt,pdt,pq,pdq, & 1370 1402 zdtrain,zdqrain,zdqsrain,zdqssnow,cloudfrac) 1371 1403 1372 pdq(1:ngrid,1:nlayer mx,igcm_h2o_vap) = pdq(1:ngrid,1:nlayermx,igcm_h2o_vap) &1373 +zdqrain(1:ngrid,1:nlayer mx,igcm_h2o_vap)1374 pdq(1:ngrid,1:nlayer mx,igcm_h2o_ice) = pdq(1:ngrid,1:nlayermx,igcm_h2o_ice) &1375 +zdqrain(1:ngrid,1:nlayer mx,igcm_h2o_ice)1376 pdt(1:ngrid,1:nlayer mx) = pdt(1:ngrid,1:nlayermx)+zdtrain(1:ngrid,1:nlayermx)1404 pdq(1:ngrid,1:nlayer,igcm_h2o_vap) = pdq(1:ngrid,1:nlayer,igcm_h2o_vap) & 1405 +zdqrain(1:ngrid,1:nlayer,igcm_h2o_vap) 1406 pdq(1:ngrid,1:nlayer,igcm_h2o_ice) = pdq(1:ngrid,1:nlayer,igcm_h2o_ice) & 1407 +zdqrain(1:ngrid,1:nlayer,igcm_h2o_ice) 1408 pdt(1:ngrid,1:nlayer) = pdt(1:ngrid,1:nlayer)+zdtrain(1:ngrid,1:nlayer) 1377 1409 dqsurf(1:ngrid,igcm_h2o_vap) = dqsurf(1:ngrid,igcm_h2o_vap)+zdqsrain(1:ngrid) ! a bug was here 1378 1410 dqsurf(1:ngrid,igcm_h2o_ice) = dqsurf(1:ngrid,igcm_h2o_ice)+zdqssnow(1:ngrid) … … 1419 1451 ! ------------- 1420 1452 if (sedimentation) then 1421 zdqsed(1:ngrid,1:nlayer mx,1:nq) = 0.01453 zdqsed(1:ngrid,1:nlayer,1:nq) = 0.0 1422 1454 zdqssed(1:ngrid,1:nq) = 0.0 1423 1455 … … 1455 1487 ! and as in rain.F90, whether it falls as rain or snow depends 1456 1488 ! only on the surface temperature 1457 pdq(1:ngrid,1:nlayer mx,1:nq) = pdq(1:ngrid,1:nlayermx,1:nq) + zdqsed(1:ngrid,1:nlayermx,1:nq)1489 pdq(1:ngrid,1:nlayer,1:nq) = pdq(1:ngrid,1:nlayer,1:nq) + zdqsed(1:ngrid,1:nlayer,1:nq) 1458 1490 dqsurf(1:ngrid,1:nq) = dqsurf(1:ngrid,1:nq) + zdqssed(1:ngrid,1:nq) 1459 1491 … … 1483 1515 if(mass_redistrib) then 1484 1516 1485 zdmassmr(1:ngrid,1:nlayer mx) = mass(1:ngrid,1:nlayermx) * &1486 ( zdqevap(1:ngrid,1:nlayer mx) &1487 + zdqrain(1:ngrid,1:nlayer mx,igcm_h2o_vap) &1488 + dqmoist(1:ngrid,1:nlayer mx,igcm_h2o_vap) &1489 + dqvaplscale(1:ngrid,1:nlayer mx) )1517 zdmassmr(1:ngrid,1:nlayer) = mass(1:ngrid,1:nlayer) * & 1518 ( zdqevap(1:ngrid,1:nlayer) & 1519 + zdqrain(1:ngrid,1:nlayer,igcm_h2o_vap) & 1520 + dqmoist(1:ngrid,1:nlayer,igcm_h2o_vap) & 1521 + dqvaplscale(1:ngrid,1:nlayer) ) 1490 1522 1491 1523 do ig = 1, ngrid 1492 zdmassmr_col(ig)=SUM(zdmassmr(ig,1:nlayer mx))1524 zdmassmr_col(ig)=SUM(zdmassmr(ig,1:nlayer)) 1493 1525 enddo 1494 1526 … … 1503 1535 1504 1536 1505 pdq(1:ngrid,1:nlayer mx,1:nq) = pdq(1:ngrid,1:nlayermx,1:nq) + zdqmr(1:ngrid,1:nlayermx,1:nq)1537 pdq(1:ngrid,1:nlayer,1:nq) = pdq(1:ngrid,1:nlayer,1:nq) + zdqmr(1:ngrid,1:nlayer,1:nq) 1506 1538 dqsurf(1:ngrid,1:nq) = dqsurf(1:ngrid,1:nq) + zdqsurfmr(1:ngrid,1:nq) 1507 pdt(1:ngrid,1:nlayer mx) = pdt(1:ngrid,1:nlayermx) + zdtmr(1:ngrid,1:nlayermx)1508 pdu(1:ngrid,1:nlayer mx) = pdu(1:ngrid,1:nlayermx) + zdumr(1:ngrid,1:nlayermx)1509 pdv(1:ngrid,1:nlayer mx) = pdv(1:ngrid,1:nlayermx) + zdvmr(1:ngrid,1:nlayermx)1539 pdt(1:ngrid,1:nlayer) = pdt(1:ngrid,1:nlayer) + zdtmr(1:ngrid,1:nlayer) 1540 pdu(1:ngrid,1:nlayer) = pdu(1:ngrid,1:nlayer) + zdumr(1:ngrid,1:nlayer) 1541 pdv(1:ngrid,1:nlayer) = pdv(1:ngrid,1:nlayer) + zdvmr(1:ngrid,1:nlayer) 1510 1542 pdpsrf(1:ngrid) = pdpsrf(1:ngrid) + zdpsrfmr(1:ngrid) 1511 1543 zdtsurf(1:ngrid) = zdtsurf(1:ngrid) + zdtsurfmr(1:ngrid) 1512 1544 1513 ! print*,'after mass redistrib, q=',pq(211,1:nlayer mx,igcm_h2o_vap)+ptimestep*pdq(211,1:nlayermx,igcm_h2o_vap)1545 ! print*,'after mass redistrib, q=',pq(211,1:nlayer,igcm_h2o_vap)+ptimestep*pdq(211,1:nlayer,igcm_h2o_vap) 1514 1546 endif 1515 1547 … … 1672 1704 1673 1705 ! temperature, zonal and meridional wind 1674 zt(1:ngrid,1:nlayer mx) = pt(1:ngrid,1:nlayermx) + pdt(1:ngrid,1:nlayermx)*ptimestep1675 zu(1:ngrid,1:nlayer mx) = pu(1:ngrid,1:nlayermx) + pdu(1:ngrid,1:nlayermx)*ptimestep1676 zv(1:ngrid,1:nlayer mx) = pv(1:ngrid,1:nlayermx) + pdv(1:ngrid,1:nlayermx)*ptimestep1706 zt(1:ngrid,1:nlayer) = pt(1:ngrid,1:nlayer) + pdt(1:ngrid,1:nlayer)*ptimestep 1707 zu(1:ngrid,1:nlayer) = pu(1:ngrid,1:nlayer) + pdu(1:ngrid,1:nlayer)*ptimestep 1708 zv(1:ngrid,1:nlayer) = pv(1:ngrid,1:nlayer) + pdv(1:ngrid,1:nlayer)*ptimestep 1677 1709 1678 1710 ! diagnostic 1679 zdtdyn(1:ngrid,1:nlayer mx) = pt(1:ngrid,1:nlayermx)-ztprevious(1:ngrid,1:nlayermx)1680 ztprevious(1:ngrid,1:nlayer mx) = zt(1:ngrid,1:nlayermx)1711 zdtdyn(1:ngrid,1:nlayer) = pt(1:ngrid,1:nlayer)-ztprevious(1:ngrid,1:nlayer) 1712 ztprevious(1:ngrid,1:nlayer) = zt(1:ngrid,1:nlayer) 1681 1713 1682 1714 if(firstcall)then 1683 zdtdyn(1:ngrid,1:nlayer mx)=0.01715 zdtdyn(1:ngrid,1:nlayer)=0.0 1684 1716 endif 1685 1717 … … 1690 1722 1691 1723 ! tracers 1692 zq(1:ngrid,1:nlayer mx,1:nq) = pq(1:ngrid,1:nlayermx,1:nq) + pdq(1:ngrid,1:nlayermx,1:nq)*ptimestep1724 zq(1:ngrid,1:nlayer,1:nq) = pq(1:ngrid,1:nlayer,1:nq) + pdq(1:ngrid,1:nlayer,1:nq)*ptimestep 1693 1725 1694 1726 ! surface pressure … … 1798 1830 do iq=1,nq 1799 1831 do ig=1,ngrid 1800 qcol(ig,iq) = SUM( zq(ig,1:nlayer mx,iq) * mass(ig,1:nlayermx))1832 qcol(ig,iq) = SUM( zq(ig,1:nlayer,iq) * mass(ig,1:nlayer)) 1801 1833 enddo 1802 1834 enddo … … 1805 1837 reffcol(1:ngrid,1:naerkind)=0.0 1806 1838 if(co2cond.and.(iaero_co2.ne.0))then 1807 call co2_reffrad(ngrid,n q,zq,reffrad(1,1,iaero_co2))1839 call co2_reffrad(ngrid,nlayer,nq,zq,reffrad(1,1,iaero_co2)) 1808 1840 do ig=1,ngrid 1809 reffcol(ig,iaero_co2) = SUM(zq(ig,1:nlayer mx,igcm_co2_ice)*reffrad(ig,1:nlayermx,iaero_co2)*mass(ig,1:nlayermx))1841 reffcol(ig,iaero_co2) = SUM(zq(ig,1:nlayer,igcm_co2_ice)*reffrad(ig,1:nlayer,iaero_co2)*mass(ig,1:nlayer)) 1810 1842 enddo 1811 1843 endif 1812 1844 if(water.and.(iaero_h2o.ne.0))then 1813 call h2o_reffrad(ngrid, zq(1,1,igcm_h2o_ice),zt, &1845 call h2o_reffrad(ngrid,nlayer,zq(1,1,igcm_h2o_ice),zt, & 1814 1846 reffrad(1,1,iaero_h2o),nueffrad(1,1,iaero_h2o)) 1815 1847 do ig=1,ngrid 1816 reffcol(ig,iaero_h2o) = SUM(zq(ig,1:nlayer mx,igcm_h2o_ice)*reffrad(ig,1:nlayermx,iaero_h2o)*mass(ig,1:nlayermx))1848 reffcol(ig,iaero_h2o) = SUM(zq(ig,1:nlayer,igcm_h2o_ice)*reffrad(ig,1:nlayer,iaero_h2o)*mass(ig,1:nlayer)) 1817 1849 enddo 1818 1850 endif … … 1989 2021 end do 1990 2022 if (water) then 1991 vmr=zq(1:ngrid,1:nlayer mx,igcm_h2o_vap)*mugaz/mmol(igcm_h2o_vap)2023 vmr=zq(1:ngrid,1:nlayer,igcm_h2o_vap)*mugaz/mmol(igcm_h2o_vap) 1992 2024 call wstats(ngrid,"vmr_h2ovapor", & 1993 2025 "H2O vapour volume mixing ratio","mol/mol", & … … 2183 2215 2184 2216 ! deallocate gas variables 2217 !$OMP BARRIER 2218 !$OMP MASTER 2185 2219 IF ( ALLOCATED( gnom ) ) DEALLOCATE( gnom ) 2186 2220 IF ( ALLOCATED( gfrac ) ) DEALLOCATE( gfrac ) ! both allocated in su_gases.F90 2221 !$OMP END MASTER 2222 !$OMP BARRIER 2187 2223 2188 2224 ! deallocate saved arrays … … 2203 2239 IF ( ALLOCATED(q2)) DEALLOCATE(q2) 2204 2240 IF ( ALLOCATED(ztprevious)) DEALLOCATE(ztprevious) 2241 IF ( ALLOCATED(hice)) DEALLOCATE(hice) 2205 2242 IF ( ALLOCATED(cloudfrac)) DEALLOCATE(cloudfrac) 2206 2243 IF ( ALLOCATED(totcloudfrac)) DEALLOCATE(totcloudfrac) … … 2223 2260 IF ( ALLOCATED(zdtsw)) DEALLOCATE(zdtsw) 2224 2261 IF ( ALLOCATED(tau_col)) DEALLOCATE(tau_col) 2262 IF ( ALLOCATED(pctsrf_sic)) DEALLOCATE(pctsrf_sic) 2263 IF ( ALLOCATED(tslab)) DEALLOCATE(tslab) 2264 IF ( ALLOCATED(tsea_ice)) DEALLOCATE(tsea_ice) 2265 IF ( ALLOCATED(sea_ice)) DEALLOCATE(sea_ice) 2266 IF ( ALLOCATED(zmasq)) DEALLOCATE(zmasq) 2267 IF ( ALLOCATED(knindex)) DEALLOCATE(knindex) 2225 2268 2226 2269 !! this is defined in comsaison_h -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/radcommon_h.F90
r222 r227 61 61 ! gIR : mean assymetry factor 62 62 63 REAL*8 BWNI(L_NSPECTI+1), WNOI(L_NSPECTI), DWNI(L_NSPECTI), WAVEI(L_NSPECTI) 64 REAL*8 BWNV(L_NSPECTV+1), WNOV(L_NSPECTV), DWNV(L_NSPECTV), WAVEV(L_NSPECTV) 63 REAL*8 BWNI(L_NSPECTI+1), WNOI(L_NSPECTI), DWNI(L_NSPECTI), WAVEI(L_NSPECTI) !BWNI read by master in setspi 64 REAL*8 BWNV(L_NSPECTV+1), WNOV(L_NSPECTV), DWNV(L_NSPECTV), WAVEV(L_NSPECTV) !BWNV read by master in setspv 65 65 REAL*8 STELLARF(L_NSPECTV), TAURAY(L_NSPECTV), TAURAYVAR(L_NSPECTV) 66 !$OMP THREADPRIVATE(WNOI,DWNI,WAVEI,& 67 !$OMP WNOV,DWNV,WAVEV,& 68 !$OMP STELLARF,TAURAY,TAURAYVAR) 66 69 67 70 REAL*8 blami(L_NSPECTI+1) 68 71 REAL*8 blamv(L_NSPECTV+1) ! these are needed by suaer.F90 72 !$OMP THREADPRIVATE(blami,blamv) 69 73 70 74 !! AS: introduced to avoid doing same computations again for continuum 71 75 INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: indi 72 76 INTEGER, DIMENSION(:,:,:), ALLOCATABLE :: indv 77 !$OMP THREADPRIVATE(indi,indv) 73 78 74 79 !!! ALLOCATABLE STUFF SO THAT DIMENSIONS ARE READ in *.dat FILES -- AS 12/2011 … … 79 84 real*8 pgasmin, pgasmax 80 85 real*8 tgasmin, tgasmax 86 !$OMP THREADPRIVATE(gasi,gasv,& !wrefvar,pgasref,tgasref,pfgasref read by master in sugas_corrk 87 !$OMP FZEROI,FZEROV) !pgasmin,pgasmax,tgasmin,tgasmax read by master in sugas_corrk 81 88 82 89 real QVISsQREF(L_NSPECTV,naerkind,nsizemax) … … 86 93 real omegair(L_NSPECTI,naerkind,nsizemax) 87 94 real gir(L_NSPECTI,naerkind,nsizemax) 95 !$OMP THREADPRIVATE(QVISsQREF,omegavis,gvis,QIRsQREF,omegair,gir) 88 96 89 97 … … 102 110 103 111 DOUBLE PRECISION :: radiustab(naerkind,2,nsizemax) 112 !$OMP THREADPRIVATE(lamrefir,lamrefvis,radiustab) !nsize read by suaer_corrk 104 113 105 114 ! Extinction coefficient at reference wavelengths; … … 121 130 122 131 real*8 gweight(L_NGAUSS) 132 !$OMP THREADPRIVATE(QREFvis,QREFir,omegaREFvis,omegaREFir,& ! gweight read by master in sugas_corrk 133 !$OMP tstellar,planckir,PTOP,TAUREF) 123 134 124 135 ! If the gas optical depth (top to the surface) is less than … … 137 148 real*8 glat_ig 138 149 save glat_ig 150 !$OMP THREADPRIVATE(Cmk,glat_ig) 139 151 140 152 ! extinction of incoming sunlight (Saturn's rings, eclipses, etc...) … … 143 155 !Latitude-dependent gravity 144 156 REAL, DIMENSION(:), ALLOCATABLE :: glat 157 !$OMP THREADPRIVATE(glat,eclipse) 145 158 146 159 end module radcommon_h -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/radii_mod.F90
r222 r227 12 12 real, save :: Nmix_h2o 13 13 real, save :: Nmix_h2o_ice 14 !$OMP THREADPRIVATE(rad_h2o,rad_h2o_ice,Nmix_h2o,Nmix_h2o_ice) 14 15 real, parameter :: coef_chaud=0.13 15 16 real, parameter :: coef_froid=0.09 … … 20 21 21 22 !================================================================== 22 subroutine su_aer_radii(ngrid, reffrad,nueffrad)23 subroutine su_aer_radii(ngrid,nlayer,reffrad,nueffrad) 23 24 !================================================================== 24 25 ! Purpose … … 32 33 !================================================================== 33 34 ! to use 'getin' 34 use ioipsl_getincom 35 ! use ioipsl_getincom 36 use ioipsl_getincom_p 35 37 use radinc_h, only: naerkind 36 38 use aerosol_mod … … 39 41 40 42 include "callkeys.h" 41 include "dimensions.h" 42 include "dimphys.h" 43 44 integer,intent(in) :: ngrid 45 46 real, intent(out) :: reffrad(ngrid,nlayermx,naerkind) !aerosols radii (K) 47 real, intent(out) :: nueffrad(ngrid,nlayermx,naerkind) !variance 43 ! include "dimensions.h" 44 ! include "dimphys.h" 45 46 integer,intent(in) :: ngrid 47 integer,intent(in) :: nlayer 48 49 real, intent(out) :: reffrad(ngrid,nlayer,naerkind) !aerosols radii (K) 50 real, intent(out) :: nueffrad(ngrid,nlayer,naerkind) !variance 48 51 49 52 logical, save :: firstcall=.true. 53 !$OMP THREADPRIVATE(firstcall) 50 54 integer :: iaer 51 55 … … 58 62 59 63 if(iaer.eq.iaero_co2)then ! CO2 ice 60 reffrad(1:ngrid,1:nlayer mx,iaer) = 1.e-461 nueffrad(1:ngrid,1:nlayer mx,iaer) = 0.164 reffrad(1:ngrid,1:nlayer,iaer) = 1.e-4 65 nueffrad(1:ngrid,1:nlayer,iaer) = 0.1 62 66 endif 63 67 64 68 if(iaer.eq.iaero_h2o)then ! H2O ice 65 reffrad(1:ngrid,1:nlayer mx,iaer) = 1.e-566 nueffrad(1:ngrid,1:nlayer mx,iaer) = 0.169 reffrad(1:ngrid,1:nlayer,iaer) = 1.e-5 70 nueffrad(1:ngrid,1:nlayer,iaer) = 0.1 67 71 endif 68 72 69 73 if(iaer.eq.iaero_dust)then ! dust 70 reffrad(1:ngrid,1:nlayer mx,iaer) = 1.e-571 nueffrad(1:ngrid,1:nlayer mx,iaer) = 0.174 reffrad(1:ngrid,1:nlayer,iaer) = 1.e-5 75 nueffrad(1:ngrid,1:nlayer,iaer) = 0.1 72 76 endif 73 77 74 78 if(iaer.eq.iaero_h2so4)then ! H2O ice 75 reffrad(1:ngrid,1:nlayer mx,iaer) = 1.e-676 nueffrad(1:ngrid,1:nlayer mx,iaer) = 0.179 reffrad(1:ngrid,1:nlayer,iaer) = 1.e-6 80 nueffrad(1:ngrid,1:nlayer,iaer) = 0.1 77 81 endif 78 82 79 83 if(iaer.eq.iaero_back2lay)then ! Two-layer aerosols 80 reffrad(1:ngrid,1:nlayer mx,iaer) = 2.e-681 nueffrad(1:ngrid,1:nlayer mx,iaer) = 0.184 reffrad(1:ngrid,1:nlayer,iaer) = 2.e-6 85 nueffrad(1:ngrid,1:nlayer,iaer) = 0.1 82 86 endif 83 87 … … 98 102 write(*,*)"radius of H2O water particles:" 99 103 rad_h2o=13. ! default value 100 call getin ("rad_h2o",rad_h2o)104 call getin_p("rad_h2o",rad_h2o) 101 105 write(*,*)" rad_h2o = ",rad_h2o 102 106 103 107 write(*,*)"radius of H2O ice particles:" 104 108 rad_h2o_ice=35. ! default value 105 call getin ("rad_h2o_ice",rad_h2o_ice)109 call getin_p("rad_h2o_ice",rad_h2o_ice) 106 110 write(*,*)" rad_h2o_ice = ",rad_h2o_ice 107 111 … … 110 114 write(*,*)"Number mixing ratio of H2O water particles:" 111 115 Nmix_h2o=1.e6 ! default value 112 call getin ("Nmix_h2o",Nmix_h2o)116 call getin_p("Nmix_h2o",Nmix_h2o) 113 117 write(*,*)" Nmix_h2o = ",Nmix_h2o 114 118 115 119 write(*,*)"Number mixing ratio of H2O ice particles:" 116 120 Nmix_h2o_ice=Nmix_h2o ! default value 117 call getin ("Nmix_h2o_ice",Nmix_h2o_ice)121 call getin_p("Nmix_h2o_ice",Nmix_h2o_ice) 118 122 write(*,*)" Nmix_h2o_ice = ",Nmix_h2o_ice 119 123 endif … … 126 130 127 131 !================================================================== 128 subroutine h2o_reffrad(ngrid, pq,pt,reffrad,nueffrad)132 subroutine h2o_reffrad(ngrid,nlayer,pq,pt,reffrad,nueffrad) 129 133 !================================================================== 130 134 ! Purpose … … 141 145 142 146 include "callkeys.h" 143 include "dimensions.h"144 include "dimphys.h"147 ! include "dimensions.h" 148 ! include "dimphys.h" 145 149 include "comcstfi.h" 146 150 147 151 integer,intent(in) :: ngrid 148 149 real, intent(in) :: pq(ngrid,nlayermx) !water ice mixing ratios (kg/kg) 150 real, intent(in) :: pt(ngrid,nlayermx) !temperature (K) 151 real, intent(out) :: reffrad(ngrid,nlayermx) !aerosol radii 152 real, intent(out) :: nueffrad(ngrid,nlayermx) ! dispersion 152 integer,intent(in) :: nlayer 153 154 real, intent(in) :: pq(ngrid,nlayer) !water ice mixing ratios (kg/kg) 155 real, intent(in) :: pt(ngrid,nlayer) !temperature (K) 156 real, intent(out) :: reffrad(ngrid,nlayer) !aerosol radii 157 real, intent(out) :: nueffrad(ngrid,nlayer) ! dispersion 153 158 154 159 integer :: ig,l … … 158 163 159 164 if (radfixed) then 160 do l=1,nlayer mx165 do l=1,nlayer 161 166 do ig=1,ngrid 162 167 zfice = 1.0 - (pt(ig,l)-T_h2O_ice_clouds) / (T_h2O_ice_liq-T_h2O_ice_clouds) … … 167 172 enddo 168 173 else 169 do l=1,nlayer mx174 do l=1,nlayer 170 175 do ig=1,ngrid 171 176 zfice = 1.0 - (pt(ig,l)-T_h2O_ice_clouds) / (T_h2O_ice_liq-T_h2O_ice_clouds) … … 186 191 187 192 !================================================================== 188 subroutine h2o_cloudrad(ngrid, pql,reffliq,reffice)193 subroutine h2o_cloudrad(ngrid,nlayer,pql,reffliq,reffice) 189 194 !================================================================== 190 195 ! Purpose … … 201 206 202 207 include "callkeys.h" 203 include "dimensions.h"204 include "dimphys.h"208 ! include "dimensions.h" 209 ! include "dimphys.h" 205 210 include "comcstfi.h" 206 211 207 212 integer,intent(in) :: ngrid 208 209 real, intent(in) :: pql(ngrid,nlayermx) !condensed water mixing ratios (kg/kg) 210 real, intent(out) :: reffliq(ngrid,nlayermx),reffice(ngrid,nlayermx) !liquid and ice water particle radii (m) 213 integer,intent(in) :: nlayer 214 215 real, intent(in) :: pql(ngrid,nlayer) !condensed water mixing ratios (kg/kg) 216 real, intent(out) :: reffliq(ngrid,nlayer),reffice(ngrid,nlayer) !liquid and ice water particle radii (m) 211 217 212 218 real,external :: CBRT … … 214 220 215 221 if (radfixed) then 216 reffliq(1:ngrid,1:nlayer mx)= rad_h2o217 reffice(1:ngrid,1:nlayer mx)= rad_h2o_ice222 reffliq(1:ngrid,1:nlayer)= rad_h2o 223 reffice(1:ngrid,1:nlayer)= rad_h2o_ice 218 224 else 219 do k=1,nlayer mx225 do k=1,nlayer 220 226 do i=1,ngrid 221 227 reffliq(i,k) = CBRT(3*pql(i,k)/(4*Nmix_h2o*pi*rhowater)) … … 234 240 235 241 !================================================================== 236 subroutine co2_reffrad(ngrid,n q,pq,reffrad)242 subroutine co2_reffrad(ngrid,nlayer,nq,pq,reffrad) 237 243 !================================================================== 238 244 ! Purpose … … 249 255 250 256 include "callkeys.h" 251 include "dimensions.h"252 include "dimphys.h"257 ! include "dimensions.h" 258 ! include "dimphys.h" 253 259 include "comcstfi.h" 254 260 255 integer,intent(in) :: ngrid,n q256 257 real, intent(in) :: pq(ngrid,nlayer mx,nq) !tracer mixing ratios (kg/kg)258 real, intent(out) :: reffrad(ngrid,nlayer mx) !co2 ice particles radii (K)261 integer,intent(in) :: ngrid,nlayer,nq 262 263 real, intent(in) :: pq(ngrid,nlayer,nq) !tracer mixing ratios (kg/kg) 264 real, intent(out) :: reffrad(ngrid,nlayer) !co2 ice particles radii (m) 259 265 260 266 integer :: ig,l … … 265 271 266 272 if (radfixed) then 267 reffrad(1:ngrid,1:nlayer mx) = 5.e-5 ! CO2 ice273 reffrad(1:ngrid,1:nlayer) = 5.e-5 ! CO2 ice 268 274 else 269 do l=1,nlayer mx275 do l=1,nlayer 270 276 do ig=1,ngrid 271 277 zrad = CBRT( 3*pq(ig,l,igcm_co2_ice)/(4*Nmix_co2*pi*rho_co2) ) … … 281 287 282 288 !================================================================== 283 subroutine dust_reffrad(ngrid, reffrad)289 subroutine dust_reffrad(ngrid,nlayer,reffrad) 284 290 !================================================================== 285 291 ! Purpose … … 294 300 Implicit none 295 301 296 include "callkeys.h" 297 include "dimensions.h" 298 include "dimphys.h" 299 300 integer,intent(in) :: ngrid 301 302 real, intent(out) :: reffrad(ngrid,nlayermx) !dust particles radii (K) 302 ! include "callkeys.h" 303 ! include "dimensions.h" 304 ! include "dimphys.h" 305 306 integer,intent(in) :: ngrid 307 integer,intent(in) :: nlayer 308 309 real, intent(out) :: reffrad(ngrid,nlayer) !dust particles radii (m) 303 310 304 reffrad(1:ngrid,1:nlayer mx) = 2.e-6 ! dust311 reffrad(1:ngrid,1:nlayer) = 2.e-6 ! dust 305 312 306 313 end subroutine dust_reffrad … … 309 316 310 317 !================================================================== 311 subroutine h2so4_reffrad(ngrid, reffrad)318 subroutine h2so4_reffrad(ngrid,nlayer,reffrad) 312 319 !================================================================== 313 320 ! Purpose … … 322 329 Implicit none 323 330 324 include "callkeys.h" 325 include "dimensions.h" 326 include "dimphys.h" 327 328 integer,intent(in) :: ngrid 329 330 real, intent(out) :: reffrad(ngrid,nlayermx) !h2so4 particle radii (K) 331 ! include "callkeys.h" 332 ! include "dimensions.h" 333 ! include "dimphys.h" 334 335 integer,intent(in) :: ngrid 336 integer,intent(in) :: nlayer 337 338 real, intent(out) :: reffrad(ngrid,nlayer) !h2so4 particle radii (m) 331 339 332 reffrad(1:ngrid,1:nlayer mx) = 1.e-6 ! h2so4340 reffrad(1:ngrid,1:nlayer) = 1.e-6 ! h2so4 333 341 334 342 end subroutine h2so4_reffrad … … 352 360 353 361 include "callkeys.h" 354 include "dimensions.h"355 include "dimphys.h"356 357 integer,intent(in) :: ngrid 358 359 real, intent(out) :: reffrad(ngrid,nlayer mx) ! particle radii362 ! include "dimensions.h" 363 ! include "dimphys.h" 364 365 integer,intent(in) :: ngrid 366 367 real, intent(out) :: reffrad(ngrid,nlayer) ! particle radii (m) 360 368 REAL,INTENT(IN) :: pplev(ngrid,nlayer+1) ! inter-layer pressure (Pa) 361 369 INTEGER,INTENT(IN) :: nlayer ! number of atmospheric layers -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/radinc_h.F90
r222 r227 63 63 ! These are set in sugas_corrk 64 64 ! [uses allocatable arrays] -- AS 12/2011 65 integer :: L_NPREF, L_NTREF, L_REFVAR, L_PINT 65 integer :: L_NPREF, L_NTREF, L_REFVAR, L_PINT !L_NPREF, L_NTREF, L_REFVAR, L_PINT read by master in sugas_corrk 66 66 67 67 integer, parameter :: L_NGAUSS = 17 … … 92 92 character (len=100) :: corrkdir 93 93 save corrkdir 94 !$OMP THREADPRIVATE(corrkdir) 94 95 95 96 character (len=100) :: banddir 96 97 save banddir 98 !$OMP THREADPRIVATE(banddir) 97 99 98 100 end module radinc_h -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/rain.F90
r222 r227 1 subroutine rain(ngrid,nq,ptimestep,pplev,pplay,t,pdt,pq,pdq,d_t,dqrain,dqsrain,dqssnow,rneb) 2 3 4 use ioipsl_getincom, only: getin 1 subroutine rain(ngrid,nlayer,nq,ptimestep,pplev,pplay,t,pdt,pq,pdq,d_t,dqrain,dqsrain,dqssnow,rneb) 2 3 4 ! to use 'getin' 5 ! use ioipsl_getincom 6 use ioipsl_getincom_p 5 7 use watercommon_h, only: T_h2O_ice_liq,T_h2O_ice_clouds, RLVTT, RCPD, RCPV, RV, RVTMP2,Psat_water,Tsat_water,rhowater 6 8 use radii_mod, only: h2o_cloudrad … … 22 24 !================================================================== 23 25 24 include "dimensions.h"25 include "dimphys.h"26 ! include "dimensions.h" 27 ! include "dimphys.h" 26 28 include "comcstfi.h" 27 29 include "callkeys.h" 28 30 29 31 ! Arguments 30 integer,intent(in) :: ngrid ! number of atmospherci columns 32 integer,intent(in) :: ngrid ! number of atmospheric columns 33 integer,intent(in) :: nlayer ! number of atmospheric layers 31 34 integer,intent(in) :: nq ! number of tracers 32 35 real,intent(in) :: ptimestep ! time interval 33 real,intent(in) :: pplev(ngrid,nlayer mx+1) ! inter-layer pressure (Pa)34 real,intent(in) :: pplay(ngrid,nlayer mx) ! mid-layer pressure (Pa)35 real,intent(in) :: t(ngrid,nlayer mx) ! input temperature (K)36 real,intent(in) :: pdt(ngrid,nlayer mx) ! input tendency on temperature (K/s)37 real,intent(in) :: pq(ngrid,nlayer mx,nq) ! tracers (kg/kg)38 real,intent(in) :: pdq(ngrid,nlayer mx,nq) ! input tendency on tracers39 real,intent(out) :: d_t(ngrid,nlayer mx) ! temperature tendency (K/s)40 real,intent(out) :: dqrain(ngrid,nlayer mx,nq) ! tendency of H2O precipitation (kg/kg.s-1)36 real,intent(in) :: pplev(ngrid,nlayer+1) ! inter-layer pressure (Pa) 37 real,intent(in) :: pplay(ngrid,nlayer) ! mid-layer pressure (Pa) 38 real,intent(in) :: t(ngrid,nlayer) ! input temperature (K) 39 real,intent(in) :: pdt(ngrid,nlayer) ! input tendency on temperature (K/s) 40 real,intent(in) :: pq(ngrid,nlayer,nq) ! tracers (kg/kg) 41 real,intent(in) :: pdq(ngrid,nlayer,nq) ! input tendency on tracers 42 real,intent(out) :: d_t(ngrid,nlayer) ! temperature tendency (K/s) 43 real,intent(out) :: dqrain(ngrid,nlayer,nq) ! tendency of H2O precipitation (kg/kg.s-1) 41 44 real,intent(out) :: dqsrain(ngrid) ! rain flux at the surface (kg.m-2.s-1) 42 45 real,intent(out) :: dqssnow(ngrid) ! snow flux at the surface (kg.m-2.s-1) 43 real,intent(in) :: rneb(ngrid,nlayer mx) ! cloud fraction44 45 REAL zt(ngrid,nlayer mx) ! working temperature (K)46 REAL ql(ngrid,nlayer mx) ! liquid water (Kg/Kg)47 REAL q(ngrid,nlayer mx) ! specific humidity (Kg/Kg)48 REAL d_q(ngrid,nlayer mx) ! water vapor increment49 REAL d_ql(ngrid,nlayer mx) ! liquid water / ice increment46 real,intent(in) :: rneb(ngrid,nlayer) ! cloud fraction 47 48 REAL zt(ngrid,nlayer) ! working temperature (K) 49 REAL ql(ngrid,nlayer) ! liquid water (Kg/Kg) 50 REAL q(ngrid,nlayer) ! specific humidity (Kg/Kg) 51 REAL d_q(ngrid,nlayer) ! water vapor increment 52 REAL d_ql(ngrid,nlayer) ! liquid water / ice increment 50 53 51 54 ! Subroutine options … … 62 65 REAL,PARAMETER :: Kboucher=1.19E8 63 66 REAL,SAVE :: c1 67 !$OMP THREADPRIVATE(precip_scheme,rainthreshold,cloud_sat,precip_timescale,Cboucher,c1) 64 68 65 69 INTEGER,PARAMETER :: ninter=5 66 70 67 71 logical,save :: evap_prec ! Does the rain evaporate? 72 !$OMP THREADPRIVATE(evap_prec) 68 73 69 74 ! for simple scheme … … 73 78 ! Local variables 74 79 INTEGER i, k, n 75 REAL zqs(ngrid,nlayer mx),Tsat(ngrid,nlayermx), zdelta, zcor80 REAL zqs(ngrid,nlayer),Tsat(ngrid,nlayer), zdelta, zcor 76 81 REAL zrfl(ngrid), zrfln(ngrid), zqev, zqevt 77 82 … … 80 85 REAL zchau(ngrid),zfroi(ngrid),zfrac(ngrid),zneb(ngrid) 81 86 82 real reffh2oliq(ngrid,nlayer mx),reffh2oice(ngrid,nlayermx)87 real reffh2oliq(ngrid,nlayer),reffh2oice(ngrid,nlayer) 83 88 84 89 real ttemp, ptemp, psat_tmp 85 real tnext(ngrid,nlayer mx)86 87 real l2c(ngrid,nlayer mx)90 real tnext(ngrid,nlayer) 91 92 real l2c(ngrid,nlayer) 88 93 real dWtot 89 94 … … 92 97 INTEGER, SAVE :: i_vap=0 ! water vapour 93 98 INTEGER, SAVE :: i_ice=0 ! water ice 99 !$OMP THREADPRIVATE(i_vap,i_ice) 94 100 95 101 LOGICAL,SAVE :: firstcall=.true. 102 !$OMP THREADPRIVATE(firstcall) 96 103 97 104 ! Online functions … … 114 121 write(*,*) "Precipitation scheme to use?" 115 122 precip_scheme=1 ! default value 116 call getin ("precip_scheme",precip_scheme)123 call getin_p("precip_scheme",precip_scheme) 117 124 write(*,*) " precip_scheme = ",precip_scheme 118 125 … … 120 127 write(*,*) "rainthreshold in simple scheme?" 121 128 rainthreshold=0. ! default value 122 call getin ("rainthreshold",rainthreshold)129 call getin_p("rainthreshold",rainthreshold) 123 130 write(*,*) " rainthreshold = ",rainthreshold 124 131 … … 126 133 write(*,*) "cloud water saturation level in non simple scheme?" 127 134 cloud_sat=2.6e-4 ! default value 128 call getin ("cloud_sat",cloud_sat)135 call getin_p("cloud_sat",cloud_sat) 129 136 write(*,*) " cloud_sat = ",cloud_sat 130 137 write(*,*) "precipitation timescale in non simple scheme?" 131 138 precip_timescale=3600. ! default value 132 call getin ("precip_timescale",precip_timescale)139 call getin_p("precip_timescale",precip_timescale) 133 140 write(*,*) " precip_timescale = ",precip_timescale 134 141 … … 136 143 write(*,*) "multiplicative constant in Boucher 95 precip scheme" 137 144 Cboucher=1. ! default value 138 call getin ("Cboucher",Cboucher)145 call getin_p("Cboucher",Cboucher) 139 146 write(*,*) " Cboucher = ",Cboucher 140 147 c1=1.00*1.097/rhowater*Cboucher*Kboucher … … 144 151 write(*,*) "re-evaporate precipitations?" 145 152 evap_prec=.true. ! default value 146 call getin ("evap_prec",evap_prec)153 call getin_p("evap_prec",evap_prec) 147 154 write(*,*) " evap_prec = ",evap_prec 148 155 … … 151 158 152 159 ! GCM -----> subroutine variables 153 DO k = 1, nlayer mx160 DO k = 1, nlayer 154 161 DO i = 1, ngrid 155 162 … … 172 179 173 180 ! Initialise the outputs 174 d_t(1:ngrid,1:nlayer mx) = 0.0175 d_q(1:ngrid,1:nlayer mx) = 0.0176 d_ql(1:ngrid,1:nlayer mx) = 0.0181 d_t(1:ngrid,1:nlayer) = 0.0 182 d_q(1:ngrid,1:nlayer) = 0.0 183 d_ql(1:ngrid,1:nlayer) = 0.0 177 184 zrfl(1:ngrid) = 0.0 178 185 zrfln(1:ngrid) = 0.0 179 186 180 187 ! calculate saturation mixing ratio 181 DO k = 1, nlayer mx188 DO k = 1, nlayer 182 189 DO i = 1, ngrid 183 190 ttemp = zt(i,k) … … 190 197 191 198 ! get column / layer conversion factor 192 DO k = 1, nlayer mx199 DO k = 1, nlayer 193 200 DO i = 1, ngrid 194 201 l2c(i,k)=(pplev(i,k)-pplev(i,k+1))/g … … 199 206 ! We carry the rain with us and calculate that added by warm/cold precipitation 200 207 ! processes and that subtracted by evaporation at each level. 201 DO k = nlayer mx, 1, -1208 DO k = nlayer, 1, -1 202 209 203 210 IF (evap_prec) THEN ! note no rneb dependence! … … 273 280 274 281 !recalculate liquid water particle radii 275 call h2o_cloudrad(ngrid, ql,reffh2oliq,reffh2oice)282 call h2o_cloudrad(ngrid,nlayer,ql,reffh2oliq,reffh2oice) 276 283 277 284 SELECT CASE(precip_scheme) … … 357 364 endif ! if precip_scheme=1 358 365 359 ENDDO ! of DO k = nlayer mx, 1, -1366 ENDDO ! of DO k = nlayer, 1, -1 360 367 361 368 ! Rain or snow on the ground … … 376 383 ! now subroutine -----> GCM variables 377 384 if (evap_prec) then 378 dqrain(1:ngrid,1:nlayer mx,i_vap)=d_q(1:ngrid,1:nlayermx)/ptimestep379 d_t(1:ngrid,1:nlayer mx)=d_t(1:ngrid,1:nlayermx)/ptimestep385 dqrain(1:ngrid,1:nlayer,i_vap)=d_q(1:ngrid,1:nlayer)/ptimestep 386 d_t(1:ngrid,1:nlayer)=d_t(1:ngrid,1:nlayer)/ptimestep 380 387 else 381 dqrain(1:ngrid,1:nlayer mx,i_vap)=0.0382 d_t(1:ngrid,1:nlayer mx)=0.0388 dqrain(1:ngrid,1:nlayer,i_vap)=0.0 389 d_t(1:ngrid,1:nlayer)=0.0 383 390 endif 384 dqrain(1:ngrid,1:nlayer mx,i_ice) = d_ql(1:ngrid,1:nlayermx)/ptimestep391 dqrain(1:ngrid,1:nlayer,i_ice) = d_ql(1:ngrid,1:nlayer)/ptimestep 385 392 386 393 end subroutine rain -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/setspi.F90
r222 r227 92 92 endif 93 93 94 !$OMP MASTER 94 95 nb=0 95 96 ierr=0 … … 121 122 BWNI(L_NSPECTI) =lastband(1) 122 123 BWNI(L_NSPECTI+1)=lastband(2) 124 !$OMP END MASTER 125 !$OMP BARRIER 123 126 124 127 print*,'' -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/setspv.F90
r222 r227 69 69 call abort 70 70 endif 71 71 72 !$OMP MASTER 72 73 nb=0 73 74 ierr=0 … … 98 99 BWNV(L_NSPECTV) =lastband(1) 99 100 BWNV(L_NSPECTV+1)=lastband(2) 100 101 !$OMP END MASTER 102 !$OMP BARRIER 101 103 102 104 print*,'setspv: VI band limits:' -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/soil.F
r222 r227 17 17 !----------------------------------------------------------------------- 18 18 19 include "dimensions.h"20 include "dimphys.h"19 ! include "dimensions.h" 20 ! include "dimphys.h" 21 21 22 22 … … 45 45 real,dimension(:,:),save,allocatable :: beta ! beta_k coefficients 46 46 real,save :: mu 47 !$OMP THREADPRIVATE(mthermdiff,thermdiff,coefq,coefd,alph,beta,mu) 47 48 48 49 ! local variables: -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/start2archive.F
r222 r227 27 27 ! to use 'getin' 28 28 USE ioipsl_getincom 29 USE planete_mod 29 30 30 31 implicit none 31 32 32 33 #include "dimensions.h" 34 integer, parameter :: ngridmx = (2+(jjm-1)*iim - 1/jjm) 33 35 #include "paramet.h" 34 36 #include "comconst.h" … … 41 43 #include "ener.h" 42 44 43 #include "dimphys.h"44 #include "planete.h"45 !#include "dimphys.h" 46 !#include "planete.h" 45 47 !#include"advtrac.h" 46 48 #include "netcdf.inc" … … 69 71 REAL tsoil(ngridmx,nsoilmx) ! Soil temperature 70 72 REAL co2ice(ngridmx) ! CO2 ice layer 71 REAL q2(ngridmx, nlayermx+1)73 REAL q2(ngridmx,llm+1) 72 74 REAL,ALLOCATABLE :: qsurf(:,:) 73 75 REAL emis(ngridmx) … … 80 82 ! added by FF for cloud fraction setup 81 83 REAL hice(ngridmx) 82 REAL cloudfrac(ngridmx, nlayermx),totalcloudfrac(ngridmx)84 REAL cloudfrac(ngridmx,llm),totalcloudfrac(ngridmx) 83 85 84 86 ! added by BC for slab ocean … … 199 201 200 202 201 CALL phyetat0 (ngridmx, fichnom,0,Lmodif,nsoilmx,nqtot,day_ini_fi,202 . timefi,203 CALL phyetat0 (ngridmx,llm,fichnom,0,Lmodif,nsoilmx,nqtot, 204 . day_ini_fi,timefi, 203 205 . tsurf,tsoil,emis,q2,qsurf, 204 206 ! change FF 05/2011 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/stellarlong.F
r222 r227 1 1 SUBROUTINE stellarlong(pday,pstellong) 2 3 USE planete_mod, ONLY: year_day, peri_day, e_elips, timeperi 2 4 IMPLICIT NONE 3 5 … … 40 42 c ------------- 41 43 42 #include "planete.h"44 !#include "planete.h" 43 45 #include "comcstfi.h" 44 46 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/stokes.F90
r222 r227 33 33 real a,b,molrad,visc 34 34 save a,b 35 !$OMP THREADPRIVATE(a,b) 35 36 36 37 LOGICAL firstcall 37 38 SAVE firstcall 38 39 DATA firstcall/.true./ 40 !$OMP THREADPRIVATE(firstcall) 39 41 40 42 if (firstcall) then -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/su_gases.F90
r222 r227 20 20 !================================================================== 21 21 22 !$OMP MASTER 22 23 ! load gas names from file 'gases.def' 23 24 open(90,file='gases.def',status='old',form='formatted',iostat=ierr) … … 122 123 endif 123 124 close(90) 125 !$OMP END MASTER 126 !$OMP BARRIER 124 127 125 128 end subroutine su_gases -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/suaer_corrk.F90
r222 r227 46 46 #include "callkeys.h" 47 47 ! Optical properties (read in external ASCII files) 48 INTEGER 49 ! the domain (VIS or IR) 48 INTEGER,SAVE :: nwvl ! Number of wavelengths in 49 ! the domain (VIS or IR), read by master 50 50 51 51 ! REAL :: solsir ! visible to infrared ratio … … 53 53 54 54 REAL, DIMENSION(:),& 55 ALLOCATABLE, SAVE :: wvl ! Wavelength axis 55 ALLOCATABLE, SAVE :: wvl ! Wavelength axis, read by master 56 56 REAL, DIMENSION(:),& 57 ALLOCATABLE, SAVE :: radiusdyn ! Particle size axis 57 ALLOCATABLE, SAVE :: radiusdyn ! Particle size axis, read by master 58 58 59 59 REAL, DIMENSION(:,:),& 60 ALLOCATABLE, SAVE :: ep,& ! Extinction coefficient Qext 61 omeg,& ! Single Scattering Albedo 62 gfactor ! Assymetry Factor 60 ALLOCATABLE, SAVE :: ep,& ! Extinction coefficient Qext, read by master 61 omeg,& ! Single Scattering Albedo, read by master 62 gfactor ! Assymetry Factor, read by master 63 63 64 64 ! Local variables: … … 98 98 99 99 CHARACTER(LEN=30), DIMENSION(naerkind,2), SAVE :: file_id 100 !$OMP THREADPRIVATE(file_id) 100 101 !---- Please indicate the names of the optical property files below 101 102 ! Please also choose the reference wavelengths of each aerosol … … 203 204 ! 1.1 Open the ASCII file 204 205 206 !$OMP MASTER 205 207 INQUIRE(FILE=TRIM(datadir)//& 206 208 '/'//TRIM(file_id(iaer,idomain)),& … … 339 341 endif 340 342 343 !$OMP END MASTER 344 !$OMP BARRIER 341 345 342 346 … … 432 436 !======================================================================== 433 437 434 DEALLOCATE(wvl) ! wvl 435 DEALLOCATE(radiusdyn) ! radiusdyn 436 DEALLOCATE(ep) ! ep 437 DEALLOCATE(omeg) ! omeg 438 DEALLOCATE(gfactor) ! g 438 !$OMP BARRIER 439 !$OMP MASTER 440 IF (ALLOCATED(wvl)) DEALLOCATE(wvl) ! wvl 441 IF (ALLOCATED(radiusdyn)) DEALLOCATE(radiusdyn) ! radiusdyn 442 IF (ALLOCATED(ep)) DEALLOCATE(ep) ! ep 443 IF (ALLOCATED(omeg)) DEALLOCATE(omeg) ! omeg 444 IF (ALLOCATED(gfactor)) DEALLOCATE(gfactor) ! g 445 !$OMP END MASTER 446 !$OMP BARRIER 439 447 440 448 END DO ! Loop on iaer -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/sugas_corrk.F90
r222 r227 28 28 29 29 use gases_h 30 use ioipsl_getincom 30 ! use ioipsl_getincom 31 use ioipsl_getincom_p 31 32 implicit none 32 33 … … 45 46 46 47 ! ALLOCATABLE ARRAYS -- AS 12/2011 47 REAL*8, DIMENSION(:,:,:,:,:), ALLOCATABLE :: gasi8, gasv848 character*20,allocatable,DIMENSION(:) :: gastype ! for check with gnom48 REAL*8, DIMENSION(:,:,:,:,:), ALLOCATABLE,SAVE :: gasi8, gasv8 !read by master 49 character*20,allocatable,DIMENSION(:),SAVE :: gastype ! for check with gnom, read by master 49 50 50 51 real*8 x, xi(4), yi(4), ans, p … … 76 77 endif 77 78 79 !$OMP MASTER 78 80 ! check that database matches varactive toggle 79 81 open(111,file=TRIM(file_path),form='formatted') … … 261 263 tgasmin = tgasref(1) 262 264 tgasmax = tgasref(L_NTREF) 265 !$OMP END MASTER 266 !$OMP BARRIER 263 267 264 268 !----------------------------------------------------------------------- … … 302 306 write(*,*)"graybody: constant absorption coefficient in visible:" 303 307 kappa_VI=-100000. 304 call getin ("kappa_VI",kappa_VI)308 call getin_p("kappa_VI",kappa_VI) 305 309 write(*,*)" kappa_VI = ",kappa_VI 306 310 kappa_VI=kappa_VI*1.e4* mugaz * 1.672621e-27 ! conversion from m^2/kg to cm^2/molecule … … 309 313 write(*,*)"graybody: constant absorption coefficient in InfraRed:" 310 314 kappa_IR=-100000. 311 call getin ("kappa_IR",kappa_IR)315 call getin_p("kappa_IR",kappa_IR) 312 316 write(*,*)" kappa_IR = ",kappa_IR 313 317 kappa_IR=kappa_IR*1.e4* mugaz * 1.672621e-27 ! conversion from m^2/kg to cm^2/molecule … … 320 324 End if 321 325 326 !$OMP MASTER 322 327 ! print*,corrkdir(1:4) 323 328 ! VISIBLE … … 361 366 gasv8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTV,1:L_NGAUSS)=0.0 362 367 endif 368 !$OMP END MASTER 369 !$OMP BARRIER 363 370 364 371 ! INFRA-RED 365 372 if ((corrkdir(1:4).eq.'null'))then !.or.(TRIM(corrkdir).eq.'null_LowTeffStar')) then 366 373 print*,'Infrared corrk gaseous absorption is set to zero if graybody=F' 374 !$OMP MASTER 367 375 gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTI,1:L_NGAUSS)=0.0 368 else 376 !$OMP END MASTER 377 !$OMP BARRIER 378 else 369 379 file_id='/corrk_data/'//trim(adjustl(banddir))//'/corrk_gcm_IR.dat' 370 380 file_path=TRIM(datadir)//TRIM(file_id) … … 379 389 endif 380 390 391 !$OMP MASTER 381 392 open(111,file=TRIM(file_path),form='formatted') 382 393 read(111,*) gasi8 383 394 close(111) 395 !$OMP END MASTER 396 !$OMP BARRIER 384 397 385 398 ! 'fzero' is a currently unused feature that allows optimisation … … 423 436 endif 424 437 438 !$OMP MASTER 425 439 if(nIR_limit.eq.0) then 426 440 gasi8(1:L_NTREF,1:L_NPREF,1:L_REFVAR,1:L_NSPECTI,1:L_NGAUSS)= & … … 465 479 end do 466 480 end do 481 !$OMP END MASTER 482 !$OMP BARRIER 467 483 468 484 ! Interpolate the values: first the longwave … … 653 669 654 670 ! Deallocate local arrays 671 !$OMP BARRIER 672 !$OMP MASTER 655 673 IF( ALLOCATED( gasi8 ) ) DEALLOCATE( gasi8 ) 656 674 IF( ALLOCATED( gasv8 ) ) DEALLOCATE( gasv8 ) 657 675 IF( ALLOCATED( pgasref ) ) DEALLOCATE( pgasref ) 658 676 IF( ALLOCATED( gastype ) ) DEALLOCATE( gastype ) 677 !$OMP END MASTER 678 !$OMP BARRIER 659 679 660 680 return -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/surf_heat_transp_mod.F90
r222 r227 7 7 CONTAINS 8 8 9 SUBROUTINE divgrad_phy(n levs,temp,delta)9 SUBROUTINE divgrad_phy(ngrid,nlevs,temp,delta) 10 10 11 11 … … 14 14 15 15 #include "dimensions.h" 16 #include "dimphys.h"16 !#include "dimphys.h" 17 17 #include "paramet.h" 18 18 #include "comcstfi.h" … … 20 20 #include "comhdiff.h" 21 21 22 INTEGER nlevs, ll 23 REAL temp(ngridmx,nlevs) 22 INTEGER,INTENT(IN) :: ngrid, nlevs 23 REAL,INTENT(IN) :: temp(ngrid,nlevs) 24 REAL,INTENT(OUT) :: delta(ngrid,nlevs) 24 25 REAL delta_2d(ip1jmp1,nlevs) 25 REAL delta(ngridmx,nlevs)26 INTEGER :: ll 26 27 REAL ghx(ip1jmp1,nlevs), ghy(ip1jm,nlevs) 27 28 28 29 29 CALL gr_fi_dyn(nlevs,ngrid mx,iip1,jjp1,temp,delta_2d)30 CALL gr_fi_dyn(nlevs,ngrid,iip1,jjp1,temp,delta_2d) 30 31 CALL grad(nlevs,delta_2d,ghx,ghy) 31 32 DO ll=1,nlevs … … 36 37 END DO 37 38 CALL diverg(nlevs,ghx,ghy,delta_2d) 38 CALL gr_dyn_fi(nlevs,iip1,jjp1,ngrid mx,delta_2d,delta)39 CALL gr_dyn_fi(nlevs,iip1,jjp1,ngrid,delta_2d,delta) 39 40 40 41 … … 43 44 44 45 45 SUBROUTINE init_masquv( zmasq)46 SUBROUTINE init_masquv(ngrid,zmasq) 46 47 47 48 IMPLICIT NONE 48 49 49 50 #include "dimensions.h" 50 #include "dimphys.h"51 !#include "dimphys.h" 51 52 #include "paramet.h" 52 53 #include "comcstfi.h" … … 55 56 56 57 57 REAL zmasq(ngridmx) 58 INTEGER,INTENT(IN) :: ngrid 59 REAL zmasq(ngrid) 58 60 REAL zmasq_2d(ip1jmp1) 59 61 REAL ff(ip1jm) … … 66 68 zmasqv=1. 67 69 68 CALL gr_fi_dyn(1,ngrid mx,iip1,jjp1,zmasq,zmasq_2d)70 CALL gr_fi_dyn(1,ngrid,iip1,jjp1,zmasq,zmasq_2d) 69 71 70 72 DO i=1,ip1jmp1-1 … … 104 106 105 107 106 SUBROUTINE slab_ekman2( tx_phy,ty_phy,ts_phy,dt_phy)108 SUBROUTINE slab_ekman2(ngrid,tx_phy,ty_phy,ts_phy,dt_phy) 107 109 108 110 use slab_ice_h … … 111 113 112 114 #include "dimensions.h" 113 #include "dimphys.h"115 !#include "dimphys.h" 114 116 #include "paramet.h" 115 117 #include "comcstfi.h" … … 118 120 #include "comhdiff.h" 119 121 122 INTEGER,INTENT(IN) :: ngrid 120 123 INTEGER ij 121 124 REAL txv(ip1jm),fluxm(ip1jm),tyv(ip1jm) … … 123 126 REAL tyu(ip1jmp1),txu(ip1jmp1),fluxz(ip1jmp1),fluxv(ip1jmp1) 124 127 REAL dt(ip1jmp1,noceanmx),ts(ip1jmp1,noceanmx) 125 REAL tx_phy(ngrid mx),ty_phy(ngridmx)126 REAL dt_phy(ngrid mx,noceanmx),ts_phy(ngridmx,noceanmx)128 REAL tx_phy(ngrid),ty_phy(ngrid) 129 REAL dt_phy(ngrid,noceanmx),ts_phy(ngrid,noceanmx) 127 130 128 131 … … 130 133 131 134 ! Passage taux,y sur grilles 2d 132 CALL gr_fi_dyn(1,ngrid mx,iip1,jjp1,tx_phy,txu)133 CALL gr_fi_dyn(1,ngrid mx,iip1,jjp1,ty_phy,tyu)135 CALL gr_fi_dyn(1,ngrid,iip1,jjp1,tx_phy,txu) 136 CALL gr_fi_dyn(1,ngrid,iip1,jjp1,ty_phy,tyu) 134 137 ! Passage grille u,v 135 138 ! Multiplication par f ou eps. … … 144 147 145 148 ! Calcul de T, Tdeep 146 CALL gr_fi_dyn(2,ngrid mx,iip1,jjp1,ts_phy,ts)149 CALL gr_fi_dyn(2,ngrid,iip1,jjp1,ts_phy,ts) 147 150 148 151 ! Flux de masse … … 222 225 223 226 ! Retour grille physique 224 CALL gr_dyn_fi(2,iip1,jjp1,ngrid mx,dt,dt_phy)227 CALL gr_dyn_fi(2,iip1,jjp1,ngrid,dt,dt_phy) 225 228 226 229 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/surface_nature.F
r222 r227 33 33 !================================================================== 34 34 35 #include "dimensions.h"36 #include "dimphys.h"35 !#include "dimensions.h" 36 !#include "dimphys.h" 37 37 #include "comcstfi.h" 38 38 #include "callkeys.h" -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/surfdat_h.F90
r222 r227 5 5 6 6 real,allocatable,dimension(:) :: albedodat ! albedo of bare ground 7 !$OMP THREADPRIVATE(albedodat) 7 8 ! Ehouarn: moved inertiedat to comsoil.h 8 9 ! real inertiedat, ! thermal inertia 9 10 real,allocatable,dimension(:) :: phisfi ! geopotential at ground level 11 !$OMP THREADPRIVATE(phisfi) 10 12 real,dimension(2) :: albedice 11 13 real,dimension(2) :: emisice ! ice emissivity; 1:Northern hemisphere 2:Southern hemisphere 12 14 real emissiv 13 15 real,dimension(2) :: iceradius, dtemisice 16 !$OMP THREADPRIVATE(albedice,emisice,emissiv,iceradius,dtemisice) 14 17 real,allocatable,dimension(:) :: zmea,zstd,zsig,zgam,zthe 18 !$OMP THREADPRIVATE(zmea,zstd,zsig,zgam,zthe) 15 19 16 20 real,allocatable,dimension(:) :: dryness !"Dryness coefficient" for grnd water ice sublimation … … 18 22 19 23 logical,allocatable,dimension(:) :: watercaptag !! was in watercap.h 24 !$OMP THREADPRIVATE(dryness,watercaptag) 20 25 21 26 end module surfdat_h -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/surfini.F
r222 r227 15 15 c Declarations: 16 16 c ------------- 17 #include "dimensions.h"18 #include "dimphys.h"17 !#include "dimensions.h" 18 !#include "dimphys.h" 19 19 #include "callkeys.h" 20 20 c -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/tabfi.F
r222 r227 34 34 c comparer avec le day_ini dynamique) 35 35 c 36 c - lmax: tab_cntrl(tab0+2) (pour test avec nlayer mx)36 c - lmax: tab_cntrl(tab0+2) (pour test avec nlayer) 37 37 c 38 38 c - p_rad … … 51 51 use iostart, only: get_var 52 52 use mod_phys_lmdz_para, only: is_parallel 53 use planete_mod, only: year_day, periastr, apoastr, peri_day, 54 & obliquit, z0, lmixmin, emin_turb 53 55 54 56 implicit none 55 57 56 58 #include "comcstfi.h" 57 #include "planete.h"59 !#include "planete.h" 58 60 #include "netcdf.inc" 59 61 #include "callkeys.h" -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/totalcloudfrac.F90
r222 r227 1 subroutine totalcloudfrac(ngrid,n q,rneb,totalrneb,pplev,pq,tau)1 subroutine totalcloudfrac(ngrid,nlayer,nq,rneb,totalrneb,pplev,pq,tau) 2 2 3 3 use watercommon_h … … 19 19 !================================================================== 20 20 21 #include "dimensions.h"22 #include "dimphys.h"21 !#include "dimensions.h" 22 !#include "dimphys.h" 23 23 #include "comcstfi.h" 24 24 #include "callkeys.h" 25 25 26 26 integer,intent(in) :: ngrid ! number of atmospheric columns 27 integer,intent(in) :: nlayer ! number of atmospheric layers 27 28 integer,intent(in) :: nq ! number of tracers 28 real,intent(in) :: rneb(ngrid,nlayer mx) ! cloud fraction29 real,intent(in) :: rneb(ngrid,nlayer) ! cloud fraction 29 30 real,intent(out) :: totalrneb(ngrid) ! total cloud fraction 30 real,intent(in) :: pplev(ngrid,nlayer mx+1) ! inter-layer pressure (Pa)31 real,intent(in) :: pq(ngrid,nlayer mx,nq) ! tracers (.../kg_of_air)32 real,intent(in) :: tau(ngrid,nlayer mx)31 real,intent(in) :: pplev(ngrid,nlayer+1) ! inter-layer pressure (Pa) 32 real,intent(in) :: pq(ngrid,nlayer,nq) ! tracers (.../kg_of_air) 33 real,intent(in) :: tau(ngrid,nlayer) 33 34 34 real, dimension(nlayer mx+1) :: masse35 real, dimension(nlayer+1) :: masse 35 36 integer, parameter :: recovery=7 36 37 integer ltau_max … … 47 48 real clear,tau_min 48 49 real, parameter :: tau_c=0.1 !threshold of optical depth for the calculation of total cloud fraction 49 real rneb2(nlayer mx)50 real rneb2(nlayer) 50 51 51 52 … … 55 56 if (recovery.eq.1) then 56 57 clear = (1.-rneb(ig,1)) 57 do l=2,nlayer mx58 do l=2,nlayer 58 59 clear = clear*(1.-rneb(ig,l)) 59 60 enddo … … 62 63 elseif (recovery.eq.2) then 63 64 totalrneb(ig) = rneb(ig,1) 64 do l=2,14 !nlayer mx65 do l=2,14 !nlayer 65 66 totalrneb(ig) = max(rneb(ig,l),totalrneb(ig)) 66 67 enddo … … 68 69 elseif (recovery.eq.3) then 69 70 totalrneb(ig) = rneb(ig,1) 70 do l=2,nlayer mx71 do l=2,nlayer 71 72 totalrneb(ig) = min(rneb(ig,l),totalrneb(ig)) 72 73 enddo … … 77 78 elseif (recovery.eq.5) then 78 79 totalrneb(ig) = rneb(ig,1) 79 do l=1,nlayer mx80 do l=1,nlayer 80 81 masse(l)=pq(ig,l,igcm_h2o_ice)*(pplev(ig,l)-pplev(ig,l+1)) 81 82 enddo … … 85 86 elseif (recovery.eq.6) then 86 87 totalrneb(ig) = 0. 87 do l=1,nlayer mx88 do l=1,nlayer 88 89 masse(l)=pq(ig,l,igcm_h2o_ice)*(pplev(ig,l)-pplev(ig,l+1)) 89 90 masse(l)=max(masse(l),0.) 90 91 enddo 91 92 massetot=sum(masse,dim=1) 92 do l=1,nlayer mx93 do l=1,nlayer 93 94 totalrneb(ig) = totalrneb(ig)+rneb(ig,l)*masse(l)/massetot 94 95 enddo … … 96 97 elseif (recovery.eq.7) then 97 98 98 rneb2(:)=rneb(ig,1:nlayer mx)99 tau_min=MIN(tau_c,MAXVAL(tau(ig,1:nlayer mx))/2.)100 do l=1,nlayer mx99 rneb2(:)=rneb(ig,1:nlayer) 100 tau_min=MIN(tau_c,MAXVAL(tau(ig,1:nlayer))/2.) 101 do l=1,nlayer 101 102 if(tau(ig,l)<tau_min) rneb2(l)=0. 102 103 enddo 103 totalrneb(ig)=maxval(rneb2(1:nlayer mx))104 totalrneb(ig)=maxval(rneb2(1:nlayer)) 104 105 105 106 endif ! (recovery=) -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/tracer_h.F90
r222 r227 19 19 real rho_co2 ! CO2 ice density (kg.m-3) 20 20 real ref_r0 ! for computing reff=ref_r0*r0 (in log.n. distribution) 21 !$OMP THREADPRIVATE(noms,mmol,radius,rho_q,qext,alpha_lift,alpha_devil,qextrhor, & 22 !$OMP varian,r3n_q,rho_dust,rho_ice,rho_co2,ref_r0) 21 23 22 24 ! tracer indexes: these are initialized in initracer and should be 0 if the … … 47 49 integer :: igcm_ar_n2 ! for simulations using co2 +neutral gaz 48 50 integer :: igcm_co2_ice ! CO2 ice 51 !$OMP THREADPRIVATE(igcm_dustbin,igcm_dust_mass,igcm_dust_number,igcm_h2o_vap,igcm_h2o_ice, & 52 !$OMP igcm_co2,igcm_co,igcm_o,igcm_o1d,igcm_o2,igcm_o3,igcm_h,igcm_h2,igcm_oh, & 53 !$OMP igcm_ho2,igcm_h2o2,igcm_n2,igcm_ar,igcm_ar_n2,igcm_co2_ice) 49 54 50 55 end module tracer_h -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/turbdiff.F90
r222 r227 41 41 ! ------------ 42 42 43 include "dimensions.h"44 include "dimphys.h"43 ! include "dimensions.h" 44 ! include "dimphys.h" 45 45 include "comcstfi.h" 46 46 include "callkeys.h" … … 48 48 ! arguments 49 49 ! --------- 50 INTEGER,INTENT(IN) :: ngrid,nlay 50 INTEGER,INTENT(IN) :: ngrid 51 INTEGER,INTENT(IN) :: nlay 51 52 REAL,INTENT(IN) :: ptimestep 52 53 REAL,INTENT(IN) :: pplay(ngrid,nlay),pplev(ngrid,nlay+1) … … 86 87 87 88 REAL z4st,zdplanck(ngrid) 88 REAL zkv(ngrid,nlay ermx+1),zkh(ngrid,nlayermx+1)89 REAL zkv(ngrid,nlay+1),zkh(ngrid,nlay+1) 89 90 REAL zcdv(ngrid),zcdh(ngrid) 90 91 REAL zcdv_true(ngrid),zcdh_true(ngrid) 91 REAL zu(ngrid,nlay ermx),zv(ngrid,nlayermx)92 REAL zh(ngrid,nlay ermx),zt(ngrid,nlayermx)92 REAL zu(ngrid,nlay),zv(ngrid,nlay) 93 REAL zh(ngrid,nlay),zt(ngrid,nlay) 93 94 REAL ztsrf(ngrid) 94 95 REAL z1(ngrid),z2(ngrid) 95 REAL zmass(ngrid,nlay ermx)96 REAL zfluxv(ngrid,nlay ermx),zfluxt(ngrid,nlayermx),zfluxq(ngrid,nlayermx)97 REAL zb0(ngrid,nlay ermx)98 REAL zExner(ngrid,nlay ermx),zovExner(ngrid,nlayermx)99 REAL zcv(ngrid,nlay ermx),zdv(ngrid,nlayermx) !inversion coefficient for winds100 REAL zct(ngrid,nlay ermx),zdt(ngrid,nlayermx) !inversion coefficient for temperature101 REAL zcq(ngrid,nlay ermx),zdq(ngrid,nlayermx) !inversion coefficient for tracers96 REAL zmass(ngrid,nlay) 97 REAL zfluxv(ngrid,nlay),zfluxt(ngrid,nlay),zfluxq(ngrid,nlay) 98 REAL zb0(ngrid,nlay) 99 REAL zExner(ngrid,nlay),zovExner(ngrid,nlay) 100 REAL zcv(ngrid,nlay),zdv(ngrid,nlay) !inversion coefficient for winds 101 REAL zct(ngrid,nlay),zdt(ngrid,nlay) !inversion coefficient for temperature 102 REAL zcq(ngrid,nlay),zdq(ngrid,nlay) !inversion coefficient for tracers 102 103 REAL zcst1 103 104 REAL zu2!, a … … 106 107 107 108 LOGICAL,SAVE :: firstcall=.true. 109 !$OMP THREADPRIVATE(firstcall) 108 110 109 111 ! Tracers 110 112 ! ------- 111 113 INTEGER iq 112 REAL zq(ngrid,nlay ermx,nq)113 REAL zqnoevap(ngrid,nlay ermx) !special for water case to compute where evaporated water goes.114 REAL pdqevap(ngrid,nlay ermx) !special for water case to compute where evaporated water goes.114 REAL zq(ngrid,nlay,nq) 115 REAL zqnoevap(ngrid,nlay) !special for water case to compute where evaporated water goes. 116 REAL pdqevap(ngrid,nlay) !special for water case to compute where evaporated water goes. 115 117 REAL zdmassevap(ngrid) 116 118 REAL rho(ngrid) ! near-surface air density … … 123 125 124 126 integer, save :: ivap, iliq, iliq_surf,iice_surf ! also make liq for clarity on surface... 127 !$OMP THREADPRIVATE(ivap,iliq,iliq_surf,iice_surf) 125 128 126 129 real, parameter :: karman=0.4 … … 241 244 ! ------------------------------------------------------ 242 245 243 call vdif_kc(ngrid, ptimestep,g,pzlev,pzlay,pu,pv,zh,zcdv_true,pq2,zkv,zkh) !JL12 why not call vdif_kc with updated winds and temperature246 call vdif_kc(ngrid,nlay,ptimestep,g,pzlev,pzlay,pu,pv,zh,zcdv_true,pq2,zkv,zkh) !JL12 why not call vdif_kc with updated winds and temperature 244 247 245 248 ! Adding eddy mixing to mimic 3D general circulation in 1D … … 262 265 263 266 !JL12 calculate the flux coefficients (tables multiplied elements by elements) 264 zfluxv(1:ngrid,1:nlay ermx)=zkv(1:ngrid,1:nlayermx)*zb0(1:ngrid,1:nlayermx)267 zfluxv(1:ngrid,1:nlay)=zkv(1:ngrid,1:nlay)*zb0(1:ngrid,1:nlay) 265 268 266 269 !----------------------------------------------------------------------- … … 364 367 ! JL12 calculate the flux coefficients (tables multiplied elements by elements) 365 368 ! --------------------------------------------------------------- 366 zfluxq(1:ngrid,1:nlay ermx)=zkh(1:ngrid,1:nlayermx)*zb0(1:ngrid,1:nlayermx) !JL12 we save zfluxq which doesn't need the Exner factor367 zfluxt(1:ngrid,1:nlay ermx)=zfluxq(1:ngrid,1:nlayermx)*zExner(1:ngrid,1:nlayermx)369 zfluxq(1:ngrid,1:nlay)=zkh(1:ngrid,1:nlay)*zb0(1:ngrid,1:nlay) !JL12 we save zfluxq which doesn't need the Exner factor 370 zfluxt(1:ngrid,1:nlay)=zfluxq(1:ngrid,1:nlay)*zExner(1:ngrid,1:nlay) 368 371 369 372 DO ig=1,ngrid -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/vdif_cd.F
r222 r227 55 55 DATA firstcal/.true./ 56 56 SAVE b,c,d,karman,c2b,c3bc,c3b,firstcal,umin2 57 !$OMP THREADPRIVATE(b,c,d,karman,c2b,c3bc,c3b,firstcal,umin2) 57 58 58 59 c----------------------------------------------------------------------- -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/vdif_kc.F
r222 r227 1 SUBROUTINE vdif_kc(ngrid, dt,g,zlev,zlay,u,v,teta,cd,q2,km,kn)1 SUBROUTINE vdif_kc(ngrid,nlay,dt,g,zlev,zlay,u,v,teta,cd,q2,km,kn) 2 2 IMPLICIT NONE 3 3 c....................................................................... 4 #include "dimensions.h"5 #include "dimphys.h"4 !#include "dimensions.h" 5 !#include "dimphys.h" 6 6 c....................................................................... 7 7 c … … 28 28 c....................................................................... 29 29 INTEGER,INTENT(IN) :: ngrid 30 INTEGER,INTENT(IN) :: nlay 30 31 REAL,INTENT(IN) :: dt,g 31 REAL,INTENT(IN) :: zlev(ngrid,nlay ermx+1)32 REAL,INTENT(IN) :: zlay(ngrid,nlay ermx)33 REAL,INTENT(IN) :: u(ngrid,nlay ermx)34 REAL,INTENT(IN) :: v(ngrid,nlay ermx)35 REAL,INTENT(IN) :: teta(ngrid,nlay ermx)32 REAL,INTENT(IN) :: zlev(ngrid,nlay+1) 33 REAL,INTENT(IN) :: zlay(ngrid,nlay) 34 REAL,INTENT(IN) :: u(ngrid,nlay) 35 REAL,INTENT(IN) :: v(ngrid,nlay) 36 REAL,INTENT(IN) :: teta(ngrid,nlay) 36 37 REAL,INTENT(IN) :: cd(ngrid) 37 REAL,INTENT(INOUT) :: q2(ngrid,nlay ermx+1)38 REAL,INTENT(OUT) :: km(ngrid,nlay ermx+1)39 REAL,INTENT(OUT) :: kn(ngrid,nlay ermx+1)38 REAL,INTENT(INOUT) :: q2(ngrid,nlay+1) 39 REAL,INTENT(OUT) :: km(ngrid,nlay+1) 40 REAL,INTENT(OUT) :: kn(ngrid,nlay+1) 40 41 c....................................................................... 41 42 c … … 50 51 c 51 52 c....................................................................... 52 INTEGER,PARAMETER :: nlay=nlayermx 53 INTEGER,PARAMETER :: nlev=nlayermx+1 54 REAL unsdz(ngrid,nlayermx) 55 REAL unsdzdec(ngrid,nlayermx+1) 56 REAL q(ngrid,nlayermx+1) 53 INTEGER :: nlev 54 REAL unsdz(ngrid,nlay) 55 REAL unsdzdec(ngrid,nlay+1) 56 REAL q(ngrid,nlay+1) 57 57 c....................................................................... 58 58 c … … 64 64 c 65 65 c....................................................................... 66 REAL kmpre(ngrid,nlay ermx+1)66 REAL kmpre(ngrid,nlay+1) 67 67 REAL qcstat 68 68 REAL q2cstat … … 72 72 c 73 73 c....................................................................... 74 REAL long(ngrid,nlay ermx+1)74 REAL long(ngrid,nlay+1) 75 75 c....................................................................... 76 76 c … … 95 95 REAL mcstat 96 96 REAL m2cstat 97 REAL m(ngrid,nlay ermx+1)98 REAL mpre(ngrid,nlay ermx+1)99 REAL m2(ngrid,nlay ermx+1)100 REAL n2(ngrid,nlay ermx+1)97 REAL m(ngrid,nlay+1) 98 REAL mpre(ngrid,nlay+1) 99 REAL m2(ngrid,nlay+1) 100 REAL n2(ngrid,nlay+1) 101 101 c....................................................................... 102 102 c … … 120 120 LOGICAL gnsup 121 121 REAL gm 122 c REAL ri(ngrid,nlaye rmx+1)123 REAL sn(ngrid,nlay ermx+1)124 REAL snq2(ngrid,nlay ermx+1)125 REAL sm(ngrid,nlay ermx+1)126 REAL smq2(ngrid,nlay ermx+1)122 c REAL ri(ngrid,nlaye+1) 123 REAL sn(ngrid,nlay+1) 124 REAL snq2(ngrid,nlay+1) 125 REAL sm(ngrid,nlay+1) 126 REAL smq2(ngrid,nlay+1) 127 127 c....................................................................... 128 128 c … … 179 179 180 180 ! initialization of local variables: 181 nlev=nlay+1 181 182 long(:,:)=0. 182 183 n2(:,:)=0. -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/vdifc.F
r222 r227 38 38 ! ------------ 39 39 40 #include "dimensions.h"41 #include "dimphys.h"40 !#include "dimensions.h" 41 !#include "dimphys.h" 42 42 #include "comcstfi.h" 43 43 #include "callkeys.h" … … 77 77 78 78 REAL z4st,zdplanck(ngrid) 79 REAL zkv(ngrid,nlay ermx+1),zkh(ngrid,nlayermx+1)79 REAL zkv(ngrid,nlay+1),zkh(ngrid,nlay+1) 80 80 REAL zcdv(ngrid),zcdh(ngrid) 81 81 REAL zcdv_true(ngrid),zcdh_true(ngrid) 82 REAL zu(ngrid,nlay ermx),zv(ngrid,nlayermx)83 REAL zh(ngrid,nlay ermx)82 REAL zu(ngrid,nlay),zv(ngrid,nlay) 83 REAL zh(ngrid,nlay) 84 84 REAL ztsrf2(ngrid) 85 85 REAL z1(ngrid),z2(ngrid) 86 REAL za(ngrid,nlay ermx),zb(ngrid,nlayermx)87 REAL zb0(ngrid,nlay ermx)88 REAL zc(ngrid,nlay ermx),zd(ngrid,nlayermx)86 REAL za(ngrid,nlay),zb(ngrid,nlay) 87 REAL zb0(ngrid,nlay) 88 REAL zc(ngrid,nlay),zd(ngrid,nlay) 89 89 REAL zcst1 90 90 REAL zu2!, a 91 REAL zcq(ngrid,nlay ermx),zdq(ngrid,nlayermx)91 REAL zcq(ngrid,nlay),zdq(ngrid,nlay) 92 92 REAL evap(ngrid) 93 93 REAL zcq0(ngrid),zdq0(ngrid) … … 96 96 LOGICAL firstcall 97 97 SAVE firstcall 98 !$OMP THREADPRIVATE(firstcall) 98 99 99 100 LOGICAL lastcall … … 101 102 ! variables added for CO2 condensation 102 103 ! ------------------------------------ 103 REAL hh !, zhcond(ngrid,nlay ermx)104 REAL hh !, zhcond(ngrid,nlay) 104 105 ! REAL latcond,tcond1mb 105 106 ! REAL acond,bcond 106 107 ! SAVE acond,bcond 108 !!$OMP THREADPRIVATE(acond,bcond) 107 109 ! DATA latcond,tcond1mb/5.9e5,136.27/ 108 110 … … 110 112 ! ------- 111 113 INTEGER iq 112 REAL zq(ngrid,nlay ermx,nq)114 REAL zq(ngrid,nlay,nq) 113 115 REAL zq1temp(ngrid) 114 116 REAL rho(ngrid) ! near-surface air density … … 123 125 real z1_T(ngrid), z2_T(ngrid) 124 126 real zb_T(ngrid) 125 real zc_T(ngrid,nlay ermx)126 real zd_T(ngrid,nlay ermx)127 real zc_T(ngrid,nlay) 128 real zd_T(ngrid,nlay) 127 129 real lat1(ngrid), lat2(ngrid) 128 130 real surfh2otot … … 132 134 integer ivap, iice ! also make liq for clarity on surface... 133 135 save ivap, iice 136 !$OMP THREADPRIVATE(ivap,iice) 134 137 135 138 real, parameter :: karman=0.4 … … 249 252 ! ------------------------------------------------------ 250 253 251 call vdif_kc(ngrid, ptimestep,g,pzlev,pzlay254 call vdif_kc(ngrid,nlay,ptimestep,g,pzlev,pzlay 252 255 & ,pu,pv,ph,zcdv_true 253 256 & ,pq2,zkv,zkh) -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/vlz_fi.F
r222 r227 1 SUBROUTINE vlz_fi(ngrid, q,pente_max,masse,w,wq)1 SUBROUTINE vlz_fi(ngrid,nlayer,q,pente_max,masse,w,wq) 2 2 c 3 3 c Auteurs: P.Le Van, F.Hourdin, F.Forget … … 16 16 IMPLICIT NONE 17 17 c 18 #include "dimensions.h"19 #include "dimphys.h"18 !#include "dimensions.h" 19 !#include "dimphys.h" 20 20 21 21 c … … 23 23 c Arguments: 24 24 c ---------- 25 integer ngrid26 real masse(ngrid,llm),pente_max27 REAL q(ngrid,llm)28 REAL w(ngrid,llm)29 REAL wq(ngrid,llm+1)25 integer,intent(in) :: ngrid, nlayer 26 real,intent(in) :: masse(ngrid,nlayer),pente_max 27 REAL,INTENT(INOUT) :: q(ngrid,nlayer) 28 REAL,INTENT(INOUT) :: w(ngrid,nlayer) 29 REAL,INTENT(OUT) :: wq(ngrid,nlayer+1) 30 30 c 31 31 c Local … … 35 35 c 36 36 37 real dzq(ngrid,llm),dzqw(ngrid,llm),adzqw(ngrid,llm),dzqmax 37 real dzq(ngrid,nlayer),dzqw(ngrid,nlayer),adzqw(ngrid,nlayer) 38 real dzqmax 38 39 real newmasse 39 40 real sigw, Mtot, MQtot … … 47 48 c sens de W 48 49 49 do l=2, llm50 do l=2,nlayer 50 51 do ij=1,ngrid 51 52 dzqw(ij,l)=q(ij,l-1)-q(ij,l) … … 54 55 enddo 55 56 56 do l=2, llm-157 do l=2,nlayer-1 57 58 do ij=1,ngrid 58 59 #ifdef CRAY … … 73 74 do ij=1,ngrid 74 75 dzq(ij,1)=0. 75 dzq(ij, llm)=0.76 dzq(ij,nlayer)=0. 76 77 enddo 77 78 c --------------------------------------------------------------- … … 83 84 c No flux at the model top: 84 85 do ij=1,ngrid 85 wq(ij, llm+1)=0.86 wq(ij,nlayer+1)=0. 86 87 enddo 87 88 … … 89 90 c =============================== 90 91 91 do l = 1, llm! loop different than when w<092 do l = 1,nlayer ! loop different than when w<0 92 93 do ij=1,ngrid 93 94 … … 107 108 Mtot = masse(ij,m) 108 109 MQtot = masse(ij,m)*q(ij,m) 109 if(m.ge. llm)goto 88110 if(m.ge.nlayer)goto 88 110 111 do while(w(ij,l).gt.(Mtot+masse(ij,m+1))) 111 112 m=m+1 112 113 Mtot = Mtot + masse(ij,m) 113 114 MQtot = MQtot + masse(ij,m)*q(ij,m) 114 if(m.ge. llm)goto 88115 if(m.ge.nlayer)goto 88 115 116 end do 116 117 88 continue 117 if (m.lt. llm) then118 if (m.lt.nlayer) then 118 119 sigw=(w(ij,l)-Mtot)/masse(ij,m+1) 119 120 wq(ij,l)=(MQtot + (w(ij,l)-Mtot)* … … 137 138 end do 138 139 139 do l = 1, llm-1 ! loop different than when w>0140 do l = 1,nlayer-1 ! loop different than when w>0 140 141 do ij=1,ngrid 141 142 if(w(ij,l+1).le.0)then … … 176 177 99 continue 177 178 178 do l=1, llm179 do l=1,nlayer 179 180 do ij=1,ngrid 180 181 … … 191 192 192 193 193 194 return195 194 end -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/watercommon_h.F90
r222 r227 19 19 20 20 real, save :: epsi, RCPD, RCPV, RV, RVTMP2 21 !$OMP THREADPRIVATE(epsi,RCPD,RCPV,RV,RVTMP2) 21 22 22 23 contains -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/writediagfi.F
r222 r227 68 68 69 69 real*4,save :: date 70 !$OMP THREADPRIVATE(date) 70 71 71 72 REAL phis(ip1jmp1) … … 78 79 integer,save :: zitau=0 79 80 character(len=20),save :: firstnom='1234567890' 81 !$OMP THREADPRIVATE(zitau,firstnom) 80 82 81 83 ! Ajouts 82 84 integer, save :: ntime=0 85 !$OMP THREADPRIVATE(ntime) 83 86 integer :: idim,varid 84 87 integer :: nid … … 95 98 character(len=120),save :: nom_def(n_nom_def_max) 96 99 logical,save :: firstcall=.true. 100 !$OMP THREADPRIVATE(firstcall) !diagfi_def,n_nom_def,nom_def read in diagfi.def 97 101 98 102 #ifndef MESOSCALE … … 128 132 firstcall=.false. 129 133 134 !$OMP MASTER 130 135 ! Open diagfi.def definition file if there is one: 131 136 open(99,file="diagfi.def",status='old',form='formatted', … … 151 156 diagfi_def=.false. 152 157 endif 158 !$OMP END MASTER 159 !$OMP BARRIER 153 160 END IF ! of IF (firstcall) 154 161 -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/writediagsoil.F90
r222 r227 46 46 character(len=20),save :: firstname="1234567890" 47 47 integer,save :: zitau=0 48 !$OMP THREADPRIVATE(date,isample,ntime,firstname,zitau) 48 49 49 50 character(len=30) :: filename="diagsoil.nc" -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/writediagspecIR.F
r222 r227 54 54 ! Commons 55 55 #include "dimensions.h" 56 #include "dimphys.h"56 !#include "dimphys.h" 57 57 #include "paramet.h" 58 58 !#include "control.h" … … 91 91 data firstnom /'1234567890'/ 92 92 data zitau /0/ 93 !$OMP THREADPRIVATE(firstnom,zitau,date) 93 94 94 95 ! Ajouts 95 96 integer, save :: ntime=0 97 !$OMP THREADPRIVATE(ntime) 96 98 integer :: idim,varid 97 99 integer :: nid -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/writediagspecVI.F
r222 r227 54 54 ! Commons 55 55 #include "dimensions.h" 56 #include "dimphys.h"56 !#include "dimphys.h" 57 57 #include "paramet.h" 58 58 !#include "control.h" … … 91 91 data firstnom /'1234567890'/ 92 92 data zitau /0/ 93 !$OMP THREADPRIVATE(firstnom,zitau,date) 93 94 94 95 ! Ajouts 95 96 integer, save :: ntime=0 97 !$OMP THREADPRIVATE(ntime) 96 98 integer :: idim,varid 97 99 integer :: nid -
codes/icosagcm/branches/SATURN_DYNAMICO/LMDZ.COMMON/libf/phystd/wstats.F90
r222 r227 7 7 8 8 #include "dimensions.h" 9 #include "dimphys.h"9 !#include "dimphys.h" 10 10 #include "comconst.h" 11 11 #include "statto.h" … … 22 22 character (len=50) :: namebis 23 23 character (len=50), save :: firstvar 24 !$OMP THREADPRIVATE(firstvar) 24 25 integer :: ierr,varid,nbdim,nid 25 26 integer :: meanid,sdid … … 30 31 31 32 integer, save :: step=0 33 !$OMP THREADPRIVATE(firstcall,indx,step) 32 34 33 35 ! Added to work in parallel mode … … 293 295 294 296 include "dimensions.h" 295 include "dimphys.h"297 !include "dimphys.h" 296 298 include "netcdf.inc" 297 299
Note: See TracChangeset
for help on using the changeset viewer.