9 SUBROUTINE reorder_mesh(communicator,nb_proc,mesh,mesh_loc,list_of_interfaces)
17 INTEGER,
DIMENSION(2) :: me_loc, mes_loc, np_loc
18 INTEGER,
DIMENSION(:),
POINTER,
OPTIONAL :: list_of_interfaces
19 INTEGER,
DIMENSION(mesh%me+1) :: xind
20 INTEGER,
DIMENSION(mesh%me) :: vwgt, vsize, part, old_m_to_new
21 LOGICAL,
DIMENSION(mesh%np) :: virgin
22 INTEGER,
DIMENSION(mesh%np) :: old_j_to_new
23 INTEGER,
DIMENSION(mesh%mes) :: old_ms_to_new, parts
24 INTEGER,
DIMENSION(2,mesh%mes) :: inter_news
25 INTEGER,
DIMENSION(SIZE(mesh%jjs,1)) :: i_loc
26 LOGICAL,
DIMENSION(mesh%mes) :: virgins
27 REAL(KIND=8),
DIMENSION(mesh%np) :: r_old_j_to_new
28 INTEGER,
DIMENSION(5) :: opts
29 INTEGER,
DIMENSION(:),
ALLOCATABLE :: xadj, adjwgt, nblmt_per_proc, start, displ
30 INTEGER,
DIMENSION(:,:),
ALLOCATABLE :: inter
31 INTEGER :: dim, nb_neigh, edge, m, ms, n, nb, numflag, p, wgtflag, new_dof, j, &
32 news, ns, nws, msop, nsop, proc, iop, ncon
33 LOGICAL :: test=.false.
36 REAL(KIND=4),
DIMENSION(:),
ALLOCATABLE :: tpwgts
37 REAL(KIND=4),
DIMENSION(1) :: ubvec
38 INTEGER,
DIMENSION(METIS_NOPTIONS) :: metis_opt
40 #include "petsc/finclude/petsc.h" 41 petscerrorcode :: ierr
42 petscmpiint :: rank , nb_proc
43 mpi_comm :: communicator
44 CALL mpi_comm_rank(communicator,rank,ierr)
48 me_loc = (/ 1, mesh%me /)
49 mes_loc = (/ 1, mesh%mes /)
50 np_loc = (/ 1, mesh%np /)
56 nb_neigh =
SIZE(mesh%neigh,1)
61 IF (mesh%neigh(n,m)==0) cycle
64 xind(m+1) = xind(m) + nb
66 ALLOCATE(xadj(xind(mesh%me+1)-1))
70 IF (mesh%neigh(n,m)==0) cycle
72 xadj(p) = mesh%neigh(n,m)
75 IF (p/=xind(mesh%me+1)-1)
THEN 80 ALLOCATE(adjwgt(
SIZE(xadj)))
88 ALLOCATE(tpwgts(nb_proc))
90 CALL metis_setdefaultoptions(metis_opt)
93 CALL metis_partgraphrecursive(mesh%me, ncon, xind(1:), xadj(1:), vwgt(1:), vsize(1:), adjwgt(1:), &
94 nb_proc, tpwgts(1:), ubvec(1:), metis_opt(1:), edge, part(1:))
100 ALLOCATE(nblmt_per_proc(nb_proc),start(nb_proc),displ(nb_proc))
103 nblmt_per_proc(part(m)) = nblmt_per_proc(part(m)) + 1
107 start(n) = start(n-1) + nblmt_per_proc(n-1)
109 me_loc(1) = start(rank+1) + 1
110 me_loc(2) = start(rank+1) + nblmt_per_proc(rank+1)
116 start(part(m)) = start(part(m)) + 1
117 old_m_to_new(m) = start(part(m))
123 nws =
SIZE( mesh%jjs,1)
126 parts = part(mesh%neighs)
127 IF (
PRESENT(list_of_interfaces))
THEN 128 IF (
SIZE(list_of_interfaces)/=0)
THEN 132 IF (.NOT.virgins(ms)) cycle
133 IF (minval(abs(mesh%sides(ms)-list_of_interfaces))/=0) cycle
134 i_loc = mesh%jjs(:,ms)
135 DO msop = 1, mesh%mes
136 IF (msop==ms .OR. .NOT.virgins(msop)) cycle
137 IF (minval(abs(mesh%sides(msop)-list_of_interfaces))/=0) cycle
141 iop = mesh%jjs(nsop,msop)
142 IF (maxval(abs(mesh%rr(:,i_loc(ns))-mesh%rr(:,iop))).LT.
epsilon)
THEN 154 CALL error_petsc(.NOT.
'BUG in create_local_mesh, test ')
156 IF (part(mesh%neighs(ms)) == part(mesh%neighs(msop))) cycle
157 proc = min(part(mesh%neighs(ms)),part(mesh%neighs(msop)))
160 virgins(ms) = .false.
161 virgins(msop) = .false.
162 IF (proc /= rank+1) cycle
164 inter_news(1,news)=ms
165 inter_news(2,news)=msop
172 ALLOCATE(inter(
SIZE(mesh%jj,1),mesh%me))
174 inter(:,old_m_to_new(m)) = mesh%jj(:,m)
181 DO n = 1,
SIZE(mesh%jj,1)
183 IF (.NOT.virgin(j)) cycle
184 new_dof = new_dof + 1
186 old_j_to_new(j) = new_dof
190 inter(:,m) = old_j_to_new(mesh%jj(:,m))
198 np_loc(1) = maxval(mesh%jj(:,displ(rank)+1:displ(rank)+nblmt_per_proc(rank))) + 1
200 np_loc(2) = maxval(mesh%jj(:,me_loc(1):me_loc(2)))
204 DO n = 1,
SIZE(mesh%rr,1)
205 r_old_j_to_new(old_j_to_new) = mesh%rr(n,:)
206 mesh%rr(n,:) = r_old_j_to_new(:)
211 ALLOCATE(inter(
SIZE(mesh%neigh,1),mesh%me))
212 dim =
SIZE(mesh%rr,1)
215 IF (mesh%neigh(n,m) /=0)
THEN 216 inter(n,old_m_to_new(m)) = old_m_to_new(mesh%neigh(n,m))
218 inter(n,old_m_to_new(m)) = 0
227 DEALLOCATE(xadj);
ALLOCATE(xadj(mesh%me))
228 xadj(old_m_to_new) = mesh%i_d
233 DEALLOCATE(xadj);
ALLOCATE(xadj(mesh%mes))
239 nblmt_per_proc(n) = nblmt_per_proc(n) + 1
243 start(n) = start(n-1) + nblmt_per_proc(n-1)
245 mes_loc(1) = start(rank+1) + 1
246 mes_loc(2) = start(rank+1) + nblmt_per_proc(rank+1)
250 start(n) = start(n) + 1
251 old_ms_to_new(ms) = start(n)
253 xadj(old_ms_to_new) = mesh%neighs
255 xadj = old_m_to_new(mesh%neighs)
260 xadj(1:news) = old_ms_to_new(inter_news(1,1:news))
261 inter_news(1,1:news) = xadj(1:news)
262 xadj(1:news) = old_ms_to_new(inter_news(2,1:news))
263 inter_news(2,1:news) = xadj(1:news)
267 xadj(old_ms_to_new) = mesh%sides
272 DO n = 1,
SIZE(mesh%jjs,1)
273 xadj(old_ms_to_new) = old_j_to_new(mesh%jjs(n,:))
279 CALL create_local_mesh(mesh, mesh_loc, me_loc, mes_loc, np_loc, news, inter_news(:,1:news))
280 CALL mpi_comm_rank(mpi_comm_world,rank,ierr)
282 DEALLOCATE(xadj,adjwgt,nblmt_per_proc,start,displ)
285 SUBROUTINE create_local_mesh(mesh, mesh_loc, me_loc, mes_loc, np_loc, news, inter_news)
290 INTEGER,
DIMENSION(2),
INTENT(IN) :: me_loc, mes_loc, np_loc
291 INTEGER,
DIMENSION(:,:),
INTENT(IN),
OPTIONAL :: inter_news
292 INTEGER,
OPTIONAL :: news
293 INTEGER,
DIMENSION(mesh%me) :: m_glob_to_loc,m_loc_to_glob
294 INTEGER,
DIMENSION(mesh%np) :: glob_to_loc,loc_to_glob
295 LOGICAL,
DIMENSION(mesh%np) :: virgin
296 INTEGER :: dim, nws, nw, m, ms, mop, msop, ns, msup, minf, dof, &
297 dom_me, nwc, dom_mes, dom_np, n, i
299 dim =
SIZE(mesh%rr,1)
300 nws =
SIZE(mesh%jjs,1)
302 nwc =
SIZE(mesh%neigh,1)
304 IF (me_loc(2) - me_loc(1) + 1==mesh%me)
THEN 305 mesh_loc%me = mesh%me
306 mesh_loc%np = mesh%np
307 mesh_loc%mes = mesh%mes
308 mesh_loc%dom_me = mesh%me
309 mesh_loc%dom_np = mesh%np
310 mesh_loc%dom_mes = mesh%mes
311 ALLOCATE(mesh_loc%jj(nw,mesh%me))
312 mesh_loc%jj = mesh%jj
313 ALLOCATE(mesh_loc%loc_to_glob(mesh%np))
315 mesh_loc%loc_to_glob(n) = n
317 ALLOCATE(mesh_loc%rr(dim,mesh%np))
319 ALLOCATE(mesh_loc%neigh(nwc,mesh%me))
320 mesh_loc%neigh = mesh%neigh
321 ALLOCATE(mesh_loc%i_d(mesh%me))
322 mesh_loc%i_d = mesh%i_d
323 ALLOCATE(mesh_loc%neighs(mesh_loc%mes))
324 mesh_loc%neighs = mesh%neighs
325 ALLOCATE(mesh_loc%sides(mesh_loc%mes))
326 mesh_loc%sides = mesh%sides
327 ALLOCATE(mesh_loc%jjs(nws,mesh_loc%mes))
328 mesh_loc%jjs = mesh%jjs
333 IF (.NOT.
PRESENT(news) .OR. .NOT.
PRESENT( inter_news))
THEN 334 CALL error_petsc(.NOT..OR..NOT.
'BUG in create_local_mesh present(news) present( inter_news)')
338 dom_me = me_loc(2) - me_loc(1) + 1
339 dom_mes = mes_loc(2) - mes_loc(1) + 1
340 dom_np = np_loc(2) - np_loc(1) + 1
341 mesh_loc%me = dom_me + news
342 mesh_loc%mes = dom_mes
343 mesh_loc%dom_me = dom_me
344 mesh_loc%dom_np = dom_np
345 mesh_loc%dom_mes = dom_mes
350 DO m = me_loc(1), me_loc(2)
353 IF(.NOT.virgin(i) .OR. i.GE.np_loc(1)) cycle
358 ALLOCATE(mesh_loc%jj(nw,mesh_loc%me))
363 DO m = me_loc(1), me_loc(2)
368 IF (i .LT. np_loc(1))
THEN 373 glob_to_loc(i) = i-np_loc(1) + 1
374 loc_to_glob(i-np_loc(1) + 1) = i
378 m_loc_to_glob(m-me_loc(1)+1) = m
379 m_glob_to_loc(m) = m-me_loc(1)+1
383 mesh_loc%jj(n,1:dom_me) = glob_to_loc(mesh%jj(n,me_loc(1):me_loc(2)))
387 msop=inter_news(2,ns)
388 IF (mesh%neighs(ms) < me_loc(1) .OR. mesh%neighs(ms) > me_loc(2))
THEN 391 m = mesh%neighs(msop)
398 IF (i.GE.np_loc(1) .AND. i.LE.np_loc(2))
THEN 399 CALL error_petsc(.GE..AND..LE.
'BUG in create_local_mesh, inp_loc(1) inp_loc(2)')
406 mesh_loc%jj(:,dom_me+ns) = glob_to_loc(mesh%jj(:,m))
407 m_loc_to_glob(dom_me+ns) = m
408 m_glob_to_loc(m) = dom_me+ns
413 IF (maxval(mesh_loc%jj)/=dof)
THEN 414 CALL error_petsc(
'BUG in create_local_mesh, mesh_loc%jj)/=dof')
417 ALLOCATE(mesh_loc%loc_to_glob(mesh_loc%np))
418 mesh_loc%loc_to_glob = loc_to_glob(1:mesh_loc%np)
422 ALLOCATE(mesh_loc%rr(dim,mesh_loc%np))
423 DO n = 1, mesh_loc%np
424 mesh_loc%rr(:,n) = mesh%rr(:,mesh_loc%loc_to_glob(n))
429 ALLOCATE(mesh_loc%neigh(nwc,mesh_loc%me))
430 msup = maxval(m_loc_to_glob)
431 minf = minval(m_loc_to_glob)
432 DO m = 1, mesh_loc%me
434 mop = mesh%neigh(n,m_loc_to_glob(m))
436 mesh_loc%neigh(n,m) = 0
437 ELSE IF(mop<minf .OR. mop>msup)
THEN 438 mesh_loc%neigh(n,m) = -mop
440 mesh_loc%neigh(n,m) = m_glob_to_loc(mop)
447 ALLOCATE(mesh_loc%i_d(mesh_loc%me))
448 mesh_loc%i_d = mesh%i_d(m_loc_to_glob(1:mesh_loc%me))
452 ALLOCATE(mesh_loc%neighs(mesh_loc%mes))
453 mesh_loc%neighs = m_glob_to_loc(mesh%neighs(mes_loc(1):mes_loc(2)))
458 ALLOCATE(mesh_loc%sides(mesh_loc%mes))
459 mesh_loc%sides = mesh%sides(mes_loc(1):mes_loc(2))
463 ALLOCATE(mesh_loc%jjs(nws,mesh_loc%mes))
465 mesh_loc%jjs(ns,:) = glob_to_loc(mesh%jjs(ns,mes_loc(1):mes_loc(2)))
480 DEALLOCATE(mesh%neigh)
481 DEALLOCATE(mesh%sides)
482 DEALLOCATE(mesh%neighs)
485 NULLIFY(mesh%loc_to_glob)
490 IF (mesh%edge_stab)
THEN 493 DEALLOCATE(mesh%jjsi)
494 DEALLOCATE(mesh%neighi)
504 mesh%edge_stab = .false.
subroutine create_local_mesh(mesh, mesh_loc, me_loc, mes_loc, np_loc, news, inter_news)
subroutine error_petsc(string)
subroutine, public free_mesh_after(mesh)
subroutine plot_const_p1_label(jj, rr, uu, file_name)
integer metis_option_numbering
subroutine, public reorder_mesh(communicator, nb_proc, mesh, mesh_loc, list_of_interfaces)