Neko 1.99.1
A portable framework for high-order spectral element flow simulations
Loading...
Searching...
No Matches
parmetis.F90
Go to the documentation of this file.
1! Copyright (c) 2021, The Neko Authors
2! All rights reserved.
3!
4! Redistribution and use in source and binary forms, with or without
5! modification, are permitted provided that the following conditions
6! are met:
7!
8! * Redistributions of source code must retain the above copyright
9! notice, this list of conditions and the following disclaimer.
10!
11! * Redistributions in binary form must reproduce the above
12! copyright notice, this list of conditions and the following
13! disclaimer in the documentation and/or other materials provided
14! with the distribution.
15!
16! * Neither the name of the authors nor the names of its
17! contributors may be used to endorse or promote products derived
18! from this software without specific prior written permission.
19!
20! THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
21! "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
22! LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
23! FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
24! COPYRIGHT OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
25! INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
26! BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
27! LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
28! CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
29! LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
30! ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
31! POSSIBILITY OF SUCH DAMAGE.
32!
35 use comm, only : pe_size, pe_rank
36 use mpi_f08, only : mpi_allgather, mpi_integer
37 use point, only : point_t
38 use utils, only : neko_error
39 use num_types, only : rp
40 use mesh_field, only : mesh_fld_t
41 use mesh, only : mesh_t
42 use, intrinsic :: iso_c_binding
43 implicit none
44 private
45
46 integer, private, parameter :: metis_ok = 1, metis_error = -4
47
49
50 interface
51 integer (c_int) function parmetis_v3_partmeshkway &
52 (elmdist, eptr, eind, elmwgt, wgtflag, numflag, ncon, &
53 ncommonnodes, nparts, tpwgts, ubvec, options, edgecut, part) &
54 bind(c, name='ParMETIS_V3_PartMeshKway_wrapper')
55 use, intrinsic :: iso_c_binding
56 implicit none
57 ! idx_t variables
58 type(c_ptr), value :: elmdist, eptr, eind, elmwgt, wgtflag, &
59 numflag, ncon, ncommonnodes, nparts, options, edgecut, part
60 ! real_t variables
61 type(c_ptr), value :: tpwgts, ubvec
62 end function parmetis_v3_partmeshkway
63 end interface
64
65 interface
66 integer (c_int) function parmetis_v3_partgeom &
67 (vtxdist, ndims, xyz, part) &
68 bind(c, name='ParMETIS_V3_PartGeom_wrapper')
69 use, intrinsic :: iso_c_binding
70 implicit none
71 ! idx_t variables
72 type(c_ptr), value :: vtxdist, ndims, part
73 ! real_t variables
74 type(c_ptr), value :: xyz
75 end function parmetis_v3_partgeom
76 end interface
77
78!
79! Define data types depending on ParMETIS' configuration
80!
81! parmetis_real converts between NEKO and ParMETIS' real representation
82!
83! parmetis_idx converts between NEKO and ParMETIS' integer representation
84! neko_idx converts between ParMETIS and NEKO's integer representation
85!
86
87#ifdef HAVE_PARMETIS
88#ifdef HAVE_PARMETIS_REAL64
89#define M_REAL c_double
90#define parmetis_real(i) real((i), i8)
91#else
92#define M_REAL c_float
93#define parmetis_real(i) real((i), 4)
94#endif
95#ifdef HAVE_PARMETIS_INT64
96#define M_INT c_int64_t
97#define parmetis_idx(i) int((i), i8)
98#define neko_idx(i) int((i), 4)
99#else
100#define M_INT c_int32_t
101#define parmetis_idx(i) (i)
102#define neko_idx(i) (i)
103#endif
104#endif
105
106contains
107
108#ifdef HAVE_PARMETIS
109
111 subroutine parmetis_partmeshkway(msh, parts, weights, nprts)
112 type(mesh_t), intent(inout) :: msh
113 type(mesh_fld_t), intent(inout) :: parts
114 type(mesh_fld_t), intent(in), optional :: weights
115 integer, intent(in), optional :: nprts
116 integer(kind=M_INT), target :: wgtflag, numflag, ncon, ncommonnodes
117 integer(kind=M_INT), target :: nparts, options(3), edgecut, rcode
118 real(kind=m_real), allocatable, target, dimension(:) :: tpwgts, ubvec
119 integer(kind=M_INT), allocatable, target, dimension(:) :: &
120 elmdist, eptr, eind, elmwgt, part
121 integer :: i, j, k, ierr
122
126 numflag = 0
127 ncon = 1
128 ncommonnodes = 2**(msh%gdim - 1)
129 options(1) = 1
130 options(2) = 1
131 options(3) = 15
132 wgtflag = 2
133
134 if (present(nprts)) then
135 nparts = nprts
136 else
137 nparts = pe_size
138 end if
139
140 allocate(elmdist(0:pe_size), eptr(0:msh%nelv))
141 allocate(eind(0:(msh%nelv * msh%npts)), part(msh%nelv))
142 allocate(elmwgt(msh%nelv), tpwgts(ncon * nparts), ubvec(ncon))
143
144 call parmetis_dist(elmdist, msh%nelv)
145
146 if (present(weights)) then
147 call parmetis_wgt(msh, elmwgt, tpwgts, ubvec, nparts, ncon, weights)
148 else
149 call parmetis_wgt(msh, elmwgt, tpwgts, ubvec, nparts, ncon)
150 end if
151
152 eptr(0) = 0
153 do i = 1, msh%nelv
154 eptr(i) = parmetis_idx(eptr(i - 1) + msh%npts)
155 end do
156
157 k = 0
158 do i = 1, msh%nelv
159 do j = 1, msh%npts
160 eind(k) = parmetis_idx(msh%elements(i)%e%pts(j)%p%id() - 1)
161 k = k + 1
162 end do
163 end do
164
165 rcode = parmetis_v3_partmeshkway(c_loc(elmdist), c_loc(eptr), c_loc(eind), &
166 c_loc(elmwgt), c_loc(wgtflag), c_loc(numflag), c_loc(ncon), &
167 c_loc(ncommonnodes), c_loc(nparts), c_loc(tpwgts), c_loc(ubvec),&
168 c_loc(options), c_loc(edgecut), c_loc(part))
169
170 if (rcode .eq. metis_ok) then
171 call parmetis_mark_parts(parts, msh, part)
172 else
173 call neko_error(rcode)
174 end if
175
176 deallocate(elmdist, eptr, eind, part, elmwgt, tpwgts, ubvec)
177
178 end subroutine parmetis_partmeshkway
179
182 subroutine parmetis_partgeom(msh, parts)
183 type(mesh_t), intent(inout) :: msh
184 type(mesh_fld_t), intent(inout) :: parts
185 integer(kind=M_INT), target :: ndims
186 real(kind=m_real), allocatable, target, dimension(:) :: xyz
187 integer(kind=M_INT), allocatable, target, dimension(:) :: vtxdist, part
188 type(point_t) :: c
189 integer :: i, j, ierr, rcode
190
191 ndims = msh%gdim
192
193 allocate(part(msh%nelv), xyz(ndims * msh%nelv))
194 allocate(vtxdist(0:pe_size))
195
196 call parmetis_dist(vtxdist, msh%nelv)
197
198 i = 1
199 do j = 1, msh%nelv
200 c = msh%elements(j)%e%centroid()
201 xyz(i) = parmetis_real(c%x(1))
202 xyz(i + 1) = parmetis_real(c%x(2))
203 xyz(i + 2) = parmetis_real(c%x(3))
204 i = i + 3
205 end do
206
207 rcode = parmetis_v3_partgeom(c_loc(vtxdist), c_loc(ndims), &
208 c_loc(xyz), c_loc(part))
209
210 if (rcode .eq. metis_ok) then
211 call parmetis_mark_parts(parts, msh, part)
212 else
213 call neko_error(rcode)
214 end if
215
216 deallocate(part, xyz, vtxdist)
217
218 end subroutine parmetis_partgeom
219
221 subroutine parmetis_mark_parts(parts, msh, part)
222 type(mesh_fld_t), intent(inout) :: parts
223 type(mesh_t), intent(in) :: msh
224 integer(kind=M_INT), allocatable, intent(in) :: part(:)
225 integer :: i
226
227 call mesh_field_init(parts, msh, 'partitions')
228
229 do i = 1, msh%nelv
230 parts%data(i) = neko_idx(part(i))
231 end do
232
233 end subroutine parmetis_mark_parts
234
236 subroutine parmetis_wgt(msh, wgt, tpwgts, ubvec, nparts, ncon, weight)
237 type(mesh_t), intent(in) :: msh
238 integer(kind=M_INT), allocatable, intent(inout) :: wgt(:)
239 real(kind=m_real), allocatable, intent(inout) :: tpwgts(:)
240 real(kind=m_real), allocatable, intent(inout) :: ubvec(:)
241 integer, intent(in) :: nparts, ncon
242 type(mesh_fld_t), intent(in), optional :: weight
243 integer :: i
244
245 if (present(weight)) then
246 do i = 1, msh%nelv
247 wgt(i) = parmetis_idx(weight%data(i))
248 end do
249 else
250 wgt = parmetis_idx(1)
251 end if
252
253 do i = 1, (ncon * nparts)
254 tpwgts(i) = parmetis_real(1) / parmetis_real(nparts)
255 end do
256
257 do i = 1, ncon
258 ubvec(i) = parmetis_real(1.05d0)
259 end do
260
261 end subroutine parmetis_wgt
262
264 subroutine parmetis_dist(dist, nelv)
265 integer(kind=M_INT), intent(inout) :: dist(0:pe_size)
266 integer, intent(in) :: nelv
267 integer(kind=M_INT) :: tmp, sum
268 integer :: i, ierr
269
270 dist(pe_rank) = parmetis_idx(nelv)
271
272 call mpi_allgather(nelv, 1, mpi_integer, dist, 1, &
273 mpi_integer, neko_comm, ierr)
274
275 sum = dist(0)
276 dist(0) = 0
277 do i = 1, pe_size
278 tmp = dist(i)
279 dist(i) = sum
280 sum = tmp + sum
281 end do
282
283 end subroutine parmetis_dist
284
285#else
286
288 subroutine parmetis_partmeshkway(msh, parts, weights, nprts)
289 type(mesh_t), intent(inout) :: msh
290 type(mesh_fld_t), intent(inout) :: parts
291 type(mesh_fld_t), intent(in), optional :: weights
292 integer, intent(in), optional :: nprts
293 call neko_error('NEKO needs to be built with ParMETIS support')
294 end subroutine parmetis_partmeshkway
295
298 subroutine parmetis_partgeom(msh, parts)
299 type(mesh_t), intent(inout) :: msh
300 type(mesh_fld_t), intent(inout) :: parts
301 call neko_error('NEKO needs to be built with ParMETIS support')
302 end subroutine parmetis_partgeom
303
304#endif
305
306end module parmetis
Definition comm.F90:1
integer, public pe_size
MPI size of communicator.
Definition comm.F90:58
integer, public pe_rank
MPI rank.
Definition comm.F90:55
type(mpi_comm), public neko_comm
MPI communicator.
Definition comm.F90:42
Defines a mesh field.
subroutine, public mesh_field_init(fld, msh, fld_name)
Defines a mesh.
Definition mesh.f90:34
integer, parameter, public rp
Global precision used in computations.
Definition num_types.f90:12
Interface to ParMETIS.
Definition parmetis.F90:34
subroutine parmetis_mark_parts(parts, msh, part)
Fill mesh field according to new partitions.
Definition parmetis.F90:222
integer, parameter, private metis_error
Definition parmetis.F90:46
subroutine parmetis_dist(dist, nelv)
Compute the (parallel) vertex distribution of the dual graph.
Definition parmetis.F90:265
subroutine, public parmetis_partmeshkway(msh, parts, weights, nprts)
Compute a k-way partitioning of a mesh msh.
Definition parmetis.F90:112
integer, parameter, private metis_ok
Definition parmetis.F90:46
subroutine parmetis_wgt(msh, wgt, tpwgts, ubvec, nparts, ncon, weight)
Setup weights and balance constraints for the dual graph.
Definition parmetis.F90:237
subroutine, public parmetis_partgeom(msh, parts)
Compute a k-way partitioning of a mesh msh using a coordinated-based space-filing curves method.
Definition parmetis.F90:183
Implements a point.
Definition point.f90:35
Utilities.
Definition utils.f90:35
A point in with coordinates .
Definition point.f90:43