SPEED
MAKE_PARTITION_AND_MPI_FILES.f90
Go to the documentation of this file.
1! Copyright (C) 2012 The SPEED FOUNDATION
2! Author: Ilario Mazzieri
3!
4! This file is part of SPEED.
5!
6! SPEED is free software; you can redistribute it and/or modify it
7! under the terms of the GNU Affero General Public License as
8! published by the Free Software Foundation, either version 3 of the
9! License, or (at your option) any later version.
10!
11! SPEED is distributed in the hope that it will be useful, but
12! WITHOUT ANY WARRANTY; without even the implied warranty of
13! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
14! Affero General Public License for more details.
15!
16! You should have received a copy of the GNU Affero General Public License
17! along with SPEED. If not, see <http://www.gnu.org/licenses/>.
18
23
26
27 subroutine make_partition_and_mpi_files()
28
29 use max_var
30 use speed_par
31
32
33 implicit none
34
35 include 'SPEED.MPI'
36
37!***************************************************************************************************************
38! Partitioning
39
40 if (mpi_id.eq.0) write(*,'(A)') '---------------------Partitioning----------------------'
41
42
43
44 if (mpi_id.eq.0) then
45 if(len_trim(mpi_file) .ne. 70) then
46 file_part = mpi_file(1:len_trim(mpi_file)) // '/elemdomain.mpi'
47 else
48 file_part = 'elemdomain.mpi'
49 endif
50
51 inquire(file=file_part,exist=filefound)
52
53! write(*,*) filefound, file_part, mpi_id
54! read(*,*)
55
56 if(mpi_np .gt. 1) then
57
58 if(filefound .eqv. .true.) then
59
60 if (mpi_id.eq.0) write(*,'(A)') 'Reading existing partitioning...'
61 unit_part = 400
62 open(unit_part,file=file_part)
63 read(unit_part,*) trash
64 do i = 1, nelem
66 enddo
67 close(unit_part)
68 if (mpi_id.eq.0) write(*,'(A)') 'Read.'
69
70 else
71
72 write(*,'(A)') 'Making new partitioning...'
74! call MPI_BARRIER(mpi_comm, mpi_ierr)
75 !sustitute this part adding output elem_domain in mesh partitioning
76 unit_part = 400
77 open(unit_part, file = file_part)
78 read(unit_part,*) trash
79 do i = 1, nelem
81 enddo
82 close(unit_part)
83 write(*,'(A)') 'Made.'
84
85 endif
86 endif
87
88 endif!if (mpi_id .eq. 0)
89
90 call mpi_barrier(mpi_comm, mpi_ierr)
91 call mpi_bcast(elem_domain,nelem,speed_integer,0,mpi_comm_world,mpi_ierr)
92
93!**************************************************************************************************************
94
95 allocate(elem_index(nelem))
96
97
98 nelem_dom = 0
99 do ie = 1,nelem
100 elem_index(ie) = 0
101 if (elem_domain(ie).eq.mpi_id) then
104 endif
105 enddo
106
107
108
109!***************************************************************************************************************
110! SPECTRAL CONNECTIVITY --- ONLY FOR MASTER PROCESS ---
111
112 if (mpi_id.eq.0) then
113 write(*,'(A)')
114 write(*,'(A)')'--------------Making Spectral connectivity-------------'
115 endif
116
117 call mpi_barrier(mpi_comm, mpi_ierr)
118
119 if(mpi_id .eq. 0) then
120
121 file_mpi = 'nodedomain.mpi';
122
123 if(len_trim(mpi_file) .ne. 70) then
124 file_mpi_new = mpi_file(1:len_trim(mpi_file)) // '/' // file_mpi
125 else
127 endif
128
129 inquire(file=file_mpi_new, exist=filefound)
130
131 if(filefound .eqv. .false.) then
132
133 allocate(node_weight(nnod_macro))
134
136
137 allocate(node_pointer(0:nnz_node_weight))
138
140
141 deallocate(node_weight)
142
143 ! Why are we not counting nodes in block nmat_nle
144 con_nnz = nelem +1
145 do ie = 1,nelem
146 do j = 1,nmat
147 if (tag_mat(j).eq.con(ie,1)) nn = sdeg_mat(j) +1
148 enddo
149 con_nnz = con_nnz + nn*nn*nn +1
150 enddo
151
152 allocate(con_spx(0:con_nnz))
153
156
157
158 deallocate(node_pointer); allocate(node_weight(nnod))
159
161
162 allocate(node_pointer(0:nnz_node_weight))
163
165
166 deallocate(node_weight)
167
168 write(*,'(A)') 'Made.'
169
170! SPECTRAL CONNECTIVITY END --- ONLY FOR MASTER PROCESS ---
171!***************************************************************************************************************
172!***************************************************************************************************************
173! NODE DOMAIN FOR MPI
174
175 allocate(node_domain(nnod))
176
177 ! Will this replace partition of repetiting node with partition where the last element occurs
178 do ie = 1, nelem
179 do i = con_spx(ie -1) +1, con_spx(ie) -1
181 enddo
182 enddo
183
184 file_mpi = 'nodedomain.mpi'; unit_mpi = 40
185
186 if(len_trim(mpi_file) .ne. 70) then
187 file_mpi_new = mpi_file(1:len_trim(mpi_file)) // '/' // file_mpi
188 else
190 endif
191
192 open(unit_mpi,file=file_mpi_new)
193 write(unit_mpi,*) nnod
194 do i = 1, nnod
195 write(unit_mpi,*) node_domain(i)
196 enddo
197 close(unit_mpi)
198
199 deallocate(node_domain)
200
201! NODE DOMAIN FOR MPI
202!***************************************************************************************************************
203!***************************************************************************************************************
204! SPECTRAL LOCAL CONNECTIVITY
205
206 write(*,'(A)')
207 write(*,'(A)') '---------------Making local connectivity---------------'
208
209 do ip = 0, mpi_np-1
210
211 nelem_loc = 0
212 do in = 1,nelem
213 if(elem_domain(in) .eq. ip) nelem_loc = nelem_loc +1
214 enddo
215
217 do ie = 1, nelem
218 do j = 1, nmat
219 if ((tag_mat(j).eq. con_spx(con_spx(ie-1))) .and. (elem_domain(ie) .eq. ip)) then
220 nn = sdeg_mat(j) +1
222 endif
223 enddo
224 enddo
225
226 allocate(con_spx_loc(0:con_nnz_loc))
227
231
232
233
234
235 file_mpi = 'cons000000.mpi'
236 unit_mpi = 40
237 if (ip.lt. 10) then
238 write(file_mpi(10:10),'(i1)') ip
239 else if (ip .lt. 100) then
240 write(file_mpi(9:10),'(i2)') ip
241 else if (ip .lt. 1000) then
242 write(file_mpi(8:10),'(i3)') ip
243 else if (ip .lt. 10000) then
244 write(file_mpi(7:10),'(i4)') ip
245 else if (ip .lt. 100000) then
246 write(file_mpi(6:10),'(i5)') ip
247 else if (ip .lt. 1000000) then
248 write(file_mpi(5:10),'(i6)') ip
249 endif
250
251 if(len_trim(mpi_file) .ne. 70) then
252 file_mpi_new = mpi_file(1:len_trim(mpi_file)) // '/' // file_mpi
253 else
255 endif
256
257 open(unit_mpi,file=file_mpi_new)
258 write(unit_mpi,*) con_nnz_loc
259 do i = 0, con_nnz_loc
260 write(unit_mpi,*) con_spx_loc(i)
261 enddo
262 close(unit_mpi)
263
264 ! We already have nelem_loc...do we need this again?
265 ic = 0
266 do in = 1,nelem
267 if (elem_domain(in) .eq. ip) ic = ic +1
268 enddo
269 nelem_loc = ic
270
271 allocate(local_el_num(nelem_loc))
272
273 ic = 0
274 do in = 1,nelem
275 if (elem_domain(in) .eq. ip) then
276 ic = ic +1
278 endif
279 enddo
280
282 nface_loc = 0
283
284 ic = 0
285 if (nface.gt.0) then
286
287 do i = 1,nface
288
291
292 if (ic .ne. 0) then !I found global index in local_el_num
293
294 do j = 1,nmat
295
296 if (tag_mat(j) .eq. con_spx_loc(con_spx_loc(ic-1))) then
297 nn = sdeg_mat(j) +1
299 nface_loc = nface_loc + 1
300 endif
301
302 enddo
303 endif
304
305 enddo
306
308
309 allocate(con_spx_bc_loc(0:con_nnz_bc_loc))
311
315 nface_loc)
316
317 endif
318
319
320 file_mpi = 'conb000000.mpi'
321 unit_mpi = 40
322 if (ip.lt. 10) then
323 write(file_mpi(10:10),'(i1)') ip
324 else if (ip .lt. 100) then
325 write(file_mpi(9:10),'(i2)') ip
326 else if (ip .lt. 1000) then
327 write(file_mpi(8:10),'(i3)') ip
328 else if (ip .lt. 10000) then
329 write(file_mpi(7:10),'(i4)') ip
330 else if (ip .lt. 100000) then
331 write(file_mpi(6:10),'(i5)') ip
332 else if (ip .lt. 1000000) then
333 write(file_mpi(5:10),'(i6)') ip
334 endif
335
336 if(len_trim(mpi_file) .ne. 70) then
337 file_mpi_new = mpi_file(1:len_trim(mpi_file)) // '/' // file_mpi
338 else
340 endif
341
342 open(unit_mpi,file=file_mpi_new)
343 write(unit_mpi,*) con_nnz_bc_loc
344 do i = 0, con_nnz_bc_loc
345 write(unit_mpi,*) con_spx_bc_loc(i)
346 enddo
347 close(unit_mpi)
348
349
351
352 enddo
353
354 deallocate(con_spx,node_pointer)
355
356 endif
357
358 endif
359
360 if(mpi_id.eq.0) write(*,'(A)') 'Made.'
361 call mpi_barrier(mpi_comm, mpi_ierr)
362
363!***************************************************************************************************************
364! LOCAL NUMERATION FOR ELEMENTS
365
366 ! This part runs in all processors.
367 ic = 0
368 do in = 1,nelem
369
370 if (elem_domain(in) .eq. mpi_id) then
371 ic = ic +1
372 endif
373 enddo
374 nelem_loc = ic
375
376 allocate(local_el_num(nelem_loc))
377
378 ic = 0
379 do in = 1,nelem
380
381 if (elem_domain(in) .eq. mpi_id) then
382 ic = ic +1
384 endif
385 enddo
386
387 call mpi_barrier(mpi_comm, mpi_ierr)
388
389! LOCAL NUMERATION FOR ELEMENTS END
390!***************************************************************************************************************
391
392 file_mpi = 'cons000000.mpi'; unit_mpi = 40
393
394 if (mpi_id .lt. 10) then
395 write(file_mpi(10:10),'(i1)') mpi_id
396 else if (mpi_id .lt. 100) then
397 write(file_mpi(9:10),'(i2)') mpi_id
398 else if (mpi_id .lt. 1000) then
399 write(file_mpi(8:10),'(i3)') mpi_id
400 else if (mpi_id .lt. 10000) then
401 write(file_mpi(7:10),'(i4)') mpi_id
402 else if (mpi_id .lt. 100000) then
403 write(file_mpi(6:10),'(i5)') mpi_id
404 else if (mpi_id .lt. 1000000) then
405 write(file_mpi(5:10),'(i6)') mpi_id
406 endif
407
408 if(len_trim(mpi_file) .ne. 70) then
409 file_mpi_new = mpi_file(1:len_trim(mpi_file)) // '/' // file_mpi
410 else
412 endif
413
414
415 open(unit_mpi,file=file_mpi_new)
416 read(unit_mpi,*) con_nnz_loc
417 allocate(con_spx_loc(0:con_nnz_loc))
418
419
420 do i = 0, con_nnz_loc
421 read(unit_mpi,*) con_spx_loc(i)
422 enddo
423 close(unit_mpi)
424
425
426 file_mpi = 'conb000000.mpi'; unit_mpi = 40
427
428 if (mpi_id .lt. 10) then
429 write(file_mpi(10:10),'(i1)') mpi_id
430 else if (mpi_id .lt. 100) then
431 write(file_mpi(9:10),'(i2)') mpi_id
432 else if (mpi_id .lt. 1000) then
433 write(file_mpi(8:10),'(i3)') mpi_id
434 else if (mpi_id .lt. 10000) then
435 write(file_mpi(7:10),'(i4)') mpi_id
436 else if (mpi_id .lt. 100000) then
437 write(file_mpi(6:10),'(i5)') mpi_id
438 else if (mpi_id .lt. 1000000) then
439 write(file_mpi(5:10),'(i6)') mpi_id
440 endif
441
442 if(len_trim(mpi_file) .ne. 70) then
443 file_mpi_new = mpi_file(1:len_trim(mpi_file)) // '/' // file_mpi
444 else
446 endif
447
448 open(unit_mpi,file=file_mpi_new)
450 allocate(con_spx_bc_loc(0:con_nnz_bc_loc))
451
452
453 do i = 0, con_nnz_bc_loc
454 read(unit_mpi,*) con_spx_bc_loc(i)
455 enddo
456 close(unit_mpi)
457
458
459! SPECTRAL LOCAL CONNECTIVITY
460!***************************************************************************************************************
461!***************************************************************************************************************
462! COUNTING LOCAL NODES & LOCAL NODE NUMBERING
463 if(mpi_id .eq.0) write(*,'(A)')
464 if(mpi_id .eq.0) write(*,'(A)') '------------------Counting local nodes-----------------'
465
467
468 if(mpi_id .eq.0) write(*,'(A)') 'Made.'
469
470 allocate(local_node_num(nnod_loc))
472
473
474 if(mpi_id .eq.0) write(*,'(A)')
475 if(mpi_id .eq.0) write(*,'(A)') '--------------Make local nodes numbering---------------'
476
477
480
481! ic = 0
482!! do im = 1,nmat
483!! nn = sdeg_mat(im) +1
484! do ie = 1,nelem_loc
485! im = con_spx_loc(con_spx_loc(ie -1) +0); nn = sdeg_mat(im) +1
486!
487!! if (con_spx_loc(con_spx_loc(ie -1) +0) .eq. tag_mat(im)) then
488! do i = 1, nn*nn*nn
489! if (find(nnod_loc, local_node_num, con_spx_loc(con_spx_loc(ie -1) + i)) .ne. 1) then
490! ic = ic + 1
491! local_node_num(ic) = con_spx_loc(con_spx_loc(ie -1) + i)
492! endif
493! enddo
494!! endif
495! enddo
496!! enddo
497
498
499 if(mpi_id .eq.0) write(*,'(A)') 'Made.'
500 if(mpi_id .eq.0) write(*,'(A)')
501
502 file_mpi = 'nodedomain.mpi'; unit_mpi = 40
503 if(len_trim(mpi_file) .ne. 70) then
504 file_mpi_new = mpi_file(1:len_trim(mpi_file)) // '/' // file_mpi
505 else
507 endif
508
509 if(mpi_id .eq.0) then
510 !start1=MPI_WTIME()
511 open(unit_mpi,file=file_mpi_new)
512 read(unit_mpi,*) nnod
513 endif
514 call mpi_bcast(nnod,1,speed_integer,0,mpi_comm_world,mpi_ierr)
515
516 allocate(node_domain(nnod))
517 if(mpi_id .eq.0) then
518 do i = 1, nnod
519 read(unit_mpi,*) node_domain(i)
520 enddo
521 close(unit_mpi)
522 !write (*,*) "File nodedomain.mpi read, broadcast to slaves"
523 endif
524 call mpi_bcast(node_domain,nnod,speed_integer,0,mpi_comm_world,mpi_ierr)
525
526 !write (*,*) "Vector received by process",mpi_id
527
528
529 ic = 0
530 do in = 1, nnod
531 if (node_domain(in).eq. mpi_id) ic = ic +1
532 enddo
533 nnode_dom = ic
534
535 ! write(*,'(A,I3,I8)')'Nodes active on proc : ', mpi_id, ic
536 write(*,'(A,I6,I15)') 'Total nodes on Proc. : ', mpi_id, nnode_dom
537
538
539 call mpi_barrier(mpi_comm, mpi_ierr)
540
541 !if(mpi_id .eq.0) then
542 !start2=MPI_WTIME()
543 !write (*,*) "Time to send nodedomain.mpi to all &
544 !& proc",start2-start1,"sec"
545 !endif
546
547! COUNTING LOCAL NODES & LOCAL NODE NUMBERING
548!***********************************************************************************
549
550 end subroutine make_partition_and_mpi_files
subroutine get_elem_from_face(nb_nz_el, list_el, v1, v2, v3, ie)
Find element index from verteces number.
subroutine get_indloc_from_indglo(local_el, nel_loc, ie, ic)
Returns local id from global id.
subroutine get_loc_node_num(nnz_loc, cs_loc, nn_loc, id)
Computes number of local nodes.
subroutine make_grid_nodes(nb_node, nb_elem, con_matrix, node_wgt, nnz
Makes pointer for connectivity of mesh nodes.
subroutine make_ren_loc_node(nnz_loc, cs_loc, nn_loc, id, loc_n_nu
Renumerates the grid nodes.
subroutine make_spx_con(nelem, con_mac, nmat, tag_mat, sdeg, nnz_pntr, node_pntr, con_nnz, c
Makes spectral connectivity vector.
subroutine make_spx_con_loc(nel_loc, nel, el_dom, nnz, cs, nnz_loc, cs_loc, nm, tag_mat, sdeg_mat, mpi_id)
Makes local connectivity vector.
subroutine make_spx_con_loc_bound(cs_nnz, cs, ne_bc, cm_bc, nm, tm, sd, ennz, ebin, cs_nnz_bc, cs_bc, loc_el
Makes spectral connectivity vector for the boundary.
subroutine make_spx_nodes(nb_node, spx_con_nnz, spx_con, node_wgt, nnz
Makes a pointer for spectral connectivity vector.
subroutine make_wgt_grid_nodes(nb_node, nb_elem, con_matrix, node_wgt
Computes multeplicity for mesh nodes.
subroutine make_wgt_spx_nodes(nb_node, spx_con_nnz, spx_con, node_wgt
Computes multeplicity for spectral nodes.
subroutine mesh_partitioning(mpi_file, nelem, nnode, nparts, conn, w_yn
Makes mesh partitioning using METIS.
Set maximal bounds.
Definition MODULES.f90:54
Contains SPEED PARAMETERS used in (SPEED, READ_INPUT_FILES, MAKE_PARTION_AND_MPI_FILES,...
Definition MODULES.f90:207
integer *4, dimension(:), allocatable elem_index
Definition MODULES.f90:322
integer *4, dimension(:), allocatable con_spx
Definition MODULES.f90:322
integer *4 nmat
Definition MODULES.f90:269
integer *4 unit_mpi
Definition MODULES.f90:308
integer *4, dimension(:), allocatable local_el_num
Definition MODULES.f90:322
integer *4 nnod
Definition MODULES.f90:269
integer *4 con_nnz_bc_loc
Definition MODULES.f90:269
integer *4 nface_loc
Definition MODULES.f90:269
integer *4 ip
Definition MODULES.f90:263
character *70 file_mpi_new
Definition MODULES.f90:243
integer *4 j
Definition MODULES.f90:263
integer *4 trash
Definition MODULES.f90:308
integer *4 nelem_dom
Definition MODULES.f90:269
integer *4 ie
Definition MODULES.f90:263
integer *4 con_nnz_loc
Definition MODULES.f90:269
integer *4 in
Definition MODULES.f90:263
integer *4, dimension(:), allocatable node_domain
Definition MODULES.f90:322
character *70 mpi_file
Definition MODULES.f90:243
integer *4, dimension(:,:), allocatable con
Definition MODULES.f90:383
integer *4 nelem
Definition MODULES.f90:269
integer *4 i
Definition MODULES.f90:263
integer *4 ic
Definition MODULES.f90:263
integer *4 nface
Definition MODULES.f90:269
integer *4, dimension(:), allocatable con_spx_bc_loc
Definition MODULES.f90:322
integer *4 mpi_comm
Definition MODULES.f90:308
character *70 file_mpi
Definition MODULES.f90:243
integer *4, dimension(:), allocatable node_weight
Definition MODULES.f90:322
integer *4, dimension(:), allocatable local_node_num
Definition MODULES.f90:322
integer *4 mpi_ierr
Definition MODULES.f90:298
integer *4 nnz_node_weight
Definition MODULES.f90:269
integer *4, dimension(:), allocatable sdeg_mat
Definition MODULES.f90:335
integer *4 nelem_loc
Definition MODULES.f90:269
integer *4, dimension(:), allocatable con_spx_loc
Definition MODULES.f90:322
logical filefound
Definition MODULES.f90:217
character *70 file_part
Definition MODULES.f90:243
integer *4 nnod_macro
Definition MODULES.f90:269
integer *4, dimension(:), allocatable tag_mat
Definition MODULES.f90:335
integer *4 nnode_dom
Definition MODULES.f90:269
integer *4 mpi_np
Definition MODULES.f90:308
integer *4 mpi_id
Definition MODULES.f90:308
integer *4, dimension(:,:), allocatable con_bc
Definition MODULES.f90:383
integer *4 con_nnz
Definition MODULES.f90:269
integer *4 nnod_loc
Definition MODULES.f90:269
integer *4, dimension(:), allocatable node_pointer
Definition MODULES.f90:322
integer *4 nn
Definition MODULES.f90:269
integer *4, dimension(:), allocatable elem_domain
Definition MODULES.f90:322
integer *4 unit_part
Definition MODULES.f90:308