Index: LMDZ5/branches/testing/libf/dyn3d/calfis.F =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/calfis.F (revision 1669) +++ LMDZ5/branches/testing/libf/dyn3d/calfis.F (revision 1707) @@ -507,5 +507,5 @@ . debut_split, !! firstcall . lafin_split, !! lastcall - . float(day_ini), !! pday <-- day_ini (dans temps.h) + . jD_cur, !! pday. see leapfrog . jH_cur_split, !! ptime "fraction of day" . zdt_split, !! ptimestep Index: LMDZ5/branches/testing/libf/dyn3d/comconst.h =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/comconst.h (revision 1669) +++ LMDZ5/branches/testing/libf/dyn3d/comconst.h (revision 1707) @@ -21,6 +21,7 @@ REAL dtdiss ! (s) time step for the dissipation REAL rad ! (m) radius of the planet - REAL r ! Gas constant R=8.31 J.K-1.mol-1 - REAL cpp ! Cp + REAL r ! Reduced Gas constant r=R/mu + ! with R=8.31.. J.K-1.mol-1, mu: mol mass of atmosphere (kg/mol) + REAL cpp ! Specific heat Cp (J.kg-1.K-1) REAL kappa ! kappa=R/Cp REAL cotot Index: LMDZ5/branches/testing/libf/dyn3d/comdissnew.h =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/comdissnew.h (revision 1669) +++ LMDZ5/branches/testing/libf/dyn3d/comdissnew.h (revision 1707) @@ -12,8 +12,14 @@ COMMON/comdissnew/ lstardis,nitergdiv,nitergrot,niterh,tetagdiv, & - & tetagrot,tetatemp,coefdis + & tetagrot,tetatemp,coefdis, vert_prof_dissip LOGICAL lstardis INTEGER nitergdiv, nitergrot, niterh + + integer vert_prof_dissip ! vertical profile of horizontal dissipation +! Allowed values: +! 0: rational fraction, function of pressure +! 1: tanh of altitude + REAL tetagdiv, tetagrot, tetatemp, coefdis Index: LMDZ5/branches/testing/libf/dyn3d/conf_gcm.F =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/conf_gcm.F (revision 1669) +++ LMDZ5/branches/testing/libf/dyn3d/conf_gcm.F (revision 1707) @@ -14,4 +14,6 @@ #endif USE infotrac, ONLY : type_trac + use assert_m, only: assert + IMPLICIT NONE c----------------------------------------------------------------------- @@ -93,5 +95,6 @@ CALL getin('lunout', lunout) IF (lunout /= 5 .and. lunout /= 6) THEN - OPEN(lunout,FILE='lmdz.out') + OPEN(UNIT=lunout,FILE='lmdz.out',ACTION='write', + & STATUS='unknown',FORM='formatted') ENDIF @@ -173,8 +176,4 @@ !Config Key = nsplit_phys -!Config Desc = nombre de pas par jour -!Config Def = 1 -!Config Help = nombre de pas par jour (multiple de iperiod) ( -!Config ici pour dt = 1 min ) nsplit_phys = 1 CALL getin('nsplit_phys',nsplit_phys) @@ -625,5 +624,4 @@ CALL getin('ok_dyn_ave',ok_dyn_ave) - write(lunout,*)' #########################################' write(lunout,*)' Configuration des parametres du gcm: ' @@ -635,4 +633,5 @@ write(lunout,*)' day_step = ', day_step write(lunout,*)' iperiod = ', iperiod + write(lunout,*)' nsplit_phys = ', nsplit_phys write(lunout,*)' iconser = ', iconser write(lunout,*)' iecri = ', iecri @@ -805,8 +804,8 @@ !Config Desc = sortie des transports zonaux dans la dynamique !Config Def = n -!Config Help = +!Config Help = Permet de mettre en route le calcul des transports !Config - ok_dynzon = .FALSE. - CALL getin('ok_dynzon',ok_dynzon) + ok_dynzon = .FALSE. + CALL getin('ok_dynzon',ok_dynzon) !Config Key = ok_dyn_ins @@ -838,5 +837,5 @@ write(lunout,*)'STOP !!!' write(lunout,*)'use_filtre_fft n est pas implemente dans dyn3d' - STOP + STOP 1 ENDIF @@ -848,4 +847,9 @@ ok_strato=.FALSE. CALL getin('ok_strato',ok_strato) + + vert_prof_dissip = merge(1, 0, ok_strato .and. llm==39) + CALL getin('vert_prof_dissip', vert_prof_dissip) + call assert(vert_prof_dissip == 0 .or. vert_prof_dissip == 1, + $ "bad value for vert_prof_dissip") !Config Key = ok_gradsfile Index: LMDZ5/branches/testing/libf/dyn3d/etat0_netcdf.F90 =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/etat0_netcdf.F90 (revision 1669) +++ (revision ) @@ -1,530 +1,0 @@ -! -! $Id$ -! -!------------------------------------------------------------------------------- -! -SUBROUTINE etat0_netcdf(ib, masque, phis, letat0) -! -!------------------------------------------------------------------------------- -! Purpose: Creates initial states -!------------------------------------------------------------------------------- -! Note: This routine is designed to work for Earth -!------------------------------------------------------------------------------- - USE control_mod -#ifdef CPP_EARTH - USE startvar - USE ioipsl - USE dimphy - USE infotrac - USE fonte_neige_mod - USE pbl_surface_mod - USE phys_state_var_mod - USE filtreg_mod - USE regr_lat_time_climoz_m, ONLY: regr_lat_time_climoz - USE conf_phys_m, ONLY: conf_phys -! For parameterization of ozone chemistry: - use regr_lat_time_coefoz_m, only: regr_lat_time_coefoz - use press_coefoz_m, only: press_coefoz - use regr_pr_o3_m, only: regr_pr_o3 - USE netcdf, ONLY : NF90_OPEN, NF90_NOWRITE, NF90_CLOSE, NF90_NOERR -#endif - IMPLICIT NONE -!------------------------------------------------------------------------------- -! Arguments: -#include "dimensions.h" -#include "paramet.h" -#include "iniprint.h" - LOGICAL, INTENT(IN) :: ib ! barycentric interpolat. - REAL, DIMENSION(iip1,jjp1), INTENT(INOUT) :: masque ! land mask - REAL, DIMENSION(iip1,jjp1), INTENT(OUT) :: phis ! geopotentiel au sol - LOGICAL, INTENT(IN) :: letat0 ! F: masque only required -#ifndef CPP_EARTH - WRITE(lunout,*)'limit_netcdf: Earth-specific routine, needs Earth physics' -#else -!------------------------------------------------------------------------------- -! Local variables: -#include "comgeom2.h" -#include "comvert.h" -#include "comconst.h" -#include "indicesol.h" -#include "dimsoil.h" -#include "temps.h" - REAL, DIMENSION(klon) :: tsol, qsol - REAL, DIMENSION(klon) :: sn, rugmer, run_off_lic_0 - REAL, DIMENSION(iip1,jjp1) :: orog, rugo, psol - REAL, DIMENSION(iip1,jjp1,llm+1) :: p3d - REAL, DIMENSION(iip1,jjp1,llm) :: uvent, t3d, tpot, qsat, qd - REAL, DIMENSION(iip1,jjm ,llm) :: vvent - REAL, DIMENSION(:,:,:,:), ALLOCATABLE :: q3d - REAL, DIMENSION(klon,nbsrf) :: qsolsrf, snsrf, evap - REAL, DIMENSION(klon,nbsrf) :: frugs, agesno - REAL, DIMENSION(klon,nsoilmx,nbsrf) :: tsoil - -!--- Local variables for sea-ice reading: - INTEGER :: iml_lic, jml_lic, llm_tmp - INTEGER :: ttm_tmp, iret, fid - INTEGER, DIMENSION(1) :: itaul - REAL, DIMENSION(1) :: lev - REAL :: date - REAL, DIMENSION(:,:), ALLOCATABLE :: lon_lic, lat_lic, fraclic - REAL, DIMENSION(:), ALLOCATABLE :: dlon_lic, dlat_lic - REAL, DIMENSION(iip1,jjp1) :: flic_tmp - -!--- Misc - CHARACTER(LEN=80) :: x, fmt - INTEGER :: i, j, l, ji - REAL, DIMENSION(iip1,jjp1,llm) :: alpha, beta, pk, pls, y - REAL, DIMENSION(ip1jmp1) :: pks - -#include "comdissnew.h" -#include "serre.h" -#include "clesphys.h" - - REAL, DIMENSION(iip1,jjp1,llm) :: masse - INTEGER :: itau, iday - REAL :: xpn, xps, time, phystep - REAL, DIMENSION(iim) :: xppn, xpps - REAL, DIMENSION(ip1jmp1,llm) :: pbaru, phi, w - REAL, DIMENSION(ip1jm ,llm) :: pbarv - REAL, DIMENSION(klon) :: fder - -!--- Local variables for ocean mask reading: - INTEGER :: nid_o2a, iml_omask, jml_omask - LOGICAL :: couple=.FALSE. - REAL, DIMENSION(:,:), ALLOCATABLE :: lon_omask, lat_omask, ocemask, ocetmp - REAL, DIMENSION(:), ALLOCATABLE :: dlon_omask,dlat_omask - REAL, DIMENSION(klon) :: ocemask_fi - INTEGER, DIMENSION(klon-2) :: isst - REAL, DIMENSION(iim,jjp1) :: zx_tmp_2d - REAL :: dummy - LOGICAL :: ok_newmicro, ok_journe, ok_mensuel, ok_instan, ok_hf - LOGICAL :: ok_LES, ok_ade, ok_aie, aerosol_couple, new_aod, callstats - INTEGER :: iflag_radia, flag_aerosol - REAL :: bl95_b0, bl95_b1, fact_cldcon, facttemps, ratqsbas, ratqshaut - REAL :: tau_ratqs - INTEGER :: iflag_cldcon, iflag_ratqs, iflag_coupl, iflag_clos, iflag_wake - INTEGER :: iflag_thermals, nsplit_thermals - INTEGER :: iflag_thermals_ed, iflag_thermals_optflux - REAL :: tau_thermals, solarlong0, seuil_inversion - INTEGER :: read_climoz ! read ozone climatology -! Allowed values are 0, 1 and 2 -! 0: do not read an ozone climatology -! 1: read a single ozone climatology that will be used day and night -! 2: read two ozone climatologies, the average day and night -! climatology and the daylight climatology -!------------------------------------------------------------------------------- - REAL :: alp_offset - logical found - -!--- Constants - pi = 4. * ATAN(1.) - rad = 6371229. - daysec = 86400. - omeg = 2.*pi/daysec - g = 9.8 - kappa = 0.2857143 - cpp = 1004.70885 - preff = 101325. - pa = 50000. - jmp1 = jjm + 1 - -!--- CONSTRUCT A GRID - CALL conf_phys( ok_journe, ok_mensuel, ok_instan, ok_hf, ok_LES, & - callstats, & - solarlong0,seuil_inversion, & - fact_cldcon, facttemps,ok_newmicro,iflag_radia, & - iflag_cldcon, & - iflag_ratqs,ratqsbas,ratqshaut,tau_ratqs, & - ok_ade, ok_aie, aerosol_couple, & - flag_aerosol, new_aod, & - bl95_b0, bl95_b1, & - read_climoz, & - alp_offset) - -! co2_ppm0 : initial value of atmospheric CO2 from .def file (co2_ppm value) - co2_ppm0 = co2_ppm - - dtvr = daysec/FLOAT(day_step) - WRITE(lunout,*)'dtvr',dtvr - - CALL iniconst() - CALL inigeom() - -!--- Initializations for tracers - CALL infotrac_init - ALLOCATE(q3d(iip1,jjp1,llm,nqtot)) - - CALL inifilr() - CALL phys_state_var_init(read_climoz) - - rlat(1) = ASIN(1.0) - DO j=2,jjm; rlat((j-2)*iim+2:(j-1)*iim+1)=rlatu(j); END DO - rlat(klon) = - ASIN(1.0) - rlat(:)=rlat(:)*(180.0/pi) - - rlon(1) = 0.0 - DO j=2,jjm; rlon((j-2)*iim+2:(j-1)*iim+1)=rlonv(1:iim); END DO - rlon(klon) = 0.0 - rlon(:)=rlon(:)*(180.0/pi) - -! For a coupled simulation, the ocean mask from ocean model is used to compute -! the weights an to insure ocean fractions are the same for atmosphere and ocean -! Otherwise, mask is created using Relief file. - - WRITE(lunout,*)'Essai de lecture masque ocean' - iret = NF90_OPEN("o2a.nc", NF90_NOWRITE, nid_o2a) - IF(iret/=NF90_NOERR) THEN - WRITE(lunout,*)'ATTENTION!! pas de fichier o2a.nc trouve' - WRITE(lunout,*)'Run force' - x='masque' - masque(:,:)=0.0 - CALL startget_phys2d(x, iip1, jjp1, rlonv, rlatu, masque, 0.0, jjm, & - & rlonu, rlatv, ib) - WRITE(lunout,*)'MASQUE construit : Masque' - WRITE(lunout,'(97I1)') nINT(masque) - CALL gr_dyn_fi(1, iip1, jjp1, klon, masque, zmasq) - WHERE( zmasq(:)=1) CALL regr_lat_time_climoz(read_climoz) - - x = 'tsol'; tsol(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,tsol,0.0,jjm,rlonu,rlatv,ib) - - x = 'qsol'; qsol(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,qsol,0.0,jjm,rlonu,rlatv,ib) - - x = 'snow'; sn(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,sn,0.0,jjm,rlonu,rlatv,ib) - - x = 'rads'; radsol(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,radsol,0.0,jjm,rlonu,rlatv,ib) - - x = 'rugmer'; rugmer(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,rugmer,0.0,jjm,rlonu,rlatv,ib) - - x = 'zmea'; zmea(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,zmea,0.0,jjm,rlonu,rlatv,ib) - - x = 'zstd'; zstd(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,zstd,0.0,jjm,rlonu,rlatv,ib) - - x = 'zsig'; zsig(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,zsig,0.0,jjm,rlonu,rlatv,ib) - - x = 'zgam'; zgam(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,zgam,0.0,jjm,rlonu,rlatv,ib) - - x = 'zthe'; zthe(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,zthe,0.0,jjm,rlonu,rlatv,ib) - - x = 'zpic'; zpic(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,zpic,0.0,jjm,rlonu,rlatv,ib) - - x = 'zval'; zval(:) = 0.0 - CALL startget_phys1d(x,iip1,jjp1,rlonv,rlatu,klon,zval,0.0,jjm,rlonu,rlatv,ib) - -! WRITE(lunout,'(48I3)') 'TSOL :', INT(tsol(2:klon)-273) - -! Soil ice file reading for soil fraction and soil ice fraction -!******************************************************************************* - CALL flininfo("landiceref.nc", iml_lic, jml_lic, llm_tmp, ttm_tmp, fid) - ALLOCATE( lat_lic(iml_lic,jml_lic),lon_lic(iml_lic, jml_lic)) - ALLOCATE(dlat_lic(jml_lic), dlon_lic(iml_lic)) - ALLOCATE( fraclic(iml_lic,jml_lic)) - CALL flinopen("landiceref.nc", .FALSE., iml_lic, jml_lic, llm_tmp, & - & lon_lic, lat_lic, lev, ttm_tmp, itaul, date, dt, fid) - CALL flinget(fid, 'landice', iml_lic, jml_lic, llm_tmp, ttm_tmp, 1,1, fraclic) - CALL flinclo(fid) - -! Interpolation on model T-grid -!******************************************************************************* - WRITE(lunout,*)'dimensions de landice iml_lic, jml_lic : ',iml_lic,jml_lic -! conversion if coordinates are in degrees - IF(MAXVAL(lon_lic)>pi) lon_lic=lon_lic*pi/180. - IF(MAXVAL(lat_lic)>pi) lat_lic=lat_lic*pi/180. - dlon_lic(:)=lon_lic(:,1) - dlat_lic(:)=lat_lic(1,:) - CALL grille_m( iml_lic, jml_lic, dlon_lic, dlat_lic, fraclic, iim,jjp1, & - & rlonv, rlatu, flic_tmp(1:iim,:) ) - flic_tmp(iip1,:)=flic_tmp(1,:) - -!--- To the physical grid - CALL gr_dyn_fi(1, iip1, jjp1, klon, flic_tmp, pctsrf(:,is_lic)) - -!--- Adequation with soil/sea mask - WHERE(pctsrf(:,is_lic)EPSFRA) THEN - IF(pctsrf(ji,is_lic)>=zmasq(ji)) THEN - pctsrf(ji,is_lic)=zmasq(ji) - pctsrf(ji,is_ter)=0. - ELSE - pctsrf(ji,is_ter)=zmasq(ji)-pctsrf(ji,is_lic) - IF(pctsrf(ji,is_ter)0.) isst=1 - -! It is checked that the sub-surfaces sum is equal to 1 -!******************************************************************************* - ji=COUNT((ABS(SUM(pctsrf(:,:),dim=2))-1.0)>EPSFRA) - IF(ji/=0) WRITE(lunout,*) 'pb repartition sous maille pour ',ji,' points' - CALL gr_fi_ecrit(1, klon, iim, jjp1, zmasq, zx_tmp_2d) -! WRITE(fmt,"(i3,')')")iim; fmt='(i'//ADJUSTL(fmt) -! WRITE(lunout,*)'zmasq = ' -! WRITE(lunout,TRIM(fmt))NINT(zx_tmp_2d) - CALL gr_fi_dyn(1, klon, iip1, jjp1, zmasq, masque) - WRITE(fmt,"(i4,'i1)')")iip1 ; fmt='('//ADJUSTL(fmt) - WRITE(lunout,*) 'MASQUE construit : Masque' - WRITE(lunout,TRIM(fmt))NINT(masque(:,:)) - -! Intermediate computation -!******************************************************************************* - CALL massdair(p3d,masse) - WRITE(lunout,*)' ALPHAX ',alphax - DO l=1,llm - xppn(:)=aire(1:iim,1 )*masse(1:iim,1 ,l) - xpps(:)=aire(1:iim,jjp1)*masse(1:iim,jjp1,l) - xpn=SUM(xppn)/apoln - xps=SUM(xpps)/apols - masse(:,1 ,l)=xpn - masse(:,jjp1,l)=xps - END DO - q3d(iip1,:,:,:)=q3d(1,:,:,:) - phis(iip1,:) = phis(1,:) - - IF(.NOT.letat0) RETURN - -! Writing -!******************************************************************************* - CALL inidissip(lstardis,nitergdiv,nitergrot,niterh,tetagdiv,tetagrot,tetatemp) - WRITE(lunout,*)'sortie inidissip' - itau=0 - itau_dyn=0 - itau_phy=0 - iday=dayref+itau/day_step - time=FLOAT(itau-(iday-dayref)*day_step)/day_step - IF(time>1.) THEN - time=time-1 - iday=iday+1 - END IF - day_ref=dayref - annee_ref=anneeref - - CALL geopot( ip1jmp1, tpot, pk, pks, phis, phi ) - WRITE(lunout,*)'sortie geopot' - - CALL caldyn0( itau, uvent, vvent, tpot, psol, masse, pk, phis, & - phi, w, pbaru, pbarv, time+iday-dayref) - WRITE(lunout,*)'sortie caldyn0' - CALL dynredem0( "start.nc", dayref, phis) - WRITE(lunout,*)'sortie dynredem0' - CALL dynredem1( "start.nc", 0.0, vvent, uvent, tpot, q3d, masse, psol) - WRITE(lunout,*)'sortie dynredem1' - -! Physical initial state writting -!******************************************************************************* - WRITE(lunout,*)'phystep ',dtvr,iphysiq,nbapp_rad - phystep = dtvr * FLOAT(iphysiq) - radpas = NINT (86400./phystep/ FLOAT(nbapp_rad) ) - WRITE(lunout,*)'phystep =', phystep, radpas - -! Init: tsol, qsol, sn, evap, tsoil, rain_fall, snow_fall, solsw, sollw, frugs -!******************************************************************************* - DO i=1,nbsrf; ftsol(:,i) = tsol; END DO - DO i=1,nbsrf; snsrf(:,i) = sn; END DO - falb1(:,is_ter) = 0.08; falb1(:,is_lic) = 0.6 - falb1(:,is_oce) = 0.5; falb1(:,is_sic) = 0.6 - falb2 = falb1 - evap(:,:) = 0. - DO i=1,nbsrf; qsolsrf(:,i)=150.; END DO - DO i=1,nbsrf; DO j=1,nsoilmx; tsoil(:,j,i) = tsol; END DO; END DO - rain_fall = 0.; snow_fall = 0. - solsw = 165.; sollw = -53. - t_ancien = 273.15 - q_ancien = 0. - agesno = 0. - frugs(:,is_oce) = rugmer(:) - frugs(:,is_ter) = MAX(1.0e-05,zstd(:)*zsig(:)/2.0) - frugs(:,is_lic) = MAX(1.0e-05,zstd(:)*zsig(:)/2.0) - frugs(:,is_sic) = 0.001 - fder = 0.0 - clwcon = 0.0 - rnebcon = 0.0 - ratqs = 0.0 - run_off_lic_0 = 0.0 - rugoro = 0.0 - -! Before phyredem calling, surface modules and values to be saved in startphy.nc -! are initialized -!******************************************************************************* - dummy = 1.0 - pbl_tke(:,:,:) = 1.e-8 - zmax0(:) = 40. - f0(:) = 1.e-5 - ema_work1(:,:) = 0. - ema_work2(:,:) = 0. - wake_deltat(:,:) = 0. - wake_deltaq(:,:) = 0. - wake_s(:) = 0. - wake_cstar(:) = 0. - wake_fip(:) = 0. - wake_pe = 0. - fm_therm = 0. - entr_therm = 0. - detr_therm = 0. - - CALL fonte_neige_init(run_off_lic_0) - CALL pbl_surface_init( qsol, fder, snsrf, qsolsrf, evap, frugs, agesno, tsoil ) - CALL phyredem( "startphy.nc" ) - -! WRITE(lunout,*)'CCCCCCCCCCCCCCCCCC REACTIVER SORTIE VISU DANS ETAT0' -! WRITE(lunout,*)'entree histclo' - CALL histclo() - -#endif -!#endif of #ifdef CPP_EARTH - RETURN - -END SUBROUTINE etat0_netcdf -! -!------------------------------------------------------------------------------- Index: LMDZ5/branches/testing/libf/dyn3d/fxhyp.F =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/fxhyp.F (revision 1669) +++ LMDZ5/branches/testing/libf/dyn3d/fxhyp.F (revision 1707) @@ -68,4 +68,25 @@ xzoom = xzoomdeg * pi/180. c + if (iim==1) then + + rlonm025(1)=-pi/2. + rlonv(1)=0. + rlonu(1)=pi + rlonp025(1)=pi/2. + rlonm025(2)=rlonm025(1)+depi + rlonv(2)=rlonv(1)+depi + rlonu(2)=rlonu(1)+depi + rlonp025(2)=rlonp025(1)+depi + + xprimm025(:)=1. + xprimv(:)=1. + xprimu(:)=1. + xprimp025(:)=1. + champmin=depi + champmax=depi + return + + endif + decalx = .75 IF( grossism.EQ.1..AND.scal180 ) THEN @@ -286,8 +307,10 @@ + IF(ik.EQ.1.and.grossism.EQ.1.) THEN xvrai(1) = xvrai(iip1)-depi xxprim(1) = xxprim(iip1) ENDIF + DO i = 1 , iim xlon(i) = xvrai(i) Index: LMDZ5/branches/testing/libf/dyn3d/gcm.F =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/gcm.F (revision 1669) +++ LMDZ5/branches/testing/libf/dyn3d/gcm.F (revision 1707) @@ -405,5 +405,5 @@ CALL inidissip( lstardis, nitergdiv, nitergrot, niterh , - * tetagdiv, tetagrot , tetatemp ) + * tetagdiv, tetagrot , tetatemp, vert_prof_dissip) c----------------------------------------------------------------------- @@ -433,6 +433,7 @@ ! Physics: #ifdef CPP_PHYS - CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys/nsplit_phys , - , latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp ) + CALL iniphysiq(ngridmx,llm,daysec,day_ini,dtphys/nsplit_phys, + & latfi,lonfi,airefi,zcufi,zcvfi,rad,g,r,cpp, + & iflag_phys) #endif call_iniphys=.false. @@ -457,12 +458,4 @@ 301 FORMAT('1'/,15x,'run du ', i2,'/',i2,'/',i4) 302 FORMAT('1'/,15x,' au ', i2,'/',i2,'/',i4) -#endif - -#ifdef CPP_PHYS -! Create start file (startphy.nc) and boundary conditions (limit.nc) -! for the Earth verstion - if (iflag_phys>=100) then - call iniaqua(ngridmx,latfi,lonfi,iflag_phys) - endif #endif Index: LMDZ5/branches/testing/libf/dyn3d/groupe.F =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/groupe.F (revision 1669) +++ LMDZ5/branches/testing/libf/dyn3d/groupe.F (revision 1707) @@ -38,13 +38,21 @@ integer i,j,l - logical firstcall - save firstcall + logical firstcall,groupe_ok + save firstcall,groupe_ok data firstcall/.true./ + data groupe_ok/.true./ + + if (iim==1) then + groupe_ok=.false. + endif if (firstcall) then - if(mod(iim,2**ngroup).ne.0) stop'probleme du nombre ede point' + if (groupe_ok) then + if(mod(iim,2**ngroup).ne.0) stop'probleme du nombre de point' + endif firstcall=.false. endif + c Champs 1D @@ -52,14 +60,12 @@ call convflu(pbaru,pbarv,llm,zconvm) -c call scopy(ijp1llm,zconvm,1,zconvmm,1) call scopy(ijmllm,pbarv,1,pbarvm,1) -c + if (groupe_ok) then call groupeun(jjp1,llm,zconvmm) call groupeun(jjm,llm,pbarvm) c Champs 3D - do l=1,llm do j=2,jjm @@ -74,4 +80,9 @@ enddo enddo + + else + pbarum(:,:,:)=pbaru(:,:,:) + pbarvm(:,:,:)=pbarv(:,:,:) + endif c integration de la convergence de masse de haut en bas ...... Index: LMDZ5/branches/testing/libf/dyn3d/inidissip.F90 =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/inidissip.F90 (revision 1669) +++ LMDZ5/branches/testing/libf/dyn3d/inidissip.F90 (revision 1707) @@ -3,5 +3,5 @@ ! SUBROUTINE inidissip ( lstardis,nitergdiv,nitergrot,niterh , & - tetagdiv,tetagrot,tetatemp ) + tetagdiv,tetagrot,tetatemp, vert_prof_dissip) !======================================================================= ! initialisation de la dissipation horizontale @@ -25,4 +25,10 @@ INTEGER,INTENT(in) :: nitergdiv,nitergrot,niterh REAL,INTENT(in) :: tetagdiv,tetagrot,tetatemp + + integer, INTENT(in):: vert_prof_dissip + ! Vertical profile of horizontal dissipation + ! Allowed values: + ! 0: rational fraction, function of pressure + ! 1: tanh of altitude ! Local variables: @@ -167,5 +173,5 @@ ! -------------------------------------------------- - if (ok_strato .and. llm==39) then + if (vert_prof_dissip == 1) then do l=1,llm pseudoz=8.*log(preff/presnivs(l)) Index: LMDZ5/branches/testing/libf/dyn3d/leapfrog.F =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/leapfrog.F (revision 1669) +++ LMDZ5/branches/testing/libf/dyn3d/leapfrog.F (revision 1707) @@ -383,4 +383,10 @@ jD_cur = jD_ref + day_ini - day_ref + & & itau/day_step + + IF (planet_type .eq."generic") THEN + ! AS: we make jD_cur to be pday + jD_cur = int(day_ini + itau/day_step) + ENDIF + jH_cur = jH_ref + start_time + & & mod(itau,day_step)/float(day_step) Index: LMDZ5/branches/testing/libf/dyn3d/limit_netcdf.F90 =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/limit_netcdf.F90 (revision 1669) +++ (revision ) @@ -1,642 +1,0 @@ -! -! $Id$ -!------------------------------------------------------------------------------- -! -SUBROUTINE limit_netcdf(interbar, extrap, oldice, masque) -! -!------------------------------------------------------------------------------- -! Author : L. Fairhead, 27/01/94 -!------------------------------------------------------------------------------- -! Purpose: Boundary conditions files building for new model using climatologies. -! Both grids have to be regular. -!------------------------------------------------------------------------------- -! Note: This routine is designed to work for Earth -!------------------------------------------------------------------------------- -! Modification history: -! * 23/03/1994: Z. X. Li -! * 09/1999: L. Fairhead (netcdf reading in LMDZ.3.3) -! * 07/2001: P. Le Van -! * 11/2009: L. Guez (ozone day & night climatos, see etat0_netcdf.F90) -! * 12/2009: D. Cugnet (f77->f90, calendars, files from coupled runs) -!------------------------------------------------------------------------------- - USE control_mod -#ifdef CPP_EARTH - USE dimphy - USE ioipsl, ONLY : ioget_year_len - USE phys_state_var_mod, ONLY : pctsrf - USE netcdf, ONLY : NF90_OPEN, NF90_CREATE, NF90_CLOSE, & - NF90_DEF_DIM, NF90_DEF_VAR, NF90_PUT_VAR, NF90_PUT_ATT, & - NF90_NOERR, NF90_NOWRITE, NF90_DOUBLE, NF90_GLOBAL, & - NF90_CLOBBER, NF90_ENDDEF, NF90_UNLIMITED, NF90_FLOAT - USE inter_barxy_m, only: inter_barxy -#endif - IMPLICIT NONE -!------------------------------------------------------------------------------- -! Arguments: -#include "dimensions.h" -#include "paramet.h" -#include "iniprint.h" - LOGICAL, INTENT(IN) :: interbar ! barycentric interpolation - LOGICAL, INTENT(IN) :: extrap ! SST extrapolation flag - LOGICAL, INTENT(IN) :: oldice ! old way ice computation - REAL, DIMENSION(iip1,jjp1), INTENT(IN) :: masque ! land mask -#ifndef CPP_EARTH - CALL abort_gcm('limit_netcdf','Earth-specific routine, needs Earth physics',1) -#else -!------------------------------------------------------------------------------- -! Local variables: -#include "logic.h" -#include "comvert.h" -#include "comgeom2.h" -#include "comconst.h" -#include "indicesol.h" - -!--- INPUT NETCDF FILES NAMES -------------------------------------------------- - CHARACTER(LEN=25) :: icefile, sstfile, dumstr - CHARACTER(LEN=25), PARAMETER :: famipsst='amipbc_sst_1x1.nc ', & - famipsic='amipbc_sic_1x1.nc ', & - fcpldsst='cpl_atm_sst.nc ', & - fcpldsic='cpl_atm_sic.nc ', & - fhistsst='histmth_sst.nc ', & - fhistsic='histmth_sic.nc ', & - frugo ='Rugos.nc ', & - falbe ='Albedo.nc ' - CHARACTER(LEN=10) :: varname -!--- OUTPUT VARIABLES FOR NETCDF FILE ------------------------------------------ - REAL, DIMENSION(klon) :: fi_ice, verif - REAL, DIMENSION(:,:), POINTER :: phy_rug=>NULL(), phy_ice=>NULL() - REAL, DIMENSION(:,:), POINTER :: phy_sst=>NULL(), phy_alb=>NULL() - REAL, DIMENSION(:,:), ALLOCATABLE :: phy_bil - REAL, DIMENSION(:,:,:), ALLOCATABLE :: pctsrf_t - INTEGER :: nbad - -!--- VARIABLES FOR OUTPUT FILE WRITING ----------------------------------------- - INTEGER :: ierr, nid, ndim, ntim, k - INTEGER, DIMENSION(2) :: dims - INTEGER :: id_tim, id_SST, id_BILS, id_RUG, id_ALB - INTEGER :: id_FOCE, id_FSIC, id_FTER, id_FLIC - INTEGER :: NF90_FORMAT - INTEGER :: ndays !--- Depending on the output calendar - -!--- INITIALIZATIONS ----------------------------------------------------------- -#ifdef NC_DOUBLE - NF90_FORMAT=NF90_DOUBLE -#else - NF90_FORMAT=NF90_FLOAT -#endif - - pi = 4.*ATAN(1.) - rad = 6371229. - daysec= 86400. - omeg = 2.*pi/daysec - g = 9.8 - kappa = 0.2857143 - cpp = 1004.70885 - dtvr = daysec/REAL(day_step) - CALL inigeom - -!--- Beware: anneeref (from gcm.def) is used to determine output time sampling - ndays=ioget_year_len(anneeref) - -!--- RUGOSITY TREATMENT -------------------------------------------------------- - IF (prt_level>1) WRITE(lunout,*) 'Traitement de la rugosite' - varname='RUGOS' - CALL get_2Dfield(frugo,varname,'RUG',interbar,ndays,phy_rug,mask=masque(1:iim,:)) - -!--- OCEAN TREATMENT ----------------------------------------------------------- - IF (prt_level>1) WRITE(lunout,*) 'Traitement de la glace oceanique' - -! Input SIC file selection -! Open file only to test if available - IF ( NF90_OPEN(TRIM(famipsic),NF90_NOWRITE,nid)==NF90_NOERR ) THEN - icefile=TRIM(famipsic) - varname='sicbcs' - ELSE IF( NF90_OPEN(TRIM(fcpldsic),NF90_NOWRITE,nid)==NF90_NOERR ) THEN - icefile=TRIM(fcpldsic) - varname='SIICECOV' - ELSE IF ( NF90_OPEN(TRIM(fhistsic),NF90_NOWRITE,nid)==NF90_NOERR ) THEN - icefile=TRIM(fhistsic) - varname='pourc_sic' - ELSE - WRITE(lunout,*) 'ERROR! No sea-ice input file was found.' - WRITE(lunout,*) 'One of following files must be availible : ',trim(famipsic),', ',trim(fcpldsic),', ',trim(fhistsic) - CALL abort_gcm('limit_netcdf','No sea-ice file was found',1) - END IF - ierr=NF90_CLOSE(nid) - IF (prt_level>=0) WRITE(lunout,*)'Pour la glace de mer a ete choisi le fichier '//TRIM(icefile) - - CALL get_2Dfield(icefile,varname, 'SIC',interbar,ndays,phy_ice,flag=oldice) - - ALLOCATE(pctsrf_t(klon,nbsrf,ndays)) - DO k=1,ndays - fi_ice=phy_ice(:,k) - WHERE(fi_ice>=1.0 ) fi_ice=1.0 - WHERE(fi_ice=1.0-zmasq) - pctsrf_t(:,is_sic,k)=1.0-zmasq - pctsrf_t(:,is_oce,k)=0.0 - ELSEWHERE - pctsrf_t(:,is_oce,k)=1.0-zmasq-pctsrf_t(:,is_sic,k) - WHERE(pctsrf_t(:,is_oce,k)0) WRITE(lunout,*) 'pb sous maille pour nb point = ',nbad - nbad=COUNT(abs(sum(pctsrf_t(:,:,k),dim=2)-1.0)>EPSFRA) - IF(nbad>0) WRITE(lunout,*) 'pb sous surface pour nb points = ',nbad - END DO - DEALLOCATE(phy_ice) - -!--- SST TREATMENT ------------------------------------------------------------- - IF (prt_level>1) WRITE(lunout,*) 'Traitement de la sst' - -! Input SST file selection -! Open file only to test if available - IF ( NF90_OPEN(TRIM(famipsst),NF90_NOWRITE,nid)==NF90_NOERR ) THEN - sstfile=TRIM(famipsst) - varname='tosbcs' - ELSE IF ( NF90_OPEN(TRIM(fcpldsst),NF90_NOWRITE,nid)==NF90_NOERR ) THEN - sstfile=TRIM(fcpldsst) - varname='SISUTESW' - ELSE IF ( NF90_OPEN(TRIM(fhistsst),NF90_NOWRITE,nid)==NF90_NOERR ) THEN - sstfile=TRIM(fhistsst) - varname='tsol_oce' - ELSE - WRITE(lunout,*) 'ERROR! No sst input file was found.' - WRITE(lunout,*) 'One of following files must be availible : ',trim(famipsst),trim(fcpldsst),trim(fhistsst) - CALL abort_gcm('limit_netcdf','No sst file was found',1) - END IF - ierr=NF90_CLOSE(nid) - IF (prt_level>=0) WRITE(lunout,*)'Pour la temperature de mer a ete choisi un fichier '//TRIM(sstfile) - - CALL get_2Dfield(sstfile,varname,'SST',interbar,ndays,phy_sst,flag=extrap) - -!--- ALBEDO TREATMENT ---------------------------------------------------------- - IF (prt_level>1) WRITE(lunout,*) 'Traitement de l albedo' - varname='ALBEDO' - CALL get_2Dfield(falbe,varname,'ALB',interbar,ndays,phy_alb) - -!--- REFERENCE GROUND HEAT FLUX TREATMENT -------------------------------------- - ALLOCATE(phy_bil(klon,ndays)); phy_bil=0.0 - -!--- OUTPUT FILE WRITING ------------------------------------------------------- - IF (prt_level>5) WRITE(lunout,*) 'Ecriture du fichier limit : debut' - - !--- File creation - ierr=NF90_CREATE("limit.nc",NF90_CLOBBER,nid) - ierr=NF90_PUT_ATT(nid,NF90_GLOBAL,"title","Fichier conditions aux limites") - - !--- Dimensions creation - ierr=NF90_DEF_DIM(nid,"points_physiques",klon,ndim) - ierr=NF90_DEF_DIM(nid,"time",NF90_UNLIMITED,ntim) - - dims=(/ndim,ntim/) - - !--- Variables creation - ierr=NF90_DEF_VAR(nid,"TEMPS",NF90_FORMAT,(/ntim/),id_tim) - ierr=NF90_DEF_VAR(nid,"FOCE", NF90_FORMAT,dims,id_FOCE) - ierr=NF90_DEF_VAR(nid,"FSIC", NF90_FORMAT,dims,id_FSIC) - ierr=NF90_DEF_VAR(nid,"FTER", NF90_FORMAT,dims,id_FTER) - ierr=NF90_DEF_VAR(nid,"FLIC", NF90_FORMAT,dims,id_FLIC) - ierr=NF90_DEF_VAR(nid,"SST", NF90_FORMAT,dims,id_SST) - ierr=NF90_DEF_VAR(nid,"BILS", NF90_FORMAT,dims,id_BILS) - ierr=NF90_DEF_VAR(nid,"ALB", NF90_FORMAT,dims,id_ALB) - ierr=NF90_DEF_VAR(nid,"RUG", NF90_FORMAT,dims,id_RUG) - - !--- Attributes creation - ierr=NF90_PUT_ATT(nid,id_tim, "title","Jour dans l annee") - ierr=NF90_PUT_ATT(nid,id_FOCE,"title","Fraction ocean") - ierr=NF90_PUT_ATT(nid,id_FSIC,"title","Fraction glace de mer") - ierr=NF90_PUT_ATT(nid,id_FTER,"title","Fraction terre") - ierr=NF90_PUT_ATT(nid,id_FLIC,"title","Fraction land ice") - ierr=NF90_PUT_ATT(nid,id_SST ,"title","Temperature superficielle de la mer") - ierr=NF90_PUT_ATT(nid,id_BILS,"title","Reference flux de chaleur au sol") - ierr=NF90_PUT_ATT(nid,id_ALB, "title","Albedo a la surface") - ierr=NF90_PUT_ATT(nid,id_RUG, "title","Rugosite") - - ierr=NF90_ENDDEF(nid) - - !--- Variables saving - ierr=NF90_PUT_VAR(nid,id_tim,(/(REAL(k),k=1,ndays)/)) - ierr=NF90_PUT_VAR(nid,id_FOCE,pctsrf_t(:,is_oce,:),(/1,1/),(/klon,ndays/)) - ierr=NF90_PUT_VAR(nid,id_FSIC,pctsrf_t(:,is_sic,:),(/1,1/),(/klon,ndays/)) - ierr=NF90_PUT_VAR(nid,id_FTER,pctsrf_t(:,is_ter,:),(/1,1/),(/klon,ndays/)) - ierr=NF90_PUT_VAR(nid,id_FLIC,pctsrf_t(:,is_lic,:),(/1,1/),(/klon,ndays/)) - ierr=NF90_PUT_VAR(nid,id_SST ,phy_sst(:,:),(/1,1/),(/klon,ndays/)) - ierr=NF90_PUT_VAR(nid,id_BILS,phy_bil(:,:),(/1,1/),(/klon,ndays/)) - ierr=NF90_PUT_VAR(nid,id_ALB ,phy_alb(:,:),(/1,1/),(/klon,ndays/)) - ierr=NF90_PUT_VAR(nid,id_RUG ,phy_rug(:,:),(/1,1/),(/klon,ndays/)) - - ierr=NF90_CLOSE(nid) - - IF (prt_level>5) WRITE(lunout,*) 'Ecriture du fichier limit : fin' - - DEALLOCATE(pctsrf_t,phy_sst,phy_bil,phy_alb,phy_rug) - - -!=============================================================================== -! - CONTAINS -! -!=============================================================================== - - -!------------------------------------------------------------------------------- -! -SUBROUTINE get_2Dfield(fnam, varname, mode, ibar, ndays, champo, flag, mask) -! -!----------------------------------------------------------------------------- -! Comments: -! There are two assumptions concerning the NetCDF files, that are satisfied -! with files that are conforming NC convention: -! 1) The last dimension of the variables used is the time record. -! 2) Dimensional variables have the same names as corresponding dimensions. -!----------------------------------------------------------------------------- - USE netcdf, ONLY: NF90_OPEN, NF90_INQ_VARID, NF90_INQUIRE_VARIABLE, & - NF90_CLOSE, NF90_INQ_DIMID, NF90_INQUIRE_DIMENSION, NF90_GET_VAR, & - NF90_GET_ATT - USE dimphy, ONLY : klon - USE phys_state_var_mod, ONLY : pctsrf - USE control_mod - use pchsp_95_m, only: pchsp_95 - use pchfe_95_m, only: pchfe_95 - use arth_m, only: arth - - IMPLICIT NONE -#include "dimensions.h" -#include "paramet.h" -#include "comgeom2.h" -#include "indicesol.h" -#include "iniprint.h" -!----------------------------------------------------------------------------- -! Arguments: - CHARACTER(LEN=*), INTENT(IN) :: fnam ! NetCDF file name - CHARACTER(LEN=10), INTENT(IN) :: varname ! NetCDF variable name - CHARACTER(LEN=3), INTENT(IN) :: mode ! RUG, SIC, SST or ALB - LOGICAL, INTENT(IN) :: ibar ! interp on pressure levels - INTEGER, INTENT(IN) :: ndays ! current year number of days - REAL, POINTER, DIMENSION(:, :) :: champo ! output field = f(t) - LOGICAL, OPTIONAL, INTENT(IN) :: flag ! extrapol. (SST) old ice (SIC) - REAL, OPTIONAL, DIMENSION(iim, jjp1), INTENT(IN) :: mask -!------------------------------------------------------------------------------ -! Local variables: -!--- NetCDF - INTEGER :: ncid, varid ! NetCDF identifiers - CHARACTER(LEN=30) :: dnam ! dimension name -!--- dimensions - INTEGER, DIMENSION(4) :: dids ! NetCDF dimensions identifiers - REAL, ALLOCATABLE, DIMENSION(:) :: dlon_ini ! initial longitudes vector - REAL, ALLOCATABLE, DIMENSION(:) :: dlat_ini ! initial latitudes vector - REAL, POINTER, DIMENSION(:) :: dlon, dlat ! reordered lon/lat vectors -!--- fields - INTEGER :: imdep, jmdep, lmdep ! dimensions of 'champ' - REAL, ALLOCATABLE, DIMENSION(:, :) :: champ ! wanted field on initial grid - REAL, ALLOCATABLE, DIMENSION(:) :: yder, timeyear - REAL, DIMENSION(iim, jjp1) :: champint ! interpolated field - REAL, ALLOCATABLE, DIMENSION(:, :, :) :: champtime - REAL, ALLOCATABLE, DIMENSION(:, :, :) :: champan -!--- input files - CHARACTER(LEN=20) :: cal_in ! calendar - CHARACTER(LEN=20) :: unit_sic ! attribute unit in sea-ice file - INTEGER :: ndays_in ! number of days -!--- misc - INTEGER :: i, j, k, l ! loop counters - REAL, ALLOCATABLE, DIMENSION(:, :) :: work ! used for extrapolation - CHARACTER(LEN=25) :: title ! for messages - LOGICAL :: extrp ! flag for extrapolation - LOGICAL :: oldice ! flag for old way ice computation - REAL :: chmin, chmax - INTEGER ierr - integer n_extrap ! number of extrapolated points - logical skip - -!------------------------------------------------------------------------------ -!---Variables depending on keyword 'mode' ------------------------------------- - NULLIFY(champo) - - SELECT CASE(mode) - CASE('RUG'); title='Rugosite' - CASE('SIC'); title='Sea-ice' - CASE('SST'); title='SST' - CASE('ALB'); title='Albedo' - END SELECT - - - extrp=.FALSE. - oldice=.FALSE. - IF ( PRESENT(flag) ) THEN - IF ( flag .AND. mode=='SST' ) extrp=.TRUE. - IF ( flag .AND. mode=='SIC' ) oldice=.TRUE. - END IF - -!--- GETTING SOME DIMENSIONAL VARIABLES FROM FILE ----------------------------- - IF (prt_level>5) WRITE(lunout,*) ' Now reading file : ',fnam - ierr=NF90_OPEN(fnam, NF90_NOWRITE, ncid); CALL ncerr(ierr, fnam) - ierr=NF90_INQ_VARID(ncid, trim(varname), varid); CALL ncerr(ierr, fnam) - ierr=NF90_INQUIRE_VARIABLE(ncid, varid, dimids=dids); CALL ncerr(ierr, fnam) - -!--- Read unit for sea-ice variable only - IF (mode=='SIC') THEN - ierr=NF90_GET_ATT(ncid, varid, 'units', unit_sic) - IF(ierr/=NF90_NOERR) THEN - IF (prt_level>5) WRITE(lunout,*) 'No unit was given in sea-ice file. Take percentage as default value' - unit_sic='X' - ELSE - IF (prt_level>5) WRITE(lunout,*) ' Sea-ice cover has unit=',unit_sic - END IF - END IF - -!--- Longitude - ierr=NF90_INQUIRE_DIMENSION(ncid, dids(1), name=dnam, len=imdep) - CALL ncerr(ierr, fnam); ALLOCATE(dlon_ini(imdep), dlon(imdep)) - ierr=NF90_INQ_VARID(ncid, dnam, varid); CALL ncerr(ierr, fnam) - ierr=NF90_GET_VAR(ncid, varid, dlon_ini); CALL ncerr(ierr, fnam) - IF (prt_level>5) WRITE(lunout, *) 'variable ', dnam, 'dimension ', imdep - -!--- Latitude - ierr=NF90_INQUIRE_DIMENSION(ncid, dids(2), name=dnam, len=jmdep) - CALL ncerr(ierr, fnam); ALLOCATE(dlat_ini(jmdep), dlat(jmdep)) - ierr=NF90_INQ_VARID(ncid, dnam, varid); CALL ncerr(ierr, fnam) - ierr=NF90_GET_VAR(ncid, varid, dlat_ini); CALL ncerr(ierr, fnam) - IF (prt_level>5) WRITE(lunout, *) 'variable ', dnam, 'dimension ', jmdep - -!--- Time (variable is not needed - it is rebuilt - but calendar is) - ierr=NF90_INQUIRE_DIMENSION(ncid, dids(3), name=dnam, len=lmdep) - CALL ncerr(ierr, fnam); ALLOCATE(timeyear(lmdep)) - ierr=NF90_INQ_VARID(ncid, dnam, varid); CALL ncerr(ierr, fnam) - cal_in=' ' - ierr=NF90_GET_ATT(ncid, varid, 'calendar', cal_in) - IF(ierr/=NF90_NOERR) THEN - SELECT CASE(mode) - CASE('RUG', 'ALB'); cal_in='360d' - CASE('SIC', 'SST'); cal_in='gregorian' - END SELECT - IF (prt_level>5) WRITE(lunout, *)'ATTENTION: variable "time" sans attribut "calendrier" ' & - // 'dans '//TRIM(fnam)//'. On choisit la valeur par defaut.' - END IF - IF (prt_level>5) WRITE(lunout, *) 'variable ', dnam, 'dimension ', lmdep, 'calendrier ', & - cal_in - - -!--- CONSTRUCTING THE INPUT TIME VECTOR FOR INTERPOLATION -------------------- - !--- Determining input file number of days, depending on calendar - ndays_in=year_len(anneeref, cal_in) - -!--- Time vector reconstruction (time vector from file is not trusted) -!--- If input records are not monthly, time sampling has to be constant ! - timeyear=mid_months(anneeref, cal_in, lmdep) - IF (lmdep /= 12) WRITE(lunout,*) 'Note : les fichiers de ', TRIM(mode), & - ' ne comportent pas 12, mais ', lmdep, ' enregistrements.' - -!--- GETTING THE FIELD AND INTERPOLATING IT ---------------------------------- - ALLOCATE(champ(imdep, jmdep), champtime(iim, jjp1, lmdep)) - IF(extrp) ALLOCATE(work(imdep, jmdep)) - - IF (prt_level>5) WRITE(lunout, *) - IF (prt_level>5) WRITE(lunout,*)'LECTURE ET INTERPOLATION HORIZ. DE ', lmdep, ' CHAMPS.' - ierr=NF90_INQ_VARID(ncid, varname, varid); CALL ncerr(ierr, fnam) - DO l=1, lmdep - ierr=NF90_GET_VAR(ncid, varid, champ, (/1, 1, l/), (/imdep, jmdep, 1/)) - CALL ncerr(ierr, fnam) - CALL conf_dat2d(title, imdep, jmdep, dlon_ini, dlat_ini, dlon, dlat, & - champ, ibar) - - IF (extrp) CALL extrapol(champ, imdep, jmdep, 999999., .TRUE., .TRUE., 2, & - work) - - IF(ibar .AND. .NOT.oldice) THEN - IF(l==1 .AND. prt_level>5) THEN - WRITE(lunout, *) '-------------------------------------------------------------------------' - WRITE(lunout, *) 'Utilisation de l''interpolation barycentrique pour ',TRIM(title),' $$$' - WRITE(lunout, *) '-------------------------------------------------------------------------' - END IF - IF(mode=='RUG') champ=LOG(champ) - CALL inter_barxy(dlon, dlat(:jmdep-1), champ, rlonu(:iim), & - rlatv, champint) - IF(mode=='RUG') THEN - champint=EXP(champint) - WHERE(NINT(mask)/=1) champint=0.001 - END IF - ELSE - SELECT CASE(mode) - CASE('RUG'); CALL rugosite(imdep, jmdep, dlon, dlat, champ, & - iim, jjp1, rlonv, rlatu, champint, mask) - CASE('SIC'); CALL sea_ice (imdep, jmdep, dlon, dlat, champ, & - iim, jjp1, rlonv, rlatu, champint) - CASE('SST', 'ALB'); CALL grille_m(imdep, jmdep, dlon, dlat, champ, & - iim, jjp1, rlonv, rlatu, champint) - END SELECT - END IF - champtime(:, :, l)=champint - END DO - ierr=NF90_CLOSE(ncid); CALL ncerr(ierr, fnam) - - DEALLOCATE(dlon_ini, dlat_ini, dlon, dlat, champ) - IF(extrp) DEALLOCATE(work) - -!--- TIME INTERPOLATION ------------------------------------------------------ - IF (prt_level>5) THEN - WRITE(lunout, *) - WRITE(lunout, *)'INTERPOLATION TEMPORELLE.' - WRITE(lunout, *)' Vecteur temps en entree: ', timeyear - WRITE(lunout, *)' Vecteur temps en sortie de 0 a ', ndays - END IF - - ALLOCATE(yder(lmdep), champan(iip1, jjp1, ndays)) - skip = .false. - n_extrap = 0 - DO j=1, jjp1 - DO i=1, iim - yder = pchsp_95(timeyear, champtime(i, j, :), ibeg=2, iend=2, & - vc_beg=0., vc_end=0.) - CALL pchfe_95(timeyear, champtime(i, j, :), yder, skip, & - arth(0., real(ndays_in) / ndays, ndays), champan(i, j, :), ierr) - if (ierr < 0) stop 1 - n_extrap = n_extrap + ierr - END DO - END DO - if (n_extrap /= 0) then - WRITE(lunout,*) "get_2Dfield pchfe_95: n_extrap = ", n_extrap - end if - champan(iip1, :, :)=champan(1, :, :) - DEALLOCATE(yder, champtime, timeyear) - -!--- Checking the result - DO j=1, jjp1 - CALL minmax(iip1, champan(1, j, 10), chmin, chmax) - IF (prt_level>5) WRITE(lunout, *)' ',TRIM(title),' au temps 10 ', chmin, chmax, j - END DO - -!--- SPECIAL FILTER FOR SST: SST>271.38 -------------------------------------- - IF(mode=='SST') THEN - IF (prt_level>5) WRITE(lunout, *) 'Filtrage de la SST: SST >= 271.38' - WHERE(champan<271.38) champan=271.38 - END IF - -!--- SPECIAL FILTER FOR SIC: 0.05) WRITE(lunout, *) 'Filtrage de la SIC: 0.0 < Sea-ice < 1.0' - - IF (unit_sic=='1') THEN - ! Nothing to be done for sea-ice field is already in fraction of 1 - ! This is the case for sea-ice in file cpl_atm_sic.nc - IF (prt_level>5) WRITE(lunout,*) 'Sea-ice field already in fraction of 1' - ELSE - ! Convert sea ice from percentage to fraction of 1 - IF (prt_level>5) WRITE(lunout,*) 'Transformt sea-ice field from percentage to fraction of 1.' - champan(:, :, :)=champan(:, :, :)/100. - END IF - - champan(iip1, :, :)=champan(1, :, :) - WHERE(champan>1.0) champan=1.0 - WHERE(champan<0.0) champan=0.0 - END IF - -!--- DYNAMICAL TO PHYSICAL GRID ---------------------------------------------- - ALLOCATE(champo(klon, ndays)) - DO k=1, ndays - CALL gr_dyn_fi(1, iip1, jjp1, klon, champan(1, 1, k), champo(1, k)) - END DO - DEALLOCATE(champan) - -END SUBROUTINE get_2Dfield -! -!------------------------------------------------------------------------------- - - - -!------------------------------------------------------------------------------- -! -FUNCTION year_len(y,cal_in) -! -!------------------------------------------------------------------------------- - USE ioipsl, ONLY : ioget_calendar,ioconf_calendar,lock_calendar,ioget_year_len - IMPLICIT NONE -!------------------------------------------------------------------------------- -! Arguments: - INTEGER :: year_len - INTEGER, INTENT(IN) :: y - CHARACTER(LEN=*), INTENT(IN) :: cal_in -!------------------------------------------------------------------------------- -! Local variables: - CHARACTER(LEN=20) :: cal_out ! calendar (for outputs) -!------------------------------------------------------------------------------- -!--- Getting the input calendar to reset at the end of the function - CALL ioget_calendar(cal_out) - -!--- Unlocking calendar and setting it to wanted one - CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_in)) - -!--- Getting the number of days in this year - year_len=ioget_year_len(y) - -!--- Back to original calendar - CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_out)) - -END FUNCTION year_len -! -!------------------------------------------------------------------------------- - - -!------------------------------------------------------------------------------- -! -FUNCTION mid_months(y,cal_in,nm) -! -!------------------------------------------------------------------------------- - USE ioipsl, ONLY : ioget_calendar,ioconf_calendar,lock_calendar,ioget_mon_len - IMPLICIT NONE -!------------------------------------------------------------------------------- -! Arguments: - INTEGER, INTENT(IN) :: y ! year - CHARACTER(LEN=*), INTENT(IN) :: cal_in ! calendar - INTEGER, INTENT(IN) :: nm ! months/year number - REAL, DIMENSION(nm) :: mid_months ! mid-month times -!------------------------------------------------------------------------------- -! Local variables: - CHARACTER(LEN=99) :: mess ! error message - CHARACTER(LEN=20) :: cal_out ! calendar (for outputs) - INTEGER, DIMENSION(nm) :: mnth ! months lengths (days) - INTEGER :: m ! months counter - INTEGER :: nd ! number of days -!------------------------------------------------------------------------------- - nd=year_len(y,cal_in) - - IF(nm==12) THEN - - !--- Getting the input calendar to reset at the end of the function - CALL ioget_calendar(cal_out) - - !--- Unlocking calendar and setting it to wanted one - CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_in)) - - !--- Getting the length of each month - DO m=1,nm; mnth(m)=ioget_mon_len(y,m); END DO - - !--- Back to original calendar - CALL lock_calendar(.FALSE.); CALL ioconf_calendar(TRIM(cal_out)) - - ELSE IF(MODULO(nd,nm)/=0) THEN - WRITE(mess,'(a,i3,a,i3,a)')'Unconsistent calendar: ',nd,' days/year, but ',& - nm,' months/year. Months number should divide days number.' - CALL abort_gcm('mid_months',TRIM(mess),1) - - ELSE - mnth=(/(m,m=1,nm,nd/nm)/) - END IF - -!--- Mid-months times - mid_months(1)=0.5*REAL(mnth(1)) - DO k=2,nm - mid_months(k)=mid_months(k-1)+0.5*REAL(mnth(k-1)+mnth(k)) - END DO - -END FUNCTION mid_months -! -!------------------------------------------------------------------------------- - - -!------------------------------------------------------------------------------- -! -SUBROUTINE ncerr(ncres,fnam) -! -!------------------------------------------------------------------------------- -! Purpose: NetCDF errors handling. -!------------------------------------------------------------------------------- - USE netcdf, ONLY : NF90_NOERR, NF90_STRERROR - IMPLICIT NONE -!------------------------------------------------------------------------------- -! Arguments: - INTEGER, INTENT(IN) :: ncres - CHARACTER(LEN=*), INTENT(IN) :: fnam -!------------------------------------------------------------------------------- -#include "iniprint.h" - IF(ncres/=NF90_NOERR) THEN - WRITE(lunout,*)'Problem with file '//TRIM(fnam)//' in routine limit_netcdf.' - CALL abort_gcm('limit_netcdf',NF90_STRERROR(ncres),1) - END IF - -END SUBROUTINE ncerr -! -!------------------------------------------------------------------------------- - -#endif -! of #ifdef CPP_EARTH - -END SUBROUTINE limit_netcdf Index: LMDZ5/branches/testing/libf/dyn3d/mod_const_mpi.F90 =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/mod_const_mpi.F90 (revision 1707) +++ LMDZ5/branches/testing/libf/dyn3d/mod_const_mpi.F90 (revision 1707) @@ -0,0 +1,16 @@ +MODULE mod_const_mpi + + INTEGER :: COMM_LMDZ + INTEGER :: MPI_REAL_LMDZ + + +CONTAINS + + SUBROUTINE Init_const_mpi + IMPLICIT NONE + + COMM_LMDZ=0 + MPI_REAL_LMDZ=0 + END SUBROUTINE Init_const_mpi + +END MODULE mod_const_mpi Index: LMDZ5/branches/testing/libf/dyn3d/mod_const_para.F90 =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/mod_const_para.F90 (revision 1669) +++ (revision ) @@ -1,16 +1,0 @@ -MODULE mod_const_mpi - - INTEGER :: COMM_LMDZ - INTEGER :: MPI_REAL_LMDZ - - -CONTAINS - - SUBROUTINE Init_const_mpi - IMPLICIT NONE - - COMM_LMDZ=0 - MPI_REAL_LMDZ=0 - END SUBROUTINE Init_const_mpi - -END MODULE mod_const_mpi Index: LMDZ5/branches/testing/libf/dyn3d/paramet.h =================================================================== --- LMDZ5/branches/testing/libf/dyn3d/paramet.h (revision 1669) +++ LMDZ5/branches/testing/libf/dyn3d/paramet.h (revision 1707) @@ -17,5 +17,5 @@ INTEGER jcfil,jcfllm - PARAMETER( iip1= iim+1-1/iim,iip2=iim+2,iip3=iim+3 & + PARAMETER( iip1= iim+1,iip2=iim+2,iip3=iim+3 & & ,jjp1=jjm+1-1/jjm) PARAMETER( llmp1 = llm+1, llmp2 = llm+2, llmm1 = llm-1 )