FrontISTR 5.2.0
Large-scale structural analysis program with finit element method
Loading...
Searching...
No Matches
hecmw_ML_wrapper.c
Go to the documentation of this file.
1/*****************************************************************************
2 * Copyright (c) 2019 FrontISTR Commons
3 * This software is released under the MIT License, see LICENSE.txt
4 *****************************************************************************/
5
6#include <stdio.h>
7#include <stdlib.h>
8#include <errno.h>
9#include "hecmw_util.h"
10
11#ifdef HECMW_WITH_ML
12
13#include "Trilinos_version.h"
14#include "ml_include.h"
15#include "hecmw_ML_helper.h"
16#include "hecmw_ML_helper_33.h"
17#include "hecmw_ML_helper_nn.h"
18
19/*
20 * static variable
21 */
22
23struct ml_info {
24 ML *ml_object;
25 ML_Aggregate *agg_object;
26 int ndof;
27};
28
29#define MAX_MI 8
30
31static struct ml_info MLInfo[MAX_MI];
32
33
34/*
35 * Settings
36 */
37
38/* Whether coasest level is solved with direct solver
39 * If not, a few sweeps of smoother is applied.
40 * (Note: Trilinos must be built with Amesos package enabled)
41 */
42static int FlgDirectSolveCoarsest = 1;
43
44/* Direct solver for the coarsest level
45 * available types: KLU, MUMPS
46 * (KLU is a serial direct solver that comes with Trilinos/Amesos)
47 */
48enum direct_solver {KLU, MUMPS};
49static enum direct_solver DirectSolver = KLU;
50
51/* Smoother type
52 * available types: Jacobi, GaussSeidel, BlockGaussSeidel, SymGaussSeidel,
53 * SymBlockGaussSeidel, Cheby, Amesos, etc.
54 * However, the following three types are currently available from this interface
55 */
56enum smoother_type {Jacobi, SymBlockGaussSeidel, Cheby};
57static enum smoother_type SmootherType = Cheby;
58
59/* Whether HEC-MW smoother is used at finest level when SmootherType is SymBlockGaussSeidel
60 */
61static int FlgUseHECMWSmoother = 1;
62
63/* Solver cycle
64 * available types: ML_MGV (V-cycle), ML_MGW (W-cycle), ML_MGFULLV (Full V-Cycle), etc.
65 */
66static int MGType = ML_MGW;
67
68/* Max num of levels
69 */
70static int MaxLevels = 10;
71
72/* Coarsening scheme
73 * available types: Coupled, Uncoupled, MIS, UncoupledCoupled, UncoupledMIS,
74 * DD, METIS, ParMETIS, Zoltan, User
75 * However, the following three types are currently available from this interface
76 */
77enum coarsen_scheme {UncoupledMIS, METIS, ParMETIS, Zoltan, DD};
78static enum coarsen_scheme CoarsenScheme = UncoupledMIS;
79
80/* Num of smoother sweeps
81 */
82static int NumSweeps = 2;
83
84
85/*
86 * public functions
87 */
88
89void hecmw_ML_wrapper_setup(int *id, int *sym, int *Ndof, int *ierr) {
90 int loglevel, myrank;
91 int N_grids, N_levels;
92 int nlocal, nlocal_allcolumns;
93 ML *ml_object;
94 ML_Aggregate *agg_object;
95
96 if (*id <= 0 && MAX_MI < *id) {
97 *ierr = HECMW_ERROR;
98 return;
99 }
100
101 hecmw_ml_get_loglevel_(id, &loglevel);
102 ML_Set_PrintLevel(loglevel);
103
105
106 /* Get options */
107 {
108 int opt1, opt2, opt3, opt4, opt5, opt6;
109
110 hecmw_ml_get_opt1_(id, &opt1, ierr);
111 if (*ierr != HECMW_SUCCESS) return;
112 switch (opt1) {
113 case 0: /* default */
114 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML using default coarse solver\n");
115 case 1:
116 FlgDirectSolveCoarsest = 0;
117 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML coarse solver is smoother\n");
118 break;
119 case 2:
120 FlgDirectSolveCoarsest = 1;
121 DirectSolver = KLU;
122 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML coarse solver is KLU\n");
123 break;
124 case 3:
125 FlgDirectSolveCoarsest = 1;
126 DirectSolver = MUMPS;
127 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML coarse solver is MUMPS\n");
128 break;
129 default:
130 if (myrank == 0) fprintf(stderr, "WARNING: invalid solver_opt1=%d (ignored)\n", opt1);
131 }
132
133 hecmw_ml_get_opt2_(id, &opt2, ierr);
134 if (*ierr != HECMW_SUCCESS) return;
135 switch (opt2) {
136 case 0: /* default */
137 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML using default smoother\n");
138 case 1:
139 SmootherType = Cheby;
140 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML smoother is Cheby\n");
141 break;
142 case 2:
143 SmootherType = SymBlockGaussSeidel;
144 FlgUseHECMWSmoother = 1;
145 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML smoother is SymBlockGaussSeidel\n");
146 break;
147 case 3:
148 SmootherType = Jacobi;
149 FlgUseHECMWSmoother = 1;
150 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML smoother is Jacobi\n");
151 break;
152 default:
153 if (myrank == 0) fprintf(stderr, "WARNING: invalid solver_opt2=%d (ignored)\n", opt2);
154 }
155
156 hecmw_ml_get_opt3_(id, &opt3, ierr);
157 if (*ierr != HECMW_SUCCESS) return;
158 switch (opt3) {
159 case 0: /* default */
160 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML using default multigrid type\n");
161 case 1:
162 MGType = ML_MGV;
163 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML multigrid type is V-cycle\n");
164 break;
165 case 2:
166 MGType = ML_MGW;
167 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML multigrid type is W-cycle\n");
168 break;
169 case 3:
170 MGType = ML_MGFULLV;
171 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML multigrid type is Full-V-cycle\n");
172 break;
173 default:
174 if (myrank == 0) fprintf(stderr, "WARNING: invalid solver_opt3=%d (ignored)\n", opt3);
175 }
176
177 hecmw_ml_get_opt4_(id, &opt4, ierr);
178 if (*ierr != HECMW_SUCCESS) return;
179 if (opt4 > 0) {
180 MaxLevels = opt4;
181 } else {
182 if (opt4 < 0) {
183 if (myrank == 0) fprintf(stderr, "WARNING: invalid solver_opt4=%d (ignored)\n", opt4);
184 }
185 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML using default MaxLevels\n");
186 MaxLevels = 10; /* default */
187 }
188 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML num of max levels is %d\n", MaxLevels);
189
190 hecmw_ml_get_opt5_(id, &opt5, ierr);
191 if (*ierr != HECMW_SUCCESS) return;
192 switch (opt5) {
193 case 0: /* default */
194 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML using default coarsening scheme\n");
195 case 1:
196 CoarsenScheme = UncoupledMIS;
197 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML coarsening scheme is UncoupledMIS\n");
198 break;
199 case 2:
200 CoarsenScheme = METIS;
201 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML coarsening scheme is METIS\n");
202 break;
203 case 3:
204 CoarsenScheme = ParMETIS;
205 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML coarsening scheme is ParMETIS\n");
206 break;
207 case 4:
208 CoarsenScheme = Zoltan;
209 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML coarsening scheme is Zoltan\n");
210 break;
211 case 5:
212 CoarsenScheme = DD;
213 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML coarsening scheme is DD\n");
214 break;
215 default:
216 if (myrank == 0) fprintf(stderr, "WARNING: invalid solver_opt5=%d (ignored)\n", opt5);
217 }
218
219 hecmw_ml_get_opt6_(id, &opt6, ierr);
220 if (*ierr != HECMW_SUCCESS) return;
221 if (opt6 > 0) {
222 NumSweeps = opt6;
223 } else {
224 if (opt6 < 0) {
225 if (myrank == 0) fprintf(stderr, "WARNING: invalid solver_opt6=%d (ignored)\n", opt6);
226 }
227 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML using default num of sweeps\n");
228 NumSweeps = 2; /* default */
229 hecmw_ml_set_opt6_(id, &NumSweeps, ierr);
230 if (*ierr != HECMW_SUCCESS) return;
231 }
232 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML num of smoother sweeps is %d\n", NumSweeps);
233
234 }
235
236 /* ML object */
237 N_grids = MaxLevels;
238 ML_Create(&ml_object, N_grids);
239 hecmw_ml_get_nlocal_(id, &nlocal, &nlocal_allcolumns, ierr);
240 if (*ierr != HECMW_SUCCESS) return;
241 ML_Init_Amatrix(ml_object, 0, nlocal, nlocal, id);
242 if (*Ndof == 3) {
243 ML_Set_Amatrix_Getrow(ml_object, 0, hecmw_ML_getrow_33, hecmw_ML_comm_33,
244 nlocal_allcolumns);
245 ML_Set_Amatrix_Matvec(ml_object, 0, hecmw_ML_matvec_33);
246 } else {
247 ML_Set_Amatrix_Getrow(ml_object, 0, hecmw_ML_getrow_nn, hecmw_ML_comm_nn,
248 nlocal_allcolumns);
249 ML_Set_Amatrix_Matvec(ml_object, 0, hecmw_ML_matvec_nn);
250 }
251
252 /* if (!(*sym)) ML_Set_Symmetrize(ml_object, ML_YES); */
253
254 /* Aggregate */
255 ML_Aggregate_Create(&agg_object);
256
257 /* Null Space (Rigid Body Mode) */
258 {
259 int num_PDE_eqns = *Ndof;
260 int null_dim;
261 double *null_vect;
262 int leng = nlocal;
263 if (*Ndof == 1) {
264 null_dim = 1;
265 } else if (*Ndof == 2) {
266 null_dim = 3;
267 } else {
268 null_dim = 6;
269 }
270 null_vect = (double *)HECMW_malloc(sizeof(double) * null_dim * leng);
271 if (!null_vect) {
272 HECMW_set_error(errno, "");
273 abort();
274 }
275 hecmw_ml_get_rbm_(id, null_vect, ierr);
276 if (*ierr != HECMW_SUCCESS) return;
277 ML_Aggregate_Set_NullSpace(agg_object, num_PDE_eqns, null_dim, null_vect,
278 leng);
279 HECMW_free(null_vect);
280 }
281 /* ML_Aggregate_Set_MaxCoarseSize(agg_object, 128); */ /* default: 128 */
282
283 /* options */
284 /* CoarsenScheme */
285 {
286 if (CoarsenScheme == UncoupledMIS) {
287 ML_Aggregate_Set_CoarsenScheme_UncoupledMIS(agg_object);
288 } else if (CoarsenScheme == METIS) {
289 ML_Aggregate_Set_CoarsenScheme_METIS(agg_object);
290 } else if (CoarsenScheme == ParMETIS) {
291 ML_Aggregate_Set_CoarsenScheme_ParMETIS(agg_object);
292 } else if (CoarsenScheme == Zoltan) {
293 ML_Aggregate_Set_CoarsenScheme_Zoltan(agg_object);
294 } else if (CoarsenScheme == DD) {
295 ML_Aggregate_Set_CoarsenScheme_DD(agg_object);
296 }
297 /*
298 if (MaxLevels == 2) {
299 ML_Aggregate_Set_LocalNumber(ml_object, agg_object, ML_ALL_LEVELS, 1);
300 } else if (MaxLevels == 3) {
301 ML_Aggregate_Set_NodesPerAggr(ml_object, agg_object, ML_ALL_LEVELS, 512);
302 ML_Aggregate_Set_ReqLocalCoarseSize(ml_object->ML_num_levels, agg_object, ML_ALL_LEVELS, 128);
303 }
304 */
305 }
306 /* ML_Aggregate_Set_Threshold(agg_object, threshold); */
307 /* ML_Aggregate_Set_DampingFactor(agg_object, dampingfactor); */
308
309 /* eigen-analysis */
310 /* ML_Set_SpectralNormScheme_PowerMethod(ml_object); */ /* power-method (default) */
311 if (*sym) {
312 ML_Set_SpectralNormScheme_Calc(ml_object); /* cg */
313 }
314 /* ML_Set_SpectralNorm_Iterations(ml_object, 10); */ /* default: 10 */
315
316 /* repartitioning */
317 /* ML_Repartition_Activate(ml_object); */
318 /* ML_Repartition_Set_LargestMinMaxRatio(ml_object, 1.3); */ /* default: 1.3 */
319 /* ML_Repartition_Set_MinPerProc(ml_object, 512); */ /* default: 512 */
320 /* ML_Repartition_Set_PutOnSingleProc(ml_object, i); */
321 /* ML_Repartition_Set_StartLevel(ml_object, 1); */ /* default: 1 */
322 /* ML_Repartition_Set_Partitioner(ml_object, ML_USEPARMETIS); */ /* default: ML_USEZOLTAN */
323
324 ML_Aggregate_Set_Dimensions(agg_object, *Ndof);
325
326 /* Generate MultiGrid */
327 /* N_levels = ML_Gen_MGHierarchy_UsingAggregation( */
328 /* ml_object, 0, ML_INCREASING, agg_object); */
329 N_levels = ML_Gen_MultiLevelHierarchy_UsingAggregation(
330 ml_object, 0, ML_INCREASING, agg_object);
331 if (loglevel > 0 && myrank == 0) fprintf(stderr, "INFO: ML generated num of levels is %d\n", N_levels);
332
333 /* Smoother */
334 /*
335 * Set pre- and post-smoother for each level
336 * level : num in (0, N_levels-1) or ML_ALL_LEVELS
337 * pre-or-post: ML_PRESMOOTHER, ML_POSTSMOOTHER or ML_BOTH
338 * omega : damping factor for Jacobi, GaussSeidel, etc. (ML_DEFAULT=1.0)
339 */
340 {
341 int level;
342 int coarsest_level = N_levels - 1;
343 /*
344 * levels other than the coarsest level
345 */
346 if (SmootherType == Jacobi) {
347 level = 0;
348 if (FlgUseHECMWSmoother) {
349 /* use HEC-MW's Block-Diag preconditioner at the finest level */
350 if (*Ndof == 3) {
351 hecmw_ml_smoother_diag_setup_33_(id, ierr);
352 if (*ierr != HECMW_SUCCESS) return;
353 ML_Set_Smoother(ml_object, 0, ML_BOTH, id, hecmw_ML_smoother_diag_apply_33, "HEC-MW");
354 } else {
355 hecmw_ml_smoother_diag_setup_nn_(id, ierr);
356 if (*ierr != HECMW_SUCCESS) return;
357 ML_Set_Smoother(ml_object, 0, ML_BOTH, id, hecmw_ML_smoother_diag_apply_nn, "HEC-MW");
358 }
359 level++;
360 }
361 /* use ML's smoother at other levels */
362 for (; level < coarsest_level; level++) {
363 ML_Gen_Smoother_Jacobi(ml_object, level,
364 ML_BOTH, NumSweeps, ML_DEFAULT);
365 }
366 } else if (SmootherType == SymBlockGaussSeidel) {
367 level = 0;
368 if (FlgUseHECMWSmoother) {
369 /* use HEC-MW's Block-SSOR preconditioner at the finest level */
370 if (*Ndof == 3) {
371 hecmw_ml_smoother_ssor_setup_33_(id, ierr);
372 if (*ierr != HECMW_SUCCESS) return;
373 ML_Set_Smoother(ml_object, 0, ML_BOTH, id, hecmw_ML_smoother_ssor_apply_33, "HEC-MW");
374 } else {
375 hecmw_ml_smoother_ssor_setup_nn_(id, ierr);
376 if (*ierr != HECMW_SUCCESS) return;
377 ML_Set_Smoother(ml_object, 0, ML_BOTH, id, hecmw_ML_smoother_ssor_apply_nn, "HEC-MW");
378 }
379 level++;
380 }
381 /* use ML's smoother at other levels */
382 for (; level < coarsest_level; level++) {
383 ML_Gen_Smoother_SymBlockGaussSeidel(ml_object, level,
384 ML_BOTH, NumSweeps, ML_DEFAULT, *Ndof);
385 }
386 } else /* if (SmootherType == Cheby) */ {
387 for (level = 0; level < coarsest_level; level++) {
388 ML_Gen_Smoother_Cheby(ml_object, level, ML_BOTH, 20.0, NumSweeps);
389 }
390 }
391 /*
392 * coarsest level
393 */
394 if (FlgDirectSolveCoarsest) {
395 if (DirectSolver == MUMPS) {
396#if TRILINOS_MAJOR_VERSION < 13
397 ML_Gen_Smoother_Amesos(ml_object, coarsest_level, ML_AMESOS_MUMPS, 1, 0.0);
398#else
399 ML_Gen_Smoother_Amesos(ml_object, coarsest_level, ML_AMESOS_MUMPS, 1, 0.0, 1);
400#endif
401 } else /* if (DirectSolver == KLU) */ {
402#if TRILINOS_MAJOR_VERSION < 13
403 ML_Gen_Smoother_Amesos(ml_object, coarsest_level, ML_AMESOS_KLU, 1, 0.0);
404#else
405 ML_Gen_Smoother_Amesos(ml_object, coarsest_level, ML_AMESOS_KLU, 1, 0.0, 1);
406#endif
407 }
408 } else {
409 if (SmootherType == Jacobi) {
410 ML_Gen_Smoother_Jacobi(ml_object, coarsest_level,
411 ML_BOTH, 3, ML_DEFAULT);
412 } else if (SmootherType == SymBlockGaussSeidel) {
413 ML_Gen_Smoother_SymBlockGaussSeidel(ml_object, coarsest_level,
414 ML_BOTH, 3, ML_DEFAULT, *Ndof);
415 } else /* if (SmootherType == Cheby) */ {
416 ML_Gen_Smoother_Cheby(ml_object, coarsest_level, ML_BOTH, 20.0, 2);
417 }
418 }
419 }
420
421 /* Solver */
422 ML_Gen_Solver(ml_object, MGType, 0, N_levels - 1);
423
424 /* Save objects */
425 MLInfo[*id - 1].ml_object = ml_object;
426 MLInfo[*id - 1].agg_object = agg_object;
427 MLInfo[*id - 1].ndof = *Ndof;
428}
429
430void hecmw_ML_wrapper_apply(int *id, double rhs[], int *ierr) {
431 int nlocal, nlocal_allcolumns;
432 double *sol;
433 int i;
434 ML *ml_object;
435 if (*id <= 0 && MAX_MI < *id) {
436 *ierr = HECMW_ERROR;
437 return;
438 }
439 ml_object = MLInfo[*id - 1].ml_object;
440 hecmw_ml_get_nlocal_(id, &nlocal, &nlocal_allcolumns, ierr);
441 if (*ierr != HECMW_SUCCESS) return;
442 sol = (double *)HECMW_malloc(sizeof(double) * nlocal_allcolumns);
443 if (!sol) {
444 HECMW_set_error(errno, "");
445 abort();
446 }
447 /* MultiGrid V-cycle */
448 ML_Solve_MGV(ml_object, rhs, sol);
449 for (i = 0; i < nlocal; i++) {
450 rhs[i] = sol[i];
451 }
452 HECMW_free(sol);
453}
454
455void hecmw_ML_wrapper_clear(int *id, int *ierr) {
456 if (*id <= 0 && MAX_MI < *id) {
457 *ierr = HECMW_ERROR;
458 return;
459 }
460 ML_Aggregate_Destroy(&(MLInfo[*id - 1].agg_object));
461 ML_Destroy(&(MLInfo[*id - 1].ml_object));
462
463 if (FlgUseHECMWSmoother) {
464 if (SmootherType == Jacobi) {
465 if (MLInfo[*id - 1].ndof == 3) {
466 hecmw_ml_smoother_diag_clear_33_(id, ierr);
467 } else {
468 hecmw_ml_smoother_diag_clear_nn_(id, ierr);
469 }
470 } else if (SmootherType == SymBlockGaussSeidel) {
471 if (MLInfo[*id - 1].ndof == 3) {
472 hecmw_ml_smoother_ssor_clear_33_(id, ierr);
473 } else {
474 hecmw_ml_smoother_ssor_clear_nn_(id, ierr);
475 }
476 }
477 }
478}
479
480#else /* WITH_ML */
481
482void hecmw_ML_wrapper_setup(int *id, int *sym, int *Ndof, int *ierr) {
483 fprintf(stderr, "ERROR: ML not enabled\n");
484 *ierr = HECMW_ERROR;
485}
486void hecmw_ML_wrapper_apply(int *id, double rhs[], int *ierr) {
487 fprintf(stderr, "ERROR: ML not enabled\n");
488 *ierr = HECMW_ERROR;
489}
490void hecmw_ML_wrapper_clear(int *id, int *ierr) {
491 fprintf(stderr, "ERROR: ML not enabled\n");
492 *ierr = HECMW_ERROR;
493}
494
495#endif /* WITH_ML */
496/* Fortran interface */
497
498void hecmw_ml_wrapper_setup_(int *id, int *sym, int *ndof, int *ierr) {
499 hecmw_ML_wrapper_setup(id, sym, ndof, ierr);
500}
501void hecmw_ml_wrapper_setup__(int *id, int *sym, int *ndof, int *ierr) {
502 hecmw_ML_wrapper_setup(id, sym, ndof, ierr);
503}
504void HECMW_ML_WRAPPER_SETUP(int *id, int *sym, int *ndof, int *ierr) {
505 hecmw_ML_wrapper_setup(id, sym, ndof, ierr);
506}
507
508void hecmw_ml_wrapper_apply_(int *id, double rhs[], int *ierr) {
509 hecmw_ML_wrapper_apply(id, rhs, ierr);
510}
511void hecmw_ml_wrapper_apply__(int *id, double rhs[], int *ierr) {
512 hecmw_ML_wrapper_apply(id, rhs, ierr);
513}
514void HECMW_ML_WRAPPER_APPLY(int *id, double rhs[], int *ierr) {
515 hecmw_ML_wrapper_apply(id, rhs, ierr);
516}
517
518void hecmw_ml_wrapper_clear_(int *id, int *ierr) {
519 hecmw_ML_wrapper_clear(id, ierr);
520}
521void hecmw_ml_wrapper_clear__(int *id, int *ierr) {
522 hecmw_ML_wrapper_clear(id, ierr);
523}
524void HECMW_ML_WRAPPER_CLEAR(int *id, int *ierr) {
525 hecmw_ML_wrapper_clear(id, ierr);
526}
void hecmw_ml_get_loglevel_(int *id, int *level)
void hecmw_ml_get_opt4_(int *id, int *opt3, int *ierr)
void hecmw_ml_get_nlocal_(int *id, int *nlocal, int *nlocal_allcolumns, int *ierr)
void hecmw_ml_get_opt6_(int *id, int *opt3, int *ierr)
void hecmw_ml_get_opt3_(int *id, int *opt3, int *ierr)
void hecmw_ml_set_opt6_(int *id, int *opt3, int *ierr)
void hecmw_ml_get_opt2_(int *id, int *opt2, int *ierr)
void hecmw_ml_get_opt5_(int *id, int *opt3, int *ierr)
void hecmw_ml_get_rbm_(int *id, double rbm[], int *ierr)
void hecmw_ml_get_opt1_(int *id, int *opt1, int *ierr)
void hecmw_ml_wrapper_clear_(int *id, int *ierr)
void hecmw_ML_wrapper_apply(int *id, double rhs[], int *ierr)
void hecmw_ml_wrapper_apply__(int *id, double rhs[], int *ierr)
void hecmw_ml_wrapper_setup__(int *id, int *sym, int *ndof, int *ierr)
void hecmw_ml_wrapper_setup_(int *id, int *sym, int *ndof, int *ierr)
void hecmw_ML_wrapper_setup(int *id, int *sym, int *Ndof, int *ierr)
void hecmw_ML_wrapper_clear(int *id, int *ierr)
void HECMW_ML_WRAPPER_CLEAR(int *id, int *ierr)
void HECMW_ML_WRAPPER_SETUP(int *id, int *sym, int *ndof, int *ierr)
void HECMW_ML_WRAPPER_APPLY(int *id, double rhs[], int *ierr)
void hecmw_ml_wrapper_apply_(int *id, double rhs[], int *ierr)
void hecmw_ml_wrapper_clear__(int *id, int *ierr)
int HECMW_Comm_rank(HECMW_Comm comm, int *rank)
Definition: hecmw_comm.c:18
HECMW_Comm HECMW_comm_get_comm(void)
Definition: hecmw_comm.c:699
#define HECMW_ERROR
Definition: hecmw_config.h:66
#define HECMW_SUCCESS
Definition: hecmw_config.h:64
int HECMW_set_error(int errorno, const char *fmt,...)
Definition: hecmw_error.c:37
#define HECMW_free(ptr)
Definition: hecmw_malloc.h:24
#define HECMW_malloc(size)
Definition: hecmw_malloc.h:20
int * level
integer(kind=kint) myrank
PARALLEL EXECUTION.
Definition: m_fstr.f90:80