Commit 06c3008b authored by MEURDESOIF Yann's avatar MEURDESOIF Yann
Browse files

Introduce limited area metric and adapt nudging for forcing field at frontier.

YM
parent f5517719
Pipeline #143541 failed with stages
in 1 minute and 27 seconds
MODULE grid_param
INTEGER :: iim_glo=40
INTEGER :: jjm_glo
INTEGER,PARAMETER :: nb_face=10
INTEGER :: llm=19
INTEGER :: nqtot ! number of tracers handled by advection scheme
INTEGER :: nqdyn ! number of dynamical tracers : 1 if dry, more if moist
......
......@@ -46,12 +46,14 @@ CONTAINS
SUBROUTINE allocate_time_field(field, field_type, data_type, dim3, dim4, name, pname)
USE field_mod, ONLY : allocate_field
USE domain_mod, ONLY : ndomain, assigned_domain
! value of argument pname must reflect the contents of XIOS XML configuration files
TYPE(t_time_field), INTENT(OUT) :: field
INTEGER, INTENT(IN) :: field_type
INTEGER, INTENT(IN) :: data_type
INTEGER, OPTIONAL, INTENT(IN) :: dim3, dim4
CHARACTER(*), OPTIONAL, INTENT(IN) :: name, pname
INTEGER :: ind
field%name = TRIM(name) ! note that name is in fact a compulsory argument
field%is_ps = PRESENT(pname)
......@@ -60,6 +62,14 @@ CONTAINS
CALL allocate_field(field%f_t, field_type, data_type, dim3, dim4, name=TRIM(name)//'_t')
CALL allocate_field(field%f_t0, field_type, data_type, dim3, dim4, name=TRIM(name)//'_t0')
CALL allocate_field(field%f_t1, field_type, data_type, dim3, dim4, name=TRIM(name)//'_t1')
DO ind=1, ndomain
IF(.NOT. assigned_domain(ind)) CYCLE
field%f_t(ind)%rval4d = 0.
field%f_t0(ind)%rval4d = 0.
field%f_t1(ind)%rval4d = 0.
END DO
END SUBROUTINE allocate_time_field
SUBROUTINE update_time_interpolation(time, interp)
......
File mode changed from 100755 to 100644
......@@ -100,6 +100,9 @@ CONTAINS
INTEGER :: l
CHARACTER(len=255) :: rayleigh_friction_key
REAL(rstd) :: mintau
INTEGER :: ip1, im1, jp1, jm1, count
LOGICAL :: flag
REAL(rstd) :: min_u, max_u, min_theta, max_theta
CHARACTER(len=50) :: vert_prof_dissip_type
! New variables added for dissipation vertical profile (SF, 19/09/18)
......@@ -193,8 +196,7 @@ CONTAINS
! set random seed to get reproductibility when using a different number of process
CALL RANDOM_SEED(size=M)
CALL RANDOM_SEED(put=(/(i,i=1,M)/))
! This cannot be ported on GPU due to compiler limitations
u=0
DO j=jj_begin,jj_end
DO i=ii_begin,ii_end
ij=(j-1)*iim+i
......@@ -205,7 +207,8 @@ CONTAINS
CALL RANDOM_NUMBER(r)
u(ij+u_ldown)=r-0.5
ENDDO
ENDDO
ENDDO
! This cannot be ported on GPU due to compiler limitations
ENDDO
!$OMP END MASTER
!$OMP BARRIER
......@@ -214,9 +217,59 @@ CONTAINS
dumax=0
DO iter=1,nitergdiv
DO ind=1,ndomain
IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE
CALL swap_dimensions(ind)
CALL swap_geometry(ind)
u=f_u(ind)
count=0
DO j=jj_begin-1,jj_end+1
DO i=ii_begin-1,ii_end+1
ij=(j-1)*iim+i
ip1=i+1 ; IF (ip1>iim) ip1=iim
im1=i-1 ; IF (im1<1) im1=1
jp1=j+1 ; IF (jp1>jjm) jp1=jjm
jm1=j-1 ; IF (jm1<1) jm1=1
flag = domain(ind)%outside(ip1,j) .OR. domain(ind)%outside(ip1,jm1) .OR. domain(ind)%outside(i,j) &
.OR. domain(ind)%outside(i,jm1) .OR. domain(ind)%outside(i,jp1) &
.OR. domain(ind)%outside(im1,j) .OR. domain(ind)%outside(im1,jp1)
IF (flag) THEN
IF (count<0) PRINT*,flag
count=count+1
u(ij+u_right)=0.
u(ij+u_rup)=0.
u(ij+u_lup)=0.
u(ij+u_left)=0.
u(ij+u_rdown)=0.
u(ij+u_ldown)=0.
ENDIF
ENDDO
ENDDO
ENDDO
CALL update_device_field(f_u)
CALL transfert_request(f_u,req_e1_vect)
DO ind=1,ndomain
IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE
CALL swap_dimensions(ind)
CALL swap_geometry(ind)
u=f_u(ind)
min_u=HUGE(min_u)
max_u=0
DO j=jj_begin,jj_end
DO i=ii_begin,ii_end
ij=(j-1)*iim+i
min_u=MIN(min_u,ABS(u(ij+u_right)),ABS(u(ij+u_rup)),ABS(u(ij+u_lup)),&
ABS(u(ij+u_left)),ABS(u(ij+u_rdown)),ABS(u(ij+u_ldown)))
max_u=MAX(max_u,ABS(u(ij+u_right)),ABS(u(ij+u_rup)),ABS(u(ij+u_lup)),&
ABS(u(ij+u_left)),ABS(u(ij+u_rdown)),ABS(u(ij+u_ldown)))
ENDDO
ENDDO
! PRINT*,"Min U ", min_u, " Max U ", max_u
ENDDO
DO ind=1,ndomain
IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE
CALL swap_dimensions(ind)
......@@ -228,6 +281,7 @@ CONTAINS
!$acc update host(u(:)) wait
du=u
ENDDO
ENDDO
CALL update_device_field(f_du)
......@@ -239,17 +293,38 @@ CONTAINS
CALL swap_dimensions(ind)
CALL swap_geometry(ind)
du=f_du(ind)
max_u=0.
min_u=HUGE(min_u)
! Not ported on GPU because all the other kernels cannot be ported
DO j=jj_begin,jj_end
DO i=ii_begin,ii_end
ij=(j-1)*iim+i
if (le(ij+u_right)>1e-100) dumax=MAX(dumax,ABS(du(ij+u_right)))
if (le(ij+u_lup)>1e-100) dumax=MAX(dumax,ABS(du(ij+u_lup)))
if (le(ij+u_ldown)>1e-100) dumax=MAX(dumax,ABS(du(ij+u_ldown)))
ij=(j-1)*iim+i
ip1=MIN(i+1,iim)
im1=MAX(i-1,1)
jp1=MIN(j+1,jjm)
jm1=MAX(j-1,1)
IF ( domain(ind)%outside(ip1,j) .OR. domain(ind)%outside(ip1,jm1) .OR. domain(ind)%outside(i,j) &
.OR. domain(ind)%outside(i,jm1) .OR. domain(ind)%outside(i,jp1) &
.OR. domain(ind)%outside(im1,j) .OR. domain(ind)%outside(im1,jp1)) THEN
du(ij+u_right)=0
du(ij+u_rup)=0
du(ij+u_lup)=0
du(ij+u_left)=0
du(ij+u_rdown)=0
du(ij+u_ldown)=0
ENDIF
IF (le(ij+u_right)>1e-100) max_u=MAX(max_u,ABS(du(ij+u_right)))
IF (le(ij+u_right)>1e-100) min_u=MIN(min_u,ABS(du(ij+u_right)))
IF (le(ij+u_lup)>1e-100) max_u=MAX(max_u,ABS(du(ij+u_lup)))
IF (le(ij+u_lup)>1e-100) min_u=MIN(min_u,ABS(du(ij+u_lup)))
IF (le(ij+u_ldown)>1e-100) max_u=MAX(max_u,ABS(du(ij+u_ldown)))
IF (le(ij+u_ldown)>1e-100) min_u=MIN(min_u,ABS(du(ij+u_ldown)))
ENDDO
ENDDO
ENDDO
! PRINT*,"Min dU ", min_u, " Max dU ", max_u
dumax=MAX(dumax,max_u)
ENDDO
IF (using_mpi) THEN
CALL reduce_max_omp(dumax,dumax1)
......@@ -289,7 +364,7 @@ CONTAINS
CALL swap_dimensions(ind)
CALL swap_geometry(ind)
u=f_u(ind)
u=0
! set random seed to get reproductibility when using a different number of process
CALL RANDOM_SEED(size=M)
CALL RANDOM_SEED(put=(/(i,i=1,M)/))
......@@ -314,8 +389,57 @@ CONTAINS
dumax=0
DO iter=1,nitergrot
DO ind=1,ndomain
IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE
CALL swap_dimensions(ind)
CALL swap_geometry(ind)
u=f_u(ind)
DO j=jj_begin-1,jj_end+1
DO i=ii_begin-1,ii_end+1
ij=(j-1)*iim+i
ip1=MIN(i+1,iim)
im1=MAX(i-1,1)
jp1=MIN(j+1,jjm)
jm1=MAX(j-1,1)
IF ( domain(ind)%outside(ip1,j) .OR. domain(ind)%outside(ip1,jm1) .OR. domain(ind)%outside(i,j) &
.OR. domain(ind)%outside(i,jm1) .OR. domain(ind)%outside(i,jp1) &
.OR. domain(ind)%outside(im1,j) .OR. domain(ind)%outside(im1,jp1)) THEN
u(ij+u_right)=0
u(ij+u_rup)=0
u(ij+u_lup)=0
u(ij+u_left)=0
u(ij+u_rdown)=0
u(ij+u_ldown)=0
ENDIF
ENDDO
ENDDO
ENDDO
CALL update_device_field(f_u)
CALL transfert_request(f_u,req_e1_vect)
DO ind=1,ndomain
IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE
CALL swap_dimensions(ind)
CALL swap_geometry(ind)
u=f_u(ind)
min_u=HUGE(min_u)
max_u=0
DO j=jj_begin,jj_end
DO i=ii_begin,ii_end
ij=(j-1)*iim+i
min_u=MIN(min_u,ABS(u(ij+u_right)),ABS(u(ij+u_rup)),ABS(u(ij+u_lup)),&
ABS(u(ij+u_left)),ABS(u(ij+u_rdown)),ABS(u(ij+u_ldown)))
max_u=MAX(max_u,ABS(u(ij+u_right)),ABS(u(ij+u_rup)),ABS(u(ij+u_lup)),&
ABS(u(ij+u_left)),ABS(u(ij+u_rdown)),ABS(u(ij+u_ldown)))
ENDDO
ENDDO
! PRINT*,"Min U ", min_u, " Max U ", max_u
ENDDO
DO ind=1,ndomain
IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE
CALL swap_dimensions(ind)
......@@ -340,14 +464,39 @@ CONTAINS
du=f_du(ind)
! Not ported on GPU because all the other kernels cannot be ported
max_u=0.
min_u=HUGE(min_u)
DO j=jj_begin,jj_end
DO i=ii_begin,ii_end
ij=(j-1)*iim+i
if (le(ij+u_right)>1e-100) dumax=MAX(dumax,ABS(du(ij+u_right)))
if (le(ij+u_lup)>1e-100) dumax=MAX(dumax,ABS(du(ij+u_lup)))
if (le(ij+u_ldown)>1e-100) dumax=MAX(dumax,ABS(du(ij+u_ldown)))
ij=(j-1)*iim+i
ip1=MIN(i+1,iim)
im1=MAX(i-1,1)
jp1=MIN(j+1,jjm)
jm1=MAX(j-1,1)
IF ( domain(ind)%outside(ip1,j) .OR. domain(ind)%outside(ip1,jm1) .OR. domain(ind)%outside(i,j) &
.OR. domain(ind)%outside(i,jm1) .OR. domain(ind)%outside(i,jp1) &
.OR. domain(ind)%outside(im1,j) .OR. domain(ind)%outside(im1,jp1)) THEN
du(ij+u_right)=0
du(ij+u_rup)=0
du(ij+u_lup)=0
du(ij+u_left)=0
du(ij+u_rdown)=0
du(ij+u_ldown)=0
ENDIF
IF (le(ij+u_right)>1e-100) max_u=MAX(max_u,ABS(du(ij+u_right)))
IF (le(ij+u_right)>1e-100) min_u=MIN(min_u,ABS(du(ij+u_right)))
IF (le(ij+u_lup)>1e-100) max_u=MAX(max_u,ABS(du(ij+u_lup)))
IF (le(ij+u_lup)>1e-100) min_u=MIN(min_u,ABS(du(ij+u_lup)))
IF (le(ij+u_ldown)>1e-100) max_u=MAX(max_u,ABS(du(ij+u_ldown)))
IF (le(ij+u_ldown)>1e-100) min_u=MIN(min_u,ABS(du(ij+u_ldown)))
ENDDO
ENDDO
! PRINT*,"Min dU ", min_u, " Max dU ", max_u
dumax=MAX(dumax,max_u)
ENDDO
IF (using_mpi) THEN
......@@ -388,7 +537,7 @@ CONTAINS
CALL swap_dimensions(ind)
CALL swap_geometry(ind)
theta=f_theta(ind)
theta=0
! set random seed to get reproductibility when using a different number of process
CALL RANDOM_SEED(size=M)
CALL RANDOM_SEED(put=(/(i,i=1,M)/))
......@@ -409,8 +558,51 @@ CONTAINS
dthetamax=0
DO iter=1,niterdivgrad
DO ind=1,ndomain
IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE
CALL swap_dimensions(ind)
CALL swap_geometry(ind)
theta=f_theta(ind)
DO j=jj_begin-1,jj_end+1
DO i=ii_begin-1,ii_end+1
ij=(j-1)*iim+i
ip1=MIN(i+1,iim)
im1=MAX(i-1,1)
jp1=MIN(j+1,jjm)
jm1=MAX(j-1,1)
IF ( domain(ind)%outside(ip1,j) .OR. domain(ind)%outside(ip1,jm1) .OR. domain(ind)%outside(i,j) &
.OR. domain(ind)%outside(i,jm1) .OR. domain(ind)%outside(i,jp1) &
.OR. domain(ind)%outside(im1,j) .OR. domain(ind)%outside(im1,jp1)) THEN
theta(ij)=0
ENDIF
ENDDO
ENDDO
ENDDO
CALL update_device_field(f_theta)
CALL transfert_request(f_theta,req_i1)
DO ind=1,ndomain
IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE
CALL swap_dimensions(ind)
CALL swap_geometry(ind)
theta=f_theta(ind)
min_theta=HUGE(min_theta)
max_theta=0
DO j=jj_begin,jj_end
DO i=ii_begin,ii_end
ij=(j-1)*iim+i
min_theta=MIN(min_theta,theta(ij))
max_theta=MIN(max_theta,theta(ij))
ENDDO
ENDDO
! PRINT*,"Min Theta ", min_theta, " Max Theta ", max_theta
ENDDO
DO ind=1,ndomain
IF (.NOT. assigned_domain(ind) .OR. .NOT. is_omp_level_master) CYCLE
CALL swap_dimensions(ind)
......@@ -435,12 +627,27 @@ CONTAINS
dtheta=f_dtheta(ind)
! Not ported on GPU because all the other kernels cannot be ported
min_theta=HUGE(min_theta)
max_theta=0
DO j=jj_begin,jj_end
DO i=ii_begin,ii_end
ij=(j-1)*iim+i
dthetamax=MAX(dthetamax,ABS(dtheta(ij)))
ip1=MIN(i+1,iim)
im1=MAX(i-1,1)
jp1=MIN(j+1,jjm)
jm1=MAX(j-1,1)
IF ( domain(ind)%outside(ip1,j) .OR. domain(ind)%outside(ip1,jm1) .OR. domain(ind)%outside(i,j) &
.OR. domain(ind)%outside(i,jm1) .OR. domain(ind)%outside(i,jp1) &
.OR. domain(ind)%outside(im1,j) .OR. domain(ind)%outside(im1,jp1)) THEN
dtheta(ij)=0
ENDIF
min_theta=MIN(min_theta,ABS(dtheta(ij)))
max_theta=MAX(max_theta,ABS(dtheta(ij)))
ENDDO
ENDDO
! PRINT*,"Min dTheta ", min_theta, " Max dTheta ", max_theta
dthetamax=MAX(dthetamax,max_theta)
ENDDO
IF (using_mpi) THEN
......
......@@ -32,13 +32,14 @@ CONTAINS
END SUBROUTINE init_guided
SUBROUTINE guided(tt, f_ps, f_theta_rhodz, f_u, f_q)
SUBROUTINE guided(tt, f_ps, f_mass, f_theta_rhodz, f_u, f_q)
USE icosa
USE guided_ncar_mod, ONLY : guided_ncar => guided
USE nudging_mod, ONLY : guided_nudging => guided
IMPLICIT NONE
REAL(rstd), INTENT(IN):: tt
TYPE(t_field),POINTER :: f_ps(:)
TYPE(t_field),POINTER :: f_mass(:)
TYPE(t_field),POINTER :: f_theta_rhodz(:)
TYPE(t_field),POINTER :: f_u(:)
TYPE(t_field),POINTER :: f_q(:)
......@@ -48,13 +49,30 @@ CONTAINS
CASE ('dcmip1')
CALL guided_ncar(tt, f_ps, f_theta_rhodz, f_u, f_q)
CASE ('nudging')
CALL guided_nudging(tt, f_ps, f_theta_rhodz, f_u, f_q)
CALL guided_nudging(tt, f_ps, f_mass, f_theta_rhodz, f_u, f_q)
CASE DEFAULT
PRINT*,"Bad selector for varaible guided_type >",TRIM(guided_type),"> option are <none>, <dcmip1>"
STOP
END SELECT
END SUBROUTINE guided
SUBROUTINE update_border(tt, f_ps, f_mass, f_theta_rhodz, f_u, f_q)
USE icosa
USE metric
USE nudging_mod, ONLY : up_border => update_border
IMPLICIT NONE
REAL(rstd), INTENT(IN):: tt
TYPE(t_field),POINTER :: f_ps(:)
TYPE(t_field),POINTER :: f_mass(:)
TYPE(t_field),POINTER :: f_theta_rhodz(:)
TYPE(t_field),POINTER :: f_u(:)
TYPE(t_field),POINTER :: f_q(:)
IF (metric_type == metric_icosa_area) CALL up_border(tt, f_ps, f_mass, f_theta_rhodz, f_u, f_q)
END SUBROUTINE update_border
END MODULE guided_mod
This diff is collapsed.
......@@ -16,7 +16,7 @@ CONTAINS
USE earth_const, ONLY : init_earth_const
USE grid_param, ONLY : init_grid_param
USE omp_para, ONLY : init_omp_para, switch_omp_no_distrib_level
USE metric, ONLY : compute_metric
USE metric, ONLY : init_metric
USE domain_mod, ONLY : compute_domain
USE transfert_mod, ONLY : init_transfert
USE write_field_mod, ONLY : init_writefield, writefield, close_files
......@@ -41,7 +41,7 @@ CONTAINS
CALL init_earth_const
CALL init_grid_param(is_mpi_master)
CALL init_omp_para(is_mpi_master)
CALL compute_metric
CALL init_metric
CALL compute_domain ! allocates domain_glo(ndomain_glo), domain(ndomain)
......
File mode changed from 100755 to 100644
......@@ -344,6 +344,7 @@ CONTAINS
SUBROUTINE write_restart_field(field,fieldId,ncid)
USE prec
USE metric
USE grid_param
USE field_mod
USE domain_mod
USE netcdf_mod
......@@ -512,6 +513,7 @@ CONTAINS
field10,field11,field12,field13,field14,field15,field16,field17,field18,field19 )
USE prec
USE metric
USE grid_param
USE field_mod
USE domain_mod
USE netcdf_mod
......@@ -678,6 +680,7 @@ CONTAINS
SUBROUTINE read_start_field(field,fieldId,ncid)
USE prec
USE metric
USE grid_param
USE field_mod
USE domain_mod
USE netcdf_mod
......
......@@ -46,18 +46,21 @@ CONTAINS
USE geometry
USE mpi_mod
USE time_mod
USE metric, ONLY : vup,vdown, cell_glo
USE metric, ONLY : vup,vdown, cell_glo, ncell_glo_tot=>ncell_glo
USE tracer_icosa_mod
IMPLICIT NONE
TYPE(xios_context) :: ctx_hdl
TYPE(xios_duration) :: dtime
REAL(rstd) :: lev_value(llm)
REAL(rstd) :: lev_valuep1(llm+1)
REAL(rstd) :: nq_value(nqtot)
INTEGER :: ncell, ncell_tot, ncell_glo(0:mpi_size-1), displ
INTEGER :: ncell, ncell_tot, ncell_mpi(0:mpi_size-1), displ(0:mpi_size-1), ncell_glo
INTEGER :: ind, i,j,k,l,ij
REAL(rstd),ALLOCATABLE :: lon(:), lat(:), bounds_lon(:,:), bounds_lat(:,:)
INTEGER, ALLOCATABLE :: ind_glo(:)
TYPE(t_domain),POINTER :: d
INTEGER,ALLOCATABLE :: cell_glo_tot(:)
INTEGER,ALLOCATABLE :: ind_glo_tot(:)
CALL xios_set_context
!$OMP BARRIER
......@@ -84,14 +87,14 @@ CONTAINS
ENDDO
ncell_i=ncell
CALL MPI_ALLGATHER(ncell,1,MPI_INTEGER,ncell_glo,1,MPI_INTEGER,comm_icosa,ierr)
CALL MPI_ALLGATHER(ncell,1,MPI_INTEGER,ncell_mpi,1,MPI_INTEGER,comm_icosa,ierr)
displ=0
DO i=1,mpi_rank
displ=displ+ncell_glo(i-1)
displ(0)=0
DO i=1,mpi_size-1
displ(i)=displ(i-1)+ncell_mpi(i-1)
ENDDO
ncell_tot=sum(ncell_glo(:))
ncell_tot=sum(ncell_mpi(:))
ALLOCATE(lon(ncell), lat(ncell), bounds_lon(0:5,ncell), bounds_lat(0:5,ncell), ind_glo(ncell))
......@@ -116,12 +119,31 @@ CONTAINS
ENDDO
ENDDO
ENDDO
ALLOCATE(ind_glo_tot(ncell_tot))
ALLOCATE(cell_glo_tot(0:ncell_glo_tot-1))
cell_glo_tot(:)= -1
CALL MPI_ALLGATHERV(ind_glo, ncell, MPI_INTEGER, ind_glo_tot, ncell_mpi, displ, MPI_INTEGER, comm_icosa,ierr)
ncell_glo=0
DO i=1,ncell_tot
IF (cell_glo_tot(ind_glo_tot(i))==-1) THEN
cell_glo_tot(ind_glo_tot(i))=ncell_glo
ncell_glo=ncell_glo+1
ENDIF
ENDDO
DO i=1,ncell
ind_glo(i)=cell_glo_tot(ind_glo(i))
ENDDO
CALL xios_set_domaingroup_attr("i",ni_glo=ncell_tot, ibegin=displ, ni=ncell)
CALL xios_set_domaingroup_attr("i",ni_glo=ncell_glo, ibegin=displ(mpi_rank), ni=ncell)
CALL xios_set_domaingroup_attr("i", data_dim=1, type='unstructured' , nvertex=6, i_index=ind_glo)
CALL xios_set_domaingroup_attr("i",lonvalue_1d=lon, latvalue_1d=lat, bounds_lon_1d=bounds_lon, bounds_lat_1d=bounds_lat)
DEALLOCATE(lon, lat, bounds_lon, bounds_lat,ind_glo)
DEALLOCATE(lon, lat, bounds_lon, bounds_lat,ind_glo, ind_glo_tot, cell_glo_tot)
......@@ -142,12 +164,12 @@ CONTAINS
ENDDO
ncell_e=ncell
CALL MPI_ALLGATHER(ncell_e,1,MPI_INTEGER,ncell_glo,1,MPI_INTEGER,comm_icosa,ierr)
displ=0
DO i=1,mpi_rank
displ=displ+ncell_glo(i-1)
CALL MPI_ALLGATHER(ncell_e,1,MPI_INTEGER,ncell_mpi,1,MPI_INTEGER,comm_icosa,ierr)
displ(0)=0
DO i=1,mpi_size-1
displ(i)=displ(i-1)+ncell_mpi(i-1)
ENDDO
ncell_tot=sum(ncell_glo(:))
ncell_tot=sum(ncell_mpi(:))
ALLOCATE(lon(ncell), lat(ncell), bounds_lon(0:1,ncell), bounds_lat(0:1,ncell),ind_glo(ncell))