Index: trunk/LMDZ.MARS/libf/phymars/iostart.F90 =================================================================== --- trunk/LMDZ.MARS/libf/phymars/iostart.F90 (revision 2903) +++ trunk/LMDZ.MARS/libf/phymars/iostart.F90 (revision 2907) @@ -692,9 +692,19 @@ ierr=NF90_REDEF(nid_restart) #ifdef NC_DOUBLE + if(field_name.eq. "tauscaling" .or. field_name.eq. "totcloudfrac" .or. field_name.eq. "wstar" ) then ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& (/idim2,idim7/),nvarid) -#else + else + ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& + (/idim2,idim8,idim7/),nvarid) + endif +#else + if(field_name.eq. "tauscaling" .or. field_name.eq. "totcloudfrac" .or. field_name.eq. "wstar" ) then ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_FLOAT,& (/idim2,idim7/),nvarid) + else + ierr=NF90_DEF_VAR(nid_restart,field_name,NF90_DOUBLE,& + (/idim2,idim8,idim7/),nvarid) + endif #endif if (ierr.ne.NF90_NOERR) then Index: trunk/LMDZ.MARS/libf/phymars/physiq_mod.F =================================================================== --- trunk/LMDZ.MARS/libf/phymars/physiq_mod.F (revision 2903) +++ trunk/LMDZ.MARS/libf/phymars/physiq_mod.F (revision 2907) @@ -815,4 +815,7 @@ dsotop(:,:)=0. dwatercap(:,:)=0 + + call compute_meshgridavg(ngrid,nq,albedo,emis,tsurf,qsurf, + & albedo_meshavg,emis_meshavg,tsurf_meshavg,qsurf_meshavg) ! Dust scenario conversion coefficient from IRabs to VISext @@ -1029,6 +1032,6 @@ clearsky=.false. ! part with clouds for both cases CLFvarying true/false CALL callradite(icount,ngrid,nlayer,nq,zday,zls,pq, - & albedo(:,:,1), - & emis(:,1),mu0,zplev,zplay,pt,tsurf(:,1),fract,dist_sol,igout, + & albedo_meshavg,emis_meshavg, + & mu0,zplev,zplay,pt,tsurf(:,1),fract,dist_sol,igout, & zdtlw,zdtsw,fluxsurf_lw(:,iflat),fluxsurf_dn_sw(:,:,iflat), & fluxsurf_up_sw, @@ -1037,5 +1040,5 @@ & tau,aerosol,dsodust,tauscaling,dust_rad_adjust,IRtoVIScoef, & taucloudtes,rdust,rice,nuice,riceco2,nuiceco2, - & qsurf(:,igcm_co2,1),rstormdust,rtopdust,totstormfract, + & qsurf_meshavg(:,igcm_co2),rstormdust,rtopdust,totstormfract, & clearatm,dsords,dsotop,nohmons,clearsky,totcloudfrac) @@ -1053,5 +1056,5 @@ clearsky=.true. CALL callradite(icount,ngrid,nlayer,nq,zday,zls,pq, - & albedo(:,:,1),emis(:,1),mu0,zplev,zplay,pt, + & albedo_meshavg,emis_meshavg,mu0,zplev,zplay,pt, & tsurf(:,1),fract, & dist_sol,igout,zdtlwclf,zdtswclf, @@ -1061,5 +1064,6 @@ & dsodust,tauscaling,dust_rad_adjust,IRtoVIScoef, & taucloudtesclf,rdust, - & rice,nuice,riceco2, nuiceco2,qsurf(:,igcm_co2,1), + & rice,nuice,riceco2, nuiceco2, + & qsurf_meshavg(:,igcm_co2), & rstormdust,rtopdust,totstormfract, & clearatm,dsords,dsotop, @@ -1313,5 +1317,5 @@ c for radiative transfer & clearatm,icount,zday,zls, - & tsurf,qsurf(:,igcm_co2,iflat),igout, + & tsurf,qsurf_meshavg(:,igcm_co2),igout, & totstormfract,tauscaling, & dust_rad_adjust,IRtoVIScoef, @@ -1377,5 +1381,5 @@ & zzlay,zdtsw,zdtlw, & icount,zday,zls,tsurf(:,iflat), - & qsurf(:,igcm_co2,iflat), + & qsurf_meshavg(:,igcm_co2), & igout,aerosol,tauscaling,dust_rad_adjust, & IRtoVIScoef,totstormfract,clearatm, @@ -1500,6 +1504,6 @@ dwatercap(ig,:)=dwatercap(ig,:)+dwatercap_dif(ig,:) ENDDO -c call compute_meshgridavg(ngrid,nq,albedo,emis,tsurf,zdqsdif, -c & albedo_meshavg,emis_meshavg,tsurf_meshavg,zdqsdif_meshavg_tmp) + call compute_meshgridavg(ngrid,nq,albedo,emis,tsurf,zdqsdif, + & albedo_meshavg,emis_meshavg,tsurf_meshavg,zdqsdif_meshavg_tmp) IF (.not.turb_resolved) THEN DO l=1,nlayer @@ -1539,7 +1543,7 @@ zdqdif(1:ngrid,1,1:nq)=0. zdqdif(1:ngrid,1,igcm_dust_number) = - . -zdqsdif(1:ngrid,igcm_dust_number,iflat) + . -zdqsdif_meshavg_tmp(1:ngrid,igcm_dust_number) zdqdif(1:ngrid,1,igcm_dust_mass) = - . -zdqsdif(1:ngrid,igcm_dust_mass,iflat) + . -zdqsdif_meshavg_tmp(1:ngrid,igcm_dust_mass) zdqdif(1:ngrid,2:nlayer,1:nq) = 0. DO iq=1, nq @@ -2605,8 +2609,7 @@ call physdem1("restartfi.nc",nsoilmx,ngrid,nlayer,nq, . ptimestep,ztime_fin, - . tsurf(:,iflat),tsoil(:,:,iflat),albedo(:,:,iflat), - . emis(:,iflat), - . q2,qsurf(:,:,iflat),tauscaling,totcloudfrac,wstar, - . watercap(:,iflat)) + . tsurf,tsoil,albedo, + . emis,q2,qsurf,tauscaling,totcloudfrac,wstar, + . watercap) ENDIF ! of IF (write_restart) @@ -3155,4 +3158,9 @@ call WRITEDIAGFI(ngrid,"emis","Surface emissivity","w.m-1",2, & emis(:,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,"emis_slope"//str2, + & "Surface emissivity","w.m-1",2, emis(:,islope)) + ENDDO c call WRITEDIAGFI(ngrid,"pplay","Pressure","Pa",3,zplay) c call WRITEDIAGFI(ngrid,"pplev","Pressure","Pa",3,zplev) @@ -3165,10 +3173,26 @@ call WRITEDIAGFI(ngrid,"tsurf","Surface temperature","K",2, & tsurf(:,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,"tsurf_slope"//str2, + & "Surface temperature","K",2, + & tsurf(:,islope)) + ENDDO call WRITEDIAGFI(ngrid,"ps","surface pressure","Pa",2,ps) call WRITEDIAGFI(ngrid,"co2ice","co2 ice thickness" & ,"kg.m-2",2,qsurf(:,igcm_co2,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,"co2ice_slope"//str2,"co2 ice thickness" + & ,"kg.m-2",2,qsurf(:,igcm_co2,islope)) + ENDDO call WRITEDIAGFI(ngrid,"watercap","Water ice thickness" & ,"kg.m-2",2,watercap(:,iflat)) - + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,"watercap_slope"//str2, + & "Water ice thickness" + & ,"kg.m-2",2,watercap(:,islope)) + ENDDO call WRITEDIAGFI(ngrid,"temp_layer1","temperature in layer 1", & "K",2,zt(1,1)) @@ -3177,6 +3201,18 @@ call WRITEDIAGFI(ngrid,"fluxsurf_lw","fluxsurf_lw","W.m-2",2, & fluxsurf_lw(:,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,"fluxsurf_lw_slope"//str2, + & "fluxsurf_lw","W.m-2",2, + & fluxsurf_lw(:,islope)) + ENDDO call WRITEDIAGFI(ngrid,"fluxsurf_dn_sw","fluxsurf_dn_sw", & "W.m-2",2,fluxsurf_dn_sw_tot(:,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,"fluxsurf_dn_sw_slope"//str2, + & "fluxsurf_dn_sw", + & "W.m-2",2,fluxsurf_dn_sw_tot(:,islope)) + ENDDO call WRITEDIAGFI(ngrid,"fluxtop_lw","fluxtop_lw","W.m-2",2, & fluxtop_lw) @@ -3306,4 +3342,10 @@ & 'kg.m-2',2, & qsurf(1:ngrid,igcm_h2o_ice,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'qsurf02_slope'//str2, + & 'surface tracer','kg.m-2',2, + & qsurf(1:ngrid,igcm_h2o_ice,islope)) + ENDDO #endif call WRITEDIAGFI(ngrid,'mtot', @@ -3400,9 +3442,21 @@ call WRITEDIAGFI(ngrid,'h2o_ice_s', & 'surface h2o_ice', - & 'kg.m-2',2,qsurf(1,igcm_h2o_ice,iflat)) + & 'kg.m-2',2,qsurf(1,igcm_h2o_ice,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'h2o_ice_s_slope'//str2, + & 'surface h2o_ice', + & 'kg.m-2',2,qsurf(1,igcm_h2o_ice,islope)) + ENDDO if (hdo) then call WRITEDIAGFI(ngrid,'hdo_ice_s', & 'surface hdo_ice', - & 'kg.m-2',2,qsurf(1,igcm_hdo_ice,iflat)) + & 'kg.m-2',2,qsurf(1,igcm_hdo_ice,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'hdo_ice_s_slope'//str2, + & 'surface hdo_ice', + & 'kg.m-2',2,qsurf(1,igcm_hdo_ice,islope)) + ENDDO do ig=1,ngrid @@ -3423,4 +3477,10 @@ & 'albedo', & '',2,albedo(1,1,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + CALL WRITEDIAGFI(ngrid,'albedo_slope'//str2, + & 'albedo', + & '',2,albedo(1,1,islope)) + ENDDO if (tifeedback) then call WRITEDIAGSOIL(ngrid,"soiltemp", @@ -3429,5 +3489,11 @@ call WRITEDIAGSOIL(ngrid,'soilti', & 'Soil Thermal Inertia', - & 'J.s-1/2.m-2.K-1',3,inertiesoil(:,:,iflat)) + & 'J.s-1/2.m-2.K-1',3,inertiesoil(:,:,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGSOIL(ngrid,'soilti_slope'//str2, + & 'Soil Thermal Inertia', + & 'J.s-1/2.m-2.K-1',3,inertiesoil(:,:,islope)) + ENDDO endif !A. Pottier @@ -3528,4 +3594,9 @@ call WRITEDIAGFI(ngrid,'qsurf'//str2,'qsurf', & 'kg.m-2',2,qsurf(1,iq,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'qsurf_slope'//str2,'qsurf', + & 'kg.m-2',2,qsurf(1,iq,islope)) + ENDDO end do endif ! (doubleq) @@ -3555,7 +3626,19 @@ & 'stormdust injection', & 'kg.m-2',2,qsurf(:,igcm_stormdust_mass,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'qsurf_slope'//str2, + & 'stormdust injection', + & 'kg.m-2',2,qsurf(:,igcm_stormdust_mass,islope)) + ENDDO call WRITEDIAGFI(ngrid,'pdqsurf', & 'tendancy stormdust mass at surface', - & 'kg.m-2',2,dqsurf(:,igcm_stormdust_mass,iflat)) + & 'kg.m-2',2,dqsurf(:,igcm_stormdust_mass,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'pdqsurf_slope'//str2, + & 'tendancy stormdust mass at surface', + & 'kg.m-2',2,dqsurf(:,igcm_stormdust_mass,islope)) + ENDDO call WRITEDIAGFI(ngrid,'wspeed','vertical speed stormdust', & 'm/s',3,wspeed(:,1:nlayer)) @@ -3620,6 +3703,18 @@ call WRITEDIAGFI(ngrid,'surfccnq','Surf nuclei mass mr', & 'kg.m-2',2,qsurf(1,igcm_ccn_mass,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'surfccnq_slope'//str2, + & 'Surf nuclei mass mr', + & 'kg.m-2',2,qsurf(1,igcm_ccn_mass,islope)) + ENDDO call WRITEDIAGFI(ngrid,'surfccnN','Surf nuclei number', & 'kg.m-2',2,qsurf(1,igcm_ccn_number,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'surfccnN_slope'//str2, + & 'Surf nuclei number', + & 'kg.m-2',2,qsurf(1,igcm_ccn_number,islope)) + ENDDO endif ! (scavenging) @@ -3772,4 +3867,10 @@ call writediagsoil(ngrid,"soiltemp","Soil temperature","K", & 3,tsoil(:,:,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call writediagsoil(ngrid,"soiltemp_slope"//str2, + & "Soil temperature","K", + & 3,tsoil(:,:,islope)) + ENDDO ! Write surface temperature ! call writediagsoil(ngrid,"tsurf","Surface temperature","K", @@ -3802,4 +3903,10 @@ call writediagsoil(ngrid,"soiltemp","Soil temperature","K", & 3,tsoil(:,:,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call writediagsoil(ngrid,"soiltemp_slope"//str2, + & "Soil temperature","K", + & 3,tsoil(:,:,islope)) + ENDDO ! THERMALS STUFF 1D @@ -3826,4 +3933,10 @@ call WRITEDIAGFI(ngrid,"tsurf","Surface temperature","K",0, & tsurf(:,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,"tsurf_slope"//str2, + & "Surface temperature","K",0, + & tsurf(:,islope)) + ENDDO call WRITEDIAGFI(ngrid,"u","u wind","m/s",1,zu) call WRITEDIAGFI(ngrid,"v","v wind","m/s",1,zv) @@ -3840,4 +3953,10 @@ call WRITEDIAGFI(ngrid,"co2ice","co2 ice thickness" & ,"kg.m-2",0,qsurf(:,igcm_co2,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,"co2ice_slope"//str2, + & "co2 ice thickness" + & ,"kg.m-2",0,qsurf(:,igcm_co2,islope)) + ENDDO if (igcm_co2.ne.0) then @@ -3878,6 +3997,18 @@ call WRITEDIAGFI(ngrid,'dqsdifdustq','diffusion', & 'kg.m-2.s-1',0,zdqsdif(1,igcm_dust_mass,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'dqsdifdustq_slope'//str2, + & 'diffusion', + & 'kg.m-2.s-1',0,zdqsdif(1,igcm_dust_mass,islope)) + ENDDO call WRITEDIAGFI(ngrid,'dqsdifrdsq','diffusion', & 'kg.m-2.s-1',0,zdqsdif(1,igcm_stormdust_mass,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'dqsdifrdsq_slope'//str2, + & 'diffusion', + & 'kg.m-2.s-1',0,zdqsdif(1,igcm_stormdust_mass,islope)) + ENDDO call WRITEDIAGFI(ngrid,'mstormdtot', & 'total mass of stormdust only', @@ -3903,7 +4034,19 @@ & 'stormdust at surface', & 'kg.m-2',0,qsurf(:,igcm_stormdust_mass,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'rdsqsurf_slope'//str2, + & 'stormdust at surface', + & 'kg.m-2',0,qsurf(:,igcm_stormdust_mass,islope)) + ENDDO call WRITEDIAGFI(ngrid,'qsurf', & 'dust mass at surface', & 'kg.m-2',0,qsurf(:,igcm_dust_mass,iflat)) + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + call WRITEDIAGFI(ngrid,'qsurf_slope'//str2, + & 'dust mass at surface', + & 'kg.m-2',0,qsurf(:,igcm_dust_mass,islope)) + ENDDO call WRITEDIAGFI(ngrid,'wspeed','vertical speed stormdust', & 'm/s',1,wspeed) @@ -3996,5 +4139,10 @@ & 'albedo', & '',2,albedo(1,1,iflat)) - + do islope=1,nslope + write(str2(1:2),'(i2.2)') islope + CALL WRITEDIAGFI(ngrid,'albedo_slope'//str2, + & 'albedo', + & '',2,albedo(1,1,islope)) + ENDDO IF (hdo) THEN CALL WRITEDIAGFI(ngrid,'mtotD',