Changeset 879 for codes/icosagcm/devel/src/unstructured

Timestamp:

06/04/19 18:03:54 (5 years ago)

Author:

dubos

Message:

devel : introduced derived type to store cell bounds

File:

• codes/icosagcm/devel/src/unstructured/init_unstructured.f90 (modified) (6 diffs)

codes/icosagcm/devel/src/unstructured/init_unstructured.f90

@@ -2 +2 @@
   USE mpipara, ONLY : is_mpi_master
   USE data_unstructured_mod
-  USE geometry, ONLY : de
+  USE geometry, ONLY : swap_geometry
   IMPLICIT NONE
   SAVE
@@ -128 +128 @@
 
 
+  SUBROUTINE read_field(id_nc, field)
+    USE field_mod, ONLY : t_field, type_integer, type_real
+    INTEGER :: id_nc
+    TYPE(t_field), POINTER :: field(:), fld
+    fld=>field(1)
+    SELECT CASE(fld%data_type)
+    CASE(type_integer)
+       CALL read_from_gridfile(id_nc, 'integer', TRIM(fld%name))
+       SELECT CASE(field(1)%ndim)
+       CASE(2)
+          fld%ival2d = Idata_read1
+       CASE(3)
+          fld%ival3d = Idata_read2
+       END SELECT
+    CASE(type_real)
+       CALL read_from_gridfile(id_nc, 'float', TRIM(fld%name))
+       SELECT CASE(field(1)%ndim)
+       CASE(2)
+          fld%rval2d = Ddata_read1
+       CASE(3)
+          fld%rval3d = Ddata_read2
+       END SELECT
+    END SELECT
+  END SUBROUTINE read_field
+
   SUBROUTINE read_local_mesh
     USE field_mod
-    USE geometry, ONLY : lon_i, lat_i, lon_e, lat_e, ep_e
+    USE geometry, ONLY : geom, lon_i, lat_i, lon_e, lat_e, ep_e
+    USE domain_mod, ONLY : domain, domain_glo, t_domain, swap_needed
     IMPLICIT NONE
     INTEGER :: ij
+    INTEGER, ALLOCATABLE :: cell_ij(:)
     REAL(rstd), ALLOCATABLE :: angle_e(:)
     REAL(rstd) :: clon, slon, clat, slat, & ! COS/SIN of lon/lat
          elon(3), elat(3) ! lon/lat unit vectors
     
-    !status = nf90_get_att(id_nc, nf90_global, "primal_num", primal_num)
-    !status = nf90_get_att(id_nc, nf90_global, "dual_num", dual_num)
-    !status = nf90_get_att(id_nc, nf90_global, "edge_num", edge_num)
-    !print *,primal_num,dual_num,edge_num
-
+    ! if possible data is read into fields allocated in geom(1) using read_field
+    CALL read_field(id_nc, geom%Ai)
+    CALL read_field(id_nc, geom%Av)
+    CALL read_field(id_nc, geom%fv)
+    CALL read_field(id_nc, geom%le_de)
+    CALL read_field(id_nc, geom%le)
+    CALL read_field(id_nc, geom%de)
+    CALL read_field(id_nc, geom%lon_i)
+    CALL read_field(id_nc, geom%lat_i)
+    CALL read_field(id_nc, geom%lon_e)
+    CALL read_field(id_nc, geom%lat_e)
+
+    swap_needed = .TRUE.
+    CALL swap_geometry(1)
+
+    ! cell centers and bounds are read into domain(1)%primal_own
+    domain(1)%primal_own%ncell = primal_num
+    domain_glo(1)%primal_own%ncell = primal_num
+    ALLOCATE(domain(1)%primal_own%lon,     source=lon_i)
+    ALLOCATE(domain_glo(1)%primal_own%lon, source=lon_i)
+    ALLOCATE(domain(1)%primal_own%lat,     source=lat_i)
+    ALLOCATE(domain_glo(1)%primal_own%lat, source=lat_i)
+
+    CALL read_from_gridfile(id_nc, 'float', 'primal_bounds_lon')
+    ALLOCATE(domain(1)%primal_own%bnds_lon,     source=Ddata_read2)
+    ALLOCATE(domain_glo(1)%primal_own%bnds_lon, source=Ddata_read2)
+    CALL read_from_gridfile(id_nc, 'float', 'primal_bounds_lat')
+    ALLOCATE(domain(1)%primal_own%bnds_lat,     source=Ddata_read2)
+    ALLOCATE(domain_glo(1)%primal_own%bnds_lat, source=Ddata_read2)
+    PRINT *, 'read_local_mesh : primal_num =', primal_num, domain_glo(1)%primal_own%ncell
+
+    ! other fields are defined only in the case of an unstructured mesh
+    ! those fields are read into simple arrays
     CALL read_from_gridfile(id_nc, 'integer', 'primal_deg')
     ALLOCATE(primal_deg, source = Idata_read1)
