From 0cadc749423a5db4562da8bf19ece6b8304d21d9 Mon Sep 17 00:00:00 2001
From: Jan Polcher <jan.polcher@lmd.jussieu.fr>
Date: Mon, 17 Jun 2019 16:14:55 +0200
Subject: [PATCH] Corrected some stuff and ensured that the diagnostics from
the F90 part are now written into files. This allows for finer diagnostics of
the computation there.
---
.gitignore | 1 +
Diags/plotgraphs.py | 10 +++++-----
DocumentationInterface | 4 +++-
F90subroutines/ioipsl.f90 | 2 +-
F90subroutines/routing_interface.f90 | 18 +++++++++++++++++-
Interface.py | 23 +++++++++++++++++++++--
RoutingPreProc.py | 4 +++-
7 files changed, 51 insertions(+), 11 deletions(-)
diff --git a/.gitignore b/.gitignore
index aca7377..5ad03b4 100644
--- a/.gitignore
+++ b/.gitignore
@@ -39,6 +39,7 @@ xcuserdata/
Weights/
__pycache__/
run.def
+Out*.txt
#
# Test configurations and directories
#
diff --git a/Diags/plotgraphs.py b/Diags/plotgraphs.py
index b705c2f..0875be4 100644
--- a/Diags/plotgraphs.py
+++ b/Diags/plotgraphs.py
@@ -6,9 +6,9 @@ import numpy as np
from netCDF4 import Dataset
import sys
-NCfile="../MEDCORDEX_test_graph.nc"
+NCfile="../SPortugal_test_graph_n4.nc"
GraphicFile = "test.png"
-title="Test"
+title="From file : "+NCfile
#
#####################################################################
#
@@ -101,8 +101,8 @@ for il in range(landindex.nbland) :
j,i = landindex.land2ij(il)
xt,yt=m(lon_full[j,i],lat_full[j,i])
plt.text(xt,yt,str(il))
- #
- for ih in range(nhtu) :
+ act_nhtu=np.count_nonzero(routetohtu[:,j,i] > 0)
+ for ih in range(act_nhtu) :
lo = cglon[ih,j,i]
la = cglat[ih,j,i]
xl,yl=m(lo, la)
@@ -117,7 +117,7 @@ for il in range(landindex.nbland) :
else :
print("Unforseen coding of outlof HTU : ",int(routetohtu[ih,j,i]-1))
sys.exit()
- #
+
if int(routetohtu[ih,j,i]-1) < nhtu :
ile = int(routetogrid[ih,j,i]-1)
je,ie = landindex.land2ij(ile)
diff --git a/DocumentationInterface b/DocumentationInterface
index 8770db9..d71e5bd 100644
--- a/DocumentationInterface
+++ b/DocumentationInterface
@@ -1,9 +1,11 @@
-initatmgrid(polygons_in,centre_in,area_in,contfrac_in,neighbours_in,[nbpt,nbseg])
+initatmgrid(rank_bn,nbcore,polygons_in,centre_in,area_in,contfrac_in,neighbours_in,[nbpt,nbseg])
Wrapper for ``initatmgrid``.
Parameters
----------
+rank_bn : input int
+nbcore : input int
polygons_in : input rank-3 array('f') with bounds (nbpt,2 * nbseg + 1,2)
centre_in : input rank-2 array('f') with bounds (nbpt,2)
area_in : input rank-1 array('f') with bounds (nbpt)
diff --git a/F90subroutines/ioipsl.f90 b/F90subroutines/ioipsl.f90
index e41dd69..9140288 100644
--- a/F90subroutines/ioipsl.f90
+++ b/F90subroutines/ioipsl.f90
@@ -4,7 +4,7 @@ MODULE ioipsl
PUBLIC :: ipslerr, ipslerr_p
- INTEGER, PARAMETER :: numout=6
+ INTEGER, SAVE :: numout=6
INTEGER :: ipslout=6, ilv_cur=0, ilv_max=0
LOGICAL :: ioipsl_debug=.FALSE., lact_mode=.TRUE.
diff --git a/F90subroutines/routing_interface.f90 b/F90subroutines/routing_interface.f90
index e6e4e20..591d6d2 100644
--- a/F90subroutines/routing_interface.f90
+++ b/F90subroutines/routing_interface.f90
@@ -1,11 +1,14 @@
-SUBROUTINE initatmgrid(nbpt, nbseg, polygons_in, centre_in, area_in, contfrac_in, neighbours_in)
+SUBROUTINE initatmgrid(rank, nbcore, nbpt, nbseg, polygons_in, centre_in, area_in, contfrac_in, neighbours_in)
!
USE grid
+ USE ioipsl
!
IMPLICIT NONE
!
!! INPUT VARIABLES
+ INTEGER, INTENT(in) :: rank !! rank of processor
+ INTEGER, INTENT(in) :: nbcore
INTEGER, INTENT(in) :: nbpt !! Atmospheric Domain size (unitless)
INTEGER, INTENT(in) :: nbseg !! Number of segments in the polygone
REAL, INTENT(in) :: polygons_in(nbpt,2*nbseg+1,2)
@@ -17,6 +20,7 @@ SUBROUTINE initatmgrid(nbpt, nbseg, polygons_in, centre_in, area_in, contfrac_in
!! Local
INTEGER :: i
INTEGER :: neighbours_loc(nbpt,nbseg*2)
+ CHARACTER(LEN=33) :: outfname
!
DO i=1,nbpt
! Move from C to F indexing.
@@ -25,8 +29,20 @@ SUBROUTINE initatmgrid(nbpt, nbseg, polygons_in, centre_in, area_in, contfrac_in
!
CALL grid_init(nbpt, nbseg, polygons_in, centre_in, area_in, contfrac_in, neighbours_loc)
!
+ numout=100+rank
+ WRITE(outfname,'(A15,I4.4,A6,I4.4,A4)') "Out_RoutingReg_",rank,"_from_",nbcore,".txt"
+ OPEN(unit=numout, file=outfname)
+ !
END SUBROUTINE initatmgrid
+SUBROUTINE closeinterface()
+ !
+ USE ioipsl
+ !
+ CLOSE(unit=numout)
+ !
+END SUBROUTINE closeinterface
+
SUBROUTINE gethydrogrid(nbpt, nbvmax_in, nbxmax_in, sub_pts, sub_index, sub_area, max_basins, min_topoind, &
& lon_rel, lat_rel, trip, basins, topoindex, fac, hierarchy, nbi, nbj, area_bx, trip_bx, basin_bx, &
& topoind_bx, fac_bx, hierarchy_bx, lon_bx, lat_bx, lshead_bx, nwbas)
diff --git a/Interface.py b/Interface.py
index ede1055..be95db1 100644
--- a/Interface.py
+++ b/Interface.py
@@ -52,11 +52,18 @@ if gendoc.lower() == "true" :
#
# initatmgrid : Initialises the grid.f90 module and passes the description of the atmospheric grid.
#
-def initatmgrid(nbpt, modelgrid) :
+def initatmgrid(rank, nbcore, nbpt, modelgrid) :
print("INITATMGRID corners", np.array(modelgrid.polyll).shape)
print("INITATMGRID area", np.array(modelgrid.area).shape)
print("INITATMGRID neighbours", np.array(modelgrid.neighbours).shape)
- routing_interface.initatmgrid( modelgrid.polyll, modelgrid.coordll, modelgrid.area, modelgrid.contfrac, modelgrid.neighbours)
+ routing_interface.initatmgrid(rank, nbcore, modelgrid.polyll, modelgrid.coordll, modelgrid.area, modelgrid.contfrac, modelgrid.neighbours)
+ return
+#
+#
+#
+def closeinterface(comm) :
+ comm.Barrier()
+ routing_interface.closeinterface()
return
#
#
@@ -267,7 +274,19 @@ class HydroGraph :
#
# The field route_togrid is with indices on the local grid. That needs to be converted to the global grid.
#
+ itarget=-1
+ for il in range(procgrid.nbland) :
+ lo = procgrid.lon_full[procgrid.indP[il][0],procgrid.indP[il][1]]
+ la = procgrid.lat_full[procgrid.indP[il][0],procgrid.indP[il][1]]
+ d=np.sqrt((lo-3.13)**2+(la-39.70)**2)
+ if d < 0.05 :
+ itarget = il
+
+ if itarget >+ 0 :
+ print part.rank, itarget, " Before route_togrid = ", self.route_togrid[itarget,:]
grgrid = part.l2glandindex(self.route_togrid[:,:])
+ if itarget >+ 0 :
+ print part.rank, itarget, " After route_togrid = ", self.route_togrid[itarget,:]
rgrid = procgrid.landscatter(grgrid, order='F')
rgrid = rgrid.astype(vtyp, copy=False)
rgrid[rgrid >= RPP.IntFillValue] = NCFillValue
diff --git a/RoutingPreProc.py b/RoutingPreProc.py
index b4e5cb2..63aa0a0 100644
--- a/RoutingPreProc.py
+++ b/RoutingPreProc.py
@@ -130,7 +130,7 @@ print("nbasmax : ", nbasmax)
#
hydrodata = HG.HydroData(hydrogrid.ncfile, hydrogrid.box, sub_index)
-IF.initatmgrid(nbpt, modelgrid)
+IF.initatmgrid(rank, nbcore, nbpt, modelgrid)
@@ -177,3 +177,5 @@ hsuper.fetch()
hgraph = IF.HydroGraph(nbasmax, hsuper)
hgraph.dumpnetcdf(OutGraphFile, gg, modelgrid, part)
+
+IF.closeinterface(comm)
--
GitLab