Index: trunk/LMDZ.MARS/libf/dynphy_lonlat/calfis.F =================================================================== --- trunk/LMDZ.MARS/libf/dynphy_lonlat/calfis.F (revision 2502) +++ (revision ) @@ -1,559 +1,0 @@ - SUBROUTINE calfis(nq, lafin, rdayvrai,rday_ecri, heure, - $ pucov,pvcov,pteta,pq,pmasse,pps,pp,ppk,pphis,pphi, - $ pducov,pdvcov,pdteta,pdq,pw, - $ pdufi,pdvfi,pdhfi,pdqfi,pdpsfi ) -c -c Auteur : P. Le Van, F. Hourdin -c ......... - - USE comvert_mod, ONLY: preff - USE comconst_mod, ONLY: dtphys,cpp,kappa,pi - USE physiq_mod, ONLY: physiq - - IMPLICIT NONE -c======================================================================= -c -c 1. rearrangement des tableaux et transformation -c variables dynamiques > variables physiques -c 2. calcul des termes physiques -c 3. retransformation des tendances physiques en tendances dynamiques -c -c remarques: -c ---------- -c -c - les vents sont donnes dans la physique par leurs composantes -c naturelles. -c - la variable thermodynamique de la physique est une variable -c intensive : T -c pour la dynamique on prend T * ( preff / p(l) ) **kappa -c - les deux seules variables dependant de la geometrie necessaires -c pour la physique sont la latitude pour le rayonnement et -c l'aire de la maille quand on veut integrer une grandeur -c horizontalement. -c - les points de la physique sont les points scalaires de la -c la dynamique; numerotation: -c 1 pour le pole nord -c (jjm-1)*iim pour l'interieur du domaine -c ngridmx pour le pole sud -c ---> ngridmx=2+(jjm-1)*iim -c -c Input : -c ------- -c ecritphy frequence d'ecriture (en jours)de histphy -c pucov covariant zonal velocity -c pvcov covariant meridional velocity -c pteta potential temperature -c pps surface pressure -c pmasse masse d'air dans chaque maille -c pts surface temperature (K) -c pw flux vertical (kg/s) -c -c Output : -c -------- -c pdufi tendency for the natural zonal velocity (ms-1) -c pdvfi tendency for the natural meridional velocity -c pdhfi tendency for the potential temperature -c pdtsfi tendency for the surface temperature -c -c======================================================================= -c -c----------------------------------------------------------------------- -c -c 0. Declarations : -c ------------------ - -#include "dimensions.h" -#include "paramet.h" - - INTEGER ngridmx,nq - PARAMETER( ngridmx = 2+(jjm-1)*iim - 1/jjm ) - -#include "comgeom2.h" -!#include "control.h" - -c Arguments : -c ----------- - LOGICAL lafin - REAL heure - - REAL pvcov(iip1,jjm,llm) - REAL pucov(iip1,jjp1,llm) - REAL pteta(iip1,jjp1,llm) - REAL pmasse(iip1,jjp1,llm) - REAL pq(iip1,jjp1,llm,nq) - REAL pphis(iip1,jjp1) - REAL pphi(iip1,jjp1,llm) -c - REAL pdvcov(iip1,jjm,llm) - REAL pducov(iip1,jjp1,llm) - REAL pdteta(iip1,jjp1,llm) - REAL pdq(iip1,jjp1,llm,nq) -c - REAL pw(iip1,jjp1,llm) -c - REAL pps(iip1,jjp1) - REAL pp(iip1,jjp1,llmp1) - REAL ppk(iip1,jjp1,llm) -c - REAL pdvfi(iip1,jjm,llm) - REAL pdufi(iip1,jjp1,llm) - REAL pdhfi(iip1,jjp1,llm) - REAL pdqfi(iip1,jjp1,llm,nq) - REAL pdpsfi(iip1,jjp1) - -c Local variables : -c ----------------- - - INTEGER i,j,l,ig0,ig,iq - REAL zpsrf(ngridmx) - REAL zplev(ngridmx,llm+1),zplay(ngridmx,llm) - REAL zphi(ngridmx,llm),zphis(ngridmx) -c - REAL zufi(ngridmx,llm), zvfi(ngridmx,llm) - REAL ztfi(ngridmx,llm),zqfi(ngridmx,llm,nq) -c -! REAL zvervel(ngridmx,llm) - REAL flxwfi(ngridmx,llm) ! vertical mass flux (kg/s) on physics grid -c - REAL zdufi(ngridmx,llm),zdvfi(ngridmx,llm) - REAL zdtfi(ngridmx,llm),zdqfi(ngridmx,llm,nq) - REAL zdpsrf(ngridmx) -c - REAL zsin(iim),zcos(iim),z1(iim) - REAL zsinbis(iim),zcosbis(iim),z1bis(iim) - REAL unskap, pksurcp -c - - EXTERNAL gr_dyn_fi,gr_fi_dyn - REAL SSUM - EXTERNAL SSUM - - REAL latfi(ngridmx),lonfi(ngridmx) - REAL airefi(ngridmx) - SAVE latfi, lonfi, airefi - - LOGICAL firstcal, debut - DATA firstcal/.true./ - SAVE firstcal,debut - REAL rdayvrai,rday_ecri -c -c----------------------------------------------------------------------- -c -c 1. Initialisations : -c -------------------- -c - - - IF (ngridmx.NE.2+(jjm-1)*iim) THEN - PRINT*,'STOP dans calfis' - PRINT*,'La dimension ngridmx doit etre egale a 2 + (jjm-1)*iim' - PRINT*,' ngridmx jjm iim ' - PRINT*,ngridmx,jjm,iim - STOP - ENDIF - -c----------------------------------------------------------------------- -c latitude, longitude et aires des mailles pour la physique: -c ---------------------------------------------------------- - -c - IF ( firstcal ) THEN - debut = .TRUE. - ELSE - debut = .FALSE. - ENDIF - -c -! IF (firstcal) THEN -! latfi(1)=rlatu(1) -! lonfi(1)=0. -! DO j=2,jjm -! DO i=1,iim -! latfi((j-2)*iim+1+i)= rlatu(j) -! lonfi((j-2)*iim+1+i)= rlonv(i) -! ENDDO -! ENDDO -! latfi(ngridmx)= rlatu(jjp1) -! lonfi(ngridmx)= 0. -! -! ! build airefi(), mesh area on physics grid -! CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,aire,airefi) -! ! Poles are single points on physics grid -! airefi(1)=airefi(1)*iim -! airefi(ngridmx)=airefi(ngridmx)*iim -! -! CALL inifis(ngridmx,llm,nq,day_ini,daysec,dtphys, -! . latfi,lonfi,airefi,rad,g,r,cpp) -! ENDIF - -c -c----------------------------------------------------------------------- -c 40. transformation des variables dynamiques en variables physiques: -c --------------------------------------------------------------- - -c 41. pressions au sol (en Pascals) -c ---------------------------------- - - - zpsrf(1) = pps(1,1) - - ig0 = 2 - DO j = 2,jjm - CALL SCOPY( iim,pps(1,j),1,zpsrf(ig0), 1 ) - ig0 = ig0+iim - ENDDO - - zpsrf(ngridmx) = pps(1,jjp1) - - -c 42. pression intercouches : -c -c ----------------------------------------------------------------- -c .... zplev definis aux (llm +1) interfaces des couches .... -c .... zplay definis aux ( llm ) milieux des couches .... -c ----------------------------------------------------------------- - -c ... Exner = cp * ( p(l) / preff ) ** kappa .... -c - unskap = 1./ kappa -c - DO l = 1, llmp1 - zplev( 1,l ) = pp(1,1,l) - ig0 = 2 - DO j = 2, jjm - DO i =1, iim - zplev( ig0,l ) = pp(i,j,l) - ig0 = ig0 +1 - ENDDO - ENDDO - zplev( ngridmx,l ) = pp(1,jjp1,l) - ENDDO -c -c - -c 43. temperature naturelle (en K) et pressions milieux couches . -c --------------------------------------------------------------- - - DO l=1,llm - - pksurcp = ppk(1,1,l) / cpp - zplay(1,l) = preff * pksurcp ** unskap - ztfi(1,l) = pteta(1,1,l) * pksurcp - ig0 = 2 - - DO j = 2, jjm - DO i = 1, iim - pksurcp = ppk(i,j,l) / cpp - zplay(ig0,l) = preff * pksurcp ** unskap - ztfi(ig0,l) = pteta(i,j,l) * pksurcp - ig0 = ig0 + 1 - ENDDO - ENDDO - - pksurcp = ppk(1,jjp1,l) / cpp - zplay(ig0,l) = preff * pksurcp ** unskap - ztfi (ig0,l) = pteta(1,jjp1,l) * pksurcp - - ENDDO - - DO l=1, llm - DO ig=1,ngridmx - if (ztfi(ig,l).lt.15) then - write(*,*) 'New Temperature below 15 K !!! ' - write(*,*) 'Stop in calfis.F ' - write(*,*) 'ig=', ig, ' l=', l - write(*,*) 'ztfi(ig,l)=',ztfi(ig,l) - stop - end if - ENDDO - ENDDO - - - -c 43.bis Taceurs (en kg/kg) -c -------------------------- - DO iq=1,nq - DO l=1,llm - zqfi(1,l,iq) = pq(1,1,l,iq) - ig0 = 2 - DO j=2,jjm - DO i = 1, iim - zqfi(ig0,l,iq) = pq(i,j,l,iq) - ig0 = ig0 + 1 - ENDDO - ENDDO - zqfi(ig0,l,iq) = pq(1,jjp1,l,iq) - ENDDO - ENDDO - -c Geopotentiel calcule par rapport a la surface locale: -c ----------------------------------------------------- - - CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,pphi,zphi) - CALL gr_dyn_fi(1,iip1,jjp1,ngridmx,pphis,zphis) - DO l=1,llm - DO ig=1,ngridmx - zphi(ig,l)=zphi(ig,l)-zphis(ig) - ENDDO - ENDDO - -c Calcul de la vitesse verticale (m/s) pour diagnostique -c ------------------------------- -c pw est en kg/s -c On interpole "lineairement" la temperature entre les couches(FF,10/95) - -! DO ig=1,ngridmx -! zvervel(ig,1)=0. -! END DO -! DO l=2,llm -! zvervel(1,l)=(pw(1,1,l)/apoln) -! & * r *0.5*(ztfi(1,l)+ztfi(1,l-1)) /zplev(1,l) -! ig0=2 -! DO j=2,jjm -! DO i = 1, iim -! zvervel(ig0,l) = pw(i,j,l) * unsaire(i,j) -! & * r *0.5*(ztfi(ig0,l)+ztfi(ig0,l-1)) /zplev(ig0,l) -! ig0 = ig0 + 1 -! ENDDO -! ENDDO -! zvervel(ig0,l)=(pw(1,jjp1,l)/apols) -! & * r *0.5*(ztfi(ig0,l)+ztfi(ig0,l-1)) /zplev(ig0,l) -! ENDDO - -c ......... Reindexation : calcul de zvervel au MILIEU des couches -! DO l=1,llm-1 -! DO ig=1,ngridmx -! zvervel(ig,l) = 0.5*(zvervel(ig,l)+zvervel(ig,l+1)) -! END DO -! END DO -c (dans la couche llm, on garde la valeur à la limite inférieure llm) - -! vertical mass flux - ! tranfer values from dynamics grid to physics grid: - CALL gr_dyn_fi(llm,iip1,jjp1,ngridmx,pw,flxwfi) - ! but mass flux is an extensive variable, so take the sum at the poles - DO l=1,llm - flxwfi(1,l)=sum(pw(1:iim,1,l)) - flxwfi(ngridmx,l)=sum(pw(1:iim,jjp1,l)) - ENDDO - -c 45. champ u: -c ------------ - - DO l=1,llm - DO j=2,jjm - ig0 = 1+(j-2)*iim - zufi(ig0+1,l)= 0.5 * - $ ( pucov(iim,j,l)/cu(iim,j) + pucov(1,j,l)/cu(1,j) ) - DO i=2,iim - zufi(ig0+i,l)= 0.5 * - $ ( pucov(i-1,j,l)/cu(i-1,j) + pucov(i,j,l)/cu(i,j) ) - ENDDO - ENDDO - ENDDO - - -c 46.champ v: -c ----------- - - DO l=1,llm - DO j=2,jjm - ig0=1+(j-2)*iim - DO i=1,iim - zvfi(ig0+i,l)= 0.5 * - $ ( pvcov(i,j-1,l)/cv(i,j-1) + pvcov(i,j,l)/cv(i,j) ) - ENDDO - ENDDO - ENDDO - - -c 47. champs de vents aux pole nord -c ------------------------------ -c U = 1 / pi * integrale [ v * cos(long) * d long ] -c V = 1 / pi * integrale [ v * sin(long) * d long ] - - DO l=1,llm - - z1(1) =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,1,l)/cv(1,1) - z1bis(1)=(rlonu(1)-rlonu(iim)+2.*pi)*pdvcov(1,1,l)/cv(1,1) - DO i=2,iim - z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,1,l)/cv(i,1) - z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,1,l)/cv(i,1) - ENDDO - - DO i=1,iim - zcos(i) = COS(rlonv(i))*z1(i) - zcosbis(i)= COS(rlonv(i))*z1bis(i) - zsin(i) = SIN(rlonv(i))*z1(i) - zsinbis(i)= SIN(rlonv(i))*z1bis(i) - ENDDO - - zufi(1,l) = SSUM(iim,zcos,1)/pi - zvfi(1,l) = SSUM(iim,zsin,1)/pi - - ENDDO - - -c 48. champs de vents aux pole sud: -c --------------------------------- -c U = 1 / pi * integrale [ v * cos(long) * d long ] -c V = 1 / pi * integrale [ v * sin(long) * d long ] - - DO l=1,llm - - z1(1) =(rlonu(1)-rlonu(iim)+2.*pi)*pvcov(1,jjm,l)/cv(1,jjm) - z1bis(1)=(rlonu(1)-rlonu(iim)+2.*pi)*pdvcov(1,jjm,l)/cv(1,jjm) - DO i=2,iim - z1(i) =(rlonu(i)-rlonu(i-1))*pvcov(i,jjm,l)/cv(i,jjm) - z1bis(i)=(rlonu(i)-rlonu(i-1))*pdvcov(i,jjm,l)/cv(i,jjm) - ENDDO - - DO i=1,iim - zcos(i) = COS(rlonv(i))*z1(i) - zcosbis(i) = COS(rlonv(i))*z1bis(i) - zsin(i) = SIN(rlonv(i))*z1(i) - zsinbis(i) = SIN(rlonv(i))*z1bis(i) - ENDDO - - zufi(ngridmx,l) = SSUM(iim,zcos,1)/pi - zvfi(ngridmx,l) = SSUM(iim,zsin,1)/pi - - ENDDO - -c----------------------------------------------------------------------- -c Appel de la physique: -c --------------------- - - - CALL physiq (ngridmx,llm,nq, - , debut,lafin, - , rday_ecri,heure,dtphys, - , zplev,zplay,zphi, - , zufi, zvfi,ztfi, zqfi, -! , zvervel, - , flxwfi, -C - sorties - s zdufi, zdvfi, zdtfi, zdqfi,zdpsrf) - - -c----------------------------------------------------------------------- -c transformation des tendances physiques en tendances dynamiques: -c --------------------------------------------------------------- - -c tendance sur la pression : -c ----------------------------------- - - CALL gr_fi_dyn(1,ngridmx,iip1,jjp1,zdpsrf,pdpsfi) -c -ccc CALL multipl(ip1jmp1,aire,pdpsfi,pdpsfi) - -c 62. enthalpie potentielle -c --------------------- - - DO l=1,llm - - DO i=1,iip1 - pdhfi(i,1,l) = cpp * zdtfi(1,l) / ppk(i, 1 ,l) - pdhfi(i,jjp1,l) = cpp * zdtfi(ngridmx,l)/ ppk(i,jjp1,l) - ENDDO - - DO j=2,jjm - ig0=1+(j-2)*iim - DO i=1,iim - pdhfi(i,j,l) = cpp * zdtfi(ig0+i,l) / ppk(i,j,l) - ENDDO - pdhfi(iip1,j,l) = pdhfi(1,j,l) - ENDDO - - ENDDO - - -c 62. humidite specifique -c --------------------- - - DO iq=1,nq - DO l=1,llm - DO i=1,iip1 - pdqfi(i,1,l,iq) = zdqfi(1,l,iq) - pdqfi(i,jjp1,l,iq) = zdqfi(ngridmx,l,iq) - ENDDO - DO j=2,jjm - ig0=1+(j-2)*iim - DO i=1,iim - pdqfi(i,j,l,iq) = zdqfi(ig0+i,l,iq) - ENDDO - pdqfi(iip1,j,l,iq) = pdqfi(1,j,l,iq) - ENDDO - ENDDO - ENDDO - -c 65. champ u: -c ------------ - - DO l=1,llm - - DO i=1,iip1 - pdufi(i,1,l) = 0. - pdufi(i,jjp1,l) = 0. - ENDDO - - DO j=2,jjm - ig0=1+(j-2)*iim - DO i=1,iim-1 - pdufi(i,j,l)= - $ 0.5*(zdufi(ig0+i,l)+zdufi(ig0+i+1,l))*cu(i,j) - ENDDO - pdufi(iim,j,l)= - $ 0.5*(zdufi(ig0+1,l)+zdufi(ig0+iim,l))*cu(iim,j) - pdufi(iip1,j,l)=pdufi(1,j,l) - ENDDO - - ENDDO - - -c 67. champ v: -c ------------ - - DO l=1,llm - - DO j=2,jjm-1 - ig0=1+(j-2)*iim - DO i=1,iim - pdvfi(i,j,l)= - $ 0.5*(zdvfi(ig0+i,l)+zdvfi(ig0+i+iim,l))*cv(i,j) - ENDDO - pdvfi(iip1,j,l) = pdvfi(1,j,l) - ENDDO - ENDDO - - -c 68. champ v pres des poles: -c --------------------------- -c v = U * cos(long) + V * SIN(long) - - DO l=1,llm - - DO i=1,iim - pdvfi(i,1,l)= - $ zdufi(1,l)*COS(rlonv(i))+zdvfi(1,l)*SIN(rlonv(i)) - pdvfi(i,jjm,l)=zdufi(ngridmx,l)*COS(rlonv(i)) - $ +zdvfi(ngridmx,l)*SIN(rlonv(i)) - pdvfi(i,1,l)= - $ 0.5*(pdvfi(i,1,l)+zdvfi(i+1,l))*cv(i,1) - pdvfi(i,jjm,l)= - $ 0.5*(pdvfi(i,jjm,l)+zdvfi(ngridmx-iip1+i,l))*cv(i,jjm) - ENDDO - - pdvfi(iip1,1,l) = pdvfi(1,1,l) - pdvfi(iip1,jjm,l)= pdvfi(1,jjm,l) - - ENDDO - -c----------------------------------------------------------------------- - -700 CONTINUE - - firstcal = .FALSE. - - RETURN - END Index: trunk/LMDZ.MARS/libf/dynphy_lonlat/gr_dyn_fi.F =================================================================== --- trunk/LMDZ.MARS/libf/dynphy_lonlat/gr_dyn_fi.F (revision 2502) +++ (revision ) @@ -1,37 +1,0 @@ - SUBROUTINE gr_dyn_fi(nfield,im,jm,ngrid,pdyn,pfi) - - IMPLICIT NONE -c======================================================================= -c passage d'un champ de la grille scalaire a la grille physique -c======================================================================= - -c----------------------------------------------------------------------- -c declarations: -c ------------- - - INTEGER im,jm,ngrid,nfield - REAL pdyn(im,jm,nfield) - REAL pfi(ngrid,nfield) - - INTEGER j,ifield,ig - EXTERNAL SCOPY - -c----------------------------------------------------------------------- -c calcul: -c ------- - - IF(ngrid.NE.2+(jm-2)*(im-1)) STOP 'probleme de dim' -c traitement des poles - CALL SCOPY(nfield,pdyn,im*jm,pfi,ngrid) - CALL SCOPY(nfield,pdyn(1,jm,1),im*jm,pfi(ngrid,1),ngrid) - -c traitement des point normaux - DO ifield=1,nfield - DO j=2,jm-1 - ig=2+(j-2)*(im-1) - CALL SCOPY(im-1,pdyn(1,j,ifield),1,pfi(ig,ifield),1) - ENDDO - ENDDO - - RETURN - END Index: trunk/LMDZ.MARS/libf/dynphy_lonlat/gr_fi_dyn.F =================================================================== --- trunk/LMDZ.MARS/libf/dynphy_lonlat/gr_fi_dyn.F (revision 2502) +++ (revision ) @@ -1,38 +1,0 @@ - SUBROUTINE gr_fi_dyn(nfield,ngrid,im,jm,pfi,pdyn) - IMPLICIT NONE -c======================================================================= -c passage d'un champ de la grille scalaire a la grille physique -c======================================================================= - -c----------------------------------------------------------------------- -c declarations: -c ------------- - - INTEGER im,jm,ngrid,nfield - REAL pdyn(im,jm,nfield) - REAL pfi(ngrid,nfield) - - INTEGER i,j,ifield,ig - EXTERNAL SCOPY - -c----------------------------------------------------------------------- -c calcul: -c ------- - - DO ifield=1,nfield -c traitement des poles - DO i=1,im - pdyn(i,1,ifield)=pfi(1,ifield) - pdyn(i,jm,ifield)=pfi(ngrid,ifield) - ENDDO - -c traitement des point normaux - DO j=2,jm-1 - ig=2+(j-2)*(im-1) - CALL SCOPY(im-1,pfi(ig,ifield),1,pdyn(1,j,ifield),1) - pdyn(im,j,ifield)=pdyn(1,j,ifield) - ENDDO - ENDDO - - RETURN - END Index: trunk/LMDZ.MARS/libf/dynphy_lonlat/inigeomphy_mod.F90 =================================================================== --- trunk/LMDZ.MARS/libf/dynphy_lonlat/inigeomphy_mod.F90 (revision 2502) +++ (revision ) @@ -1,230 +1,0 @@ -MODULE inigeomphy_mod - -CONTAINS - -SUBROUTINE inigeomphy(iim,jjm,nlayer, & - nbp, communicator, & - rlatu,rlatv,rlonu,rlonv,aire,cu,cv) - USE mod_grid_phy_lmdz, ONLY: klon_glo, & ! number of atmospheric columns (on full grid) - regular_lonlat, & ! regular longitude-latitude grid type - nbp_lon, nbp_lat, nbp_lev - USE mod_phys_lmdz_para, ONLY: klon_omp, & ! number of columns (on local omp grid) - klon_omp_begin, & ! start index of local omp subgrid - klon_omp_end, & ! end index of local omp subgrid - klon_mpi_begin ! start indes of columns (on local mpi grid) - USE geometry_mod, ONLY : init_geometry - USE physics_distribution_mod, ONLY : init_physics_distribution - USE regular_lonlat_mod, ONLY : init_regular_lonlat, & - east, west, north, south, & - north_east, north_west, & - south_west, south_east - USE mod_interface_dyn_phys, ONLY : init_interface_dyn_phys - USE nrtype, ONLY: pi - USE comvert_mod, ONLY: preff, ap, bp, aps, bps, presnivs, & - scaleheight, pseudoalt - USE vertical_layers_mod, ONLY: init_vertical_layers - IMPLICIT NONE - - ! ======================================================================= - ! Initialisation of the physical constants and some positional and - ! geometrical arrays for the physics - ! ======================================================================= - - include "iniprint.h" - - INTEGER, INTENT (IN) :: nlayer ! number of atmospheric layers - INTEGER, INTENT (IN) :: iim ! number of atmospheric columns along longitudes - INTEGER, INTENT (IN) :: jjm ! number of atompsheric columns along latitudes - INTEGER, INTENT(IN) :: nbp ! number of physics columns for this MPI process - INTEGER, INTENT(IN) :: communicator ! MPI communicator - REAL, INTENT (IN) :: rlatu(jjm+1) ! latitudes of the physics grid - REAL, INTENT (IN) :: rlatv(jjm) ! latitude boundaries of the physics grid - REAL, INTENT (IN) :: rlonv(iim+1) ! longitudes of the physics grid - REAL, INTENT (IN) :: rlonu(iim+1) ! longitude boundaries of the physics grid - REAL, INTENT (IN) :: aire(iim+1,jjm+1) ! area of the dynamics grid (m2) - REAL, INTENT (IN) :: cu((iim+1)*(jjm+1)) ! cu coeff. (u_covariant = cu * u) - REAL, INTENT (IN) :: cv((iim+1)*jjm) ! cv coeff. (v_covariant = cv * v) - - INTEGER :: ibegin, iend, offset - INTEGER :: i,j,k - CHARACTER (LEN=20) :: modname = 'iniphysiq' - CHARACTER (LEN=80) :: abort_message - REAL :: total_area_phy, total_area_dyn - - ! boundaries, on global grid - REAL,ALLOCATABLE :: boundslon_reg(:,:) - REAL,ALLOCATABLE :: boundslat_reg(:,:) - - ! global array, on full physics grid: - REAL,ALLOCATABLE :: latfi_glo(:) - REAL,ALLOCATABLE :: lonfi_glo(:) - REAL,ALLOCATABLE :: cufi_glo(:) - REAL,ALLOCATABLE :: cvfi_glo(:) - REAL,ALLOCATABLE :: airefi_glo(:) - REAL,ALLOCATABLE :: boundslonfi_glo(:,:) - REAL,ALLOCATABLE :: boundslatfi_glo(:,:) - - ! local arrays, on given MPI/OpenMP domain: - REAL,ALLOCATABLE,SAVE :: latfi(:) - REAL,ALLOCATABLE,SAVE :: lonfi(:) - REAL,ALLOCATABLE,SAVE :: cufi(:) - REAL,ALLOCATABLE,SAVE :: cvfi(:) - REAL,ALLOCATABLE,SAVE :: airefi(:) - REAL,ALLOCATABLE,SAVE :: boundslonfi(:,:) - REAL,ALLOCATABLE,SAVE :: boundslatfi(:,:) - INTEGER,ALLOCATABLE,SAVE :: ind_cell_glo_fi(:) -!$OMP THREADPRIVATE (latfi,lonfi,cufi,cvfi,airefi,boundslonfi,boundslatfi,ind_cell_glo_fi) - - ! Initialize Physics distibution and parameters and interface with dynamics - IF (iim*jjm>1) THEN ! general 3D case - CALL init_physics_distribution(regular_lonlat,4, & - nbp,iim,jjm+1,nlayer,communicator) - ELSE ! For 1D model - CALL init_physics_distribution(regular_lonlat,4, & - 1,1,1,nlayer,communicator) - ENDIF - CALL init_interface_dyn_phys - - ! init regular global longitude-latitude grid points and boundaries - ALLOCATE(boundslon_reg(iim,2)) - ALLOCATE(boundslat_reg(jjm+1,2)) - - DO i=1,iim - boundslon_reg(i,east)=rlonu(i) - boundslon_reg(i,west)=rlonu(i+1) - ENDDO - - boundslat_reg(1,north)= PI/2 - boundslat_reg(1,south)= rlatv(1) - DO j=2,jjm - boundslat_reg(j,north)=rlatv(j-1) - boundslat_reg(j,south)=rlatv(j) - ENDDO - boundslat_reg(jjm+1,north)= rlatv(jjm) - boundslat_reg(jjm+1,south)= -PI/2 - - ! Write values in module regular_lonlat_mod - CALL init_regular_lonlat(iim,jjm+1, rlonv(1:iim), rlatu, & - boundslon_reg, boundslat_reg) - - ! Generate global arrays on full physics grid - ALLOCATE(latfi_glo(klon_glo),lonfi_glo(klon_glo)) - ALLOCATE(cufi_glo(klon_glo),cvfi_glo(klon_glo)) - ALLOCATE(airefi_glo(klon_glo)) - ALLOCATE(boundslonfi_glo(klon_glo,4)) - ALLOCATE(boundslatfi_glo(klon_glo,4)) - - IF (klon_glo>1) THEN ! general case - ! North pole - latfi_glo(1)=rlatu(1) - lonfi_glo(1)=0. - cufi_glo(1) = cu(1) - cvfi_glo(1) = cv(1) - boundslonfi_glo(1,north_east)=0 - boundslatfi_glo(1,north_east)=PI/2 - boundslonfi_glo(1,north_west)=2*PI - boundslatfi_glo(1,north_west)=PI/2 - boundslonfi_glo(1,south_west)=2*PI - boundslatfi_glo(1,south_west)=rlatv(1) - boundslonfi_glo(1,south_east)=0 - boundslatfi_glo(1,south_east)=rlatv(1) - DO j=2,jjm - DO i=1,iim - k=(j-2)*iim+1+i - latfi_glo(k)= rlatu(j) - lonfi_glo(k)= rlonv(i) - cufi_glo(k) = cu((j-1)*(iim+1)+i) - cvfi_glo(k) = cv((j-1)*(iim+1)+i) - boundslonfi_glo(k,north_east)=rlonu(i) - boundslatfi_glo(k,north_east)=rlatv(j-1) - boundslonfi_glo(k,north_west)=rlonu(i+1) - boundslatfi_glo(k,north_west)=rlatv(j-1) - boundslonfi_glo(k,south_west)=rlonu(i+1) - boundslatfi_glo(k,south_west)=rlatv(j) - boundslonfi_glo(k,south_east)=rlonu(i) - boundslatfi_glo(k,south_east)=rlatv(j) - ENDDO - ENDDO - ! South pole - latfi_glo(klon_glo)= rlatu(jjm+1) - lonfi_glo(klon_glo)= 0. - cufi_glo(klon_glo) = cu((iim+1)*jjm+1) - cvfi_glo(klon_glo) = cv((iim+1)*jjm-iim) - boundslonfi_glo(klon_glo,north_east)= 0 - boundslatfi_glo(klon_glo,north_east)= rlatv(jjm) - boundslonfi_glo(klon_glo,north_west)= 2*PI - boundslatfi_glo(klon_glo,north_west)= rlatv(jjm) - boundslonfi_glo(klon_glo,south_west)= 2*PI - boundslatfi_glo(klon_glo,south_west)= -PI/2 - boundslonfi_glo(klon_glo,south_east)= 0 - boundslatfi_glo(klon_glo,south_east)= -Pi/2 - - ! build airefi(), mesh area on physics grid - CALL gr_dyn_fi(1,iim+1,jjm+1,klon_glo,aire,airefi_glo) - ! Poles are single points on physics grid - airefi_glo(1)=sum(aire(1:iim,1)) - airefi_glo(klon_glo)=sum(aire(1:iim,jjm+1)) - - ! Sanity check: do total planet area match between physics and dynamics? - total_area_dyn=sum(aire(1:iim,1:jjm+1)) - total_area_phy=sum(airefi_glo(1:klon_glo)) - IF (total_area_dyn/=total_area_phy) THEN - WRITE (lunout, *) 'iniphysiq: planet total surface discrepancy !!!' - WRITE (lunout, *) ' in the dynamics total_area_dyn=', total_area_dyn - WRITE (lunout, *) ' but in the physics total_area_phy=', total_area_phy - IF (abs(total_area_dyn-total_area_phy)>0.00001*total_area_dyn) THEN - ! stop here if the relative difference is more than 0.001% - abort_message = 'planet total surface discrepancy' - CALL abort_gcm(modname, abort_message, 1) - ENDIF - ENDIF - ELSE ! klon_glo==1, running the 1D model - ! just copy over input values - latfi_glo(1)=rlatu(1) - lonfi_glo(1)=rlonv(1) - cufi_glo(1)=cu(1) - cvfi_glo(1)=cv(1) - airefi_glo(1)=aire(1,1) - boundslonfi_glo(1,north_east)=rlonu(1) - boundslatfi_glo(1,north_east)=PI/2 - boundslonfi_glo(1,north_west)=rlonu(2) - boundslatfi_glo(1,north_west)=PI/2 - boundslonfi_glo(1,south_west)=rlonu(2) - boundslatfi_glo(1,south_west)=rlatv(1) - boundslonfi_glo(1,south_east)=rlonu(1) - boundslatfi_glo(1,south_east)=rlatv(1) - ENDIF ! of IF (klon_glo>1) - -!$OMP PARALLEL - ! Now generate local lon/lat/cu/cv/area/bounds arrays - ALLOCATE(latfi(klon_omp),lonfi(klon_omp),cufi(klon_omp),cvfi(klon_omp)) - ALLOCATE(airefi(klon_omp)) - ALLOCATE(boundslonfi(klon_omp,4)) - ALLOCATE(boundslatfi(klon_omp,4)) - ALLOCATE(ind_cell_glo_fi(klon_omp)) - - - offset = klon_mpi_begin - 1 - airefi(1:klon_omp) = airefi_glo(offset+klon_omp_begin:offset+klon_omp_end) - cufi(1:klon_omp) = cufi_glo(offset+klon_omp_begin:offset+klon_omp_end) - cvfi(1:klon_omp) = cvfi_glo(offset+klon_omp_begin:offset+klon_omp_end) - lonfi(1:klon_omp) = lonfi_glo(offset+klon_omp_begin:offset+klon_omp_end) - latfi(1:klon_omp) = latfi_glo(offset+klon_omp_begin:offset+klon_omp_end) - boundslonfi(1:klon_omp,:) = boundslonfi_glo(offset+klon_omp_begin:offset+klon_omp_end,:) - boundslatfi(1:klon_omp,:) = boundslatfi_glo(offset+klon_omp_begin:offset+klon_omp_end,:) - ind_cell_glo_fi(1:klon_omp)=(/ (i,i=offset+klon_omp_begin,offset+klon_omp_end) /) - - ! copy over local grid longitudes and latitudes - CALL init_geometry(klon_omp,lonfi,latfi,boundslonfi,boundslatfi, & - airefi,ind_cell_glo_fi,cufi,cvfi) - - ! copy over preff , ap(), bp(), etc - CALL init_vertical_layers(nlayer,preff,scaleheight, & - ap,bp,aps,bps,presnivs,pseudoalt) - -!$OMP END PARALLEL - - -END SUBROUTINE inigeomphy - -END MODULE inigeomphy_mod Index: trunk/LMDZ.MARS/libf/dynphy_lonlat/mod_interface_dyn_phys.F90 =================================================================== --- trunk/LMDZ.MARS/libf/dynphy_lonlat/mod_interface_dyn_phys.F90 (revision 2502) +++ (revision ) @@ -1,62 +1,0 @@ -! -! $Id: mod_interface_dyn_phys.F90 2351 2015-08-25 15:14:59Z emillour $ -! -MODULE mod_interface_dyn_phys - INTEGER,SAVE,dimension(:),allocatable :: index_i - INTEGER,SAVE,dimension(:),allocatable :: index_j - - -CONTAINS - -#ifdef CPP_PARA -! Interface with parallel physics, - SUBROUTINE Init_interface_dyn_phys - USE mod_phys_lmdz_mpi_data - IMPLICIT NONE - include 'dimensions.h' - - INTEGER :: i,j,k - - ALLOCATE(index_i(klon_mpi)) - ALLOCATE(index_j(klon_mpi)) - - k=1 - IF (is_north_pole_dyn) THEN - index_i(k)=1 - index_j(k)=1 - k=2 - ELSE - DO i=ii_begin,iim - index_i(k)=i - index_j(k)=jj_begin - k=k+1 - ENDDO - ENDIF - - DO j=jj_begin+1,jj_end-1 - DO i=1,iim - index_i(k)=i - index_j(k)=j - k=k+1 - ENDDO - ENDDO - - IF (is_south_pole_dyn) THEN - index_i(k)=1 - index_j(k)=jj_end - ELSE - DO i=1,ii_end - index_i(k)=i - index_j(k)=jj_end - k=k+1 - ENDDO - ENDIF - - END SUBROUTINE Init_interface_dyn_phys -#else - SUBROUTINE Init_interface_dyn_phys - ! dummy routine for seq case - END SUBROUTINE Init_interface_dyn_phys -#endif -! of #ifdef CPP_PARA -END MODULE mod_interface_dyn_phys Index: trunk/LMDZ.MARS/libf/phymars/dyn1d/mod_interface_dyn_phys.F90 =================================================================== --- trunk/LMDZ.MARS/libf/phymars/dyn1d/mod_interface_dyn_phys.F90 (revision 2502) +++ trunk/LMDZ.MARS/libf/phymars/dyn1d/mod_interface_dyn_phys.F90 (revision 2505) @@ -1,1 +1,1 @@ -link ../../dynphy_lonlat/mod_interface_dyn_phys.F90 +link ../../../../LMDZ.COMMON/libf/dynphy_lonlat/mod_interface_dyn_phys.F90