summaryrefslogtreecommitdiff
path: root/src/include/cctk_Reduction.h
blob: 941794b6df96b731cbacdf4e262309ddae37c846 (plain)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
 /*@@
   @header    cctk_Reduction.h
   @date
   @author    Gabrielle Allen
   @desc
              Header file for using reduction operators
   @enddesc
   @version   $Header$
 @@*/

#ifndef _CCTK_REDUCTION_H_
#define _CCTK_REDUCTION_H_ 1

#define REDUCTION_OPERATOR_REGISTER_ARGLIST  \
          const cGH *arg_GH, \
          int arg_proc, \
          int arg_num_outvals, \
          int arg_outtype, \
          void *arg_outvals, \
          int arg_num_invars, \
          const int arg_varlist []


#define REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST \
  const cGH *arg_GH, \
  int arg_proc, \
  int arg_nDims, \
  const int arg_dims [], \
  int arg_nArrays, \
  const void *const arg_inArrays [], \
  int arg_inType, \
  int arg_nOutVals, \
  void *arg_outVals, \
  int arg_outType

#define REDUCTION_LOCAL_ARRAY_OPERATOR_REGISTER_ARGLIST \
  int N_dims, \
  int operator_handle, \
  int param_table_handle, \
  int N_input_arrays, \
  const CCTK_INT input_array_dims[], \
  const CCTK_INT input_array_type_codes[], \
  const void *const input_arrays[], \
  int M_output_numbers, \
  const CCTK_INT output_number_type_codes[], \
  void *const output_numbers[]

#define REDUCTION_GRID_ARRAY_OPERATOR_REGISTER_ARGLIST \
  const cGH *GH,  \
  int dest_proc,  \
  int local_reduce_handle, \
  int param_table_handle,  \
  int N_input_arrays,  \
  const CCTK_INT input_array_variable_indices[],  \
  int M_output_values,  \
  const CCTK_INT output_value_type_codes[],  \
  void* const output_values[]

#define REDUCTION_ARRAYS_GLOBALLY_OPERATOR_REGISTER_ARGLIST \
  const cGH *GH,  \
  int dest_proc,  \
  int local_reduce_handle, \
  int param_table_handle,  \
  int N_input_arrays,  \
  const void *const input_arrays[], \
  int N_dims,\
  const CCTK_INT input_array_dims[], \
  const CCTK_INT input_array_type_codes[], \
  int M_output_values,  \
  const CCTK_INT output_value_type_codes[],  \
  void* const output_values[]
    