@@ -148 +203 @@
     CALL read_from_gridfile(id_nc, 'integer', 'primal_ne')
     ALLOCATE(primal_ne, source = Idata_read2)
+    CALL read_from_gridfile(id_nc, 'integer', 'primal_vertex') 
+    ALLOCATE(primal_vertex, source = Idata_read2)
+
     CALL read_from_gridfile(id_nc, 'integer', 'dual_deg')
     ALLOCATE(dual_deg, source = Idata_read1)
@@ -154 +212 @@
     CALL read_from_gridfile(id_nc, 'integer', 'dual_ne')
     ALLOCATE(dual_ne, source = Idata_read2)
-    CALL read_from_gridfile(id_nc, 'integer', 'primal_vertex') 
-    ALLOCATE(primal_vertex, source = Idata_read2)
+    CALL read_from_gridfile(id_nc, 'integer', 'dual_vertex') 
+    ALLOCATE(dual_vertex, source = Idata_read2)
+
     CALL read_from_gridfile(id_nc, 'integer', 'left')
     ALLOCATE(left, source = Idata_read1)
@@ -164 +223 @@
     CALL read_from_gridfile(id_nc, 'integer', 'down')
     ALLOCATE(down, source = Idata_read1)
+    CALL read_from_gridfile(id_nc, 'float', 'angle_e')
+    ALLOCATE(angle_e, source = Ddata_read1)
     CALL read_from_gridfile(id_nc, 'integer', 'trisk_deg')
     ALLOCATE(trisk_deg, source = Idata_read1)
     CALL read_from_gridfile(id_nc, 'integer', 'trisk')
     ALLOCATE(trisk, source = Idata_read2)
-    CALL read_from_gridfile(id_nc, 'float', 'Ai')
-    ALLOCATE(Ai, source = Ddata_read1)
-    CALL read_from_gridfile(id_nc, 'float', 'Av')
-    ALLOCATE(Av, source = Ddata_read1)
-    CALL read_from_gridfile(id_nc, 'float', 'le_de')
-    ALLOCATE(le_de, source = Ddata_read1)
-    CALL read_from_gridfile(id_nc, 'float', 'le')
-    ALLOCATE(le, source = Ddata_read1)
-    CALL read_from_gridfile(id_nc, 'float', 'de')
-    ALLOCATE(de, source = Ddata_read1)
     CALL read_from_gridfile(id_nc, 'float', 'Riv2')
     ALLOCATE(Riv2, source = Ddata_read2)
     CALL read_from_gridfile(id_nc, 'float', 'wee')
     ALLOCATE(wee, source = Ddata_read2)
-    CALL read_from_gridfile(id_nc, 'float', 'fv') ! 
-    ALLOCATE(fv, source = Ddata_read1)
-    CALL read_from_gridfile(id_nc, 'integer', 'dual_vertex') 
-    ALLOCATE(dual_vertex, source = Idata_read2)
-
-    CALL read_from_gridfile(id_nc, 'float', 'lon_i')
-    ALLOCATE(lon_i, source = Ddata_read1)
-    CALL read_from_gridfile(id_nc, 'float', 'lat_i')
-    ALLOCATE(lat_i, source = Ddata_read1)
-    CALL read_from_gridfile(id_nc, 'float', 'lon_e')
-    ALLOCATE(lon_e, source = Ddata_read1)
-    CALL read_from_gridfile(id_nc, 'float', 'lat_e')
-    ALLOCATE(lat_e, source = Ddata_read1)
-    CALL read_from_gridfile(id_nc, 'float', 'angle_e')
-    ALLOCATE(angle_e, source = Ddata_read1)
 
     CALL free_data_read ! free buffers after reading all data from grid file
 
-    edge_num = SIZE(le_de)
-    primal_num = SIZE(Ai)
-    dual_num = SIZE(Av)
+
+!    edge_num = SIZE(le_de)
+!    primal_num = SIZE(Ai)
+!    dual_num = SIZE(Av)
     max_primal_deg = SIZE(primal_edge,1)
     max_dual_deg = SIZE(dual_edge,1)
@@ -221 +258 @@
     
     DEALLOCATE(angle_e)
+
+    ALLOCATE(domain(1)%primal_own%ij(primal_num))
+    domain(1)%primal_own%ij(:) = (/(ij,ij=1,primal_num)/)
+    PRINT *, 'read_local_mesh :', SHAPE(domain), primal_num, &
+         SHAPE(domain(1)%primal_own%ij)
+
+    ALLOCATE(cell_ij(dual_num))
+    cell_ij(:) = (/(ij,ij=1,dual_num)/)
+    ALLOCATE(domain(1)%dual_own%ij, source=cell_ij)
+    DEALLOCATE(cell_ij)
+
+    CALL swap_geometry(1)
+    swap_needed = .FALSE.
   END SUBROUTINE read_local_mesh