Changeset 1528 for trunk/LMDZ.MARS/libf/phymars/writediagfi.F

Timestamp:

Apr 2, 2016, 4:09:43 PM (9 years ago)

Author:

emillour

Message:

Mars GCM:

• Got rid of references to "dimensions.h" from physics packages: use nbp_lon (=iim), nbp_lat (==jjp1) and nbp_lev from module mod_grid_phy_lmdz (in phy_common) instead.
• Added "ioipsl_getin_p_mod.F90" (getin_p routine) in phy_common to correctly read in parameters from *.def files in a parallel environment.

EM

File:

• trunk/LMDZ.MARS/libf/phymars/writediagfi.F (modified) (13 diffs)

trunk/LMDZ.MARS/libf/phymars/writediagfi.F

@@ -40 +40 @@
 !=================================================================
       use surfdat_h, only: phisfi
+      use comgeomphy, only: airephy
       use time_phylmdz_mod, only: ecritphy, day_step, iphysiq, day_ini
       USE mod_phys_lmdz_para, only : is_parallel, is_mpi_root,
      &                               is_master, gather
-      USE mod_grid_phy_lmdz, only : klon_glo, Grid1Dto2D_glo
+      USE mod_grid_phy_lmdz, only : klon_glo, Grid1Dto2D_glo,
+     &                              nbp_lon, nbp_lat, nbp_lev
       implicit none
 
 ! Commons
-#include "dimensions.h"
-#include "paramet.h"
-#include "comgeom.h"
-#include "netcdf.inc"
+      include "netcdf.inc"
 
 ! Arguments on input:
@@ -56 +55 @@
       character (len=*),intent(in) :: nom,titre,unite
       integer,intent(in) :: dim
-      real,intent(in) :: px(ngrid,llm)
+      real,intent(in) :: px(ngrid,nbp_lev)
 
 ! Local variables:
 
-      real*4 dx3(iip1,jjp1,llm) ! to store a 3D data set
-      real*4 dx2(iip1,jjp1)     ! to store a 2D (surface) data set
-      real*4 dx1(llm)           ! to store a 1D (column) data set
+      real*4 dx3(nbp_lon+1,nbp_lat,nbp_lev) ! to store a 3D data set
+      real*4 dx2(nbp_lon+1,nbp_lat)     ! to store a 2D (surface) data set
+      real*4 dx1(nbp_lev)           ! to store a 1D (column) data set
       real*4 dx0
 
       real*4,save :: date
 
-      REAL phis(ip1jmp1)
+      REAL phis((nbp_lon+1),nbp_lat)
+      REAL area((nbp_lon+1),nbp_lat)
 
       integer irythme
       integer ierr,ierr2
-      integer iq
-      integer i,j,l,zmax , ig0
+      integer i,j,l, ig0
 
       integer,save :: zitau=0
@@ -98 +97 @@
 #ifdef CPP_PARA
 ! Added to work in parallel mode
-      real dx3_glop(klon_glo,llm)
-      real dx3_glo(iim,jjp1,llm) ! to store a global 3D data set
+      real dx3_glop(klon_glo,nbp_lev)
+      real dx3_glo(nbp_lon,nbp_lat,nbp_lev) ! to store a global 3D data set
       real dx2_glop(klon_glo)
-      real dx2_glo(iim,jjp1)     ! to store a global 2D (surface) data set
+      real dx2_glo(nbp_lon,nbp_lat)     ! to store a global 2D (surface) data set
       real px2(ngrid)
-!      real dx1_glo(llm)          ! to store a 1D (column) data set
+!      real dx1_glo(nbp_lev)          ! to store a 1D (column) data set
 !      real dx0_glo
       real phisfi_glo(klon_glo) ! surface geopotential on global physics grid
+      real areafi_glo(klon_glo) ! mesh area on global physics grid
 #else
       real phisfi_glo(ngrid) ! surface geopotential on global physics grid
+      real areafi_glo(ngrid) ! mesh area on global physics grid
 #endif
 
@@ -181 +182 @@
           ! Gather phisfi() geopotential on physics grid
           call Gather(phisfi,phisfi_glo)
+          ! Gather airephy() mesh area on physics grid
+          call Gather(airephy,areafi_glo)
 #else
          phisfi_glo(:)=phisfi(:)
+         areafi_glo(:)=airephy(:)
 #endif
 
@@ -209 +213 @@
          ierr = NF_ENDDEF(nid)
 