#ifdef __cplusplus
extern "C"
{
#endif

/* prototype for reduction operator routine */
typedef int (*cReduceOperator) (const cGH *GH,
                                int arg_proc,
                                int arg_num_outvals,
                                int arg_outtype,
                                void *arg_outvals,
                                int arg_num_invars,
                                const int arg_varlist[]);

/* prototype for local array reduction operator routine */
typedef int (*cLocalArrayReduceOperator) (int N_dims, int operator_handle, 
                          int param_table_handle,   int N_input_arrays,
                          const CCTK_INT input_array_dims[], 
                          const CCTK_INT input_array_type_codes[],
                          const void *const input_arrays[],
                          int M_output_numbers,
                          const CCTK_INT output_number_type_codes[],
                          void *const output_numbers[]);

/* prototype for GA reduction operator routine */
typedef int (*cGridArrayReduceOperator) (const cGH *GH,
                                         int dest_proc,
                                         int local_reduce_handle,
                                         int param_table_handle,
                                         int N_input_arrays,
                                         const CCTK_INT input_array_variable_indices[],
                                         int M_output_values,
                                         const CCTK_INT output_value_type_codes[],
                                         void* const output_values[]);

/* prototype for global array reduction operator routine */
typedef int (*cReduceArraysGloballyOperator) (const cGH *GH,
                                         int dest_proc,
                                         int local_reduce_handle,
                                         int param_table_handle,
                                         int N_input_arrays,
                                         const void * const input_arrays[],
                                         int input_dims,
                                         const CCTK_INT input_array_dims[],
                                         const CCTK_INT input_array_type_codes[],
                                         int M_output_values,
                                         const CCTK_INT output_value_type_codes[],
                                         void* const output_values[]);

int CCTK_Reduce(const cGH *GH,
                int proc,
                int operation_handle,
                int num_out_vals,
                int type_out_vals,
                void *out_vals,
                int num_in_fields, ...);

int CCTK_ReductionHandle(const char *reduction);

#define CCTK_RegisterReductionOperator(a,b) \
        CCTKi_RegisterReductionOperator(CCTK_THORNSTRING,a,b)

int CCTKi_RegisterReductionOperator(const char *thorn,
                                    cReduceOperator operatorGV,
                                    const char *name);

int CCTK_ReductionArrayHandle(const char *reduction);

int CCTK_RegisterReductionArrayOperator(
         int (*function)(REDUCTION_ARRAY_OPERATOR_REGISTER_ARGLIST),
         const char *name);

const char *CCTK_ReduceOperatorImplementation(int handle);

const char *CCTK_ReduceOperator (int handle);

int CCTK_NumReduceOperators(void);

/* new local array reduction API */
int CCTK_ReduceLocalArrays(int N_dims,
                           int local_reduce_handle,
                           int param_table_handle,
                           int N_input_arrays,
                           const CCTK_INT input_array_sizes[],
                           const CCTK_INT input_array_type_codes[],
                           const void *const input_arrays[],
                           int M_output_values,
                           const CCTK_INT output_value_type_codes[],
                           void *const output_values[]);

int CCTK_LocalArrayReductionHandle(const char *reduction);

#define CCTK_RegisterLocalArrayReductionOperator(a,b) \
        CCTKi_RegisterLocalArrayReductionOperator(CCTK_THORNSTRING,a,b)

int CCTKi_RegisterLocalArrayReductionOperator(const char *thorn,
                           cLocalArrayReduceOperator operatorGV,
                                              const char *name);

int CCTK_RegisterReductionLocalArrayOperator(
         int (*function)(REDUCTION_LOCAL_ARRAY_OPERATOR_REGISTER_ARGLIST),
         const char *name);

const char *CCTK_LocalArrayReduceOperatorImplementation(int handle);

const char *CCTK_LocalArrayReduceOperator (int handle);
int CCTK_NumLocalArrayReduceOperators(void);

/* new GA reduction API */
int CCTK_ReduceGridArrays(const cGH *GH,
                          int dest_proc,
                          int local_reduce_handle,
                          int param_table_handle,
                          int N_input_arrays,
                          const CCTK_INT input_array_variable_indices[],
                          int M_output_values,
                          const CCTK_INT output_value_type_codes[],
                          void* const output_values[]);

#define CCTK_RegisterGridArrayReductionOperator(a) \
        CCTKi_RegisterGridArrayReductionOperator(CCTK_THORNSTRING,a)

int CCTKi_RegisterGridArrayReductionOperator(const char *thorn, cGridArrayReduceOperator 
        operatorGV);

const char *CCTK_GridArrayReductionOperator(void);
int CCTK_NumGridArrayReductionOperators(void);

/* new pointwise reduction API */
int CCTK_ReduceArraysGlobally(const cGH *GH,
                          int dest_proc,
                          int local_reduce_handle,
                          int param_table_handle,
                          int N_input_arrays,
                          const void * const input_arrays[],
                          int input_dims,
                          const CCTK_INT input_array_dims[],
                          const CCTK_INT input_array_type_codes[],
                          int M_output_values,
                          const CCTK_INT output_value_type_codes[],
                          void* const output_values[]);

#define CCTK_RegisterReduceArraysGloballyOperator(a) \
        CCTKi_RegisterReduceArraysGloballyOperator(CCTK_THORNSTRING,a)

int CCTKi_RegisterReduceArraysGloballyOperator(const char *thorn, cReduceArraysGloballyOperator 
        operatorGV);

const char *CCTK_ReductionArraysGloballyOperator(void);
int CCTK_NumReductionArraysGloballyOperators(void);

/* FIXME: old interface - should go */
int CCTK_ReduceLocalScalar (const cGH *GH, int proc, int operation_handle,
                            const void *inScalar, void *outScalar, int dataType);


/* FIXME: old interface - should go */
int CCTK_ReduceLocalArray1D (const cGH *GH, int proc, int operation_handle,
                            const void *in_array1d, void *out_array1d,
                            int num_in_array1d, int data_type);

int CCTK_ReduceLocScalar(const cGH *GH, int proc, int operation_handle,
                         const void *in_scalar, void *out_scalar, int data_type);

int CCTK_ReduceLocArrayToArray1D(const cGH *GH, int proc, int operation_handle,
                                 const void *in_array1d, void *out_array1d,
                                 int num_in_array1d,
                                 int data_type);

int CCTK_ReduceLocArrayToArray2D(const cGH *GH, int proc, int operation_handle,
                                 const void *in_array2d, void *out_array2d,
                                 int xsize, int ysize,
                                 int data_type);

int CCTK_ReduceLocArrayToArray3D(const cGH *GH, int proc, int operation_handle,
                                 const void *in_array3d, void *out_array3d,
                                 int xsize, int  ysize, int zsize,
                                 int data_type);


int CCTK_ReduceArray(const cGH *GH,
                     int proc,
                     int operation_handle,
                     int num_out_vals,
                     int type_out_vals,
                     void *out_vals,
                     int num_dims,
                     int num_in_arrays,
                     int type_in_arrays,
                     ... );
                          


#ifdef __cplusplus
}
#endif

#endif /* _CCTK_REDUCTION_H_ */