[12] | 1 | MODULE domain_mod |
---|
| 2 | USE domain_param |
---|
[879] | 3 | IMPLICIT NONE |
---|
| 4 | SAVE |
---|
[12] | 5 | |
---|
[879] | 6 | TYPE t_cellset |
---|
| 7 | ! derived type containing arrays for lat,lon of center and bounds of cells of primal or dual mesh |
---|
| 8 | ! this information is needed to write CF-compliant NetCDF files |
---|
| 9 | ! and by XIOS |
---|
| 10 | INTEGER :: ncell |
---|
| 11 | INTEGER, ALLOCATABLE :: ij(:), & ! ij(n) = local index of cell n |
---|
| 12 | ind_glo(:), & ! ind_glo(n) = global index of cell n |
---|
| 13 | sgn(:) ! sign to apply before writing to disk (velocities) |
---|
| 14 | REAL, ALLOCATABLE :: lon(:), lat(:), bnds_lon(:,:), bnds_lat(:,:) |
---|
| 15 | END type t_cellset |
---|
| 16 | |
---|
[12] | 17 | TYPE t_domain |
---|
| 18 | INTEGER :: face |
---|
| 19 | INTEGER :: iim |
---|
| 20 | INTEGER :: jjm |
---|
| 21 | INTEGER :: ii_begin |
---|
| 22 | INTEGER :: jj_begin |
---|
| 23 | INTEGER :: ii_end |
---|
| 24 | INTEGER :: jj_end |
---|
| 25 | INTEGER :: ii_nb |
---|
| 26 | INTEGER :: jj_nb |
---|
| 27 | INTEGER :: ii_begin_glo |
---|
| 28 | INTEGER :: jj_begin_glo |
---|
| 29 | INTEGER :: ii_end_glo |
---|
| 30 | INTEGER :: jj_end_glo |
---|
| 31 | INTEGER,POINTER :: assign_domain(:,:) |
---|
[266] | 32 | INTEGER,POINTER :: assign_cell_glo(:,:) |
---|
[12] | 33 | INTEGER,POINTER :: assign_i(:,:) |
---|
| 34 | INTEGER,POINTER :: assign_j(:,:) |
---|
[21] | 35 | INTEGER,POINTER :: edge_assign_domain(:,:,:) |
---|
| 36 | INTEGER,POINTER :: edge_assign_i(:,:,:) |
---|
| 37 | INTEGER,POINTER :: edge_assign_j(:,:,:) |
---|
| 38 | INTEGER,POINTER :: edge_assign_pos(:,:,:) |
---|
[146] | 39 | INTEGER,POINTER :: edge_assign_sign(:,:,:) |
---|
[726] | 40 | INTEGER,POINTER :: vertex_assign_domain(:,:,:) |
---|
| 41 | INTEGER,POINTER :: vertex_assign_i(:,:,:) |
---|
| 42 | INTEGER,POINTER :: vertex_assign_j(:,:,:) |
---|
| 43 | INTEGER,POINTER :: vertex_assign_pos(:,:,:) |
---|
[12] | 44 | REAL,POINTER :: xyz(:,:,:) |
---|
| 45 | REAL,POINTER :: neighbour(:,:,:,:) |
---|
| 46 | INTEGER,POINTER :: delta(:,:) |
---|
| 47 | REAL,POINTER :: vertex(:,:,:,:) |
---|
| 48 | LOGICAL,POINTER :: own(:,:) |
---|
| 49 | INTEGER,POINTER :: ne(:,:,:) |
---|
| 50 | |
---|
| 51 | INTEGER :: t_right |
---|
| 52 | INTEGER :: t_rup |
---|
| 53 | INTEGER :: t_lup |
---|
| 54 | INTEGER :: t_left |
---|
| 55 | INTEGER :: t_ldown |
---|
| 56 | INTEGER :: t_rdown |
---|
| 57 | |
---|
| 58 | INTEGER :: u_right |
---|
| 59 | INTEGER :: u_rup |
---|
| 60 | INTEGER :: u_lup |
---|
| 61 | INTEGER :: u_left |
---|
| 62 | INTEGER :: u_ldown |
---|
| 63 | INTEGER :: u_rdown |
---|
| 64 | |
---|
| 65 | INTEGER :: z_rup |
---|
| 66 | INTEGER :: z_up |
---|
| 67 | INTEGER :: z_lup |
---|
| 68 | INTEGER :: z_ldown |
---|
| 69 | INTEGER :: z_down |
---|
| 70 | INTEGER :: z_rdown |
---|
| 71 | |
---|
| 72 | INTEGER :: t_pos(6) |
---|
| 73 | INTEGER :: u_pos(6) |
---|
| 74 | INTEGER :: z_pos(6) |
---|
[879] | 75 | |
---|
| 76 | TYPE(t_cellset), POINTER :: primal_own, dual_own, edge_own, & |
---|
| 77 | primal_all, dual_all, edge_all ! halo included, used for debug by writefield |
---|
[12] | 78 | END TYPE t_domain |
---|
| 79 | |
---|
[879] | 80 | TYPE t_mesh |
---|
| 81 | INTEGER :: ndomain, primal_num, edge_num, dual_num |
---|
| 82 | TYPE(t_domain), ALLOCATABLE :: domain(:) |
---|
| 83 | TYPE(t_cellset), ALLOCATABLE :: primal_own(:), dual_own(:), edge_own(:), & |
---|
| 84 | primal_all(:), dual_all(:), edge_all(:) |
---|
| 85 | END type t_mesh |
---|
[186] | 86 | |
---|
[879] | 87 | INTEGER :: ndomain |
---|
| 88 | INTEGER :: ndomain_glo |
---|
| 89 | |
---|
[856] | 90 | LOGICAL :: swap_needed=.TRUE. ! .FALSE. if a thread always computes on the same domain |
---|
| 91 | !$OMP THREADPRIVATE(swap_needed) |
---|
| 92 | |
---|
[879] | 93 | TYPE(t_mesh), TARGET :: mesh_loc, mesh_glo |
---|
| 94 | TYPE(t_domain), POINTER :: domain(:), domain_glo(:) |
---|
[12] | 95 | |
---|
[879] | 96 | INTEGER,ALLOCATABLE :: domglo_rank(:) ! size : ndomain_glo : mpi rank assigned to a domain |
---|
| 97 | INTEGER,ALLOCATABLE :: domglo_loc_ind(:) ! size : ndomain_glo : corresponding local indice for a global domain indice |
---|
| 98 | INTEGER,ALLOCATABLE :: domloc_glo_ind(:) ! size : ndomain : corresponding global indice for a local domain indice |
---|
[186] | 99 | |
---|
[879] | 100 | LOGICAL, ALLOCATABLE :: assigned_domain(:) ! size : ndomain : is an omp task is assigned to this domain |
---|
[186] | 101 | !$OMP THREADPRIVATE(assigned_domain) |
---|
| 102 | |
---|
[12] | 103 | CONTAINS |
---|
| 104 | |
---|
[879] | 105 | SUBROUTINE allocate_mesh(ndom, mesh) |
---|
| 106 | INTEGER, INTENT(IN) :: ndom |
---|
| 107 | TYPE(t_mesh), TARGET, INTENT(OUT) :: mesh |
---|
| 108 | INTEGER :: ind |
---|
| 109 | mesh%ndomain = ndom |
---|
| 110 | ALLOCATE(mesh%domain(ndom)) |
---|
| 111 | ALLOCATE(mesh%primal_own(ndom)) |
---|
| 112 | ALLOCATE(mesh%primal_all(ndom)) |
---|
| 113 | ALLOCATE(mesh%dual_own(ndom)) |
---|
| 114 | ALLOCATE(mesh%dual_all(ndom)) |
---|
| 115 | ALLOCATE(mesh%edge_own(ndom)) |
---|
| 116 | ALLOCATE(mesh%edge_all(ndom)) |
---|
| 117 | DO ind=1, ndom |
---|
| 118 | mesh%domain(ind)%primal_own => mesh%primal_own(ind) |
---|
| 119 | mesh%domain(ind)%primal_all => mesh%primal_all(ind) |
---|
| 120 | mesh%domain(ind)%dual_own => mesh%dual_own(ind) |
---|
| 121 | mesh%domain(ind)%dual_all => mesh%dual_all(ind) |
---|
| 122 | mesh%domain(ind)%edge_own => mesh%edge_own(ind) |
---|
| 123 | mesh%domain(ind)%edge_all => mesh%edge_all(ind) |
---|
| 124 | END DO |
---|
| 125 | END SUBROUTINE allocate_mesh |
---|
| 126 | |
---|
[12] | 127 | SUBROUTINE create_domain |
---|
| 128 | USE grid_param |
---|
[151] | 129 | USE mpipara |
---|
[491] | 130 | USE ioipsl |
---|
[12] | 131 | INTEGER :: ind,nf,ni,nj,i,j |
---|
| 132 | INTEGER :: quotient, rest |
---|
[174] | 133 | INTEGER :: halo_i,halo_j |
---|
| 134 | |
---|
[12] | 135 | TYPE(t_domain),POINTER :: d |
---|
| 136 | |
---|
[26] | 137 | ndomain_glo=nsplit_i*nsplit_j*nb_face |
---|
[879] | 138 | CALL allocate_mesh(ndomain_glo, mesh_glo) |
---|
| 139 | domain_glo => mesh_glo%domain |
---|
[26] | 140 | ALLOCATE(domglo_rank(ndomain_glo)) |
---|
| 141 | ALLOCATE(domglo_loc_ind(ndomain_glo)) |
---|
[174] | 142 | |
---|
| 143 | halo_i=0 |
---|
| 144 | CALL getin("halo_i",halo_i) |
---|
| 145 | halo_j=1 |
---|
| 146 | CALL getin("halo_j",halo_j) |
---|
[26] | 147 | |
---|
[12] | 148 | ind=0 |
---|
| 149 | DO nf=1,nb_face |
---|
[21] | 150 | DO nj=1,nsplit_j |
---|
| 151 | DO ni=1,nsplit_i |
---|
[12] | 152 | ind=ind+1 |
---|
[26] | 153 | d=>domain_glo(ind) |
---|
[12] | 154 | quotient=(iim_glo/nsplit_i) |
---|
| 155 | rest=MOD(iim_glo,nsplit_i) |
---|
| 156 | IF (ni-1 < rest) THEN |
---|
| 157 | d%ii_nb=quotient+1 |
---|
| 158 | d%ii_begin_glo=(quotient+1)*(ni-1)+1 |
---|
| 159 | ELSE |
---|
| 160 | d%ii_nb=quotient |
---|
| 161 | d%ii_begin_glo=(quotient+1)*rest+(ni-1-rest) * quotient+1 |
---|
| 162 | ENDIF |
---|
| 163 | d%ii_end_glo=d%ii_begin_glo+d%ii_nb-1 |
---|
[21] | 164 | |
---|
[151] | 165 | IF (ni/=1) THEN |
---|
[21] | 166 | d%ii_nb=d%ii_nb+1 |
---|
[151] | 167 | d%ii_begin_glo=d%ii_begin_glo-1 |
---|
[21] | 168 | ENDIF |
---|
| 169 | |
---|
[12] | 170 | |
---|
| 171 | quotient=(jjm_glo/nsplit_j) |
---|
| 172 | rest=MOD(jjm_glo,nsplit_j) |
---|
| 173 | IF (nj-1 < rest) THEN |
---|
| 174 | d%jj_nb=quotient+1 |
---|
| 175 | d%jj_begin_glo=(quotient+1)*(nj-1)+1 |
---|
| 176 | ELSE |
---|
| 177 | d%jj_nb=quotient |
---|
| 178 | d%jj_begin_glo=(quotient+1)*rest+(nj-1-rest) * quotient+1 |
---|
| 179 | ENDIF |
---|
[21] | 180 | d%jj_end_glo=d%jj_begin_glo+d%jj_nb-1 |
---|
[12] | 181 | |
---|
[151] | 182 | IF (nj/=1) THEN |
---|
[21] | 183 | d%jj_nb=d%jj_nb+1 |
---|
[151] | 184 | d%jj_begin_glo=d%jj_begin_glo-1 |
---|
[21] | 185 | ENDIF |
---|
| 186 | |
---|
| 187 | |
---|
[174] | 188 | d%iim=d%ii_nb+2*halo + halo_i*2 |
---|
| 189 | d%jjm=d%jj_nb+2*halo + halo_j*2 |
---|
| 190 | d%ii_begin=halo+1 + halo_i |
---|
| 191 | d%jj_begin=halo+1 + halo_j |
---|
[12] | 192 | d%ii_end=d%ii_begin+d%ii_nb-1 |
---|
| 193 | d%jj_end=d%jj_begin+d%jj_nb-1 |
---|
| 194 | d%face=nf |
---|
| 195 | ALLOCATE(d%assign_domain(d%iim,d%jjm)) |
---|
[266] | 196 | ALLOCATE(d%assign_cell_glo(d%iim,d%jjm)) |
---|
[12] | 197 | ALLOCATE(d%assign_i(d%iim,d%jjm)) |
---|
| 198 | ALLOCATE(d%assign_j(d%iim,d%jjm)) |
---|
[21] | 199 | ALLOCATE(d%edge_assign_domain(0:5,d%iim,d%jjm)) |
---|
| 200 | ALLOCATE(d%edge_assign_i(0:5,d%iim,d%jjm)) |
---|
| 201 | ALLOCATE(d%edge_assign_j(0:5,d%iim,d%jjm)) |
---|
| 202 | ALLOCATE(d%edge_assign_pos(0:5,d%iim,d%jjm)) |
---|
[146] | 203 | ALLOCATE(d%edge_assign_sign(0:5,d%iim,d%jjm)) |
---|
[726] | 204 | ALLOCATE(d%vertex_assign_domain(0:5,d%iim,d%jjm)) |
---|
| 205 | ALLOCATE(d%vertex_assign_i(0:5,d%iim,d%jjm)) |
---|
| 206 | ALLOCATE(d%vertex_assign_j(0:5,d%iim,d%jjm)) |
---|
| 207 | ALLOCATE(d%vertex_assign_pos(0:5,d%iim,d%jjm)) |
---|
[12] | 208 | ALLOCATE(d%delta(d%iim,d%jjm)) |
---|
| 209 | ALLOCATE(d%xyz(3,d%iim,d%jjm)) |
---|
| 210 | ALLOCATE(d%vertex(3,0:5,d%iim,d%jjm)) |
---|
| 211 | ALLOCATE(d%neighbour(3,0:5,d%iim,d%jjm)) |
---|
| 212 | ALLOCATE(d%own(d%iim,d%jjm)) |
---|
| 213 | ALLOCATE(d%ne(0:5,d%iim,d%jjm)) |
---|
[151] | 214 | |
---|
| 215 | IF (is_mpi_root) PRINT *,"Domain ",ind," : size ",d%ii_nb,"x",d%jj_nb |
---|
| 216 | |
---|
[12] | 217 | END DO |
---|
| 218 | END DO |
---|
| 219 | END DO |
---|
| 220 | |
---|
| 221 | END SUBROUTINE create_domain |
---|
| 222 | |
---|
[26] | 223 | SUBROUTINE copy_domain(d1,d2) |
---|
| 224 | INTEGER :: face |
---|
| 225 | TYPE(t_domain),TARGET,INTENT(IN) :: d1 |
---|
| 226 | TYPE(t_domain), INTENT(OUT) :: d2 |
---|
[12] | 227 | |
---|
[26] | 228 | d2%iim = d1%iim |
---|
| 229 | d2%jjm = d1%jjm |
---|
| 230 | d2%ii_begin = d1%ii_begin |
---|
| 231 | d2%jj_begin = d1%jj_begin |
---|
| 232 | d2%ii_end = d1%ii_end |
---|
| 233 | d2%jj_end = d1%jj_end |
---|
| 234 | d2%ii_nb = d1%ii_nb |
---|
| 235 | d2%jj_nb = d1%jj_nb |
---|
| 236 | d2%ii_begin_glo = d1%ii_begin_glo |
---|
| 237 | d2%jj_begin_glo = d1%jj_begin_glo |
---|
| 238 | d2%ii_end_glo = d1%ii_end_glo |
---|
| 239 | d2%jj_end_glo = d1%jj_end_glo |
---|
| 240 | d2%assign_domain => d1%assign_domain |
---|
[266] | 241 | d2%assign_cell_glo => d1%assign_cell_glo |
---|
[26] | 242 | d2%assign_i => d1%assign_i |
---|
| 243 | d2%assign_j => d1%assign_j |
---|
| 244 | d2%edge_assign_domain => d1%edge_assign_domain |
---|
| 245 | d2%edge_assign_i => d1%edge_assign_i |
---|
| 246 | d2%edge_assign_j => d1%edge_assign_j |
---|
| 247 | d2%edge_assign_pos => d1%edge_assign_pos |
---|
[146] | 248 | d2%edge_assign_sign => d1%edge_assign_sign |
---|
[726] | 249 | d2%vertex_assign_domain => d1%vertex_assign_domain |
---|
| 250 | d2%vertex_assign_i => d1%vertex_assign_i |
---|
| 251 | d2%vertex_assign_j => d1%vertex_assign_j |
---|
| 252 | d2%vertex_assign_pos => d1%vertex_assign_pos |
---|
[26] | 253 | d2%xyz => d1%xyz |
---|
| 254 | d2%neighbour => d1%neighbour |
---|
| 255 | d2%delta => d1%delta |
---|
| 256 | d2%vertex => d1%vertex |
---|
| 257 | d2%own => d1%own |
---|
| 258 | d2%ne => d1%ne |
---|
| 259 | |
---|
| 260 | d2%t_right = d1%t_right |
---|
| 261 | d2%t_rup = d1%t_rup |
---|
| 262 | d2%t_lup = d1%t_lup |
---|
| 263 | d2%t_left = d1%t_left |
---|
| 264 | d2%t_ldown = d1%t_ldown |
---|
| 265 | d2%t_rdown = d1%t_rdown |
---|
| 266 | |
---|
| 267 | d2%u_right = d1%u_right |
---|
| 268 | d2%u_rup = d1%u_rup |
---|
| 269 | d2%u_lup = d1%u_lup |
---|
| 270 | d2%u_left = d1%u_left |
---|
| 271 | d2%u_ldown = d1%u_ldown |
---|
| 272 | d2%u_rdown = d1%u_rdown |
---|
| 273 | |
---|
| 274 | d2%z_rup = d1%z_rup |
---|
| 275 | d2%z_up = d1%z_up |
---|
| 276 | d2%z_lup = d1%z_lup |
---|
| 277 | d2%z_ldown = d1%z_ldown |
---|
| 278 | d2%z_down = d1%z_down |
---|
| 279 | d2%z_rdown = d1%z_rdown |
---|
| 280 | |
---|
| 281 | d2%t_pos = d1%t_pos |
---|
| 282 | d2%u_pos = d1%u_pos |
---|
| 283 | d2%z_pos = d1%z_pos |
---|
| 284 | |
---|
| 285 | END SUBROUTINE copy_domain |
---|
| 286 | |
---|
| 287 | |
---|
[12] | 288 | SUBROUTINE assign_cell |
---|
| 289 | USE metric |
---|
[726] | 290 | INTEGER :: ind_d,ind,ind2,e,v |
---|
[146] | 291 | INTEGER :: nf,nf2 |
---|
[12] | 292 | INTEGER :: i,j,k,ii,jj |
---|
| 293 | TYPE(t_domain),POINTER :: d |
---|
[21] | 294 | INTEGER :: delta |
---|
[12] | 295 | |
---|
| 296 | |
---|
[26] | 297 | DO ind_d=1,ndomain_glo |
---|
| 298 | d=>domain_glo(ind_d) |
---|
[12] | 299 | nf=d%face |
---|
| 300 | DO j=d%jj_begin,d%jj_end |
---|
| 301 | DO i=d%ii_begin,d%ii_end |
---|
| 302 | ii=d%ii_begin_glo-d%ii_begin+i |
---|
| 303 | jj=d%jj_begin_glo-d%jj_begin+j |
---|
| 304 | ind=vertex_glo(ii,jj,nf)%ind |
---|
[21] | 305 | delta=vertex_glo(ii,jj,nf)%delta |
---|
[12] | 306 | IF (cell_glo(ind)%assign_face==nf) THEN |
---|
| 307 | cell_glo(ind)%assign_domain=ind_d |
---|
| 308 | cell_glo(ind)%assign_i=i |
---|
| 309 | cell_glo(ind)%assign_j=j |
---|
| 310 | ENDIF |
---|
[21] | 311 | IF ( i+1 <= d%ii_end ) CALL assign_edge(ind_d,ind,i,j,delta,0) |
---|
| 312 | IF ( j+1 <= d%jj_end ) CALL assign_edge(ind_d,ind,i,j,delta,1) |
---|
| 313 | IF ( i-1 >= d%ii_begin .AND. j+1<=d%jj_end ) CALL assign_edge(ind_d,ind,i,j,delta,2) |
---|
| 314 | IF ( i-1 >= d%ii_begin ) CALL assign_edge(ind_d,ind,i,j,delta,3) |
---|
| 315 | IF ( j-1 >= d%jj_begin ) CALL assign_edge(ind_d,ind,i,j,delta,4) |
---|
| 316 | IF ( i+1 <= d%ii_end .AND. j-1 >=d%jj_begin ) CALL assign_edge(ind_d,ind,i,j,delta,5) |
---|
[726] | 317 | |
---|
| 318 | IF ( i+1 <= d%ii_end .AND. j+1 <= d%jj_end) CALL assign_vertex(ind_d,ind,i,j,delta,0) |
---|
| 319 | IF ( i-1 >= d%ii_begin .AND. j+1 <= d%jj_end) CALL assign_vertex(ind_d,ind,i,j,delta,1) |
---|
| 320 | IF ( i-1 >= d%ii_begin .AND. j+1 <= d%jj_end) CALL assign_vertex(ind_d,ind,i,j,delta,2) |
---|
| 321 | IF ( i-1 >= d%ii_begin .AND. j-1 >= d%jj_begin) CALL assign_vertex(ind_d,ind,i,j,delta,3) |
---|
| 322 | IF ( i+1 <= d%ii_end .AND. j-1 >= d%jj_begin) CALL assign_vertex(ind_d,ind,i,j,delta,4) |
---|
| 323 | IF ( i+1 <= d%ii_end .AND. j-1 >= d%jj_begin) CALL assign_vertex(ind_d,ind,i,j,delta,5) |
---|
[12] | 324 | ENDDO |
---|
| 325 | ENDDO |
---|
| 326 | ENDDO |
---|
| 327 | |
---|
[21] | 328 | |
---|
[26] | 329 | DO ind_d=1,ndomain_glo |
---|
| 330 | d=>domain_glo(ind_d) |
---|
[12] | 331 | nf=d%face |
---|
| 332 | DO j=d%jj_begin-1,d%jj_end+1 |
---|
| 333 | DO i=d%ii_begin-1,d%ii_end+1 |
---|
| 334 | ii=d%ii_begin_glo-d%ii_begin+i |
---|
| 335 | jj=d%jj_begin_glo-d%jj_begin+j |
---|
| 336 | ind=vertex_glo(ii,jj,nf)%ind |
---|
[266] | 337 | d%assign_cell_glo(i,j) = ind |
---|
[12] | 338 | d%assign_domain(i,j)=cell_glo(ind)%assign_domain |
---|
| 339 | d%assign_i(i,j)=cell_glo(ind)%assign_i |
---|
| 340 | d%assign_j(i,j)=cell_glo(ind)%assign_j |
---|
[21] | 341 | delta=vertex_glo(ii,jj,nf)%delta |
---|
[12] | 342 | d%delta(i,j)=vertex_glo(ii,jj,nf)%delta |
---|
| 343 | DO k=0,5 |
---|
| 344 | ind2=vertex_glo(ii,jj,nf)%neighbour(k) |
---|
| 345 | d%neighbour(:,k,i,j)=cell_glo(ind2)%xyz(:) |
---|
[146] | 346 | |
---|
| 347 | d%ne(k,i,j)=1-2*MOD(k,2) |
---|
| 348 | |
---|
[21] | 349 | e=cell_glo(ind)%edge(MOD(k+delta+6,6)) |
---|
| 350 | d%edge_assign_domain(k,i,j)=edge_glo(e)%assign_domain |
---|
| 351 | d%edge_assign_i(k,i,j)=edge_glo(e)%assign_i |
---|
| 352 | d%edge_assign_j(k,i,j)=edge_glo(e)%assign_j |
---|
| 353 | d%edge_assign_pos(k,i,j)=edge_glo(e)%assign_pos |
---|
[146] | 354 | nf2=domain_glo(edge_glo(e)%assign_domain)%face |
---|
| 355 | d%edge_assign_sign(k,i,j)=1-2*MOD(12+tab_index(nf,nf2,0),2) |
---|
[149] | 356 | IF (MOD(6+k+tab_index(nf,nf2,0),6)/=edge_glo(e)%assign_pos .AND. MOD(6+k+tab_index(nf,nf2,0),6) & |
---|
| 357 | /= MOD(edge_glo(e)%assign_pos+3,6)) THEN |
---|
[146] | 358 | d%edge_assign_sign(k,i,j)=-d%edge_assign_sign(k,i,j) |
---|
| 359 | ENDIF |
---|
[726] | 360 | |
---|
| 361 | v=cell_glo(ind)%vertex(MOD(k+delta+6,6)) |
---|
| 362 | d%vertex_assign_domain(k,i,j)=vertices_glo(v)%assign_domain |
---|
| 363 | d%vertex_assign_i(k,i,j)=vertices_glo(v)%assign_i |
---|
| 364 | d%vertex_assign_j(k,i,j)=vertices_glo(v)%assign_j |
---|
| 365 | d%vertex_assign_pos(k,i,j)=vertices_glo(v)%assign_pos |
---|
[146] | 366 | |
---|
[12] | 367 | ENDDO |
---|
| 368 | d%xyz(:,i,j)=cell_glo(ind)%xyz(:) |
---|
| 369 | IF (d%assign_domain(i,j)==ind_d) THEN |
---|
| 370 | d%own(i,j)=.TRUE. |
---|
| 371 | ELSE |
---|
| 372 | d%own(i,j)=.FALSE. |
---|
| 373 | ENDIF |
---|
| 374 | ENDDO |
---|
| 375 | ENDDO |
---|
[21] | 376 | ENDDO |
---|
| 377 | |
---|
| 378 | CONTAINS |
---|
| 379 | |
---|
| 380 | SUBROUTINE assign_edge(ind_d,ind,i,j,delta,k) |
---|
| 381 | INTEGER :: ind_d,ind,i,j,delta,k |
---|
| 382 | INTEGER :: e |
---|
| 383 | e=cell_glo(ind)%edge(MOD(k+delta+6,6)) |
---|
| 384 | edge_glo(e)%assign_domain=ind_d |
---|
| 385 | edge_glo(e)%assign_i=i |
---|
| 386 | edge_glo(e)%assign_j=j |
---|
| 387 | edge_glo(e)%assign_pos=k |
---|
[146] | 388 | edge_glo(e)%assign_delta=delta |
---|
[151] | 389 | |
---|
[726] | 390 | END SUBROUTINE assign_edge |
---|
| 391 | |
---|
| 392 | SUBROUTINE assign_vertex(ind_d,ind,i,j,delta,k) |
---|
| 393 | INTEGER :: ind_d,ind,i,j,delta,k |
---|
| 394 | INTEGER :: e |
---|
| 395 | |
---|
| 396 | e=cell_glo(ind)%vertex(MOD(k+delta+6,6)) |
---|
| 397 | vertices_glo(e)%assign_domain=ind_d |
---|
| 398 | vertices_glo(e)%assign_i=i |
---|
| 399 | vertices_glo(e)%assign_j=j |
---|
| 400 | vertices_glo(e)%assign_pos=k |
---|
| 401 | vertices_glo(e)%assign_delta=delta |
---|
| 402 | |
---|
| 403 | END SUBROUTINE assign_vertex |
---|
| 404 | |
---|
| 405 | END SUBROUTINE assign_cell |
---|
[21] | 406 | |
---|
[12] | 407 | SUBROUTINE compute_boundary |
---|
| 408 | USE spherical_geom_mod |
---|
| 409 | INTEGER :: ind_d |
---|
| 410 | INTEGER :: i,j,k |
---|
| 411 | TYPE(t_domain),POINTER :: d |
---|
| 412 | REAL(rstd) :: ng1(3),ng2(3) |
---|
| 413 | |
---|
[26] | 414 | DO ind_d=1,ndomain_glo |
---|
| 415 | d=>domain_glo(ind_d) |
---|
[12] | 416 | DO j=d%jj_begin-1,d%jj_end+1 |
---|
| 417 | DO i=d%ii_begin-1,d%ii_end+1 |
---|
| 418 | DO k=0,5 |
---|
| 419 | ng1=d%neighbour(:,MOD(k,6),i,j) |
---|
| 420 | ng2=d%neighbour(:,MOD(k+1,6),i,j) |
---|
| 421 | IF (sqrt(sum((ng1-ng2)**2))<1e-16) ng2=d%neighbour(:,MOD(k+2,6),i,j) |
---|
[15] | 422 | CALL circumcenter(d%xyz(:,i,j),ng1,ng2,d%vertex(:,k,i,j)) |
---|
[12] | 423 | ENDDO |
---|
| 424 | ENDDO |
---|
| 425 | ENDDO |
---|
| 426 | ENDDO |
---|
| 427 | END SUBROUTINE compute_boundary |
---|
| 428 | |
---|
| 429 | SUBROUTINE set_neighbour_indice |
---|
| 430 | USE metric |
---|
| 431 | INTEGER :: ind_d |
---|
| 432 | TYPE(t_domain),POINTER :: d |
---|
| 433 | |
---|
[26] | 434 | DO ind_d=1,ndomain_glo |
---|
| 435 | d=>domain_glo(ind_d) |
---|
[12] | 436 | d%t_right=1 |
---|
| 437 | d%t_left=-1 |
---|
| 438 | d%t_rup=d%iim |
---|
| 439 | d%t_lup=d%iim-1 |
---|
| 440 | d%t_ldown=-d%iim |
---|
| 441 | d%t_rdown=-d%iim+1 |
---|
| 442 | |
---|
| 443 | d%u_right=0 |
---|
| 444 | d%u_lup=d%iim*d%jjm |
---|
| 445 | d%u_ldown=2*d%iim*d%jjm |
---|
| 446 | |
---|
| 447 | d%u_rup=d%t_rup+d%u_ldown |
---|
| 448 | d%u_left=d%t_left+d%u_right |
---|
| 449 | d%u_rdown=d%t_rdown+d%u_lup |
---|
| 450 | |
---|
| 451 | d%z_up=0 |
---|
| 452 | d%z_down=d%iim*d%jjm |
---|
| 453 | d%z_rup=d%t_rup+d%z_down |
---|
| 454 | d%z_lup=d%t_lup+d%z_down |
---|
| 455 | d%z_ldown=d%t_ldown+d%z_up |
---|
| 456 | d%z_rdown=d%t_rdown+d%z_up |
---|
| 457 | |
---|
| 458 | d%t_pos(right)=d%t_right |
---|
| 459 | d%t_pos(rup)=d%t_rup |
---|
| 460 | d%t_pos(lup)=d%t_lup |
---|
| 461 | d%t_pos(left)=d%t_left |
---|
| 462 | d%t_pos(ldown)=d%t_ldown |
---|
| 463 | d%t_pos(rdown)=d%t_rdown |
---|
| 464 | |
---|
| 465 | d%u_pos(right)=d%u_right |
---|
| 466 | d%u_pos(rup)=d%u_rup |
---|
| 467 | d%u_pos(lup)=d%u_lup |
---|
| 468 | d%u_pos(left)=d%u_left |
---|
| 469 | d%u_pos(ldown)=d%u_ldown |
---|
| 470 | d%u_pos(rdown)=d%u_rdown |
---|
| 471 | |
---|
| 472 | d%z_pos(vrup)=d%z_rup |
---|
| 473 | d%z_pos(vup)=d%z_up |
---|
| 474 | d%z_pos(vlup)=d%z_lup |
---|
| 475 | d%z_pos(vldown)=d%z_ldown |
---|
| 476 | d%z_pos(vdown)=d%z_down |
---|
| 477 | d%z_pos(vrdown)=d%z_rdown |
---|
| 478 | |
---|
| 479 | ENDDO |
---|
| 480 | |
---|
| 481 | END SUBROUTINE set_neighbour_indice |
---|
| 482 | |
---|
[26] | 483 | SUBROUTINE assign_domain |
---|
| 484 | USE mpipara |
---|
[327] | 485 | USE grid_param |
---|
[26] | 486 | INTEGER :: nb_domain(0:mpi_size-1) |
---|
| 487 | INTEGER :: rank, ind,ind_glo |
---|
[327] | 488 | INTEGER :: block_j,jb,i,j,nd_glo,n,nf |
---|
| 489 | LOGICAL :: exit |
---|
[26] | 490 | |
---|
| 491 | DO rank=0,mpi_size-1 |
---|
| 492 | nb_domain(rank)=ndomain_glo/mpi_size |
---|
| 493 | IF ( rank < MOD(ndomain_glo,mpi_size) ) nb_domain(rank)=nb_domain(rank)+1 |
---|
| 494 | ENDDO |
---|
| 495 | |
---|
| 496 | ndomain=nb_domain(mpi_rank) |
---|
[879] | 497 | |
---|
| 498 | ! ALLOCATE(domain(ndomain)) |
---|
| 499 | CALL allocate_mesh(ndomain, mesh_loc) |
---|
| 500 | domain => mesh_loc%domain |
---|
| 501 | |
---|
[26] | 502 | ALLOCATE(domloc_glo_ind(ndomain)) |
---|
| 503 | |
---|
[327] | 504 | |
---|
| 505 | block_j=sqrt(nsplit_i*nsplit_j*nb_face*1./mpi_size) |
---|
| 506 | exit=.FALSE. |
---|
| 507 | jb=1 |
---|
| 508 | i=1 |
---|
| 509 | j=1 |
---|
| 510 | ind=1 |
---|
| 511 | nd_glo=0 |
---|
[26] | 512 | rank=0 |
---|
[327] | 513 | DO WHILE (.NOT. exit) |
---|
| 514 | |
---|
| 515 | IF (j==MIN(jb+block_j,nsplit_j*nb_face+1)) THEN |
---|
| 516 | j=jb |
---|
| 517 | i=i+1 |
---|
| 518 | ENDIF |
---|
| 519 | |
---|
| 520 | IF (i>nsplit_i) THEN |
---|
| 521 | i=1 |
---|
| 522 | jb=jb+block_j |
---|
| 523 | j=jb |
---|
| 524 | ENDIF |
---|
| 525 | |
---|
| 526 | IF (ind>nb_domain(rank)) THEN |
---|
| 527 | rank=rank+1 |
---|
| 528 | ind=1 |
---|
| 529 | ENDIF |
---|
| 530 | ind_glo=(j-1)*nsplit_i+i |
---|
| 531 | |
---|
| 532 | nd_glo=nd_glo+1 |
---|
| 533 | IF (nd_glo==ndomain_glo) THEN |
---|
| 534 | |
---|
| 535 | exit=.TRUE. |
---|
| 536 | IF (.NOT. (rank==mpi_size-1 .AND. ind==nb_domain(rank) )) THEN |
---|
| 537 | PRINT *, "Distribution problem in assign_domain" |
---|
| 538 | STOP |
---|
| 539 | ENDIF |
---|
| 540 | |
---|
| 541 | ENDIF |
---|
| 542 | |
---|
[26] | 543 | domglo_rank(ind_glo)=rank |
---|
| 544 | domglo_loc_ind(ind_glo)=ind |
---|
| 545 | IF (rank==mpi_rank) THEN |
---|
| 546 | CALL copy_domain(domain_glo(ind_glo),domain(ind)) |
---|
| 547 | domloc_glo_ind(ind)=ind_glo |
---|
| 548 | ENDIF |
---|
[12] | 549 | |
---|
[327] | 550 | j=j+1 |
---|
| 551 | ind=ind+1 |
---|
| 552 | |
---|
[26] | 553 | ENDDO |
---|
[186] | 554 | |
---|
[327] | 555 | IF (is_mpi_master) THEN |
---|
| 556 | |
---|
| 557 | ind_glo=0 |
---|
[358] | 558 | PRINT *,'' |
---|
[327] | 559 | PRINT*, ' MPI PROCESS DISTRIBUTION' |
---|
[358] | 560 | PRINT *,'' |
---|
[327] | 561 | |
---|
| 562 | WRITE(*,"(' ')", ADVANCE='NO') |
---|
| 563 | DO n=1,nsplit_i*7-1 |
---|
| 564 | WRITE(*,"('=')", ADVANCE='NO') |
---|
| 565 | ENDDO |
---|
[358] | 566 | PRINT *,'' |
---|
[327] | 567 | |
---|
| 568 | DO nf=1,nb_face |
---|
| 569 | DO j=1,nsplit_j |
---|
| 570 | IF (j>1) THEN |
---|
| 571 | WRITE(*,"(' ')", ADVANCE='NO') |
---|
| 572 | DO n=1,nsplit_i*7-1 |
---|
| 573 | WRITE(*,"('-')", ADVANCE='NO') |
---|
| 574 | ENDDO |
---|
[358] | 575 | PRINT *,'' |
---|
[327] | 576 | ENDIF |
---|
| 577 | |
---|
| 578 | WRITE(*,"('|')", ADVANCE='NO') |
---|
| 579 | DO i=1,nsplit_i |
---|
| 580 | WRITE(*,"(' ',' ',' |')",ADVANCE='NO') |
---|
| 581 | ENDDO |
---|
[358] | 582 | PRINT *,'' |
---|
[327] | 583 | |
---|
| 584 | WRITE(*,"('|')", ADVANCE='NO') |
---|
| 585 | DO i=1,nsplit_i |
---|
| 586 | ind_glo=ind_glo+1 |
---|
[821] | 587 | WRITE(*,"(' ',i4.4 ,' |')", ADVANCE='NO') domglo_rank(ind_glo) |
---|
[327] | 588 | END DO |
---|
[358] | 589 | PRINT *,'' |
---|
[327] | 590 | |
---|
| 591 | WRITE(*,"('|')", ADVANCE='NO') |
---|
| 592 | DO i=1,nsplit_i |
---|
[821] | 593 | WRITE(*,"(' ',' ',' |')", ADVANCE='NO') |
---|
[327] | 594 | ENDDO |
---|
[358] | 595 | PRINT *,'' |
---|
[327] | 596 | |
---|
| 597 | ENDDO |
---|
| 598 | |
---|
| 599 | WRITE(*,"(' ')", ADVANCE='NO') |
---|
| 600 | DO n=1,nsplit_i*7-1 |
---|
| 601 | WRITE(*,"('=')", ADVANCE='NO') |
---|
| 602 | ENDDO |
---|
[358] | 603 | PRINT *,'' |
---|
[327] | 604 | ENDDO |
---|
| 605 | ENDIF |
---|
| 606 | |
---|
| 607 | ! rank=0 |
---|
| 608 | ! ind=0 |
---|
| 609 | ! DO ind_glo=1,ndomain_glo |
---|
| 610 | ! ind=ind+1 |
---|
| 611 | ! domglo_rank(ind_glo)=rank |
---|
| 612 | ! domglo_loc_ind(ind_glo)=ind |
---|
| 613 | ! IF (rank==mpi_rank) THEN |
---|
| 614 | ! CALL copy_domain(domain_glo(ind_glo),domain(ind)) |
---|
| 615 | ! domloc_glo_ind(ind)=ind_glo |
---|
| 616 | ! ENDIF |
---|
| 617 | ! |
---|
| 618 | ! IF (ind==nb_domain(rank)) THEN |
---|
| 619 | ! rank=rank+1 |
---|
| 620 | ! ind=0 |
---|
| 621 | ! ENDIF |
---|
| 622 | ! ENDDO |
---|
| 623 | |
---|
[186] | 624 | !$OMP PARALLEL |
---|
| 625 | CALL assign_domain_omp |
---|
| 626 | !$OMP END PARALLEL |
---|
[813] | 627 | |
---|
[26] | 628 | END SUBROUTINE assign_domain |
---|
[186] | 629 | |
---|
| 630 | SUBROUTINE assign_domain_omp |
---|
| 631 | USE omp_para |
---|
| 632 | USE mpipara |
---|
| 633 | INTEGER :: nb_domain |
---|
| 634 | INTEGER :: rank, ind, i |
---|
| 635 | |
---|
| 636 | ALLOCATE(assigned_domain(ndomain)) |
---|
[26] | 637 | |
---|
[186] | 638 | ind=0 |
---|
[295] | 639 | DO rank=0,omp_domain_size-1 |
---|
| 640 | nb_domain=ndomain/omp_domain_size |
---|
| 641 | IF ( rank < MOD(ndomain,omp_domain_size) ) nb_domain=nb_domain+1 |
---|
[186] | 642 | |
---|
| 643 | DO i=1,nb_domain |
---|
| 644 | ind=ind+1 |
---|
[295] | 645 | IF (rank==omp_domain_rank) THEN |
---|
[186] | 646 | assigned_domain(ind)=.TRUE. |
---|
[402] | 647 | WRITE(*,'(A,I6.6,A,I4.4,A,I4.4,A,I8.8)') "Rank ",mpi_rank," task ",rank," local domain ",ind, & |
---|
| 648 | " global domain ",domloc_glo_ind(ind) |
---|
[186] | 649 | ELSE |
---|
| 650 | assigned_domain(ind)=.FALSE. |
---|
| 651 | ENDIF |
---|
| 652 | ENDDO |
---|
| 653 | |
---|
| 654 | ENDDO |
---|
| 655 | |
---|
| 656 | END SUBROUTINE assign_domain_omp |
---|
| 657 | |
---|
[867] | 658 | SUBROUTINE init_domain_unst(dom) |
---|
| 659 | USE grid_param, ONLY : primal_num |
---|
| 660 | TYPE(t_domain) :: dom |
---|
| 661 | dom%ii_begin=1 |
---|
| 662 | dom%ii_end=primal_num |
---|
| 663 | dom%iim=primal_num |
---|
| 664 | dom%jj_begin=1 |
---|
| 665 | dom%jj_end=1 |
---|
| 666 | dom%jjm=1 |
---|
| 667 | ALLOCATE(dom%assign_domain(primal_num,1)) |
---|
| 668 | dom%assign_domain(:,:)=1 ! all cells are assigned to MPI rank 0 |
---|
| 669 | END SUBROUTINE init_domain_unst |
---|
| 670 | |
---|
[12] | 671 | SUBROUTINE compute_domain |
---|
[813] | 672 | USE grid_param, ONLY : grid_type, grid_unst, grid_ico |
---|
[879] | 673 | INTEGER :: ind |
---|
[813] | 674 | SELECT CASE(grid_type) |
---|
| 675 | CASE(grid_unst) |
---|
[867] | 676 | ! FIXME temporary, sequential hack |
---|
| 677 | ndomain_glo=1 |
---|
[879] | 678 | CALL allocate_mesh(ndomain_glo, mesh_glo) |
---|
| 679 | domain_glo => mesh_glo%domain |
---|
[867] | 680 | ALLOCATE(domglo_rank(ndomain_glo)) |
---|
| 681 | ALLOCATE(domglo_loc_ind(ndomain_glo)) |
---|
| 682 | domglo_rank(:)=0 |
---|
| 683 | domglo_loc_ind(:)=1 |
---|
| 684 | |
---|
[813] | 685 | ndomain=1 |
---|
[879] | 686 | CALL allocate_mesh(ndomain, mesh_loc) |
---|
| 687 | domain => mesh_loc%domain |
---|
[867] | 688 | ALLOCATE(domloc_glo_ind(ndomain)) |
---|
| 689 | ALLOCATE(assigned_domain(ndomain)) |
---|
| 690 | domloc_glo_ind(:)=1 |
---|
| 691 | assigned_domain(:)=.TRUE. |
---|
| 692 | |
---|
| 693 | CALL init_domain_unst(domain(1)) |
---|
| 694 | CALL init_domain_unst(domain_glo(1)) |
---|
| 695 | |
---|
[813] | 696 | CASE DEFAULT |
---|
| 697 | CALL init_domain_param |
---|
| 698 | CALL create_domain |
---|
| 699 | CALL assign_cell |
---|
| 700 | CALL compute_boundary |
---|
| 701 | CALL set_neighbour_indice |
---|
| 702 | CALL assign_domain |
---|
| 703 | END SELECT |
---|
[879] | 704 | |
---|
| 705 | ! ALLOCATE(primal_cells(ndomain)) |
---|
| 706 | ! ALLOCATE(dual_cells(ndomain)) |
---|
| 707 | ! ALLOCATE(primal_cells_single(ndomain)) |
---|
| 708 | ! ALLOCATE(dual_cells_single(ndomain)) |
---|
| 709 | ! DO ind=1, ndomain |
---|
| 710 | ! domain(ind)%primal_cells => primal_cells(ind) |
---|
| 711 | ! domain(ind)%dual_cells => dual_cells(ind) |
---|
| 712 | ! domain(ind)%primal_cells_single => primal_cells_single(ind) |
---|
| 713 | ! domain(ind)%dual_cells_single => dual_cells_single(ind) |
---|
| 714 | ! END DO |
---|
| 715 | |
---|
[12] | 716 | END SUBROUTINE compute_domain |
---|
| 717 | |
---|
| 718 | END MODULE domain_mod |
---|