+         ! Build phis() and area()
+         do i=1,nbp_lon+1 ! poles
+           phis(i,1)=phisfi_glo(1)
+           phis(i,nbp_lat)=phisfi_glo(klon_glo)
+           ! for area, divide at the poles by nbp_lon
+           area(i,1)=areafi_glo(1)/nbp_lon
+           area(i,nbp_lat)=areafi_glo(klon_glo)/nbp_lon
+         enddo
+         do j=2,nbp_lat-1
+           ig0= 1+(j-2)*nbp_lon
+           do i=1,nbp_lon
+              phis(i,j)=phisfi_glo(ig0+i)
+              area(i,j)=areafi_glo(ig0+i)
+           enddo
+           ! handle redundant point in longitude
+           phis(nbp_lon+1,j)=phis(1,j)
+           area(nbp_lon+1,j)=area(1,j)
+         enddo
+         
          ! write "header" of file (longitudes, latitudes, geopotential, ...)
-         call gr_fi_dyn(1,size(phisfi_glo),iip1,jjp1,phisfi_glo,phis)
-         call iniwrite(nid,day_ini,phis)
+         call iniwrite(nid,day_ini,phis,area)
 
          endif ! of if (is_master)
@@ -290 +312 @@
             call Grid1Dto2D_glo(dx3_glop,dx3_glo)
             ! copy dx3_glo() to dx3(:) and add redundant longitude
-            dx3(1:iim,:,:)=dx3_glo(1:iim,:,:)
-            dx3(iip1,:,:)=dx3(1,:,:)
+            dx3(1:nbp_lon,:,:)=dx3_glo(1:nbp_lon,:,:)
+            dx3(nbp_lon+1,:,:)=dx3(1,:,:)
           endif
 !$OMP END MASTER
@@ -298 +320 @@
 !         Passage variable physique -->  variable dynamique
 !         recast (copy) variable from physics grid to dynamics grid
-           DO l=1,llm
-             DO i=1,iip1
+           DO l=1,nbp_lev
+             DO i=1,nbp_lon+1
                 dx3(i,1,l)=px(1,l)
-                dx3(i,jjp1,l)=px(ngrid,l)
+                dx3(i,nbp_lat,l)=px(ngrid,l)
              ENDDO
-             DO j=2,jjm
-                ig0= 1+(j-2)*iim
-                DO i=1,iim
+             DO j=2,nbp_lat-1
+                ig0= 1+(j-2)*nbp_lon
+                DO i=1,nbp_lon
                    dx3(i,j,l)=px(ig0+i,l)
                 ENDDO
-                dx3(iip1,j,l)=dx3(1,j,l)
+                dx3(nbp_lon+1,j,l)=dx3(1,j,l)
              ENDDO
            ENDDO
@@ -338 +360 @@
            corner(4)=ntime
 
-           edges(1)=iip1
-           edges(2)=jjp1
-           edges(3)=llm
+           edges(1)=nbp_lon+1
+           edges(2)=nbp_lat
+           edges(3)=nbp_lev
            edges(4)=1
 !#ifdef NC_DOUBLE
@@ -374 +396 @@
             call Grid1Dto2D_glo(dx2_glop,dx2_glo)
             ! copy dx2_glo() to dx2(:) and add redundant longitude
-            dx2(1:iim,:)=dx2_glo(1:iim,:)
-            dx2(iip1,:)=dx2(1,:)
+            dx2(1:nbp_lon,:)=dx2_glo(1:nbp_lon,:)
+            dx2(nbp_lon+1,:)=dx2(1,:)
           endif
 !$OMP END MASTER
@@ -384 +406 @@
 !         recast (copy) variable from physics grid to dynamics grid
 
-             DO i=1,iip1
+             DO i=1,nbp_lon+1
                 dx2(i,1)=px(1,1)
-                dx2(i,jjp1)=px(ngrid,1)
+                dx2(i,nbp_lat)=px(ngrid,1)
              ENDDO
-             DO j=2,jjm
-                ig0= 1+(j-2)*iim
-                DO i=1,iim
+             DO j=2,nbp_lat-1
+                ig0= 1+(j-2)*nbp_lon
+                DO i=1,nbp_lon
                    dx2(i,j)=px(ig0+i,1)
                 ENDDO
-                dx2(iip1,j)=dx2(1,j)
+                dx2(nbp_lon+1,j)=dx2(1,j)
              ENDDO
 #endif
@@ -420 +442 @@
            corner(2)=1
            corner(3)=ntime
-           edges(1)=iip1
-           edges(2)=jjp1
+           edges(1)=nbp_lon+1
+           edges(2)=nbp_lat
            edges(3)=1
 
@@ -451 +473 @@
 !         Passage variable physique -->  physique dynamique
 !         recast (copy) variable from physics grid to dynamics grid
-          do l=1,llm
+          do l=1,nbp_lev
             dx1(l)=px(1,l)
           enddo
@@ -473 +495 @@
            corner(2)=ntime
            
-           edges(1)=llm
+           edges(1)=nbp_lev
            edges(2)=1
 !#ifdef NC_DOUBLE