1 /*============================================================================
2 * Sparse Linear Equation Solvers using HYPRE
3 *============================================================================*/
4
5 /*
6 This file is part of Code_Saturne, a general-purpose CFD tool.
7
8 Copyright (C) 1998-2021 EDF S.A.
9
10 This program is free software; you can redistribute it and/or modify it under
11 the terms of the GNU General Public License as published by the Free Software
12 Foundation; either version 2 of the License, or (at your option) any later
13 version.
14
15 This program is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
18 details.
19
20 You should have received a copy of the GNU General Public License along with
21 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
22 Street, Fifth Floor, Boston, MA 02110-1301, USA.
23 */
24
25 /*----------------------------------------------------------------------------*/
26
27 #include "cs_defs.h"
28
29 /*----------------------------------------------------------------------------
30 * Standard C library headers
31 *----------------------------------------------------------------------------*/
32
33 #include <stdarg.h>
34 #include <stdio.h>
35 #include <stdlib.h>
36 #include <string.h>
37 #include <assert.h>
38 #include <math.h>
39
40 #if defined(HAVE_MPI)
41 #include <mpi.h>
42 #endif
43
44 /*----------------------------------------------------------------------------
45 * HYPRE headers
46 *----------------------------------------------------------------------------*/
47
48 /* Avoid warnings due to previous values */
49 #undef PACKAGE_BUGREPORT
50 #undef PACKAGE_NAME
51 #undef PACKAGE_STRING
52 #undef PACKAGE_TARNAME
53 #undef PACKAGE_URL
54 #undef PACKAGE_VERSION
55
56 #include <HYPRE_krylov.h>
57 #include <HYPRE_parcsr_ls.h>
58 #include <HYPRE_utilities.h>
59
60 #if !defined(HYPRE_RELEASE_NUMBER)
61 #define HYPRE_RELEASE_NUMBER 0
62 #endif
63
64 /*----------------------------------------------------------------------------
65 * Local headers
66 *----------------------------------------------------------------------------*/
67
68 #include "bft_mem.h"
69 #include "bft_error.h"
70 #include "bft_printf.h"
71
72 #include "cs_base.h"
73 #include "cs_base_accel.h"
74 #include "cs_log.h"
75 #include "cs_fp_exception.h"
76 #include "cs_halo.h"
77 #include "cs_matrix.h"
78 #include "cs_matrix_default.h"
79 #include "cs_matrix_hypre.h"
80 #include "cs_matrix_hypre_priv.h"
81 #include "cs_timer.h"
82
83 /*----------------------------------------------------------------------------
84 * Header for the current file
85 *----------------------------------------------------------------------------*/
86
87 #include "cs_sles.h"
88 #include "cs_sles_hypre.h"
89
90 /*----------------------------------------------------------------------------*/
91
92 BEGIN_C_DECLS
93
94 /*=============================================================================
95 * Additional doxygen documentation
96 *============================================================================*/
97
98 /*!
99 \file cs_sles_hypre.c
100
101 \brief handling of HYPRE-based linear solvers
102
103 \page sles_hypre HYPRE-based linear solvers.
104
105 \typedef cs_sles_hypre_setup_hook_t
106
107 \brief Function pointer for settings of a HYPRE solver setup.
108
109 This function is called during the setup stage for a HYPRE solver.
110
111 When first called, the solver argument is NULL, and must be created
112 using HYPRE functions.
113
114 Note: if the context pointer is non-NULL, it must point to valid data
115 when the selection function is called so that value or structure should
116 not be temporary (i.e. local);
117
118 \param[in] verbosity verbosity level
119 \param[in, out] context pointer to optional (untyped) value or structure
120 \param[in, out] solver handle to HYPRE solver (to be cast as HYPRE_Solver)
121 */
122
123 /*! \cond DOXYGEN_SHOULD_SKIP_THIS */
124
125 /*=============================================================================
126 * Local Macro Definitions
127 *============================================================================*/
128
129 /*=============================================================================
130 * Local Structure Definitions
131 *============================================================================*/
132
133 /* Basic per linear system options and logging */
134 /*---------------------------------------------*/
135
136 typedef struct _cs_sles_hypre_setup_t {
137
138 cs_matrix_coeffs_hypre_t *coeffs; /* HYPRE matrix and vectors */
139
140 HYPRE_Solver solver; /* Solver data */
141 HYPRE_Solver precond; /* Preconditioner data */
142
143 } cs_sles_hypre_setup_t;
144
145 typedef struct _cs_sles_hypre_t {
146
147 cs_sles_hypre_type_t solver_type; /* Solver type */
148 cs_sles_hypre_type_t precond_type; /* Preconditioner type */
149
150 int use_device; /* O for host, 1 for device */
151
152 /* Performance data */
153
154 int n_setups; /* Number of times system setup */
155 int n_solves; /* Number of times system solved */
156
157 int n_iterations_last; /* Number of iterations for last
158 system resolution */
159 int n_iterations_min; /* Minimum number of iterations
160 in system resolution history */
161 int n_iterations_max; /* Maximum number of iterations
162 in system resolution history */
163 int long long n_iterations_tot; /* Total accumulated number of
164 iterations */
165
166 cs_timer_counter_t t_setup; /* Total setup */
167 cs_timer_counter_t t_solve; /* Total time used */
168
169 /* Additional setup options */
170
171 void *hook_context; /* Optional user context */
172 cs_sles_hypre_setup_hook_t *setup_hook; /* Post setup function */
173
174 cs_sles_hypre_setup_t *setup_data;
175
176 } cs_sles_hypre_t;
177
178 /*============================================================================
179 * Global variables
180 *============================================================================*/
181
182 static int _n_hypre_systems = 0;
183
184 /*============================================================================
185 * Private function definitions
186 *============================================================================*/
187
188 /*----------------------------------------------------------------------------*/
189 /*!
190 * \brief Return name of hypre solver type.
191 *
192 * \param[in] solver_type solver type id
193 *
194 * \return name od associated solver type.
195 */
196 /*----------------------------------------------------------------------------*/
197
198 static const char *
_cs_hypre_type_name(cs_sles_hypre_type_t solver_type)199 _cs_hypre_type_name(cs_sles_hypre_type_t solver_type)
200 {
201 switch(solver_type) {
202
203 case CS_SLES_HYPRE_BOOMERAMG:
204 return N_("BoomerAMG");
205 break;
206 case CS_SLES_HYPRE_HYBRID:
207 return N_("Hybrid");
208 break;
209 case CS_SLES_HYPRE_ILU:
210 return N_("ILU");
211 break;
212 case CS_SLES_HYPRE_BICGSTAB:
213 return N_("BiCGSTAB");
214 break;
215 case CS_SLES_HYPRE_GMRES:
216 return N_("GMRES");
217 break;
218 case CS_SLES_HYPRE_FLEXGMRES:
219 return N_("Flexible GMRES");
220 break;
221 case CS_SLES_HYPRE_LGMRES:
222 return N_("LGMRES");
223 break;
224 case CS_SLES_HYPRE_PCG:
225 return N_("PCG");
226 break;
227 case CS_SLES_HYPRE_EUCLID:
228 return N_("EUCLID");
229 break;
230 case CS_SLES_HYPRE_PARASAILS:
231 return N_("ParaSails");
232 break;
233 case CS_SLES_HYPRE_NONE:
234 return N_("None");
235 break;
236
237 default:
238 return NULL;
239 }
240 }
241
242 /*----------------------------------------------------------------------------*/
243 /*!
244 * \brief Ensure MPI is initialized if present and needed
245 */
246 /*----------------------------------------------------------------------------*/
247
248 static void
_ensure_mpi_init(void)249 _ensure_mpi_init(void)
250 {
251 #if defined(HAVE_MPI) && defined(HYPRE_HAVE_MPI)
252 int flag = 0;
253 MPI_Initialized(&flag);
254 if (!flag) {
255 int mpi_threads;
256 MPI_Init_thread(NULL, NULL, MPI_THREAD_FUNNELED, &mpi_threads);
257 }
258 #endif
259 }
260
261 /*============================================================================
262 * Public function definitions
263 *============================================================================*/
264
265 /*----------------------------------------------------------------------------*/
266 /*!
267 * \brief Define and associate a HYPRE linear system solver
268 * for a given field or equation name.
269 *
270 * If this system did not previously exist, it is added to the list of
271 * "known" systems. Otherwise, its definition is replaced by the one
272 * defined here.
273 *
274 * This is a utility function: if finer control is needed, see
275 * \ref cs_sles_define and \ref cs_sles_petsc_create.
276 *
277 * The associated solver required that the matrix passed to it is a HYPRE
278 * matrix (see cs_matrix_set_type_hypre).
279 *
280 * Note that this function returns a pointer directly to the iterative solver
281 * management structure. This may be used to set further options.
282 * If needed, \ref cs_sles_find may be used to obtain a pointer to the matching
283 * \ref cs_sles_t container.
284 *
285 * \param[in] f_id associated field id, or < 0
286 * \param[in] name associated name if f_id < 0, or NULL
287 * \param[in] solver_type HYPRE solver type
288 * \param[in] precond_type HYPRE preconditioner type
289 * \param[in] setup_hook pointer to optional setup epilogue function
290 * \param[in,out] context pointer to optional (untyped) value or
291 * structure for setup_hook, or NULL
292 *
293 * \return pointer to newly created iterative solver info object.
294 */
295 /*----------------------------------------------------------------------------*/
296
297 cs_sles_hypre_t *
cs_sles_hypre_define(int f_id,const char * name,cs_sles_hypre_type_t solver_type,cs_sles_hypre_type_t precond_type,cs_sles_hypre_setup_hook_t * setup_hook,void * context)298 cs_sles_hypre_define(int f_id,
299 const char *name,
300 cs_sles_hypre_type_t solver_type,
301 cs_sles_hypre_type_t precond_type,
302 cs_sles_hypre_setup_hook_t *setup_hook,
303 void *context)
304 {
305 if (solver_type < 0 || solver_type >= CS_SLES_HYPRE_NONE)
306 bft_error(__FILE__, __LINE__, 0,
307 _("Incorrect solver type argument %d for HYPRE."),
308 (int)solver_type);
309
310 cs_sles_hypre_t *c = cs_sles_hypre_create(solver_type,
311 precond_type,
312 setup_hook,
313 context);
314
315 cs_sles_t *sc = cs_sles_define(f_id,
316 name,
317 c,
318 "cs_sles_hypre_t",
319 cs_sles_hypre_setup,
320 cs_sles_hypre_solve,
321 cs_sles_hypre_free,
322 cs_sles_hypre_log,
323 cs_sles_hypre_copy,
324 cs_sles_hypre_destroy);
325
326 cs_sles_set_error_handler(sc,
327 cs_sles_hypre_error_post_and_abort);
328
329 return c;
330 }
331
332 /*----------------------------------------------------------------------------*/
333 /*!
334 * \brief Create HYPRE linear system solver info and context.
335 *
336 * In case of rotational periodicity for a block (non-scalar) matrix,
337 * the matrix type will be forced to MATSHELL ("shell") regardless
338 * of the option used.
339 *
340 * \param[in] solver_type HYPRE solver type
341 * \param[in] precond_type HYPRE preconditioner type
342 * \param[in] setup_hook pointer to optional setup epilogue function
343 * \param[in] context pointer to optional (untyped) value or structure
344 * for setup_hook, or NULL
345 *
346 * \return pointer to newly created linear system object.
347 */
348 /*----------------------------------------------------------------------------*/
349
350 cs_sles_hypre_t *
cs_sles_hypre_create(cs_sles_hypre_type_t solver_type,cs_sles_hypre_type_t precond_type,cs_sles_hypre_setup_hook_t * setup_hook,void * context)351 cs_sles_hypre_create(cs_sles_hypre_type_t solver_type,
352 cs_sles_hypre_type_t precond_type,
353 cs_sles_hypre_setup_hook_t *setup_hook,
354 void *context)
355 {
356 cs_sles_hypre_t *c;
357
358 if (_n_hypre_systems == 0) {
359 _ensure_mpi_init();
360 HYPRE_Init(); /* Note: ideally, HYPRE should provide a function to
361 check if it is already initialized or not */
362 }
363 _n_hypre_systems += 1;
364
365
366 BFT_MALLOC(c, 1, cs_sles_hypre_t);
367
368 c->solver_type = solver_type;
369 c->precond_type = precond_type;
370
371 int device_id = cs_get_device_id();
372 c->use_device = (device_id < 0) ? 0 : 1;
373
374 c->n_setups = 0;
375 c->n_solves = 0;
376 c->n_iterations_last = 0;
377 c->n_iterations_min = 0;
378 c->n_iterations_max = 0;
379 c->n_iterations_tot = 0;
380
381 CS_TIMER_COUNTER_INIT(c->t_setup);
382 CS_TIMER_COUNTER_INIT(c->t_solve);
383
384 /* Options */
385
386 c->hook_context = context;
387 c->setup_hook = setup_hook;
388
389 /* Setup data */
390 c->setup_data = NULL;
391
392 return c;
393 }
394
395 /*----------------------------------------------------------------------------*/
396 /*!
397 * \brief Destroy iterative sparse linear system solver info and context.
398 *
399 * \param[in, out] context pointer to iterative solver info and context
400 * (actual type: cs_sles_hypre_t **)
401 */
402 /*----------------------------------------------------------------------------*/
403
404 void
cs_sles_hypre_destroy(void ** context)405 cs_sles_hypre_destroy(void **context)
406 {
407 cs_sles_hypre_t *c = (cs_sles_hypre_t *)(*context);
408 if (c != NULL) {
409 cs_sles_hypre_free(c);
410 BFT_FREE(c);
411 *context = c;
412 }
413
414 _n_hypre_systems -= 1;
415 if (_n_hypre_systems == 0) {
416 HYPRE_Finalize();
417 }
418 }
419
420 /*----------------------------------------------------------------------------*/
421 /*!
422 * \brief Create HYPRE sparse linear system solver info and context
423 * based on existing info and context.
424 *
425 * \param[in] context pointer to reference info and context
426 * (actual type: cs_sles_hypre_t *)
427 *
428 * \return pointer to newly created solver info object.
429 * (actual type: cs_sles_hypre_t *)
430 */
431 /*----------------------------------------------------------------------------*/
432
433 void *
cs_sles_hypre_copy(const void * context)434 cs_sles_hypre_copy(const void *context)
435 {
436 cs_sles_hypre_t *d = NULL;
437
438 if (context != NULL) {
439 const cs_sles_hypre_t *c = context;
440 d = cs_sles_hypre_create(c->solver_type,
441 c->precond_type,
442 c->setup_hook,
443 c->hook_context);
444 }
445
446 return d;
447 }
448
449 /*----------------------------------------------------------------------------*/
450 /*!
451 * \brief Error handler for HYPRE solver.
452 *
453 * In case of divergence or breakdown, this error handler outputs an error
454 * message
455 * It does nothing in case the maximum iteration count is reached.
456
457 * \param[in, out] sles pointer to solver object
458 * \param[in] state convergence state
459 * \param[in] a matrix
460 * \param[in] rhs right hand side
461 * \param[in, out] vx system solution
462 *
463 * \return false (do not attempt new solve)
464 */
465 /*----------------------------------------------------------------------------*/
466
467 bool
cs_sles_hypre_error_post_and_abort(cs_sles_t * sles,cs_sles_convergence_state_t state,const cs_matrix_t * a,const cs_real_t * rhs,cs_real_t * vx)468 cs_sles_hypre_error_post_and_abort(cs_sles_t *sles,
469 cs_sles_convergence_state_t state,
470 const cs_matrix_t *a,
471 const cs_real_t *rhs,
472 cs_real_t *vx)
473 {
474 CS_UNUSED(a);
475 CS_UNUSED(rhs);
476 CS_UNUSED(vx);
477
478 if (state >= CS_SLES_BREAKDOWN)
479 return false;
480
481 const char *name = cs_sles_get_name(sles);
482
483 const cs_sles_hypre_t *c = cs_sles_get_context(sles);
484 CS_UNUSED(c);
485
486 const char *error_type[] = {N_("divergence"), N_("breakdown")};
487 int err_id = (state == CS_SLES_BREAKDOWN) ? 1 : 0;
488
489 bft_error(__FILE__, __LINE__, 0,
490 _("HYPRE: error (%s) solving for %s"),
491 _(error_type[err_id]),
492 name);
493
494 return false;
495 }
496
497 /*----------------------------------------------------------------------------*/
498 /*!
499 * \brief Log sparse linear equation solver info.
500 *
501 * \param[in] context pointer to iterative solver info and context
502 * (actual type: cs_sles_hypre_t *)
503 * \param[in] log_type log type
504 */
505 /*----------------------------------------------------------------------------*/
506
507 void
cs_sles_hypre_log(const void * context,cs_log_t log_type)508 cs_sles_hypre_log(const void *context,
509 cs_log_t log_type)
510 {
511 const cs_sles_hypre_t *c = context;
512
513 if (log_type == CS_LOG_SETUP) {
514
515 cs_log_printf(log_type,
516 _(" Solver type: HYPRE (%s)\n"),
517 _cs_hypre_type_name(c->solver_type));
518 if (c->precond_type < CS_SLES_HYPRE_NONE)
519 cs_log_printf(log_type,
520 _(" Preconditioning: %s\n"),
521 _cs_hypre_type_name(c->precond_type));
522 if (c->use_device)
523 cs_log_printf(log_type,
524 _(" Accelerated device: enabled\n"));
525
526 }
527
528 else if (log_type == CS_LOG_PERFORMANCE) {
529
530 int n_calls = c->n_solves;
531 int n_it_min = c->n_iterations_min;
532 int n_it_max = c->n_iterations_max;
533 int n_it_mean = 0;
534
535 if (n_calls > 0)
536 n_it_mean = (int)( c->n_iterations_tot
537 / ((int long long)n_calls));
538
539 cs_log_printf(log_type,
540 _("\n"
541 " Solver type: HYPRE (%s)\n"),
542 _cs_hypre_type_name(c->solver_type));
543 if (c->precond_type < CS_SLES_HYPRE_NONE)
544 cs_log_printf(log_type,
545 _(" Preconditioning: %s\n"),
546 _cs_hypre_type_name(c->precond_type));
547 if (c->use_device)
548 cs_log_printf(log_type,
549 _(" Accelerated device: enabled\n"));
550 cs_log_printf(log_type,
551 _(" Number of setups: %12d\n"
552 " Number of calls: %12d\n"
553 " Minimum number of iterations: %12d\n"
554 " Maximum number of iterations: %12d\n"
555 " Mean number of iterations: %12d\n"
556 " Construction: %12.3f\n"
557 " Resolution: %12.3f\n"),
558 c->n_setups, n_calls, n_it_min, n_it_max, n_it_mean,
559 c->t_setup.nsec*1e-9,
560 c->t_solve.nsec*1e-9);
561
562 }
563 }
564
565 /*----------------------------------------------------------------------------*/
566 /*!
567 * \brief Setup iterative sparse linear equation solver.
568 *
569 * \param[in, out] context pointer to iterative solver info and context
570 * (actual type: cs_sles_hypre_t *)
571 * \param[in] name pointer to system name
572 * \param[in] a associated matrix
573 * \param[in] verbosity associated verbosity
574 */
575 /*----------------------------------------------------------------------------*/
576
577 void
cs_sles_hypre_setup(void * context,const char * name,const cs_matrix_t * a,int verbosity)578 cs_sles_hypre_setup(void *context,
579 const char *name,
580 const cs_matrix_t *a,
581 int verbosity)
582 {
583 CS_UNUSED(name);
584
585 cs_timer_t t0;
586 t0 = cs_timer_time();
587
588 cs_sles_hypre_t *c = context;
589 cs_sles_hypre_setup_t *sd = c->setup_data;
590 if (sd == NULL) {
591 BFT_MALLOC(c->setup_data, 1, cs_sles_hypre_setup_t);
592 sd = c->setup_data;
593 sd->solver = NULL;
594 sd->precond = NULL;
595 }
596
597 const char expected_matrix_type[] = "HYPRE_PARCSR";
598 if (strcmp(cs_matrix_get_type_name(a), expected_matrix_type))
599 bft_error(__FILE__, __LINE__, 0,
600 _("HYPRE [%s]:\n"
601 " expected matrix type: %s\n"
602 " provided matrix type: %s"),
603 name, expected_matrix_type, cs_matrix_get_type_name(a));
604
605 sd->coeffs = cs_matrix_hypre_get_coeffs(a);
606
607 MPI_Comm comm = cs_glob_mpi_comm;
608 if (comm == MPI_COMM_NULL)
609 comm = MPI_COMM_WORLD;
610
611 bool have_set_pc = true;
612 HYPRE_PtrToParSolverFcn solve_ftn[2] = {NULL, NULL};
613 HYPRE_PtrToParSolverFcn setup_ftn[2] = {NULL, NULL};
614
615 for (int i = 0; i < 2; i++) {
616
617 HYPRE_Solver hs = (i == 0) ? sd->precond : sd->solver;
618 cs_sles_hypre_type_t hs_type = (i == 0) ? c->precond_type : c->solver_type;
619
620 /* hs should be NULL at this point, unless we do not relly free
621 it (when cs_sles_hypre_free is called (for example to amortize setup) */
622
623 if (hs != NULL || hs_type >= CS_SLES_HYPRE_NONE)
624 continue;
625
626 switch(hs_type) {
627
628 case CS_SLES_HYPRE_BOOMERAMG:
629 {
630 HYPRE_BoomerAMGCreate(&hs);
631
632 if (verbosity > 2 ) {
633 HYPRE_BoomerAMGSetPrintLevel(hs, 1);
634 HYPRE_BoomerAMGSetPrintFileName(hs, "hypre.log");
635 }
636
637 /* Default settings for device */
638 if (c->use_device == 1) {
639
640 HYPRE_BoomerAMGSetRelaxType(hs, 6); /* 3, 4, 6, 7, 18, 11, 12 */
641 HYPRE_BoomerAMGSetRelaxOrder(hs, 0); /* must be false */
642 HYPRE_BoomerAMGSetCoarsenType(hs, 8); /* PMIS */
643 HYPRE_BoomerAMGSetInterpType(hs, 17); /* 3, 15, 6, 14, 17, 18 */
644 HYPRE_BoomerAMGSetAggInterpType(hs, 7); /* 5 or 7 */
645 HYPRE_BoomerAMGSetKeepTranspose(hs, 1); /* keep transpose to
646 avoid SpMTV */
647 HYPRE_BoomerAMGSetRAP2(hs, 0); /* RAP in two multiplications
648 (default: FALSE) */
649 }
650
651 /* Default settings for host */
652 else if (c->use_device == 0) {
653 HYPRE_BoomerAMGSetCoarsenType(hs, 10); /* HMIS */
654 HYPRE_BoomerAMGSetPMaxElmts(hs, 4);
655 HYPRE_BoomerAMGSetInterpType(hs, 7); /* extended+e */
656 HYPRE_BoomerAMGSetRelaxType(hs, 6); /* Sym G.S./Jacobi hybrid */
657 HYPRE_BoomerAMGSetRelaxOrder(hs, 0);
658 }
659
660 /* Defaults for both host and device */
661 HYPRE_BoomerAMGSetAggNumLevels(hs, 2);
662 HYPRE_BoomerAMGSetStrongThreshold(hs, 0.5); /* 2d=>0.25 3d=>0.5 */
663
664 solve_ftn[i] = HYPRE_BoomerAMGSolve;
665 setup_ftn[i] = HYPRE_BoomerAMGSetup;
666
667 if (i == 0) { /* preconditioner */
668 HYPRE_BoomerAMGSetTol(hs, 0.0);
669 HYPRE_BoomerAMGSetMaxIter(hs, 1);
670 }
671 else { /* solver */
672 have_set_pc = false;
673 HYPRE_BoomerAMGSetMaxIter(hs, 1000);
674 }
675
676 }
677 break;
678
679 case CS_SLES_HYPRE_HYBRID:
680 {
681 HYPRE_ParCSRHybridCreate(&hs);
682
683 if (verbosity > 2 ) {
684 HYPRE_ParCSRHybridSetPrintLevel(hs, 2); /* Print solve info */
685 HYPRE_ParCSRHybridSetLogging(hs, 1); /* Needed to get info later */
686 }
687
688 solve_ftn[i] = HYPRE_ParCSRHybridSolve;
689 setup_ftn[i] = HYPRE_ParCSRHybridSetup;
690
691 if (i == 1 && solve_ftn[0] != NULL) { /* solver */
692 HYPRE_ParCSRHybridSetPrecond(hs,
693 solve_ftn[0],
694 setup_ftn[0],
695 sd->precond);
696 }
697 }
698 break;
699
700 case CS_SLES_HYPRE_ILU:
701 {
702 HYPRE_ILUCreate(&hs);
703
704 solve_ftn[i] = HYPRE_ILUSolve;
705 setup_ftn[i] = HYPRE_ILUSetup;
706
707 if (i == 0) { /* preconditioner */
708 HYPRE_ILUSetTol(sd->solver, 0.);
709 }
710 else { /* solver */
711 have_set_pc = false;
712 }
713
714 if (verbosity > 2 ) {
715 HYPRE_ILUSetPrintLevel(hs, 2); /* Print solve info */
716 HYPRE_ILUSetLogging(hs, 1); /* Needed to get info later */
717 }
718 }
719 break;
720
721 case CS_SLES_HYPRE_BICGSTAB:
722 {
723 HYPRE_ParCSRBiCGSTABCreate(comm, &hs);
724
725 if (verbosity > 2 ) {
726 HYPRE_BiCGSTABSetPrintLevel(hs, 2); /* Print solve info */
727 HYPRE_BiCGSTABSetLogging(hs, 1); /* Needed to get run info later */
728 }
729
730 solve_ftn[i] = HYPRE_ParCSRBiCGSTABSolve;
731 setup_ftn[i] = HYPRE_ParCSRBiCGSTABSetup;
732
733 if (i == 0) { /* preconditioner */
734 HYPRE_ParCSRBiCGSTABSetTol(hs, 0.0);
735 HYPRE_ParCSRBiCGSTABSetMaxIter(hs, 1);
736 solve_ftn[i] = HYPRE_ParCSRBiCGSTABSolve;
737 setup_ftn[i] = HYPRE_ParCSRBiCGSTABSetup;
738 }
739 else { /* solver */
740 HYPRE_BiCGSTABSetMaxIter(hs, 1000); /* Max iterations */
741 if (solve_ftn[0] != NULL)
742 HYPRE_ParCSRBiCGSTABSetPrecond(hs,
743 solve_ftn[0],
744 setup_ftn[0],
745 sd->precond);
746 }
747 }
748 break;
749
750 case CS_SLES_HYPRE_GMRES:
751 {
752 HYPRE_ParCSRGMRESCreate(comm, &hs);
753
754 if (verbosity > 2 ) {
755 HYPRE_GMRESSetPrintLevel(hs, 2); /* Print solve info */
756 HYPRE_GMRESSetLogging(hs, 1); /* Needed to get run info later */
757 solve_ftn[i] = HYPRE_ParCSRGMRESSolve;
758 setup_ftn[i] = HYPRE_ParCSRGMRESSetup;
759 }
760
761 solve_ftn[i] = HYPRE_ParCSRGMRESSolve;
762 setup_ftn[i] = HYPRE_ParCSRGMRESSetup;
763
764 if (i == 0) { /* preconditioner */
765 HYPRE_ParCSRGMRESSetTol(hs, 0.0);
766 HYPRE_ParCSRGMRESSetMaxIter(hs, 1);
767 }
768 else { /* solver */
769 HYPRE_GMRESSetMaxIter(hs, 1000); /* Max iterations */
770 if (solve_ftn[0] != NULL)
771 HYPRE_ParCSRGMRESSetPrecond(hs,
772 solve_ftn[0],
773 setup_ftn[0],
774 sd->precond);
775 }
776 }
777 break;
778
779 case CS_SLES_HYPRE_FLEXGMRES:
780 {
781 HYPRE_ParCSRFlexGMRESCreate(comm, &hs);
782
783 if (verbosity > 2 ) {
784 HYPRE_FlexGMRESSetPrintLevel(hs, 2); /* Print solve info */
785 HYPRE_FlexGMRESSetLogging(hs, 1); /* Needed to get run info later */
786 }
787
788 solve_ftn[i] = HYPRE_ParCSRFlexGMRESSolve;
789 setup_ftn[i] = HYPRE_ParCSRFlexGMRESSetup;
790
791 if (i == 0) { /* preconditioner */
792 HYPRE_ParCSRFlexGMRESSetTol(hs, 0.0);
793 HYPRE_ParCSRFlexGMRESSetMaxIter(hs, 1);
794 }
795 else { /* solver */
796 HYPRE_FlexGMRESSetMaxIter(hs, 1000); /* Max iterations */
797 if (solve_ftn[0] != NULL)
798 HYPRE_ParCSRFlexGMRESSetPrecond(hs,
799 solve_ftn[0],
800 setup_ftn[0],
801 sd->precond);
802 }
803 }
804 break;
805
806 case CS_SLES_HYPRE_LGMRES:
807 {
808 HYPRE_ParCSRLGMRESCreate(comm, &hs);
809
810 if (verbosity > 2 ) {
811 HYPRE_LGMRESSetPrintLevel(hs, 2); /* Print solve info */
812 HYPRE_LGMRESSetLogging(hs, 1); /* Needed to get run info later */
813 }
814
815 solve_ftn[i] = HYPRE_ParCSRLGMRESSolve;
816 setup_ftn[i] = HYPRE_ParCSRLGMRESSetup;
817
818 if (i == 0) { /* preconditioner */
819 HYPRE_ParCSRLGMRESSetTol(hs, 0.0);
820 HYPRE_ParCSRLGMRESSetMaxIter(hs, 1);
821 }
822 else { /* solver */
823 HYPRE_LGMRESSetMaxIter(hs, 1000); /* Max iterations */
824 if (solve_ftn[0] != NULL)
825 HYPRE_ParCSRLGMRESSetPrecond(hs,
826 solve_ftn[0],
827 setup_ftn[0],
828 sd->precond);
829 }
830 }
831 break;
832
833 case CS_SLES_HYPRE_PCG:
834 {
835 HYPRE_ParCSRPCGCreate(comm, &hs);
836
837 if (verbosity > 2 ) {
838 HYPRE_PCGSetPrintLevel(hs, 2); /* Print solve info */
839 HYPRE_PCGSetLogging(hs, 1); /* Needed to get run info later */
840 }
841
842 solve_ftn[i] = HYPRE_ParCSRPCGSolve;
843 setup_ftn[i] = HYPRE_ParCSRPCGSetup;
844
845 if (i == 0) { /* preconditioner */
846 HYPRE_ParCSRPCGSetTol(hs, 0.0);
847 HYPRE_ParCSRPCGSetMaxIter(hs, 1);
848 }
849 else { /* solver */
850 HYPRE_PCGSetMaxIter(hs, 1000); /* Max iterations */
851 if (solve_ftn[0] != NULL)
852 HYPRE_ParCSRPCGSetPrecond(hs,
853 solve_ftn[0],
854 setup_ftn[0],
855 sd->precond);
856 }
857 }
858 break;
859
860 case CS_SLES_HYPRE_EUCLID:
861 {
862 HYPRE_EuclidCreate(comm, &hs);
863
864 solve_ftn[i] = HYPRE_EuclidSolve;
865 setup_ftn[i] = HYPRE_EuclidSetup;
866
867 if (i > 0) /* solver */
868 bft_error(__FILE__, __LINE__, 0,
869 _("HYPRE: type (%s) is a preconditioner, not a solver."),
870 _cs_hypre_type_name(c->solver_type));
871
872 }
873 break;
874
875 case CS_SLES_HYPRE_PARASAILS:
876 {
877 HYPRE_ParCSRParaSailsCreate(comm, &hs);
878
879 if (verbosity > 2 ) {
880 HYPRE_ParCSRParaSailsSetLogging(hs, 1); /* Needed to get run
881 info later */
882 }
883
884 solve_ftn[i] = HYPRE_ParCSRParaSailsSolve;
885 setup_ftn[i] = HYPRE_ParCSRParaSailsSetup;
886 if (i > 0) /* solver */
887 bft_error(__FILE__, __LINE__, 0,
888 _("HYPRE: type (%s) is a preconditioner, not a solver."),
889 _cs_hypre_type_name(c->solver_type));
890 }
891 break;
892
893 default:
894 bft_error(__FILE__, __LINE__, 0,
895 _("HYPRE: solver type (%s) not currently handled."),
896 _cs_hypre_type_name(c->solver_type));
897
898 }
899
900 if (i == 0)
901 sd->precond = hs;
902 else
903 sd->solver = hs;
904 }
905
906 if (sd->precond != NULL && have_set_pc == false)
907 bft_error(__FILE__, __LINE__, 0,
908 _("HYPRE: solver (%s) will ignore preconditioner (%s)."),
909 _cs_hypre_type_name(c->solver_type),
910 _cs_hypre_type_name(c->precond_type));
911
912 /* Call optional setup hook for user setting changes */
913
914 if (c->setup_hook != NULL)
915 c->setup_hook(verbosity, c->hook_context, &(sd->solver));
916
917 /* Now setup systems (where rhs and vx values may be different
918 when solving, but their shapes and addresses are the same) */
919
920 HYPRE_ParCSRMatrix par_a; /* Associted matrix */
921 HYPRE_ParVector p_x, p_rhs;
922
923 HYPRE_IJMatrixGetObject(sd->coeffs->hm, (void **)&par_a);
924 HYPRE_IJVectorGetObject(sd->coeffs->hx, (void **)&p_x);
925 HYPRE_IJVectorGetObject(sd->coeffs->hy, (void **)&p_rhs);
926
927 if (setup_ftn[1] != NULL)
928 setup_ftn[1](sd->solver, par_a, p_rhs, p_x);
929 else
930 bft_error(__FILE__, __LINE__, 0,
931 _("HYPRE: setup function for solver type (%s) not set."),
932 _cs_hypre_type_name(c->solver_type));
933
934 /* Update return values */
935 c->n_setups += 1;
936
937 cs_timer_t t1 = cs_timer_time();
938 cs_timer_counter_add_diff(&(c->t_setup), &t0, &t1);
939 }
940
941 /*----------------------------------------------------------------------------*/
942 /*!
943 * \brief Call HYPRE linear equation solver.
944 *
945 * \param[in, out] context pointer to iterative solver info and context
946 * (actual type: cs_sles_hypre_t *)
947 * \param[in] name pointer to system name
948 * \param[in] a matrix
949 * \param[in] verbosity associated verbosity
950 * \param[in] precision solver precision
951 * \param[in] r_norm residue normalization
952 * \param[out] n_iter number of "equivalent" iterations
953 * \param[out] residue residue
954 * \param[in] rhs right hand side
955 * \param[in, out] vx system solution
956 * \param[in] aux_size number of elements in aux_vectors (in bytes)
957 * \param aux_vectors optional working area
958 * (internal allocation if NULL)
959 *
960 * \return convergence state
961 */
962 /*----------------------------------------------------------------------------*/
963
964 cs_sles_convergence_state_t
cs_sles_hypre_solve(void * context,const char * name,const cs_matrix_t * a,int verbosity,double precision,double r_norm,int * n_iter,double * residue,const cs_real_t * rhs,cs_real_t * vx,size_t aux_size,void * aux_vectors)965 cs_sles_hypre_solve(void *context,
966 const char *name,
967 const cs_matrix_t *a,
968 int verbosity,
969 double precision,
970 double r_norm,
971 int *n_iter,
972 double *residue,
973 const cs_real_t *rhs,
974 cs_real_t *vx,
975 size_t aux_size,
976 void *aux_vectors)
977 {
978 CS_UNUSED(aux_size);
979 CS_UNUSED(aux_vectors);
980
981 cs_timer_t t0;
982 t0 = cs_timer_time();
983 cs_sles_convergence_state_t cvg = CS_SLES_ITERATING;
984 cs_sles_hypre_t *c = context;
985 cs_sles_hypre_setup_t *sd = c->setup_data;
986 cs_lnum_t n_rows = cs_matrix_get_n_rows(a);
987 HYPRE_Int its;
988 double res;
989
990 precision = precision;
991 if (sd == NULL) {
992 cs_sles_hypre_setup(c, name, a, verbosity);
993 sd = c->setup_data;
994 }
995
996 HYPRE_ParCSRMatrix par_a; /* Associted matrix */
997 HYPRE_IJMatrixGetObject(sd->coeffs->hm, (void **)&par_a);
998
999 cs_alloc_mode_t amode = CS_ALLOC_HOST;
1000 if (sd->coeffs->memory_location != HYPRE_MEMORY_HOST)
1001 amode = CS_ALLOC_HOST_DEVICE_SHARED;
1002
1003 HYPRE_Real *_t = NULL;
1004
1005 /* Set RHS and starting solution */
1006
1007 if (sizeof(cs_real_t) == sizeof(HYPRE_Real) && amode == CS_ALLOC_HOST) {
1008 HYPRE_IJVectorSetValues(sd->coeffs->hx, n_rows, NULL, vx);
1009 HYPRE_IJVectorSetValues(sd->coeffs->hy, n_rows, NULL, rhs);
1010 }
1011 else {
1012 CS_MALLOC_HD(_t, n_rows, HYPRE_Real, amode);
1013 for (HYPRE_BigInt ii = 0; ii < n_rows; ii++) {
1014 _t[ii] = vx[ii];;
1015 }
1016 HYPRE_IJVectorSetValues(sd->coeffs->hx, n_rows, NULL, _t);
1017 for (HYPRE_BigInt ii = 0; ii < n_rows; ii++) {
1018 _t[ii] = rhs[ii];;
1019 }
1020 HYPRE_IJVectorSetValues(sd->coeffs->hy, n_rows, NULL, _t);
1021 }
1022
1023 HYPRE_IJVectorAssemble(sd->coeffs->hx);
1024 HYPRE_IJVectorAssemble(sd->coeffs->hy);
1025
1026 HYPRE_ParVector p_x, p_rhs;
1027 HYPRE_IJVectorGetObject(sd->coeffs->hx, (void **)&p_x);
1028 HYPRE_IJVectorGetObject(sd->coeffs->hy, (void **)&p_rhs);
1029
1030 switch(c->solver_type) {
1031
1032 case CS_SLES_HYPRE_BOOMERAMG:
1033 {
1034 /* Finalize setup and solve; no absolute tolerance is available,
1035 so we use the available function. */
1036 HYPRE_BoomerAMGSetTol(sd->solver, precision*r_norm);
1037 HYPRE_BoomerAMGSolve(sd->solver, par_a, p_rhs, p_x);
1038
1039 /* Get solution and information */
1040 HYPRE_BoomerAMGGetFinalRelativeResidualNorm(sd->solver, &res);
1041 HYPRE_BoomerAMGGetNumIterations(sd->solver, &its);
1042 }
1043 break;
1044
1045 case CS_SLES_HYPRE_HYBRID:
1046 {
1047 /* Finalize setup and solve */
1048 HYPRE_ParCSRHybridSetAbsoluteTol(sd->solver, precision*r_norm);
1049 HYPRE_ParCSRHybridSolve(sd->solver, par_a, p_rhs, p_x);
1050
1051 /* Get solution and information */
1052 HYPRE_ParCSRHybridGetFinalRelativeResidualNorm(sd->solver, &res);
1053 HYPRE_ParCSRHybridGetNumIterations(sd->solver, &its);
1054 }
1055 break;
1056
1057 case CS_SLES_HYPRE_ILU:
1058 {
1059 /* Finalize setup and solve */
1060 HYPRE_ILUSetTol(sd->solver, precision*r_norm);
1061 HYPRE_ILUSolve(sd->solver, par_a, p_rhs, p_x);
1062
1063 /* Get solution and information */
1064 HYPRE_ILUGetFinalRelativeResidualNorm(sd->solver, &res);
1065 HYPRE_ILUGetNumIterations(sd->solver, &its);
1066 }
1067 break;
1068
1069 case CS_SLES_HYPRE_BICGSTAB:
1070 {
1071 /* Finalize setup and solve */
1072 HYPRE_BiCGSTABSetAbsoluteTol(sd->solver, precision*r_norm);
1073 HYPRE_ParCSRBiCGSTABSolve(sd->solver, par_a, p_rhs, p_x);
1074
1075 /* Get solution and information */
1076 HYPRE_ParCSRBiCGSTABGetFinalRelativeResidualNorm(sd->solver, &res);
1077 HYPRE_ParCSRBiCGSTABGetNumIterations(sd->solver, &its);
1078 }
1079 break;
1080
1081 case CS_SLES_HYPRE_GMRES:
1082 {
1083 /* Finalize setup and solve */
1084 HYPRE_GMRESSetAbsoluteTol(sd->solver, precision*r_norm);
1085 HYPRE_ParCSRGMRESSolve(sd->solver, par_a, p_rhs, p_x);
1086
1087 /* Get solution and information */
1088 HYPRE_ParCSRGMRESGetFinalRelativeResidualNorm(sd->solver, &res);
1089 HYPRE_ParCSRGMRESGetNumIterations(sd->solver, &its);
1090 }
1091 break;
1092
1093 case CS_SLES_HYPRE_FLEXGMRES:
1094 {
1095 /* Finalize setup and solve */
1096 HYPRE_FlexGMRESSetAbsoluteTol(sd->solver, precision*r_norm);
1097 HYPRE_ParCSRFlexGMRESSolve(sd->solver, par_a, p_rhs, p_x);
1098
1099 /* Get solution and information */
1100 HYPRE_ParCSRFlexGMRESGetFinalRelativeResidualNorm(sd->solver, &res);
1101 HYPRE_ParCSRFlexGMRESGetNumIterations(sd->solver, &its);
1102 }
1103 break;
1104
1105 case CS_SLES_HYPRE_LGMRES:
1106 {
1107 /* Finalize setup and solve */
1108 HYPRE_LGMRESSetAbsoluteTol(sd->solver, precision*r_norm);
1109 HYPRE_ParCSRLGMRESSolve(sd->solver, par_a, p_rhs, p_x);
1110
1111 /* Get solution and information */
1112 HYPRE_ParCSRLGMRESGetFinalRelativeResidualNorm(sd->solver, &res);
1113 HYPRE_ParCSRLGMRESGetNumIterations(sd->solver, &its);
1114 }
1115 break;
1116
1117 case CS_SLES_HYPRE_PCG:
1118 {
1119 /* Finalize setup and solve */
1120 HYPRE_PCGSetAbsoluteTol(sd->solver, precision*r_norm);
1121 HYPRE_ParCSRPCGSolve(sd->solver, par_a, p_rhs, p_x);
1122
1123 /* Get solution and information */
1124 HYPRE_ParCSRPCGGetFinalRelativeResidualNorm(sd->solver, &res);
1125 HYPRE_ParCSRPCGGetNumIterations(sd->solver, &its);
1126 }
1127 break;
1128
1129 default:
1130 bft_error(__FILE__, __LINE__, 0,
1131 _("HYPRE: solver type (%s) not handled."),
1132 _cs_hypre_type_name(c->solver_type));
1133 }
1134
1135 if (_t == NULL) {
1136 HYPRE_IJVectorGetValues(sd->coeffs->hx, n_rows, NULL, vx);
1137 }
1138 else {
1139 HYPRE_IJVectorGetValues(sd->coeffs->hx, n_rows, NULL, _t);
1140 for (HYPRE_BigInt ii = 0; ii < n_rows; ii++) {
1141 vx[ii] = _t[ii];
1142 }
1143 CS_FREE_HD(_t);
1144 }
1145
1146 *residue = res;
1147 *n_iter = its;
1148
1149 /* Update return values */
1150 if (c->n_solves == 0)
1151 c->n_iterations_min = its;
1152
1153 c->n_iterations_last = its;
1154 c->n_iterations_tot += its;
1155 if (c->n_iterations_min > its)
1156 c->n_iterations_min = its;
1157 if (c->n_iterations_max < its)
1158 c->n_iterations_max = its;
1159 c->n_solves += 1;
1160 cs_timer_t t1 = cs_timer_time();
1161 cs_timer_counter_add_diff(&(c->t_solve), &t0, &t1);
1162
1163 return cvg;
1164 }
1165
1166 /*----------------------------------------------------------------------------*/
1167 /*!
1168 * \brief Free HYPRE linear equation solver setup context.
1169 *
1170 * This function frees resolution-related data, such as
1171 * buffers and preconditioning but does not free the whole context,
1172 * as info used for logging (especially performance data) is maintained.
1173 *
1174 * \param[in, out] context pointer to iterative solver info and context
1175 * (actual type: cs_sles_hypre_t *)
1176 */
1177 /*----------------------------------------------------------------------------*/
1178
1179 void
cs_sles_hypre_free(void * context)1180 cs_sles_hypre_free(void *context)
1181 {
1182 cs_timer_t t0;
1183 t0 = cs_timer_time();
1184
1185 cs_sles_hypre_t *c = context;
1186
1187 if (c->setup_data != NULL) {
1188
1189 cs_sles_hypre_setup_t *sd = c->setup_data;
1190
1191 for (int i = 0; i < 2; i++) {
1192 HYPRE_Solver hs = (i == 0) ? sd->solver : sd->precond;
1193 cs_sles_hypre_type_t hs_type = (i == 0) ? c->solver_type : c->precond_type;
1194
1195 if (hs == NULL)
1196 continue;
1197
1198 switch(hs_type) {
1199
1200 case CS_SLES_HYPRE_BOOMERAMG:
1201 HYPRE_BoomerAMGDestroy(hs);
1202 break;
1203
1204 case CS_SLES_HYPRE_HYBRID:
1205 HYPRE_ParCSRHybridDestroy(hs);
1206 break;
1207
1208 case CS_SLES_HYPRE_ILU:
1209 HYPRE_ILUDestroy(hs);
1210 break;
1211
1212 case CS_SLES_HYPRE_BICGSTAB:
1213 HYPRE_ParCSRBiCGSTABDestroy(hs);
1214 break;
1215
1216 case CS_SLES_HYPRE_GMRES:
1217 HYPRE_ParCSRGMRESDestroy(hs);
1218 break;
1219
1220 case CS_SLES_HYPRE_FLEXGMRES:
1221 HYPRE_ParCSRFlexGMRESDestroy(hs);
1222 break;
1223
1224 case CS_SLES_HYPRE_LGMRES:
1225 HYPRE_ParCSRLGMRESDestroy(hs);
1226 break;
1227
1228 case CS_SLES_HYPRE_PCG:
1229 HYPRE_ParCSRPCGDestroy(hs);
1230 break;
1231
1232 default:
1233 bft_error(__FILE__, __LINE__, 0,
1234 _("HYPRE: solver type (%s) not handled."),
1235 _cs_hypre_type_name(c->solver_type));
1236 }
1237
1238 if (i == 0)
1239 sd->solver = NULL;
1240 else
1241 sd->precond = NULL;
1242 }
1243
1244 BFT_FREE(c->setup_data);
1245 }
1246
1247 cs_timer_t t1 = cs_timer_time();
1248 cs_timer_counter_add_diff(&(c->t_setup), &t0, &t1);
1249 }
1250
1251 /*----------------------------------------------------------------------------*/
1252 /*!
1253 * \brief Define whether the solver should run on the host or
1254 * accelerated device.
1255 *
1256 * If no device is available, this setting is ignored.
1257 *
1258 * \param[in,out] context pointer to HYPRE linear solver info
1259 * \param[in] use_device 0 for host, 1 for device (GPU)
1260 */
1261 /*----------------------------------------------------------------------------*/
1262
1263 void
cs_sles_hypre_set_host_device(cs_sles_hypre_t * context,int use_device)1264 cs_sles_hypre_set_host_device(cs_sles_hypre_t *context,
1265 int use_device)
1266 {
1267 int device_id = cs_get_device_id();
1268 if (device_id < 0)
1269 use_device = 0;
1270
1271 context->use_device = use_device;
1272 }
1273
1274 /*----------------------------------------------------------------------------*/
1275 /*!
1276 * \brief Query whether the solver should run on the host or accelerated device.
1277 *
1278 * \param[in,out] context pointer to HYPRE linear solver info
1279 *
1280 * \return 0 for host, 1 for device (GPU)
1281 */
1282 /*----------------------------------------------------------------------------*/
1283
1284 int
cs_sles_hypre_get_host_device(const cs_sles_hypre_t * context)1285 cs_sles_hypre_get_host_device(const cs_sles_hypre_t *context)
1286 {
1287 return context->use_device;
1288 }
1289
1290 /*----------------------------------------------------------------------------*/
1291 /*!
1292 * \brief Print information on hypre library.
1293 *
1294 * \param[in] log_type log type
1295 */
1296 /*----------------------------------------------------------------------------*/
1297
1298 void
cs_sles_hypre_library_info(cs_log_t log_type)1299 cs_sles_hypre_library_info(cs_log_t log_type)
1300 {
1301 char hypre_config_options[256] = "";
1302 size_t l_max = 255;
1303 size_t l = l_max;
1304
1305 snprintf(hypre_config_options, l, "%s %s (",
1306 HYPRE_RELEASE_NAME, HYPRE_RELEASE_VERSION);
1307
1308 l -= strlen(hypre_config_options);
1309
1310 #if defined(HYPRE_USING_GPU)
1311 strncat(hypre_config_options, "GPU support, ", l);
1312 l -= strlen(hypre_config_options);
1313 #endif
1314
1315 #if defined(HYPRE_SEQUENTIAL)
1316 strncat(hypre_config_options, "sequential, ", l);
1317 l -= strlen(hypre_config_options);
1318 #endif
1319
1320 #if defined(HYPRE_USING_OPENMP)
1321 strncat(hypre_config_options, "OpenMP, ", l);
1322 l -= strlen(hypre_config_options);
1323 #endif
1324
1325
1326 #if defined(HYPRE_BIGINT)
1327 strncat(hypre_config_options, "large integers, ", l);
1328 l -= strlen(hypre_config_options);
1329 #endif
1330
1331 l = strlen(hypre_config_options) - 2;
1332 if (hypre_config_options[l] == ',') {
1333 hypre_config_options[l] = ')';
1334 hypre_config_options[l+1] = '\0';
1335 }
1336 else {
1337 hypre_config_options[l] = '\0';
1338 }
1339
1340 cs_log_printf(log_type,
1341 " %s\n", hypre_config_options);
1342 }
1343
1344 /*----------------------------------------------------------------------------*/
1345
1346 END_C_DECLS
1347