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