Actual source code: tfs.h
1: #pragma once
2: #if defined(__GNUC__) || defined(__clang__)
3: #pragma GCC diagnostic ignored "-Wconversion"
4: #endif
6: /**********************************const.h*************************************
8: Author: Henry M. Tufo III
10: e-mail: hmt@cs.brown.edu
12: snail-mail:
13: Division of Applied Mathematics
14: Brown University
15: Providence, RI 02912
17: Last Modification:
18: 6.21.97
19: ***********************************const.h************************************/
21: /**********************************const.h*************************************
22: File Description:
23: -----------------
25: ***********************************const.h************************************/
26: #include <petscsys.h>
27: #include <petscblaslapack.h>
29: #define X 0
30: #define Y 1
31: #define Z 2
32: #define XY 3
33: #define XZ 4
34: #define YZ 5
36: #define THRESH 0.2
37: #define N_HALF 4096
38: #define PRIV_BUF_SZ 45
40: /*4096 8192 32768 65536 1048576 */
41: #define MAX_MSG_BUF 32768
43: #define FULL 2
44: #define PARTIAL 1
45: #define NONE 0
47: #define BYTE 8
48: #define BIT_0 0x1
49: #define BIT_1 0x2
50: #define BIT_2 0x4
51: #define BIT_3 0x8
52: #define BIT_4 0x10
53: #define BIT_5 0x20
54: #define BIT_6 0x40
55: #define BIT_7 0x80
56: #define TOP_BIT PETSC_INT_MIN
58: #define C 0
60: #define MAX_VEC 1674
61: #define FORMAT 30
62: #define MAX_COL_LEN 100
63: #define MAX_LINE FORMAT *MAX_COL_LEN
64: #define DELIM " \n \t"
65: #define LINE 12
66: #define C_LINE 80
68: #define UT 5 /* dump upper 1/2 */
69: #define LT 6 /* dump lower 1/2 */
70: #define SYMM 8 /* we assume symm and dump upper 1/2 */
71: #define NON_SYMM 9
73: #define ROW 10
74: #define COL 11
76: #define EPS 1.0e-14
77: #define EPS2 1.0e-07
79: #define MPI 1
80: #define NX 2
82: #define SWAP(a, b) \
83: do { \
84: temp = (a); \
85: (a) = (b); \
86: (b) = temp; \
87: } while (0)
88: #define P_SWAP(a, b) \
89: do { \
90: ptr = (a); \
91: (a) = (b); \
92: (b) = ptr; \
93: } while (0)
95: #define MAX_FABS(x, y) (PetscAbsScalar(x) > PetscAbsScalar(y)) ? ((PetscScalar)x) : ((PetscScalar)y)
96: #define MIN_FABS(x, y) (PetscAbsScalar(x) < PetscAbsScalar(y)) ? ((PetscScalar)x) : ((PetscScalar)y)
98: /* specer's existence ... can be done w/MAX_ABS */
99: #define EXISTS(x, y) ((x) == 0.0) ? (y) : (x)
101: #define MULT_NEG_ONE(a) (a) *= -1;
102: #define NEG(a) (a) |= BIT_31;
103: #define POS(a) (a) &= INT_MAX;
105: /**********************************types.h*************************************
107: Author: Henry M. Tufo III
109: e-mail: hmt@cs.brown.edu
111: snail-mail:
112: Division of Applied Mathematics
113: Brown University
114: Providence, RI 02912
116: Last Modification:
117: 6.21.97
118: ***********************************types.h************************************/
120: typedef PetscErrorCode (*vfp)(void *, void *, PetscInt, ...);
121: typedef PetscErrorCode (*rbfp)(PetscScalar *, PetscScalar *, PetscInt);
122: typedef PetscInt (*bfp)(void *, void *, PetscInt *, MPI_Datatype *);
124: /***********************************comm.h*************************************
126: Author: Henry M. Tufo III
128: e-mail: hmt@cs.brown.edu
130: snail-mail:
131: Division of Applied Mathematics
132: Brown University
133: Providence, RI 02912
135: Last Modification:
136: 6.21.97
137: ***********************************comm.h*************************************/
138: PETSC_INTERN PetscMPIInt PCTFS_my_id;
139: PETSC_INTERN PetscMPIInt PCTFS_num_nodes;
140: PETSC_INTERN PetscMPIInt PCTFS_floor_num_nodes;
141: PETSC_INTERN PetscMPIInt PCTFS_i_log2_num_nodes;
143: PETSC_INTERN PetscErrorCode PCTFS_giop(PetscInt *, PetscInt *, PetscInt, PetscInt *);
144: PETSC_INTERN PetscErrorCode PCTFS_grop(PetscScalar *, PetscScalar *, PetscInt, PetscInt *);
145: PETSC_INTERN PetscErrorCode PCTFS_comm_init(void);
146: PETSC_INTERN PetscErrorCode PCTFS_giop_hc(PetscInt *, PetscInt *, PetscInt, PetscInt *, PetscInt);
147: PETSC_INTERN PetscErrorCode PCTFS_grop_hc(PetscScalar *, PetscScalar *, PetscInt, PetscInt *, PetscInt);
148: PETSC_INTERN PetscErrorCode PCTFS_ssgl_radd(PetscScalar *, PetscScalar *, PetscInt, PetscInt *);
150: #define MSGTAG0 101
151: #define MSGTAG1 1001
152: #define MSGTAG2 76207
153: #define MSGTAG3 100001
154: #define MSGTAG4 163841
155: #define MSGTAG5 249439
156: #define MSGTAG6 10000001
158: #define NON_UNIFORM 0
159: #define GL_MAX 1
160: #define GL_MIN 2
161: #define GL_MULT 3
162: #define GL_ADD 4
163: #define GL_B_XOR 5
164: #define GL_B_OR 6
165: #define GL_B_AND 7
166: #define GL_L_XOR 8
167: #define GL_L_OR 9
168: #define GL_L_AND 10
169: #define GL_MAX_ABS 11
170: #define GL_MIN_ABS 12
171: #define GL_EXISTS 13
173: PETSC_INTERN PetscInt *PCTFS_ivec_copy(PetscInt *, PetscInt *, PetscInt);
174: PETSC_INTERN PetscErrorCode PCTFS_ivec_zero(PetscInt *, PetscInt);
175: PETSC_INTERN PetscErrorCode PCTFS_ivec_set(PetscInt *, PetscInt, PetscInt);
177: PETSC_INTERN PetscInt PCTFS_ivec_lb(PetscInt *, PetscInt);
178: PETSC_INTERN PetscInt PCTFS_ivec_ub(PetscInt *, PetscInt);
179: PETSC_INTERN PetscInt PCTFS_ivec_sum(PetscInt *, PetscInt);
180: PETSC_INTERN vfp PCTFS_ivec_fct_addr(PetscInt);
182: PETSC_INTERN PetscErrorCode PCTFS_ivec_non_uniform(PetscInt *, PetscInt *, PetscInt, ...);
183: PETSC_INTERN PetscErrorCode PCTFS_ivec_max(PetscInt *, PetscInt *, PetscInt);
184: PETSC_INTERN PetscErrorCode PCTFS_ivec_min(PetscInt *, PetscInt *, PetscInt);
185: PETSC_INTERN PetscErrorCode PCTFS_ivec_mult(PetscInt *, PetscInt *, PetscInt);
186: PETSC_INTERN PetscErrorCode PCTFS_ivec_add(PetscInt *, PetscInt *, PetscInt);
187: PETSC_INTERN PetscErrorCode PCTFS_ivec_xor(PetscInt *, PetscInt *, PetscInt);
188: PETSC_INTERN PetscErrorCode PCTFS_ivec_or(PetscInt *, PetscInt *, PetscInt);
189: PETSC_INTERN PetscErrorCode PCTFS_ivec_and(PetscInt *, PetscInt *, PetscInt);
190: PETSC_INTERN PetscErrorCode PCTFS_ivec_lxor(PetscInt *, PetscInt *, PetscInt);
191: PETSC_INTERN PetscErrorCode PCTFS_ivec_lor(PetscInt *, PetscInt *, PetscInt);
192: PETSC_INTERN PetscErrorCode PCTFS_ivec_land(PetscInt *, PetscInt *, PetscInt);
193: PETSC_INTERN PetscErrorCode PCTFS_ivec_and3(PetscInt *, PetscInt *, PetscInt *, PetscInt);
195: PETSC_INTERN PetscErrorCode PCTFS_ivec_sort_companion(PetscInt *, PetscInt *, PetscInt);
196: PETSC_INTERN PetscErrorCode PCTFS_ivec_sort(PetscInt *, PetscInt);
197: PETSC_INTERN PetscErrorCode PCTFS_SMI_sort(void *, void *, PetscInt, PetscInt);
198: PETSC_INTERN PetscInt PCTFS_ivec_binary_search(PetscInt, PetscInt *, PetscInt);
199: PETSC_INTERN PetscInt PCTFS_ivec_linear_search(PetscInt, PetscInt *, PetscInt);
201: PETSC_INTERN PetscErrorCode PCTFS_ivec_sort_companion_hack(PetscInt *, PetscInt **, PetscInt);
203: #define SORT_INTEGER 1
204: #define SORT_INT_PTR 2
206: PETSC_INTERN PetscErrorCode PCTFS_rvec_zero(PetscScalar *, PetscInt);
207: PETSC_INTERN PetscErrorCode PCTFS_rvec_one(PetscScalar *, PetscInt);
208: PETSC_INTERN PetscErrorCode PCTFS_rvec_set(PetscScalar *, PetscScalar, PetscInt);
209: PETSC_INTERN PetscErrorCode PCTFS_rvec_copy(PetscScalar *, PetscScalar *, PetscInt);
210: PETSC_INTERN PetscErrorCode PCTFS_rvec_scale(PetscScalar *, PetscScalar, PetscInt);
212: PETSC_INTERN vfp PCTFS_rvec_fct_addr(PetscInt);
213: PETSC_INTERN PetscErrorCode PCTFS_rvec_add(PetscScalar *, PetscScalar *, PetscInt);
214: PETSC_INTERN PetscErrorCode PCTFS_rvec_mult(PetscScalar *, PetscScalar *, PetscInt);
215: PETSC_INTERN PetscErrorCode PCTFS_rvec_max(PetscScalar *, PetscScalar *, PetscInt);
216: PETSC_INTERN PetscErrorCode PCTFS_rvec_max_abs(PetscScalar *, PetscScalar *, PetscInt);
217: PETSC_INTERN PetscErrorCode PCTFS_rvec_min(PetscScalar *, PetscScalar *, PetscInt);
218: PETSC_INTERN PetscErrorCode PCTFS_rvec_min_abs(PetscScalar *, PetscScalar *, PetscInt);
219: PETSC_INTERN PetscErrorCode PCTFS_vec_exists(PetscScalar *, PetscScalar *, PetscInt);
221: /***********************************gs.h***************************************
223: Author: Henry M. Tufo III
225: e-mail: hmt@cs.brown.edu
227: snail-mail:
228: Division of Applied Mathematics
229: Brown University
230: Providence, RI 02912
232: Last Modification:
233: 6.21.97
234: ************************************gs.h**************************************/
236: typedef struct gather_scatter_id *PCTFS_gs_ADT;
238: PETSC_INTERN PCTFS_gs_ADT PCTFS_gs_init(PetscInt *, PetscInt, PetscInt);
239: PETSC_INTERN PetscErrorCode PCTFS_gs_gop_vec(PCTFS_gs_ADT, PetscScalar *, const char *, PetscInt);
240: PETSC_INTERN PetscErrorCode PCTFS_gs_gop_hc(PCTFS_gs_ADT, PetscScalar *, const char *, PetscInt);
241: PETSC_INTERN PetscErrorCode PCTFS_gs_free(PCTFS_gs_ADT);
242: PETSC_INTERN PetscErrorCode PCTFS_gs_init_msg_buf_sz(PetscInt);
243: PETSC_INTERN PetscErrorCode PCTFS_gs_init_vec_sz(PetscInt);
245: /*************************************xxt.h************************************
246: Module Name: xxt
247: Module Info: need xxt.{c,h} gs.{c,h} comm.{c,h} ivec.{c,h} error.{c,h}
249: author: Henry M. Tufo III
250: e-mail: hmt@asci.uchicago.edu
251: contact:
252: +--------------------------------+--------------------------------+
253: |MCS Division - Building 221 |Department of Computer Science |
254: |Argonne National Laboratory |Ryerson 152 |
255: |9700 S. Cass Avenue |The University of Chicago |
256: |Argonne, IL 60439 |Chicago, IL 60637 |
257: |(630) 252-5354/5986 ph/fx |(773) 702-6019/8487 ph/fx |
258: +--------------------------------+--------------------------------+
260: Last Modification: 3.20.01
261: **************************************xxt.h***********************************/
263: typedef struct xxt_CDT *xxt_ADT;
265: /*************************************xxt.h************************************
266: Function: XXT_new()
268: Return: ADT ptr or NULL upon failure.
269: Description: This function allocates and returns an xxt handle
270: Usage: xxt_handle = xxt_new();
271: **************************************xxt.h***********************************/
272: PETSC_INTERN xxt_ADT XXT_new(void);
274: /*************************************xxt.h************************************
275: Function: XXT_free()
277: Input : pointer to ADT.
279: Description: This function frees the storage associated with an xxt handle
280: Usage: XXT_free(xxt_handle);
281: **************************************xxt.h***********************************/
282: PETSC_INTERN PetscErrorCode XXT_free(xxt_ADT);
284: /*************************************xxt.h************************************
285: Function: XXT_factor
287: Input : ADT ptr, and pointer to object
288: Return: 0 on failure, 1 on success
289: Description: This function sets the xxt solver
291: xxt assumptions: given n rows of global coarse matrix (E_loc) where
292: o global dofs N = sum_p(n), p=0,P-1
293: (i.e. row dist. with no dof replication)
294: (5.21.00 will handle dif replication case)
295: o m is the number of columns in E_loc (m>=n)
296: o local2global holds global number of column i (i=0,...,m-1)
297: o local2global holds global number of row i (i=0,...,n-1)
298: o mylocmatvec performs E_loc . x_loc where x_loc is an vector of
299: length m in 1-1 correspondence with local2global
300: (note that gs package takes care of communication).
301: (note do not zero out upper m-n entries!)
302: o mylocmatvec(void *grid_data, double *in, double *out)
304: ML beliefs/usage: move this to ML_XXT_factor routine
305: o my_ml holds address of ML struct associated w/E_loc, grid_data, grid_tag
306: o grid_tag, grid_data, my_ml used in
307: ML_Set_CSolve(my_ml, grid_tag, grid_data, ML_Do_CoarseDirect);
308: o grid_data used in
309: A_matvec(grid_data,v,u);
311: Usage:
312: **************************************xxt.h***********************************/
313: PETSC_INTERN PetscErrorCode XXT_factor(xxt_ADT, /* prev. allocated xxt handle */
314: PetscInt *, /* global column mapping */
315: PetscInt, /* local num rows */
316: PetscInt, /* local num cols */
317: PetscErrorCode (*)(void *, PetscScalar *, PetscScalar *), /* b_loc=A_local.x_loc */
318: void *); /* grid data for matvec */
320: /*************************************xxt.h************************************
321: Function: XXT_solve
323: Input : ADT ptr, b (rhs)
324: Output: x (soln)
325: Return:
326: Description: This function performs x = E^-1.b
327: Usage:
328: XXT_solve(xxt_handle, double *x, double *b)
329: XXT_solve(xxt_handle, double *x, NULL)
330: assumes x has been initialized to be b
331: **************************************xxt.h***********************************/
332: PETSC_INTERN PetscErrorCode XXT_solve(xxt_ADT, PetscScalar *, PetscScalar *);
334: /*************************************xxt.h************************************
335: Function: XXT_sp_1()
337: Input : pointer to ADT
338: Output:
339: Return:
340: Description: sets xxt parameter 1 in xxt_handle
341: Usage: implement later
343: void XXT_sp_1(xxt_handle,parameter 1 value)
344: **************************************xxt.h***********************************/
346: /*************************************xyt.h************************************
347: Module Name: xyt
348: Module Info: need xyt.{c,h} gs.{c,h} comm.{c,h} ivec.{c,h} error.{c,h}
350: author: Henry M. Tufo III
351: e-mail: hmt@asci.uchicago.edu
352: contact:
353: +--------------------------------+--------------------------------+
354: |MCS Division - Building 221 |Department of Computer Science |
355: |Argonne National Laboratory |Ryerson 152 |
356: |9700 S. Cass Avenue |The University of Chicago |
357: |Argonne, IL 60439 |Chicago, IL 60637 |
358: |(630) 252-5354/5986 ph/fx |(773) 702-6019/8487 ph/fx |
359: +--------------------------------+--------------------------------+
361: Last Modification: 3.20.01
362: **************************************xyt.h***********************************/
364: typedef struct xyt_CDT *xyt_ADT;
366: /*************************************xyt.h************************************
367: Function: XYT_new()
369: Return: ADT ptr or NULL upon failure.
370: Description: This function allocates and returns an xyt handle
371: Usage: xyt_handle = xyt_new();
372: **************************************xyt.h***********************************/
373: PETSC_INTERN xyt_ADT XYT_new(void);
375: /*************************************xyt.h************************************
376: Function: XYT_free()
378: Input : pointer to ADT.
379: Description: This function frees the storage associated with an xyt handle
380: Usage: XYT_free(xyt_handle);
381: **************************************xyt.h***********************************/
382: PETSC_INTERN PetscErrorCode XYT_free(xyt_ADT);
384: /*************************************xyt.h************************************
385: Function: XYT_factor
387: Input : ADT ptr, and pointer to object
388: Output:
389: Return: 0 on failure, 1 on success
390: Description: This function sets the xyt solver
392: xyt assumptions: given n rows of global coarse matrix (E_loc) where
393: o global dofs N = sum_p(n), p=0,P-1
394: (i.e. row dist. with no dof replication)
395: (5.21.00 will handle dif replication case)
396: o m is the number of columns in E_loc (m>=n)
397: o local2global holds global number of column i (i=0,...,m-1)
398: o local2global holds global number of row i (i=0,...,n-1)
399: o mylocmatvec performs E_loc . x_loc where x_loc is an vector of
400: length m in 1-1 correspondence with local2global
401: (note that gs package takes care of communication).
402: (note do not zero out upper m-n entries!)
403: o mylocmatvec(void *grid_data, double *in, double *out)
405: ML beliefs/usage: move this to ML_XYT_factor routine
406: o my_ml holds address of ML struct associated w/E_loc, grid_data, grid_tag
407: o grid_tag, grid_data, my_ml used in
408: ML_Set_CSolve(my_ml, grid_tag, grid_data, ML_Do_CoarseDirect);
409: o grid_data used in
410: A_matvec(grid_data,v,u);
412: Usage:
413: **************************************xyt.h***********************************/
414: PETSC_INTERN PetscErrorCode XYT_factor(xyt_ADT, /* prev. allocated xyt handle */
415: PetscInt *, /* global column mapping */
416: PetscInt, /* local num rows */
417: PetscInt, /* local num cols */
418: PetscErrorCode (*)(void *, PetscScalar *, PetscScalar *), /* b_loc=A_local.x_loc */
419: void *); /* grid data for matvec */
421: /*************************************xyt.h************************************
422: Function: XYT_solve
424: Input : ADT ptr, b (rhs)
425: Output: x (soln)
426: Return:
427: Description: This function performs x = E^-1.b
428: Usage: XYT_solve(xyt_handle, double *x, double *b)
429: **************************************xyt.h***********************************/
430: PETSC_INTERN PetscErrorCode XYT_solve(xyt_ADT, PetscScalar *, PetscScalar *);
432: /*************************************xyt.h************************************
433: Function: XYT_stats
435: Input : handle
436: **************************************xyt.h***********************************/
437: PETSC_INTERN PetscErrorCode XYT_stats(xyt_ADT);
439: /********************************bit_mask.h************************************
441: Author: Henry M. Tufo III
443: e-mail: hmt@cs.brown.edu
445: snail-mail:
446: Division of Applied Mathematics
447: Brown University
448: Providence, RI 02912
450: Last Modification:
451: 11.21.97
452: *********************************bit_mask.h***********************************/
453: PETSC_INTERN PetscInt PCTFS_div_ceil(PetscInt, PetscInt);
454: PETSC_INTERN PetscErrorCode PCTFS_set_bit_mask(PetscInt *, PetscInt, PetscInt);
455: PETSC_INTERN PetscInt PCTFS_len_bit_mask(PetscInt);
456: PETSC_INTERN PetscInt PCTFS_ct_bits(char *, PetscInt);
457: PETSC_INTERN PetscErrorCode PCTFS_bm_to_proc(char *, PetscInt, PetscInt *);
458: PETSC_INTERN PetscInt PCTFS_len_buf(PetscInt, PetscInt);