source: CIVL/examples/mpi-omp/AMG2013/parcsr_ls/par_difconv.c

main
Last change on this file was ea777aa, checked in by Alex Wilton <awilton@…>, 3 years ago

Moved examples, include, build_default.properties, common.xml, and README out from dev.civl.com into the root of the repo.

git-svn-id: ​svn://vsl.cis.udel.edu/civl/trunk@5704 fb995dde-84ed-4084-dfe6-e5aef3e2452c
  • Property mode set to 100644
File size: 9.4 KB
Line 
1/*BHEADER**********************************************************************
2 * Copyright (c) 2008, Lawrence Livermore National Security, LLC.
3 * Produced at the Lawrence Livermore National Laboratory.
4 * This file is part of HYPRE. See file COPYRIGHT for details.
5 *
6 * HYPRE is free software; you can redistribute it and/or modify it under the
7 * terms of the GNU Lesser General Public License (as published by the Free
8 * Software Foundation) version 2.1 dated February 1999.
9 *
10 * $Revision: 2.4 $
11 ***********************************************************************EHEADER*/
12
13
14
15
16#include "headers.h"
17
18/*--------------------------------------------------------------------------
19 * hypre_GenerateDifConv
20 *--------------------------------------------------------------------------*/
21
22HYPRE_ParCSRMatrix
23GenerateDifConv( MPI_Comm comm,
24 int nx,
25 int ny,
26 int nz,
27 int P,
28 int Q,
29 int R,
30 int p,
31 int q,
32 int r,
33 double *value )
34{
35 hypre_ParCSRMatrix *A;
36 hypre_CSRMatrix *diag;
37 hypre_CSRMatrix *offd;
38
39 int *diag_i;
40 int *diag_j;
41 double *diag_data;
42
43 int *offd_i;
44 int *offd_j;
45 double *offd_data;
46
47 int *global_part;
48 int ix, iy, iz;
49 int cnt, o_cnt;
50 int local_num_rows;
51 int *col_map_offd;
52 int row_index;
53 int i,j;
54
55 int nx_local, ny_local, nz_local;
56 int nx_size, ny_size, nz_size;
57 int num_cols_offd;
58 int grid_size;
59
60 int *nx_part;
61 int *ny_part;
62 int *nz_part;
63
64 int num_procs, my_id;
65 int P_busy, Q_busy, R_busy;
66
67 MPI_Comm_size(comm,&num_procs);
68 MPI_Comm_rank(comm,&my_id);
69
70 grid_size = nx*ny*nz;
71
72 hypre_GeneratePartitioning(nx,P,&nx_part);
73 hypre_GeneratePartitioning(ny,Q,&ny_part);
74 hypre_GeneratePartitioning(nz,R,&nz_part);
75
76 global_part = hypre_CTAlloc(int,P*Q*R+1);
77
78 global_part[0] = 0;
79 cnt = 1;
80 for (iz = 0; iz < R; iz++)
81 {
82 nz_size = nz_part[iz+1]-nz_part[iz];
83 for (iy = 0; iy < Q; iy++)
84 {
85 ny_size = ny_part[iy+1]-ny_part[iy];
86 for (ix = 0; ix < P; ix++)
87 {
88 nx_size = nx_part[ix+1] - nx_part[ix];
89 global_part[cnt] = global_part[cnt-1];
90 global_part[cnt++] += nx_size*ny_size*nz_size;
91 }
92 }
93 }
94
95 nx_local = nx_part[p+1] - nx_part[p];
96 ny_local = ny_part[q+1] - ny_part[q];
97 nz_local = nz_part[r+1] - nz_part[r];
98
99 my_id = r*(P*Q) + q*P + p;
100 num_procs = P*Q*R;
101
102 local_num_rows = nx_local*ny_local*nz_local;
103 diag_i = hypre_CTAlloc(int, local_num_rows+1);
104 offd_i = hypre_CTAlloc(int, local_num_rows+1);
105
106 P_busy = hypre_min(nx,P);
107 Q_busy = hypre_min(ny,Q);
108 R_busy = hypre_min(nz,R);
109
110 num_cols_offd = 0;
111 if (p) num_cols_offd += ny_local*nz_local;
112 if (p < P_busy-1) num_cols_offd += ny_local*nz_local;
113 if (q) num_cols_offd += nx_local*nz_local;
114 if (q < Q_busy-1) num_cols_offd += nx_local*nz_local;
115 if (r) num_cols_offd += nx_local*ny_local;
116 if (r < R_busy-1) num_cols_offd += nx_local*ny_local;
117
118 if (!local_num_rows) num_cols_offd = 0;
119
120 col_map_offd = hypre_CTAlloc(int, num_cols_offd);
121
122 cnt = 1;
123 o_cnt = 1;
124 diag_i[0] = 0;
125 offd_i[0] = 0;
126 for (iz = nz_part[r]; iz < nz_part[r+1]; iz++)
127 {
128 for (iy = ny_part[q]; iy < ny_part[q+1]; iy++)
129 {
130 for (ix = nx_part[p]; ix < nx_part[p+1]; ix++)
131 {
132 diag_i[cnt] = diag_i[cnt-1];
133 offd_i[o_cnt] = offd_i[o_cnt-1];
134 diag_i[cnt]++;
135 if (iz > nz_part[r])
136 diag_i[cnt]++;
137 else
138 {
139 if (iz)
140 {
141 offd_i[o_cnt]++;
142 }
143 }
144 if (iy > ny_part[q])
145 diag_i[cnt]++;
146 else
147 {
148 if (iy)
149 {
150 offd_i[o_cnt]++;
151 }
152 }
153 if (ix > nx_part[p])
154 diag_i[cnt]++;
155 else
156 {
157 if (ix)
158 {
159 offd_i[o_cnt]++;
160 }
161 }
162 if (ix+1 < nx_part[p+1])
163 diag_i[cnt]++;
164 else
165 {
166 if (ix+1 < nx)
167 {
168 offd_i[o_cnt]++;
169 }
170 }
171 if (iy+1 < ny_part[q+1])
172 diag_i[cnt]++;
173 else
174 {
175 if (iy+1 < ny)
176 {
177 offd_i[o_cnt]++;
178 }
179 }
180 if (iz+1 < nz_part[r+1])
181 diag_i[cnt]++;
182 else
183 {
184 if (iz+1 < nz)
185 {
186 offd_i[o_cnt]++;
187 }
188 }
189 cnt++;
190 o_cnt++;
191 }
192 }
193 }
194
195 diag_j = hypre_CTAlloc(int, diag_i[local_num_rows]);
196 diag_data = hypre_CTAlloc(double, diag_i[local_num_rows]);
197
198 if (num_procs > 1)
199 {
200 offd_j = hypre_CTAlloc(int, offd_i[local_num_rows]);
201 offd_data = hypre_CTAlloc(double, offd_i[local_num_rows]);
202 }
203
204 row_index = 0;
205 cnt = 0;
206 o_cnt = 0;
207 for (iz = nz_part[r]; iz < nz_part[r+1]; iz++)
208 {
209 for (iy = ny_part[q]; iy < ny_part[q+1]; iy++)
210 {
211 for (ix = nx_part[p]; ix < nx_part[p+1]; ix++)
212 {
213 diag_j[cnt] = row_index;
214 diag_data[cnt++] = value[0];
215 if (iz > nz_part[r])
216 {
217 diag_j[cnt] = row_index-nx_local*ny_local;
218 diag_data[cnt++] = value[3];
219 }
220 else
221 {
222 if (iz)
223 {
224 offd_j[o_cnt] = hypre_map(ix,iy,iz-1,p,q,r-1,P,Q,R,
225 nx_part,ny_part,nz_part,global_part, &col_map_offd[o_cnt]);
226 offd_data[o_cnt++] = value[3];
227 }
228 }
229 if (iy > ny_part[q])
230 {
231 diag_j[cnt] = row_index-nx_local;
232 diag_data[cnt++] = value[2];
233 }
234 else
235 {
236 if (iy)
237 {
238 offd_j[o_cnt] = hypre_map(ix,iy-1,iz,p,q-1,r,P,Q,R,
239 nx_part,ny_part,nz_part,global_part, &col_map_offd[o_cnt]);
240 offd_data[o_cnt++] = value[2];
241 }
242 }
243 if (ix > nx_part[p])
244 {
245 diag_j[cnt] = row_index-1;
246 diag_data[cnt++] = value[1];
247 }
248 else
249 {
250 if (ix)
251 {
252 offd_j[o_cnt] = hypre_map(ix-1,iy,iz,p-1,q,r,P,Q,R,
253 nx_part,ny_part,nz_part,global_part, &col_map_offd[o_cnt]);
254 offd_data[o_cnt++] = value[1];
255 }
256 }
257 if (ix+1 < nx_part[p+1])
258 {
259 diag_j[cnt] = row_index+1;
260 diag_data[cnt++] = value[4];
261 }
262 else
263 {
264 if (ix+1 < nx)
265 {
266 offd_j[o_cnt] = hypre_map(ix+1,iy,iz,p+1,q,r,P,Q,R,
267 nx_part,ny_part,nz_part,global_part, &col_map_offd[o_cnt]);
268 offd_data[o_cnt++] = value[4];
269 }
270 }
271 if (iy+1 < ny_part[q+1])
272 {
273 diag_j[cnt] = row_index+nx_local;
274 diag_data[cnt++] = value[5];
275 }
276 else
277 {
278 if (iy+1 < ny)
279 {
280 offd_j[o_cnt] = hypre_map(ix,iy+1,iz,p,q+1,r,P,Q,R,
281 nx_part,ny_part,nz_part,global_part, &col_map_offd[o_cnt]);
282 offd_data[o_cnt++] = value[5];
283 }
284 }
285 if (iz+1 < nz_part[r+1])
286 {
287 diag_j[cnt] = row_index+nx_local*ny_local;
288 diag_data[cnt++] = value[6];
289 }
290 else
291 {
292 if (iz+1 < nz)
293 {
294 offd_j[o_cnt] = hypre_map(ix,iy,iz+1,p,q,r+1,P,Q,R,
295 nx_part,ny_part,nz_part,global_part, &col_map_offd[o_cnt]);
296 offd_data[o_cnt++] = value[6];
297 }
298 }
299 row_index++;
300 }
301 }
302 }
303
304 if (num_procs > 1)
305 {
306 for (i=0; i < num_cols_offd; i++)
307 col_map_offd[i] = offd_j[i];
308
309 qsort0(col_map_offd, 0, num_cols_offd-1);
310
311 for (i=0; i < num_cols_offd; i++)
312 for (j=0; j < num_cols_offd; j++)
313 if (offd_j[i] == col_map_offd[j])
314 {
315 offd_j[i] = j;
316 break;
317 }
318 }
319
320 A = hypre_ParCSRMatrixCreate(comm, grid_size, grid_size,
321 global_part, global_part, num_cols_offd,
322 diag_i[local_num_rows],
323 offd_i[local_num_rows]);
324
325 hypre_ParCSRMatrixColMapOffd(A) = col_map_offd;
326
327 diag = hypre_ParCSRMatrixDiag(A);
328 hypre_CSRMatrixI(diag) = diag_i;
329 hypre_CSRMatrixJ(diag) = diag_j;
330 hypre_CSRMatrixData(diag) = diag_data;
331
332 offd = hypre_ParCSRMatrixOffd(A);
333 hypre_CSRMatrixI(offd) = offd_i;
334 if (num_cols_offd)
335 {
336 hypre_CSRMatrixJ(offd) = offd_j;
337 hypre_CSRMatrixData(offd) = offd_data;
338 }
339
340 hypre_TFree(nx_part);
341 hypre_TFree(ny_part);
342 hypre_TFree(nz_part);
343
344 return (HYPRE_ParCSRMatrix) A;
345}
346
Note: See TracBrowser for help on using the repository browser.