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