1 /******************************************************************************
2 * Copyright 1998-2019 Lawrence Livermore National Security, LLC and other
3 * HYPRE Project Developers. See the top-level COPYRIGHT file for details.
4 *
5 * SPDX-License-Identifier: (Apache-2.0 OR MIT)
6 ******************************************************************************/
7
8 /******************************************************************************
9 *
10 *
11 *****************************************************************************/
12
13 #include "_hypre_parcsr_ls.h"
14
15 /*--------------------------------------------------------------------------
16 * hypre_AMGHybridData:
17 *--------------------------------------------------------------------------*/
18
19 typedef struct
20 {
21 HYPRE_Real tol;
22 HYPRE_Real a_tol;
23 HYPRE_Real cf_tol;
24 HYPRE_Int dscg_max_its;
25 HYPRE_Int pcg_max_its;
26 HYPRE_Int two_norm;
27 HYPRE_Int stop_crit;
28 HYPRE_Int rel_change;
29 HYPRE_Int recompute_residual;
30 HYPRE_Int recompute_residual_p;
31 HYPRE_Int solver_type;
32 HYPRE_Int k_dim;
33
34 HYPRE_Int pcg_default; /* boolean */
35 HYPRE_Int (*pcg_precond_solve)(void*,void*,void*,void*);
36 HYPRE_Int (*pcg_precond_setup)(void*,void*,void*,void*);
37 void *pcg_precond;
38 void *pcg_solver;
39
40 /* log info (always logged) */
41 HYPRE_Int dscg_num_its;
42 HYPRE_Int pcg_num_its;
43 HYPRE_Real final_rel_res_norm;
44 HYPRE_Int time_index;
45
46 HYPRE_Real setup_time1;
47 HYPRE_Real setup_time2;
48 HYPRE_Real solve_time1;
49 HYPRE_Real solve_time2;
50
51 MPI_Comm comm;
52
53 /* additional information (place-holder currently used to print norms) */
54 HYPRE_Int logging;
55 HYPRE_Int print_level;
56
57 /* info for BoomerAMG */
58 HYPRE_Real strong_threshold;
59 HYPRE_Real max_row_sum;
60 HYPRE_Real trunc_factor;
61 HYPRE_Int pmax;
62 HYPRE_Int setup_type;
63 HYPRE_Int max_levels;
64 HYPRE_Int measure_type;
65 HYPRE_Int coarsen_type;
66 HYPRE_Int interp_type;
67 HYPRE_Int cycle_type;
68 HYPRE_Int relax_order;
69 HYPRE_Int keepT;
70 HYPRE_Int max_coarse_size;
71 HYPRE_Int min_coarse_size;
72 HYPRE_Int seq_threshold;
73 HYPRE_Int *num_grid_sweeps;
74 HYPRE_Int *grid_relax_type;
75 HYPRE_Int **grid_relax_points;
76 HYPRE_Real *relax_weight;
77 HYPRE_Real *omega;
78 HYPRE_Int num_paths;
79 HYPRE_Int agg_num_levels;
80 HYPRE_Int agg_interp_type;
81 HYPRE_Int num_functions;
82 HYPRE_Int nodal;
83 HYPRE_Int *dof_func;
84
85 /* data needed for non-Galerkin option */
86 HYPRE_Int nongalerk_num_tol;
87 HYPRE_Real *nongalerkin_tol;
88 } hypre_AMGHybridData;
89
90 /*--------------------------------------------------------------------------
91 * hypre_AMGHybridCreate
92 *--------------------------------------------------------------------------*/
93
94 void *
hypre_AMGHybridCreate()95 hypre_AMGHybridCreate( )
96 {
97 hypre_AMGHybridData *AMGhybrid_data;
98
99 AMGhybrid_data = hypre_CTAlloc(hypre_AMGHybridData, 1, HYPRE_MEMORY_HOST);
100
101 (AMGhybrid_data -> time_index) = hypre_InitializeTiming("AMGHybrid");
102
103 /* set defaults */
104 (AMGhybrid_data -> tol) = 1.0e-06;
105 (AMGhybrid_data -> a_tol) = 0.0;
106 (AMGhybrid_data -> cf_tol) = 0.90;
107 (AMGhybrid_data -> dscg_max_its) = 1000;
108 (AMGhybrid_data -> pcg_max_its) = 200;
109 (AMGhybrid_data -> two_norm) = 0;
110 (AMGhybrid_data -> stop_crit) = 0;
111 (AMGhybrid_data -> rel_change) = 0;
112 (AMGhybrid_data -> pcg_default) = 1;
113 (AMGhybrid_data -> solver_type) = 1;
114 (AMGhybrid_data -> pcg_precond_solve) = NULL;
115 (AMGhybrid_data -> pcg_precond_setup) = NULL;
116 (AMGhybrid_data -> pcg_precond) = NULL;
117 (AMGhybrid_data -> pcg_solver) = NULL;
118 (AMGhybrid_data -> setup_time1) = 0.0;
119 (AMGhybrid_data -> setup_time2) = 0.0;
120 (AMGhybrid_data -> solve_time1) = 0.0;
121 (AMGhybrid_data -> solve_time2) = 0.0;
122
123 /* initialize */
124 (AMGhybrid_data -> dscg_num_its) = 0;
125 (AMGhybrid_data -> pcg_num_its) = 0;
126 (AMGhybrid_data -> logging) = 0;
127 (AMGhybrid_data -> print_level) = 0;
128 (AMGhybrid_data -> k_dim) = 5;
129
130 /* BoomerAMG info */
131 (AMGhybrid_data -> setup_type) = 1;
132 (AMGhybrid_data -> strong_threshold) = 0.25;
133 (AMGhybrid_data -> max_row_sum) = 0.9;
134 (AMGhybrid_data -> trunc_factor) = 0.0;
135 (AMGhybrid_data -> pmax) = 4;
136 (AMGhybrid_data -> max_levels) = 25;
137 (AMGhybrid_data -> measure_type) = 0;
138 (AMGhybrid_data -> coarsen_type) = 10;
139 (AMGhybrid_data -> interp_type) = 6;
140 (AMGhybrid_data -> cycle_type) = 1;
141 (AMGhybrid_data -> relax_order) = 0;
142 (AMGhybrid_data -> keepT) = 0;
143 (AMGhybrid_data -> max_coarse_size) = 9;
144 (AMGhybrid_data -> min_coarse_size) = 1;
145 (AMGhybrid_data -> seq_threshold) = 0;
146 (AMGhybrid_data -> num_grid_sweeps) = NULL;
147 (AMGhybrid_data -> grid_relax_type) = NULL;
148 (AMGhybrid_data -> grid_relax_points) = NULL;
149 (AMGhybrid_data -> relax_weight) = NULL;
150 (AMGhybrid_data -> omega) = NULL;
151 (AMGhybrid_data -> agg_num_levels) = 0;
152 (AMGhybrid_data -> agg_interp_type) = 4;
153 (AMGhybrid_data -> num_paths) = 1;
154 (AMGhybrid_data -> num_functions) = 1;
155 (AMGhybrid_data -> nodal) = 0;
156 (AMGhybrid_data -> dof_func) = NULL;
157 (AMGhybrid_data -> nongalerk_num_tol) = 0;
158 (AMGhybrid_data -> nongalerkin_tol) = NULL;
159
160 return (void *) AMGhybrid_data;
161 }
162
163 /*-------------------------------------------------------------------------- *
164 hypre_AMGHybridDestroy
165 *--------------------------------------------------------------------------*/
166
167 HYPRE_Int
hypre_AMGHybridDestroy(void * AMGhybrid_vdata)168 hypre_AMGHybridDestroy( void *AMGhybrid_vdata )
169 {
170 hypre_AMGHybridData *AMGhybrid_data = (hypre_AMGHybridData *)AMGhybrid_vdata;
171 HYPRE_Int i;
172
173 if (AMGhybrid_data)
174 {
175 HYPRE_Int solver_type = (AMGhybrid_data -> solver_type);
176 /*HYPRE_Int pcg_default = (AMGhybrid_data -> pcg_default);*/
177 void *pcg_solver = (AMGhybrid_data -> pcg_solver);
178 void *pcg_precond = (AMGhybrid_data -> pcg_precond);
179
180 if (pcg_precond) hypre_BoomerAMGDestroy(pcg_precond);
181 if (solver_type == 1) hypre_PCGDestroy(pcg_solver);
182 if (solver_type == 2) hypre_GMRESDestroy(pcg_solver);
183 if (solver_type == 3) hypre_BiCGSTABDestroy(pcg_solver);
184
185 if (AMGhybrid_data -> num_grid_sweeps)
186 {
187 hypre_TFree( (AMGhybrid_data -> num_grid_sweeps) , HYPRE_MEMORY_HOST);
188 (AMGhybrid_data -> num_grid_sweeps) = NULL;
189 }
190 if (AMGhybrid_data -> grid_relax_type)
191 {
192 hypre_TFree( (AMGhybrid_data -> grid_relax_type) , HYPRE_MEMORY_HOST);
193 (AMGhybrid_data -> grid_relax_type) = NULL;
194 }
195 if (AMGhybrid_data -> grid_relax_points)
196 {
197 for (i=0; i < 4; i++)
198 {
199 hypre_TFree( (AMGhybrid_data -> grid_relax_points)[i] , HYPRE_MEMORY_HOST);
200 }
201 hypre_TFree( (AMGhybrid_data -> grid_relax_points) , HYPRE_MEMORY_HOST);
202 (AMGhybrid_data -> grid_relax_points) = NULL;
203 }
204 if (AMGhybrid_data -> relax_weight)
205 {
206 hypre_TFree( (AMGhybrid_data -> relax_weight) , HYPRE_MEMORY_HOST);
207 (AMGhybrid_data -> relax_weight) = NULL;
208 }
209 if (AMGhybrid_data -> omega)
210 {
211 hypre_TFree( (AMGhybrid_data -> omega) , HYPRE_MEMORY_HOST);
212 (AMGhybrid_data -> omega) = NULL;
213 }
214 if (AMGhybrid_data -> dof_func)
215 {
216 hypre_TFree( (AMGhybrid_data -> dof_func) , HYPRE_MEMORY_HOST);
217 (AMGhybrid_data -> dof_func) = NULL;
218 }
219 hypre_TFree(AMGhybrid_data, HYPRE_MEMORY_HOST);
220 }
221
222 return hypre_error_flag;
223 }
224
225 /*--------------------------------------------------------------------------
226 * hypre_AMGHybridSetTol
227 *--------------------------------------------------------------------------*/
228
229 HYPRE_Int
hypre_AMGHybridSetTol(void * AMGhybrid_vdata,HYPRE_Real tol)230 hypre_AMGHybridSetTol( void *AMGhybrid_vdata,
231 HYPRE_Real tol )
232 {
233 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
234
235 if (!AMGhybrid_data)
236 {
237 hypre_error_in_arg(1);
238 return hypre_error_flag;
239 }
240 if (tol < 0 || tol > 1)
241 {
242 hypre_error_in_arg(2);
243 return hypre_error_flag;
244 }
245 (AMGhybrid_data -> tol) = tol;
246
247 return hypre_error_flag;
248 }
249 /*--------------------------------------------------------------------------
250 * hypre_AMGHybridSetAbsoluteTol
251 *--------------------------------------------------------------------------*/
252
253 HYPRE_Int
hypre_AMGHybridSetAbsoluteTol(void * AMGhybrid_vdata,HYPRE_Real a_tol)254 hypre_AMGHybridSetAbsoluteTol( void *AMGhybrid_vdata,
255 HYPRE_Real a_tol )
256 {
257 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
258
259 if (!AMGhybrid_data)
260 {
261 hypre_error_in_arg(1);
262 return hypre_error_flag;
263 }
264 if (a_tol < 0 || a_tol > 1)
265 {
266 hypre_error_in_arg(2);
267 return hypre_error_flag;
268 }
269 (AMGhybrid_data -> a_tol) = a_tol;
270
271 return hypre_error_flag;
272 }
273 /*--------------------------------------------------------------------------
274 * hypre_AMGHybridSetConvergenceTol
275 *--------------------------------------------------------------------------*/
276
277 HYPRE_Int
hypre_AMGHybridSetConvergenceTol(void * AMGhybrid_vdata,HYPRE_Real cf_tol)278 hypre_AMGHybridSetConvergenceTol( void *AMGhybrid_vdata,
279 HYPRE_Real cf_tol )
280 {
281 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
282 if (!AMGhybrid_data)
283 {
284 hypre_error_in_arg(1);
285 return hypre_error_flag;
286 }
287 if (cf_tol < 0 || cf_tol > 1)
288 {
289 hypre_error_in_arg(2);
290 return hypre_error_flag;
291 }
292
293 (AMGhybrid_data -> cf_tol) = cf_tol;
294
295 return hypre_error_flag;
296 }
297
298 /*--------------------------------------------------------------------------
299 * hypre_AMGHybridSetNonGalerkinTol
300 *--------------------------------------------------------------------------*/
301
302 HYPRE_Int
hypre_AMGHybridSetNonGalerkinTol(void * AMGhybrid_vdata,HYPRE_Int nongalerk_num_tol,HYPRE_Real * nongalerkin_tol)303 hypre_AMGHybridSetNonGalerkinTol( void *AMGhybrid_vdata,
304 HYPRE_Int nongalerk_num_tol,
305 HYPRE_Real *nongalerkin_tol )
306 {
307 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
308 if (!AMGhybrid_data)
309 {
310 hypre_error_in_arg(1);
311 return hypre_error_flag;
312 }
313 if (nongalerk_num_tol < 0)
314 {
315 hypre_error_in_arg(2);
316 return hypre_error_flag;
317 }
318
319 (AMGhybrid_data -> nongalerk_num_tol) = nongalerk_num_tol;
320 (AMGhybrid_data -> nongalerkin_tol) = nongalerkin_tol;
321
322 return hypre_error_flag;
323 }
324
325 /*--------------------------------------------------------------------------
326 * hypre_AMGHybridSetDSCGMaxIter
327 *--------------------------------------------------------------------------*/
328
329 HYPRE_Int
hypre_AMGHybridSetDSCGMaxIter(void * AMGhybrid_vdata,HYPRE_Int dscg_max_its)330 hypre_AMGHybridSetDSCGMaxIter( void *AMGhybrid_vdata,
331 HYPRE_Int dscg_max_its )
332 {
333 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
334 if (!AMGhybrid_data)
335 {
336 hypre_error_in_arg(1);
337 return hypre_error_flag;
338 }
339 if (dscg_max_its < 0)
340 {
341 hypre_error_in_arg(2);
342 return hypre_error_flag;
343 }
344
345 (AMGhybrid_data -> dscg_max_its) = dscg_max_its;
346
347 return hypre_error_flag;
348 }
349
350 /*--------------------------------------------------------------------------
351 * hypre_AMGHybridSetPCGMaxIter
352 *--------------------------------------------------------------------------*/
353
354 HYPRE_Int
hypre_AMGHybridSetPCGMaxIter(void * AMGhybrid_vdata,HYPRE_Int pcg_max_its)355 hypre_AMGHybridSetPCGMaxIter( void *AMGhybrid_vdata,
356 HYPRE_Int pcg_max_its )
357 {
358 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
359 if (!AMGhybrid_data)
360 {
361 hypre_error_in_arg(1);
362 return hypre_error_flag;
363 }
364 if (pcg_max_its < 0)
365 {
366 hypre_error_in_arg(2);
367 return hypre_error_flag;
368 }
369
370 (AMGhybrid_data -> pcg_max_its) = pcg_max_its;
371
372 return hypre_error_flag;
373 }
374
375 /*--------------------------------------------------------------------------
376 * hypre_AMGHybridSetSetupType
377 *--------------------------------------------------------------------------*/
378
379 HYPRE_Int
hypre_AMGHybridSetSetupType(void * AMGhybrid_vdata,HYPRE_Int setup_type)380 hypre_AMGHybridSetSetupType( void *AMGhybrid_vdata,
381 HYPRE_Int setup_type )
382 {
383 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
384 if (!AMGhybrid_data)
385 {
386 hypre_error_in_arg(1);
387 return hypre_error_flag;
388 }
389
390 (AMGhybrid_data -> setup_type) = setup_type;
391
392 return hypre_error_flag;
393 }
394
395 /*--------------------------------------------------------------------------
396 * hypre_AMGHybridSetSolverType
397 *--------------------------------------------------------------------------*/
398
399 HYPRE_Int
hypre_AMGHybridSetSolverType(void * AMGhybrid_vdata,HYPRE_Int solver_type)400 hypre_AMGHybridSetSolverType( void *AMGhybrid_vdata,
401 HYPRE_Int solver_type )
402 {
403 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
404 if (!AMGhybrid_data)
405 {
406 hypre_error_in_arg(1);
407 return hypre_error_flag;
408 }
409
410 (AMGhybrid_data -> solver_type) = solver_type;
411
412 return hypre_error_flag;
413 }
414
415 /*--------------------------------------------------------------------------
416 * hypre_AMGHybridSetRecomputeResidual
417 *--------------------------------------------------------------------------*/
418
419 HYPRE_Int
hypre_AMGHybridSetRecomputeResidual(void * AMGhybrid_vdata,HYPRE_Int recompute_residual)420 hypre_AMGHybridSetRecomputeResidual( void *AMGhybrid_vdata,
421 HYPRE_Int recompute_residual )
422 {
423 hypre_AMGHybridData *AMGhybrid_data = (hypre_AMGHybridData *)AMGhybrid_vdata;
424 if (!AMGhybrid_data)
425 {
426 hypre_error_in_arg(1);
427 return hypre_error_flag;
428 }
429
430 (AMGhybrid_data -> recompute_residual) = recompute_residual;
431
432 return hypre_error_flag;
433 }
434
435 HYPRE_Int
hypre_AMGHybridGetRecomputeResidual(void * AMGhybrid_vdata,HYPRE_Int * recompute_residual)436 hypre_AMGHybridGetRecomputeResidual( void *AMGhybrid_vdata,
437 HYPRE_Int *recompute_residual )
438 {
439 hypre_AMGHybridData *AMGhybrid_data = (hypre_AMGHybridData *)AMGhybrid_vdata;
440 if (!AMGhybrid_data)
441 {
442 hypre_error_in_arg(1);
443 return hypre_error_flag;
444 }
445
446 *recompute_residual = (AMGhybrid_data -> recompute_residual);
447
448 return hypre_error_flag;
449 }
450
451 /*--------------------------------------------------------------------------
452 * hypre_AMGHybridSetRecomputeResidualP
453 *--------------------------------------------------------------------------*/
454
455 HYPRE_Int
hypre_AMGHybridSetRecomputeResidualP(void * AMGhybrid_vdata,HYPRE_Int recompute_residual_p)456 hypre_AMGHybridSetRecomputeResidualP( void *AMGhybrid_vdata,
457 HYPRE_Int recompute_residual_p )
458 {
459 hypre_AMGHybridData *AMGhybrid_data = (hypre_AMGHybridData *)AMGhybrid_vdata;
460 if (!AMGhybrid_data)
461 {
462 hypre_error_in_arg(1);
463 return hypre_error_flag;
464 }
465
466 (AMGhybrid_data -> recompute_residual_p) = recompute_residual_p;
467
468 return hypre_error_flag;
469 }
470
471 HYPRE_Int
hypre_AMGHybridGetRecomputeResidualP(void * AMGhybrid_vdata,HYPRE_Int * recompute_residual_p)472 hypre_AMGHybridGetRecomputeResidualP( void *AMGhybrid_vdata,
473 HYPRE_Int *recompute_residual_p )
474 {
475 hypre_AMGHybridData *AMGhybrid_data = (hypre_AMGHybridData *)AMGhybrid_vdata;
476 if (!AMGhybrid_data)
477 {
478 hypre_error_in_arg(1);
479 return hypre_error_flag;
480 }
481
482 *recompute_residual_p = (AMGhybrid_data -> recompute_residual_p);
483
484 return hypre_error_flag;
485 }
486
487 /*--------------------------------------------------------------------------
488 * hypre_AMGHybridSetKDim
489 *--------------------------------------------------------------------------*/
490
491 HYPRE_Int
hypre_AMGHybridSetKDim(void * AMGhybrid_vdata,HYPRE_Int k_dim)492 hypre_AMGHybridSetKDim( void *AMGhybrid_vdata,
493 HYPRE_Int k_dim )
494 {
495 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
496 if (!AMGhybrid_data)
497 {
498 hypre_error_in_arg(1);
499 return hypre_error_flag;
500 }
501 if (k_dim < 1)
502 {
503 hypre_error_in_arg(2);
504 return hypre_error_flag;
505 }
506
507 (AMGhybrid_data -> k_dim) = k_dim;
508
509 return hypre_error_flag;
510 }
511
512 /*--------------------------------------------------------------------------
513 * hypre_AMGHybridSetStopCrit
514 *--------------------------------------------------------------------------*/
515
516 HYPRE_Int
hypre_AMGHybridSetStopCrit(void * AMGhybrid_vdata,HYPRE_Int stop_crit)517 hypre_AMGHybridSetStopCrit( void *AMGhybrid_vdata,
518 HYPRE_Int stop_crit )
519 {
520 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
521 if (!AMGhybrid_data)
522 {
523 hypre_error_in_arg(1);
524 return hypre_error_flag;
525 }
526
527 (AMGhybrid_data -> stop_crit) = stop_crit;
528
529 return hypre_error_flag;
530 }
531
532 /*--------------------------------------------------------------------------
533 * hypre_AMGHybridSetTwoNorm
534 *--------------------------------------------------------------------------*/
535
536 HYPRE_Int
hypre_AMGHybridSetTwoNorm(void * AMGhybrid_vdata,HYPRE_Int two_norm)537 hypre_AMGHybridSetTwoNorm( void *AMGhybrid_vdata,
538 HYPRE_Int two_norm )
539 {
540 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
541 if (!AMGhybrid_data)
542 {
543 hypre_error_in_arg(1);
544 return hypre_error_flag;
545 }
546
547 (AMGhybrid_data -> two_norm) = two_norm;
548
549 return hypre_error_flag;
550 }
551
552 /*--------------------------------------------------------------------------
553 * hypre_AMGHybridSetRelChange
554 *--------------------------------------------------------------------------*/
555
556 HYPRE_Int
hypre_AMGHybridSetRelChange(void * AMGhybrid_vdata,HYPRE_Int rel_change)557 hypre_AMGHybridSetRelChange( void *AMGhybrid_vdata,
558 HYPRE_Int rel_change )
559 {
560 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
561 if (!AMGhybrid_data)
562 {
563 hypre_error_in_arg(1);
564 return hypre_error_flag;
565 }
566
567 (AMGhybrid_data -> rel_change) = rel_change;
568
569 return hypre_error_flag;
570 }
571
572 /*--------------------------------------------------------------------------
573 * hypre_AMGHybridSetPrecond
574 *--------------------------------------------------------------------------*/
575
576 HYPRE_Int
hypre_AMGHybridSetPrecond(void * pcg_vdata,HYPRE_Int (* pcg_precond_solve)(void *,void *,void *,void *),HYPRE_Int (* pcg_precond_setup)(void *,void *,void *,void *),void * pcg_precond)577 hypre_AMGHybridSetPrecond( void *pcg_vdata,
578 HYPRE_Int (*pcg_precond_solve)(void*,void*,void*,void*),
579 HYPRE_Int (*pcg_precond_setup)(void*,void*,void*,void*),
580 void *pcg_precond )
581 {
582 hypre_AMGHybridData *pcg_data =(hypre_AMGHybridData *) pcg_vdata;
583 if (!pcg_data)
584 {
585 hypre_error_in_arg(1);
586 return hypre_error_flag;
587 }
588
589 (pcg_data -> pcg_default) = 0;
590 (pcg_data -> pcg_precond_solve) = pcg_precond_solve;
591 (pcg_data -> pcg_precond_setup) = pcg_precond_setup;
592 (pcg_data -> pcg_precond) = pcg_precond;
593
594 return hypre_error_flag;
595 }
596
597 /*--------------------------------------------------------------------------
598 * hypre_AMGHybridSetLogging
599 *--------------------------------------------------------------------------*/
600
601 HYPRE_Int
hypre_AMGHybridSetLogging(void * AMGhybrid_vdata,HYPRE_Int logging)602 hypre_AMGHybridSetLogging( void *AMGhybrid_vdata,
603 HYPRE_Int logging )
604 {
605 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
606 if (!AMGhybrid_data)
607 {
608 hypre_error_in_arg(1);
609 return hypre_error_flag;
610 }
611
612 (AMGhybrid_data -> logging) = logging;
613
614 return hypre_error_flag;
615 }
616
617 /*--------------------------------------------------------------------------
618 * hypre_AMGHybridSetPrintLevel
619 *--------------------------------------------------------------------------*/
620
621 HYPRE_Int
hypre_AMGHybridSetPrintLevel(void * AMGhybrid_vdata,HYPRE_Int print_level)622 hypre_AMGHybridSetPrintLevel( void *AMGhybrid_vdata,
623 HYPRE_Int print_level )
624 {
625 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
626 if (!AMGhybrid_data)
627 {
628 hypre_error_in_arg(1);
629 return hypre_error_flag;
630 }
631
632 (AMGhybrid_data -> print_level) = print_level;
633
634 return hypre_error_flag;
635 }
636
637 /*--------------------------------------------------------------------------
638 * hypre_AMGHybridSetStrongThreshold
639 *--------------------------------------------------------------------------*/
640
641 HYPRE_Int
hypre_AMGHybridSetStrongThreshold(void * AMGhybrid_vdata,HYPRE_Real strong_threshold)642 hypre_AMGHybridSetStrongThreshold( void *AMGhybrid_vdata,
643 HYPRE_Real strong_threshold)
644 {
645 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
646 if (!AMGhybrid_data)
647 {
648 hypre_error_in_arg(1);
649 return hypre_error_flag;
650 }
651 if (strong_threshold < 0 || strong_threshold > 1)
652 {
653 hypre_error_in_arg(2);
654 return hypre_error_flag;
655 }
656
657 (AMGhybrid_data -> strong_threshold) = strong_threshold;
658
659 return hypre_error_flag;
660 }
661
662 /*--------------------------------------------------------------------------
663 * hypre_AMGHybridSetMaxRowSum
664 *--------------------------------------------------------------------------*/
665
666 HYPRE_Int
hypre_AMGHybridSetMaxRowSum(void * AMGhybrid_vdata,HYPRE_Real max_row_sum)667 hypre_AMGHybridSetMaxRowSum( void *AMGhybrid_vdata,
668 HYPRE_Real max_row_sum )
669 {
670 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
671 if (!AMGhybrid_data)
672 {
673 hypre_error_in_arg(1);
674 return hypre_error_flag;
675 }
676 if (max_row_sum < 0 || max_row_sum > 1)
677 {
678 hypre_error_in_arg(2);
679 return hypre_error_flag;
680 }
681
682 (AMGhybrid_data -> max_row_sum) = max_row_sum;
683
684 return hypre_error_flag;
685 }
686
687 /*--------------------------------------------------------------------------
688 * hypre_AMGHybridSetTruncFactor
689 *--------------------------------------------------------------------------*/
690
691 HYPRE_Int
hypre_AMGHybridSetTruncFactor(void * AMGhybrid_vdata,HYPRE_Real trunc_factor)692 hypre_AMGHybridSetTruncFactor( void *AMGhybrid_vdata,
693 HYPRE_Real trunc_factor )
694 {
695 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
696 if (!AMGhybrid_data)
697 {
698 hypre_error_in_arg(1);
699 return hypre_error_flag;
700 }
701 if (trunc_factor < 0 || trunc_factor > 1)
702 {
703 hypre_error_in_arg(2);
704 return hypre_error_flag;
705 }
706
707 (AMGhybrid_data -> trunc_factor) = trunc_factor;
708
709 return hypre_error_flag;
710 }
711
712 /*--------------------------------------------------------------------------
713 * hypre_AMGHybridSetPMaxElmts
714 *--------------------------------------------------------------------------*/
715
716 HYPRE_Int
hypre_AMGHybridSetPMaxElmts(void * AMGhybrid_vdata,HYPRE_Int P_max_elmts)717 hypre_AMGHybridSetPMaxElmts( void *AMGhybrid_vdata,
718 HYPRE_Int P_max_elmts )
719 {
720
721 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
722 if (!AMGhybrid_data)
723 {
724 hypre_error_in_arg(1);
725 return hypre_error_flag;
726 }
727
728 if (P_max_elmts < 0)
729 {
730 hypre_error_in_arg(2);
731 return hypre_error_flag;
732 }
733
734 (AMGhybrid_data -> pmax) = P_max_elmts;
735
736 return hypre_error_flag;
737 }
738
739 /*--------------------------------------------------------------------------
740 * hypre_AMGHybridSetMaxLevels
741 *--------------------------------------------------------------------------*/
742
743 HYPRE_Int
hypre_AMGHybridSetMaxLevels(void * AMGhybrid_vdata,HYPRE_Int max_levels)744 hypre_AMGHybridSetMaxLevels( void *AMGhybrid_vdata,
745 HYPRE_Int max_levels )
746 {
747 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
748 if (!AMGhybrid_data)
749 {
750 hypre_error_in_arg(1);
751 return hypre_error_flag;
752 }
753 if (max_levels < 1)
754 {
755 hypre_error_in_arg(2);
756 return hypre_error_flag;
757 }
758
759 (AMGhybrid_data -> max_levels) = max_levels;
760
761 return hypre_error_flag;
762 }
763
764 /*--------------------------------------------------------------------------
765 * hypre_AMGHybridSetMeasureType
766 *--------------------------------------------------------------------------*/
767
768 HYPRE_Int
hypre_AMGHybridSetMeasureType(void * AMGhybrid_vdata,HYPRE_Int measure_type)769 hypre_AMGHybridSetMeasureType( void *AMGhybrid_vdata,
770 HYPRE_Int measure_type )
771 {
772 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
773 if (!AMGhybrid_data)
774 {
775 hypre_error_in_arg(1);
776 return hypre_error_flag;
777 }
778
779 (AMGhybrid_data -> measure_type) = measure_type;
780
781 return hypre_error_flag;
782 }
783
784 /*--------------------------------------------------------------------------
785 * hypre_AMGHybridSetCoarsenType
786 *--------------------------------------------------------------------------*/
787
788 HYPRE_Int
hypre_AMGHybridSetCoarsenType(void * AMGhybrid_vdata,HYPRE_Int coarsen_type)789 hypre_AMGHybridSetCoarsenType( void *AMGhybrid_vdata,
790 HYPRE_Int coarsen_type )
791 {
792 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
793 if (!AMGhybrid_data)
794 {
795 hypre_error_in_arg(1);
796 return hypre_error_flag;
797 }
798
799 (AMGhybrid_data -> coarsen_type) = coarsen_type;
800
801 return hypre_error_flag;
802 }
803
804 /*--------------------------------------------------------------------------
805 * hypre_AMGHybridSetInterpType
806 *--------------------------------------------------------------------------*/
807
808 HYPRE_Int
hypre_AMGHybridSetInterpType(void * AMGhybrid_vdata,HYPRE_Int interp_type)809 hypre_AMGHybridSetInterpType( void *AMGhybrid_vdata,
810 HYPRE_Int interp_type )
811 {
812 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
813 if (!AMGhybrid_data)
814 {
815 hypre_error_in_arg(1);
816 return hypre_error_flag;
817 }
818 if (interp_type < 0)
819 {
820 hypre_error_in_arg(2);
821 return hypre_error_flag;
822 }
823
824 (AMGhybrid_data -> interp_type) = interp_type;
825
826 return hypre_error_flag;
827 }
828
829 /*--------------------------------------------------------------------------
830 * hypre_AMGHybridSetCycleType
831 *--------------------------------------------------------------------------*/
832
833 HYPRE_Int
hypre_AMGHybridSetCycleType(void * AMGhybrid_vdata,HYPRE_Int cycle_type)834 hypre_AMGHybridSetCycleType( void *AMGhybrid_vdata,
835 HYPRE_Int cycle_type )
836 {
837 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
838 if (!AMGhybrid_data)
839 {
840 hypre_error_in_arg(1);
841 return hypre_error_flag;
842 }
843 if (cycle_type < 1 || cycle_type > 2)
844 {
845 hypre_error_in_arg(2);
846 return hypre_error_flag;
847 }
848
849 (AMGhybrid_data -> cycle_type) = cycle_type;
850
851 return hypre_error_flag;
852 }
853
854 /*--------------------------------------------------------------------------
855 * hypre_AMGHybridSetNumSweeps
856 *--------------------------------------------------------------------------*/
857
858 HYPRE_Int
hypre_AMGHybridSetNumSweeps(void * AMGhybrid_vdata,HYPRE_Int num_sweeps)859 hypre_AMGHybridSetNumSweeps( void *AMGhybrid_vdata,
860 HYPRE_Int num_sweeps )
861 {
862 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
863 HYPRE_Int *num_grid_sweeps;
864 HYPRE_Int i;
865 if (!AMGhybrid_data)
866 {
867 hypre_error_in_arg(1);
868 return hypre_error_flag;
869 }
870 if (num_sweeps < 1)
871 {
872 hypre_error_in_arg(2);
873 return hypre_error_flag;
874 }
875
876 if ((AMGhybrid_data -> num_grid_sweeps) == NULL)
877 {
878 (AMGhybrid_data -> num_grid_sweeps) = hypre_CTAlloc(HYPRE_Int, 4, HYPRE_MEMORY_HOST);
879 }
880 num_grid_sweeps = (AMGhybrid_data -> num_grid_sweeps);
881 for (i=0; i < 3; i++)
882 {
883 num_grid_sweeps[i] = num_sweeps;
884 }
885 num_grid_sweeps[3] = 1;
886 return hypre_error_flag;
887 }
888
889 /*--------------------------------------------------------------------------
890 * hypre_AMGHybridSetCycleNumSweeps
891 *--------------------------------------------------------------------------*/
892
893 HYPRE_Int
hypre_AMGHybridSetCycleNumSweeps(void * AMGhybrid_vdata,HYPRE_Int num_sweeps,HYPRE_Int k)894 hypre_AMGHybridSetCycleNumSweeps( void *AMGhybrid_vdata,
895 HYPRE_Int num_sweeps,
896 HYPRE_Int k)
897 {
898 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
899 HYPRE_Int *num_grid_sweeps;
900 HYPRE_Int i;
901
902 if (!AMGhybrid_data)
903 {
904 hypre_error_in_arg(1);
905 return hypre_error_flag;
906 }
907 if (num_sweeps < 1)
908 {
909 hypre_error_in_arg(2);
910 return hypre_error_flag;
911 }
912 if (k < 1 || k > 3)
913 {
914 if (AMGhybrid_data -> print_level)
915 {
916 hypre_printf (" Warning! Invalid cycle! num_sweeps not set!\n");
917 }
918 hypre_error_in_arg(3);
919 return hypre_error_flag;
920 }
921
922 num_grid_sweeps = (AMGhybrid_data -> num_grid_sweeps);
923 if (num_grid_sweeps == NULL)
924 {
925 (AMGhybrid_data -> num_grid_sweeps) = hypre_CTAlloc(HYPRE_Int, 4, HYPRE_MEMORY_HOST);
926 num_grid_sweeps = (AMGhybrid_data -> num_grid_sweeps);
927 for (i=0; i < 4; i++)
928 {
929 num_grid_sweeps[i] = 1;
930 }
931 }
932 num_grid_sweeps[k] = num_sweeps;
933 return hypre_error_flag;
934 }
935
936 /*--------------------------------------------------------------------------
937 * hypre_AMGHybridSetRelaxType
938 *--------------------------------------------------------------------------*/
939
940 HYPRE_Int
hypre_AMGHybridSetRelaxType(void * AMGhybrid_vdata,HYPRE_Int relax_type)941 hypre_AMGHybridSetRelaxType( void *AMGhybrid_vdata,
942 HYPRE_Int relax_type )
943 {
944 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
945 HYPRE_Int *grid_relax_type;
946 HYPRE_Int i;
947 if (!AMGhybrid_data)
948 {
949 hypre_error_in_arg(1);
950 return hypre_error_flag;
951 }
952
953 if ((AMGhybrid_data -> grid_relax_type) == NULL )
954 {
955 (AMGhybrid_data -> grid_relax_type) = hypre_CTAlloc(HYPRE_Int, 4, HYPRE_MEMORY_HOST);
956 }
957 grid_relax_type = (AMGhybrid_data -> grid_relax_type);
958 for (i=0; i < 3; i++)
959 {
960 grid_relax_type[i] = relax_type;
961 }
962 grid_relax_type[3] = 9;
963
964 return hypre_error_flag;
965 }
966
967 /*--------------------------------------------------------------------------
968 * hypre_AMGHybridSetCycleRelaxType
969 *--------------------------------------------------------------------------*/
970
971 HYPRE_Int
hypre_AMGHybridSetCycleRelaxType(void * AMGhybrid_vdata,HYPRE_Int relax_type,HYPRE_Int k)972 hypre_AMGHybridSetCycleRelaxType( void *AMGhybrid_vdata,
973 HYPRE_Int relax_type,
974 HYPRE_Int k )
975 {
976 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
977 HYPRE_Int *grid_relax_type;
978 if (!AMGhybrid_data)
979 {
980 hypre_error_in_arg(1);
981 return hypre_error_flag;
982 }
983
984 if (k<1 || k > 3)
985 {
986 if (AMGhybrid_data -> print_level)
987 {
988 hypre_printf (" Warning! Invalid cycle! Relax type not set!\n");
989 }
990 hypre_error_in_arg(3);
991 return hypre_error_flag;
992 }
993
994 grid_relax_type = (AMGhybrid_data -> grid_relax_type);
995 if (grid_relax_type == NULL )
996 {
997 (AMGhybrid_data -> grid_relax_type) = hypre_CTAlloc(HYPRE_Int, 4, HYPRE_MEMORY_HOST);
998 grid_relax_type = (AMGhybrid_data -> grid_relax_type);
999
1000 grid_relax_type[1] = 13;
1001 grid_relax_type[2] = 14;
1002 grid_relax_type[3] = 9;
1003 }
1004 grid_relax_type[k] = relax_type;
1005
1006 return hypre_error_flag;
1007 }
1008
1009 /*--------------------------------------------------------------------------
1010 * hypre_AMGHybridSetRelaxOrder
1011 *--------------------------------------------------------------------------*/
1012
1013 HYPRE_Int
hypre_AMGHybridSetRelaxOrder(void * AMGhybrid_vdata,HYPRE_Int relax_order)1014 hypre_AMGHybridSetRelaxOrder( void *AMGhybrid_vdata,
1015 HYPRE_Int relax_order )
1016 {
1017 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1018 if (!AMGhybrid_data)
1019 {
1020 hypre_error_in_arg(1);
1021 return hypre_error_flag;
1022 }
1023
1024 (AMGhybrid_data -> relax_order) = relax_order;
1025
1026 return hypre_error_flag;
1027 }
1028
1029 /*--------------------------------------------------------------------------
1030 * hypre_AMGHybridSetKeepTranspose
1031 *--------------------------------------------------------------------------*/
1032
1033 HYPRE_Int
hypre_AMGHybridSetKeepTranspose(void * AMGhybrid_vdata,HYPRE_Int keepT)1034 hypre_AMGHybridSetKeepTranspose( void *AMGhybrid_vdata,
1035 HYPRE_Int keepT )
1036 {
1037 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1038 if (!AMGhybrid_data)
1039 {
1040 hypre_error_in_arg(1);
1041 return hypre_error_flag;
1042 }
1043
1044 (AMGhybrid_data -> keepT) = keepT;
1045
1046 return hypre_error_flag;
1047 }
1048
1049 /*--------------------------------------------------------------------------
1050 * hypre_AMGHybridSetMaxCoarseSize
1051 *--------------------------------------------------------------------------*/
1052
1053 HYPRE_Int
hypre_AMGHybridSetMaxCoarseSize(void * AMGhybrid_vdata,HYPRE_Int max_coarse_size)1054 hypre_AMGHybridSetMaxCoarseSize( void *AMGhybrid_vdata,
1055 HYPRE_Int max_coarse_size )
1056 {
1057 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1058 if (!AMGhybrid_data)
1059 {
1060 hypre_error_in_arg(1);
1061 return hypre_error_flag;
1062 }
1063 if (max_coarse_size < 1)
1064 {
1065 hypre_error_in_arg(2);
1066 return hypre_error_flag;
1067 }
1068
1069 (AMGhybrid_data -> max_coarse_size) = max_coarse_size;
1070
1071 return hypre_error_flag;
1072 }
1073
1074 /*--------------------------------------------------------------------------
1075 * hypre_AMGHybridSetMinCoarseSize
1076 *--------------------------------------------------------------------------*/
1077
1078 HYPRE_Int
hypre_AMGHybridSetMinCoarseSize(void * AMGhybrid_vdata,HYPRE_Int min_coarse_size)1079 hypre_AMGHybridSetMinCoarseSize( void *AMGhybrid_vdata,
1080 HYPRE_Int min_coarse_size )
1081 {
1082 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1083 if (!AMGhybrid_data)
1084 {
1085 hypre_error_in_arg(1);
1086 return hypre_error_flag;
1087 }
1088 if (min_coarse_size < 0)
1089 {
1090 hypre_error_in_arg(2);
1091 return hypre_error_flag;
1092 }
1093
1094 (AMGhybrid_data -> min_coarse_size) = min_coarse_size;
1095
1096 return hypre_error_flag;
1097 }
1098
1099 /*--------------------------------------------------------------------------
1100 * hypre_AMGHybridSetSeqThreshold
1101 *--------------------------------------------------------------------------*/
1102
1103 HYPRE_Int
hypre_AMGHybridSetSeqThreshold(void * AMGhybrid_vdata,HYPRE_Int seq_threshold)1104 hypre_AMGHybridSetSeqThreshold( void *AMGhybrid_vdata,
1105 HYPRE_Int seq_threshold )
1106 {
1107 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1108 if (!AMGhybrid_data)
1109 {
1110 hypre_error_in_arg(1);
1111 return hypre_error_flag;
1112 }
1113 if (seq_threshold < 0)
1114 {
1115 hypre_error_in_arg(2);
1116 return hypre_error_flag;
1117 }
1118
1119 (AMGhybrid_data -> seq_threshold) = seq_threshold;
1120
1121 return hypre_error_flag;
1122 }
1123
1124 /*--------------------------------------------------------------------------
1125 * hypre_AMGHybridSetNumGridSweeps
1126 *--------------------------------------------------------------------------*/
1127
1128 HYPRE_Int
hypre_AMGHybridSetNumGridSweeps(void * AMGhybrid_vdata,HYPRE_Int * num_grid_sweeps)1129 hypre_AMGHybridSetNumGridSweeps( void *AMGhybrid_vdata,
1130 HYPRE_Int *num_grid_sweeps )
1131 {
1132 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1133 if (!AMGhybrid_data)
1134 {
1135 hypre_error_in_arg(1);
1136 return hypre_error_flag;
1137 }
1138 if (!num_grid_sweeps)
1139 {
1140 hypre_error_in_arg(2);
1141 return hypre_error_flag;
1142 }
1143
1144 if ((AMGhybrid_data -> num_grid_sweeps) != NULL)
1145 {
1146 hypre_TFree((AMGhybrid_data -> num_grid_sweeps), HYPRE_MEMORY_HOST);
1147 }
1148 (AMGhybrid_data -> num_grid_sweeps) = num_grid_sweeps;
1149
1150 return hypre_error_flag;
1151 }
1152
1153 /*--------------------------------------------------------------------------
1154 * hypre_AMGHybridSetGridRelaxType
1155 *--------------------------------------------------------------------------*/
1156
1157 HYPRE_Int
hypre_AMGHybridSetGridRelaxType(void * AMGhybrid_vdata,HYPRE_Int * grid_relax_type)1158 hypre_AMGHybridSetGridRelaxType( void *AMGhybrid_vdata,
1159 HYPRE_Int *grid_relax_type )
1160 {
1161 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1162 if (!AMGhybrid_data)
1163 {
1164 hypre_error_in_arg(1);
1165 return hypre_error_flag;
1166 }
1167 if (!grid_relax_type)
1168 {
1169 hypre_error_in_arg(2);
1170 return hypre_error_flag;
1171 }
1172
1173 if ((AMGhybrid_data -> grid_relax_type) != NULL )
1174 {
1175 hypre_TFree((AMGhybrid_data -> grid_relax_type), HYPRE_MEMORY_HOST);
1176 }
1177 (AMGhybrid_data -> grid_relax_type) = grid_relax_type;
1178
1179 return hypre_error_flag;
1180 }
1181
1182 /*--------------------------------------------------------------------------
1183 * hypre_AMGHybridSetGridRelaxPoints
1184 *--------------------------------------------------------------------------*/
1185
1186 HYPRE_Int
hypre_AMGHybridSetGridRelaxPoints(void * AMGhybrid_vdata,HYPRE_Int ** grid_relax_points)1187 hypre_AMGHybridSetGridRelaxPoints( void *AMGhybrid_vdata,
1188 HYPRE_Int **grid_relax_points )
1189 {
1190 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1191 if (!AMGhybrid_data)
1192 {
1193 hypre_error_in_arg(1);
1194 return hypre_error_flag;
1195 }
1196 if (!grid_relax_points)
1197 {
1198 hypre_error_in_arg(2);
1199 return hypre_error_flag;
1200 }
1201
1202 if ((AMGhybrid_data -> grid_relax_points) != NULL )
1203 {
1204 hypre_TFree((AMGhybrid_data -> grid_relax_points), HYPRE_MEMORY_HOST);
1205 }
1206 (AMGhybrid_data -> grid_relax_points) = grid_relax_points;
1207
1208 return hypre_error_flag;
1209 }
1210
1211 /*--------------------------------------------------------------------------
1212 * hypre_AMGHybridSetRelaxWeight
1213 *--------------------------------------------------------------------------*/
1214
1215 HYPRE_Int
hypre_AMGHybridSetRelaxWeight(void * AMGhybrid_vdata,HYPRE_Real * relax_weight)1216 hypre_AMGHybridSetRelaxWeight( void *AMGhybrid_vdata,
1217 HYPRE_Real *relax_weight )
1218 {
1219 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1220 if (!AMGhybrid_data)
1221 {
1222 hypre_error_in_arg(1);
1223 return hypre_error_flag;
1224 }
1225 if (!relax_weight)
1226 {
1227 hypre_error_in_arg(2);
1228 return hypre_error_flag;
1229 }
1230
1231 if ((AMGhybrid_data -> relax_weight) != NULL )
1232 {
1233 hypre_TFree((AMGhybrid_data -> relax_weight), HYPRE_MEMORY_HOST);
1234 }
1235 (AMGhybrid_data -> relax_weight) = relax_weight;
1236
1237 return hypre_error_flag;
1238 }
1239
1240 /*--------------------------------------------------------------------------
1241 * hypre_AMGHybridSetOmega
1242 *--------------------------------------------------------------------------*/
1243
1244 HYPRE_Int
hypre_AMGHybridSetOmega(void * AMGhybrid_vdata,HYPRE_Real * omega)1245 hypre_AMGHybridSetOmega( void *AMGhybrid_vdata,
1246 HYPRE_Real *omega )
1247 {
1248 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1249 if (!AMGhybrid_data)
1250 {
1251 hypre_error_in_arg(1);
1252 return hypre_error_flag;
1253 }
1254 if (!omega)
1255 {
1256 hypre_error_in_arg(2);
1257 return hypre_error_flag;
1258 }
1259
1260 if ((AMGhybrid_data -> omega) != NULL )
1261 {
1262 hypre_TFree((AMGhybrid_data -> omega), HYPRE_MEMORY_HOST);
1263 }
1264 (AMGhybrid_data -> omega) = omega;
1265
1266 return hypre_error_flag;
1267 }
1268
1269 /*--------------------------------------------------------------------------
1270 * hypre_AMGHybridSetRelaxWt
1271 *--------------------------------------------------------------------------*/
1272
1273 HYPRE_Int
hypre_AMGHybridSetRelaxWt(void * AMGhybrid_vdata,HYPRE_Real relax_wt)1274 hypre_AMGHybridSetRelaxWt( void *AMGhybrid_vdata,
1275 HYPRE_Real relax_wt )
1276 {
1277 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1278 HYPRE_Int i , num_levels;
1279 HYPRE_Real *relax_wt_array;
1280 if (!AMGhybrid_data)
1281 {
1282 hypre_error_in_arg(1);
1283 return hypre_error_flag;
1284 }
1285
1286 num_levels = (AMGhybrid_data -> max_levels);
1287 relax_wt_array = (AMGhybrid_data -> relax_weight);
1288 if (relax_wt_array == NULL)
1289 {
1290 relax_wt_array = hypre_CTAlloc(HYPRE_Real, num_levels, HYPRE_MEMORY_HOST);
1291 (AMGhybrid_data -> relax_weight) = relax_wt_array;
1292 }
1293 for (i=0; i < num_levels; i++)
1294 {
1295 relax_wt_array[i] = relax_wt;
1296 }
1297
1298 return hypre_error_flag;
1299 }
1300
1301 /*--------------------------------------------------------------------------
1302 * hypre_AMGHybridSetLevelRelaxWt
1303 *--------------------------------------------------------------------------*/
1304
1305 HYPRE_Int
hypre_AMGHybridSetLevelRelaxWt(void * AMGhybrid_vdata,HYPRE_Real relax_wt,HYPRE_Int level)1306 hypre_AMGHybridSetLevelRelaxWt( void *AMGhybrid_vdata,
1307 HYPRE_Real relax_wt,
1308 HYPRE_Int level )
1309 {
1310 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1311 HYPRE_Int i , num_levels;
1312 HYPRE_Real *relax_wt_array;
1313 if (!AMGhybrid_data)
1314 {
1315 hypre_error_in_arg(1);
1316 return hypre_error_flag;
1317 }
1318
1319 num_levels = (AMGhybrid_data -> max_levels);
1320 if (level > num_levels-1)
1321 {
1322 if (AMGhybrid_data -> print_level)
1323 {
1324 hypre_printf (" Warning! Invalid level! Relax weight not set!\n");
1325 }
1326 hypre_error_in_arg(3);
1327 return hypre_error_flag;
1328 }
1329 relax_wt_array = (AMGhybrid_data -> relax_weight);
1330 if (relax_wt_array == NULL)
1331 {
1332 relax_wt_array = hypre_CTAlloc(HYPRE_Real, num_levels, HYPRE_MEMORY_HOST);
1333 for (i=0; i < num_levels; i++)
1334 {
1335 relax_wt_array[i] = 1.0;
1336 }
1337 (AMGhybrid_data -> relax_weight) = relax_wt_array;
1338 }
1339 relax_wt_array[level] = relax_wt;
1340
1341 return hypre_error_flag;
1342 }
1343
1344 /*--------------------------------------------------------------------------
1345 * hypre_AMGHybridSetOuterWt
1346 *--------------------------------------------------------------------------*/
1347
1348 HYPRE_Int
hypre_AMGHybridSetOuterWt(void * AMGhybrid_vdata,HYPRE_Real outer_wt)1349 hypre_AMGHybridSetOuterWt( void *AMGhybrid_vdata,
1350 HYPRE_Real outer_wt )
1351 {
1352 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1353 HYPRE_Int i , num_levels;
1354 HYPRE_Real *outer_wt_array;
1355 if (!AMGhybrid_data)
1356 {
1357 hypre_error_in_arg(1);
1358 return hypre_error_flag;
1359 }
1360
1361 num_levels = (AMGhybrid_data -> max_levels);
1362 outer_wt_array = (AMGhybrid_data -> omega);
1363 if (outer_wt_array == NULL)
1364 {
1365 outer_wt_array = hypre_CTAlloc(HYPRE_Real, num_levels, HYPRE_MEMORY_HOST);
1366 (AMGhybrid_data -> omega) = outer_wt_array;
1367 }
1368 for (i=0; i < num_levels; i++)
1369 {
1370 outer_wt_array[i] = outer_wt;
1371 }
1372
1373 return hypre_error_flag;
1374 }
1375
1376 /*--------------------------------------------------------------------------
1377 * hypre_AMGHybridSetLevelOuterWt
1378 *--------------------------------------------------------------------------*/
1379
1380 HYPRE_Int
hypre_AMGHybridSetLevelOuterWt(void * AMGhybrid_vdata,HYPRE_Real outer_wt,HYPRE_Int level)1381 hypre_AMGHybridSetLevelOuterWt( void *AMGhybrid_vdata,
1382 HYPRE_Real outer_wt,
1383 HYPRE_Int level )
1384 {
1385 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1386 HYPRE_Int i , num_levels;
1387 HYPRE_Real *outer_wt_array;
1388 if (!AMGhybrid_data)
1389 {
1390 hypre_error_in_arg(1);
1391 return hypre_error_flag;
1392 }
1393
1394 num_levels = (AMGhybrid_data -> max_levels);
1395 if (level > num_levels-1)
1396 {
1397 if (AMGhybrid_data -> print_level)
1398 {
1399 hypre_printf (" Warning! Invalid level! Outer weight not set!\n");
1400 }
1401 hypre_error_in_arg(3);
1402 return hypre_error_flag;
1403 }
1404 outer_wt_array = (AMGhybrid_data -> omega);
1405 if (outer_wt_array == NULL)
1406 {
1407 outer_wt_array = hypre_CTAlloc(HYPRE_Real, num_levels, HYPRE_MEMORY_HOST);
1408 for (i=0; i < num_levels; i++)
1409 {
1410 outer_wt_array[i] = 1.0;
1411 }
1412 (AMGhybrid_data -> omega) = outer_wt_array;
1413 }
1414 outer_wt_array[level] = outer_wt;
1415
1416 return hypre_error_flag;
1417 }
1418
1419 /*--------------------------------------------------------------------------
1420 * hypre_AMGHybridSetNumPaths
1421 *--------------------------------------------------------------------------*/
1422
1423 HYPRE_Int
hypre_AMGHybridSetNumPaths(void * AMGhybrid_vdata,HYPRE_Int num_paths)1424 hypre_AMGHybridSetNumPaths( void *AMGhybrid_vdata,
1425 HYPRE_Int num_paths )
1426 {
1427 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1428 if (!AMGhybrid_data)
1429 {
1430 hypre_error_in_arg(1);
1431 return hypre_error_flag;
1432 }
1433 if (num_paths < 1)
1434 {
1435 hypre_error_in_arg(2);
1436 return hypre_error_flag;
1437 }
1438
1439 (AMGhybrid_data -> num_paths) = num_paths;
1440
1441 return hypre_error_flag;
1442 }
1443
1444 /*--------------------------------------------------------------------------
1445 * hypre_AMGHybridSetDofFunc
1446 *--------------------------------------------------------------------------*/
1447
1448 HYPRE_Int
hypre_AMGHybridSetDofFunc(void * AMGhybrid_vdata,HYPRE_Int * dof_func)1449 hypre_AMGHybridSetDofFunc( void *AMGhybrid_vdata,
1450 HYPRE_Int *dof_func )
1451 {
1452 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1453 if (!AMGhybrid_data)
1454 {
1455 hypre_error_in_arg(1);
1456 return hypre_error_flag;
1457 }
1458 if (!dof_func)
1459 {
1460 hypre_error_in_arg(2);
1461 return hypre_error_flag;
1462 }
1463
1464 if ((AMGhybrid_data -> dof_func) != NULL )
1465 {
1466 hypre_TFree((AMGhybrid_data -> dof_func), HYPRE_MEMORY_HOST);
1467 }
1468 (AMGhybrid_data -> dof_func) = dof_func;
1469
1470 return hypre_error_flag;
1471 }
1472
1473 /*--------------------------------------------------------------------------
1474 * hypre_AMGHybridSetAggNumLevels
1475 *--------------------------------------------------------------------------*/
1476
1477 HYPRE_Int
hypre_AMGHybridSetAggNumLevels(void * AMGhybrid_vdata,HYPRE_Int agg_num_levels)1478 hypre_AMGHybridSetAggNumLevels( void *AMGhybrid_vdata,
1479 HYPRE_Int agg_num_levels )
1480 {
1481 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1482 if (!AMGhybrid_data)
1483 {
1484 hypre_error_in_arg(1);
1485 return hypre_error_flag;
1486 }
1487 if (agg_num_levels < 0)
1488 {
1489 hypre_error_in_arg(2);
1490 return hypre_error_flag;
1491 }
1492
1493 (AMGhybrid_data -> agg_num_levels) = agg_num_levels;
1494
1495 return hypre_error_flag;
1496 }
1497
1498 /*--------------------------------------------------------------------------
1499 * hypre_AMGHybridSetAggInterpType
1500 *--------------------------------------------------------------------------*/
1501
1502 HYPRE_Int
hypre_AMGHybridSetAggInterpType(void * AMGhybrid_vdata,HYPRE_Int agg_interp_type)1503 hypre_AMGHybridSetAggInterpType( void *AMGhybrid_vdata,
1504 HYPRE_Int agg_interp_type )
1505 {
1506 hypre_AMGHybridData *AMGhybrid_data = (hypre_AMGHybridData *) AMGhybrid_vdata;
1507 if (!AMGhybrid_data)
1508 {
1509 hypre_error_in_arg(1);
1510 return hypre_error_flag;
1511 }
1512
1513 (AMGhybrid_data -> agg_interp_type) = agg_interp_type;
1514
1515 return hypre_error_flag;
1516 }
1517
1518 /*--------------------------------------------------------------------------
1519 * hypre_AMGHybridSetNumFunctions
1520 *--------------------------------------------------------------------------*/
1521
1522 HYPRE_Int
hypre_AMGHybridSetNumFunctions(void * AMGhybrid_vdata,HYPRE_Int num_functions)1523 hypre_AMGHybridSetNumFunctions( void *AMGhybrid_vdata,
1524 HYPRE_Int num_functions )
1525 {
1526 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1527 if (!AMGhybrid_data)
1528 {
1529 hypre_error_in_arg(1);
1530 return hypre_error_flag;
1531 }
1532 if (num_functions < 1)
1533 {
1534 hypre_error_in_arg(2);
1535 return hypre_error_flag;
1536 }
1537
1538 (AMGhybrid_data -> num_functions) = num_functions;
1539
1540 return hypre_error_flag;
1541 }
1542
1543 /*--------------------------------------------------------------------------
1544 * hypre_AMGHybridSetNodal
1545 *--------------------------------------------------------------------------*/
1546
1547 HYPRE_Int
hypre_AMGHybridSetNodal(void * AMGhybrid_vdata,HYPRE_Int nodal)1548 hypre_AMGHybridSetNodal( void *AMGhybrid_vdata,
1549 HYPRE_Int nodal )
1550 {
1551 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1552 if (!AMGhybrid_data)
1553 {
1554 hypre_error_in_arg(1);
1555 return hypre_error_flag;
1556 }
1557
1558 (AMGhybrid_data -> nodal) = nodal;
1559
1560 return hypre_error_flag;
1561 }
1562
1563 /*--------------------------------------------------------------------------
1564 * hypre_AMGHybridGetNumIterations
1565 *--------------------------------------------------------------------------*/
1566 HYPRE_Int
hypre_AMGHybridGetSetupSolveTime(void * AMGhybrid_vdata,HYPRE_Real * time)1567 hypre_AMGHybridGetSetupSolveTime( void *AMGhybrid_vdata,
1568 HYPRE_Real *time )
1569 {
1570 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1571 if (!AMGhybrid_data)
1572 {
1573 hypre_error_in_arg(1);
1574 return hypre_error_flag;
1575 }
1576
1577 HYPRE_Real t[4];
1578 t[0] = AMGhybrid_data->setup_time1;
1579 t[1] = AMGhybrid_data->solve_time1;
1580 t[2] = AMGhybrid_data->setup_time2;
1581 t[3] = AMGhybrid_data->solve_time2;
1582
1583 MPI_Comm comm = AMGhybrid_data->comm;
1584
1585 hypre_MPI_Allreduce(t, time, 4, hypre_MPI_REAL, hypre_MPI_MAX, comm);
1586
1587 return hypre_error_flag;
1588 }
1589
1590 HYPRE_Int
hypre_AMGHybridGetNumIterations(void * AMGhybrid_vdata,HYPRE_Int * num_its)1591 hypre_AMGHybridGetNumIterations( void *AMGhybrid_vdata,
1592 HYPRE_Int *num_its )
1593 {
1594 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1595 if (!AMGhybrid_data)
1596 {
1597 hypre_error_in_arg(1);
1598 return hypre_error_flag;
1599 }
1600
1601 *num_its = (AMGhybrid_data -> dscg_num_its) + (AMGhybrid_data -> pcg_num_its);
1602
1603 return hypre_error_flag;
1604 }
1605
1606 /*--------------------------------------------------------------------------
1607 * hypre_AMGHybridGetDSCGNumIterations
1608 *--------------------------------------------------------------------------*/
1609
1610 HYPRE_Int
hypre_AMGHybridGetDSCGNumIterations(void * AMGhybrid_vdata,HYPRE_Int * dscg_num_its)1611 hypre_AMGHybridGetDSCGNumIterations( void *AMGhybrid_vdata,
1612 HYPRE_Int *dscg_num_its )
1613 {
1614 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1615 if (!AMGhybrid_data)
1616 {
1617 hypre_error_in_arg(1);
1618 return hypre_error_flag;
1619 }
1620
1621 *dscg_num_its = (AMGhybrid_data -> dscg_num_its);
1622
1623 return hypre_error_flag;
1624 }
1625
1626 /*--------------------------------------------------------------------------
1627 * hypre_AMGHybridGetPCGNumIterations
1628 *--------------------------------------------------------------------------*/
1629
1630 HYPRE_Int
hypre_AMGHybridGetPCGNumIterations(void * AMGhybrid_vdata,HYPRE_Int * pcg_num_its)1631 hypre_AMGHybridGetPCGNumIterations( void *AMGhybrid_vdata,
1632 HYPRE_Int *pcg_num_its )
1633 {
1634 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1635 if (!AMGhybrid_data)
1636 {
1637 hypre_error_in_arg(1);
1638 return hypre_error_flag;
1639 }
1640
1641 *pcg_num_its = (AMGhybrid_data -> pcg_num_its);
1642
1643 return hypre_error_flag;
1644 }
1645
1646 /*--------------------------------------------------------------------------
1647 * hypre_AMGHybridGetFinalRelativeResidualNorm
1648 *--------------------------------------------------------------------------*/
1649
1650 HYPRE_Int
hypre_AMGHybridGetFinalRelativeResidualNorm(void * AMGhybrid_vdata,HYPRE_Real * final_rel_res_norm)1651 hypre_AMGHybridGetFinalRelativeResidualNorm( void *AMGhybrid_vdata,
1652 HYPRE_Real *final_rel_res_norm )
1653 {
1654 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1655 if (!AMGhybrid_data)
1656 {
1657 hypre_error_in_arg(1);
1658 return hypre_error_flag;
1659 }
1660
1661 *final_rel_res_norm = (AMGhybrid_data -> final_rel_res_norm);
1662
1663 return hypre_error_flag;
1664 }
1665
1666 /*--------------------------------------------------------------------------
1667 * hypre_AMGHybridSetup
1668 *--------------------------------------------------------------------------*/
1669
1670 HYPRE_Int
hypre_AMGHybridSetup(void * AMGhybrid_vdata,hypre_ParCSRMatrix * A,hypre_ParVector * b,hypre_ParVector * x)1671 hypre_AMGHybridSetup( void *AMGhybrid_vdata,
1672 hypre_ParCSRMatrix *A,
1673 hypre_ParVector *b,
1674 hypre_ParVector *x )
1675 {
1676 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1677 if (!AMGhybrid_data)
1678 {
1679 hypre_error_in_arg(1);
1680 return hypre_error_flag;
1681 }
1682
1683 return hypre_error_flag;
1684 }
1685
1686 /*--------------------------------------------------------------------------
1687 * hypre_AMGHybridSolve
1688 *--------------------------------------------------------------------------
1689 *
1690 * This solver is designed to solve Ax=b using a AMGhybrid algorithm. First
1691 * the solver uses diagonally scaled conjugate gradients. If sufficient
1692 * progress is not made, the algorithm switches to preconditioned
1693 * conjugate gradients with user-specified preconditioner.
1694 *
1695 *--------------------------------------------------------------------------*/
1696
1697 HYPRE_Int
hypre_AMGHybridSolve(void * AMGhybrid_vdata,hypre_ParCSRMatrix * A,hypre_ParVector * b,hypre_ParVector * x)1698 hypre_AMGHybridSolve( void *AMGhybrid_vdata,
1699 hypre_ParCSRMatrix *A,
1700 hypre_ParVector *b,
1701 hypre_ParVector *x )
1702 {
1703 hypre_AMGHybridData *AMGhybrid_data =(hypre_AMGHybridData *) AMGhybrid_vdata;
1704
1705 HYPRE_Real tol;
1706 HYPRE_Real a_tol;
1707 HYPRE_Real cf_tol;
1708 HYPRE_Int dscg_max_its;
1709 HYPRE_Int pcg_max_its;
1710 HYPRE_Int two_norm;
1711 HYPRE_Int stop_crit;
1712 HYPRE_Int rel_change;
1713 HYPRE_Int recompute_residual;
1714 HYPRE_Int recompute_residual_p;
1715 HYPRE_Int logging;
1716 HYPRE_Int print_level;
1717 HYPRE_Int setup_type;
1718 HYPRE_Int solver_type;
1719 HYPRE_Int k_dim;
1720 /* BoomerAMG info */
1721 HYPRE_Real strong_threshold;
1722 HYPRE_Real max_row_sum;
1723 HYPRE_Real trunc_factor;
1724 HYPRE_Int pmax;
1725 HYPRE_Int max_levels;
1726 HYPRE_Int measure_type;
1727 HYPRE_Int coarsen_type;
1728 HYPRE_Int interp_type;
1729 HYPRE_Int cycle_type;
1730 HYPRE_Int num_paths;
1731 HYPRE_Int agg_num_levels;
1732 HYPRE_Int agg_interp_type;
1733 HYPRE_Int num_functions;
1734 HYPRE_Int nodal;
1735 HYPRE_Int relax_order;
1736 HYPRE_Int keepT;
1737 HYPRE_Int *num_grid_sweeps;
1738 HYPRE_Int *grid_relax_type;
1739 HYPRE_Int **grid_relax_points;
1740 HYPRE_Real *relax_weight;
1741 HYPRE_Real *omega;
1742 HYPRE_Int *dof_func;
1743
1744 HYPRE_Int *boom_ngs;
1745 HYPRE_Int *boom_grt;
1746 HYPRE_Int *boom_dof_func;
1747 HYPRE_Int **boom_grp;
1748 HYPRE_Real *boom_rlxw;
1749 HYPRE_Real *boom_omega;
1750
1751 HYPRE_Int pcg_default;
1752 HYPRE_Int (*pcg_precond_solve)(void*,void*,void*,void*);
1753 HYPRE_Int (*pcg_precond_setup)(void*,void*,void*,void*);
1754 void *pcg_precond;
1755
1756 void *pcg_solver;
1757 hypre_PCGFunctions *pcg_functions;
1758 hypre_GMRESFunctions *gmres_functions;
1759 hypre_BiCGSTABFunctions *bicgstab_functions;
1760
1761 HYPRE_Int dscg_num_its=0;
1762 HYPRE_Int pcg_num_its=0;
1763 HYPRE_Int converged=0;
1764 HYPRE_Int num_variables = hypre_VectorSize(hypre_ParVectorLocalVector(b));
1765 HYPRE_Real res_norm;
1766
1767 HYPRE_Int i, j;
1768 HYPRE_Int sol_print_level; /* print_level for solver */
1769 HYPRE_Int pre_print_level; /* print_level for preconditioner */
1770 HYPRE_Int max_coarse_size, seq_threshold;
1771 HYPRE_Int min_coarse_size;
1772 HYPRE_Int nongalerk_num_tol;
1773 HYPRE_Real *nongalerkin_tol;
1774
1775 HYPRE_Real tt1, tt2;
1776
1777 if (!AMGhybrid_data)
1778 {
1779 hypre_error_in_arg(1);
1780 return hypre_error_flag;
1781 }
1782
1783 AMGhybrid_data->setup_time1 = 0.0;
1784 AMGhybrid_data->setup_time2 = 0.0;
1785 AMGhybrid_data->solve_time1 = 0.0;
1786 AMGhybrid_data->solve_time2 = 0.0;
1787 MPI_Comm comm = hypre_ParCSRMatrixComm(A);
1788 (AMGhybrid_data -> comm) = comm;
1789 /*-----------------------------------------------------------------------
1790 * Setup diagonal scaled solver
1791 *-----------------------------------------------------------------------*/
1792 tol = (AMGhybrid_data -> tol);
1793 a_tol = (AMGhybrid_data -> a_tol);
1794 cf_tol = (AMGhybrid_data -> cf_tol);
1795 dscg_max_its = (AMGhybrid_data -> dscg_max_its);
1796 pcg_max_its = (AMGhybrid_data -> pcg_max_its);
1797 two_norm = (AMGhybrid_data -> two_norm);
1798 stop_crit = (AMGhybrid_data -> stop_crit);
1799 rel_change = (AMGhybrid_data -> rel_change);
1800 recompute_residual = (AMGhybrid_data -> recompute_residual);
1801 recompute_residual_p = (AMGhybrid_data -> recompute_residual_p);
1802 logging = (AMGhybrid_data -> logging);
1803 print_level = (AMGhybrid_data -> print_level);
1804 setup_type = (AMGhybrid_data -> setup_type);
1805 solver_type = (AMGhybrid_data -> solver_type);
1806 k_dim = (AMGhybrid_data -> k_dim);
1807 strong_threshold = (AMGhybrid_data -> strong_threshold);
1808 max_row_sum = (AMGhybrid_data -> max_row_sum);
1809 trunc_factor = (AMGhybrid_data -> trunc_factor);
1810 pmax = (AMGhybrid_data -> pmax);
1811 max_levels = (AMGhybrid_data -> max_levels);
1812 measure_type = (AMGhybrid_data -> measure_type);
1813 coarsen_type = (AMGhybrid_data -> coarsen_type);
1814 interp_type = (AMGhybrid_data -> interp_type);
1815 cycle_type = (AMGhybrid_data -> cycle_type);
1816 num_paths = (AMGhybrid_data -> num_paths);
1817 agg_num_levels = (AMGhybrid_data -> agg_num_levels);
1818 agg_interp_type = (AMGhybrid_data -> agg_interp_type);
1819 num_functions = (AMGhybrid_data -> num_functions);
1820 nodal = (AMGhybrid_data -> nodal);
1821 num_grid_sweeps = (AMGhybrid_data -> num_grid_sweeps);
1822 grid_relax_type = (AMGhybrid_data -> grid_relax_type);
1823 grid_relax_points = (AMGhybrid_data -> grid_relax_points);
1824 relax_weight = (AMGhybrid_data -> relax_weight);
1825 relax_order = (AMGhybrid_data -> relax_order);
1826 keepT = (AMGhybrid_data -> keepT);
1827 omega = (AMGhybrid_data -> omega);
1828 max_coarse_size = (AMGhybrid_data -> max_coarse_size);
1829 min_coarse_size = (AMGhybrid_data -> min_coarse_size);
1830 seq_threshold = (AMGhybrid_data -> seq_threshold);
1831 dof_func = (AMGhybrid_data -> dof_func);
1832 pcg_default = (AMGhybrid_data -> pcg_default);
1833 nongalerk_num_tol = (AMGhybrid_data -> nongalerk_num_tol);
1834 nongalerkin_tol = (AMGhybrid_data -> nongalerkin_tol);
1835 if (!b)
1836 {
1837 hypre_error_in_arg(3);
1838 return hypre_error_flag;
1839 }
1840 num_variables = hypre_VectorSize(hypre_ParVectorLocalVector(b));
1841 if (!A)
1842 {
1843 hypre_error_in_arg(2);
1844 return hypre_error_flag;
1845 }
1846 if (!x)
1847 {
1848 hypre_error_in_arg(4);
1849 return hypre_error_flag;
1850 }
1851
1852 /* print_level definitions: xy, sol_print_level = y, pre_print_level = x */
1853 pre_print_level = print_level/10;
1854 sol_print_level = print_level - pre_print_level*10;
1855
1856 pcg_solver = (AMGhybrid_data -> pcg_solver);
1857 pcg_precond = (AMGhybrid_data -> pcg_precond);
1858 (AMGhybrid_data -> dscg_num_its) = 0;
1859 (AMGhybrid_data -> pcg_num_its) = 0;
1860
1861 if (setup_type || pcg_precond == NULL)
1862 {
1863 if (pcg_precond)
1864 {
1865 hypre_BoomerAMGDestroy(pcg_precond);
1866 pcg_precond = NULL;
1867 (AMGhybrid_data -> pcg_precond) = NULL;
1868 }
1869 if (solver_type == 1)
1870 {
1871 tt1 = hypre_MPI_Wtime();
1872
1873 if (pcg_solver == NULL)
1874 {
1875 pcg_functions =
1876 hypre_PCGFunctionsCreate(
1877 hypre_ParKrylovCAlloc, hypre_ParKrylovFree,
1878 hypre_ParKrylovCommInfo,
1879 hypre_ParKrylovCreateVector,
1880 hypre_ParKrylovDestroyVector, hypre_ParKrylovMatvecCreate,
1881 hypre_ParKrylovMatvec,
1882 hypre_ParKrylovMatvecDestroy,
1883 hypre_ParKrylovInnerProd, hypre_ParKrylovCopyVector,
1884 hypre_ParKrylovClearVector,
1885 hypre_ParKrylovScaleVector, hypre_ParKrylovAxpy,
1886 hypre_ParKrylovIdentitySetup, hypre_ParKrylovIdentity );
1887 pcg_solver = hypre_PCGCreate( pcg_functions );
1888
1889 hypre_PCGSetTol(pcg_solver, tol);
1890 hypre_PCGSetAbsoluteTol(pcg_solver, a_tol);
1891 hypre_PCGSetTwoNorm(pcg_solver, two_norm);
1892 hypre_PCGSetStopCrit(pcg_solver, stop_crit);
1893 hypre_PCGSetRelChange(pcg_solver, rel_change);
1894 hypre_PCGSetRecomputeResidual(pcg_solver, recompute_residual);
1895 hypre_PCGSetRecomputeResidualP(pcg_solver, recompute_residual_p);
1896 hypre_PCGSetLogging(pcg_solver, logging);
1897 hypre_PCGSetPrintLevel(pcg_solver, sol_print_level);
1898 hypre_PCGSetHybrid(pcg_solver,-1);
1899
1900 pcg_precond = NULL;
1901 }
1902
1903 hypre_PCGSetMaxIter(pcg_solver, dscg_max_its);
1904 hypre_PCGSetConvergenceFactorTol(pcg_solver, cf_tol);
1905 hypre_PCGSetPrecond((void*) pcg_solver,
1906 (HYPRE_Int (*)(void*, void*, void*, void*)) HYPRE_ParCSRDiagScale,
1907 (HYPRE_Int (*)(void*, void*, void*, void*)) HYPRE_ParCSRDiagScaleSetup,
1908 (void*) pcg_precond);
1909
1910 hypre_PCGSetup(pcg_solver, (void*) A, (void*) b, (void*) x);
1911 (AMGhybrid_data -> pcg_solver) = pcg_solver;
1912
1913 tt2 = hypre_MPI_Wtime();
1914 AMGhybrid_data->setup_time1 = tt2 - tt1;
1915
1916 /*---------------------------------------------------------------------
1917 * Solve with DSCG.
1918 *---------------------------------------------------------------------*/
1919 tt1 = tt2;
1920
1921 hypre_PCGSolve(pcg_solver, (void*) A, (void*) b, (void*) x);
1922
1923 /*---------------------------------------------------------------------
1924 * Get information for DSCG.
1925 *---------------------------------------------------------------------*/
1926 hypre_PCGGetNumIterations(pcg_solver, &dscg_num_its);
1927 (AMGhybrid_data -> dscg_num_its) = dscg_num_its;
1928 hypre_PCGGetFinalRelativeResidualNorm(pcg_solver, &res_norm);
1929
1930 hypre_PCGGetConverged(pcg_solver, &converged);
1931
1932 tt2 = hypre_MPI_Wtime();
1933 AMGhybrid_data->solve_time1 = tt2 - tt1;
1934 }
1935 else if (solver_type == 2)
1936 {
1937 tt1 = hypre_MPI_Wtime();
1938
1939 if (pcg_solver == NULL)
1940 {
1941 gmres_functions =
1942 hypre_GMRESFunctionsCreate(
1943 hypre_ParKrylovCAlloc, hypre_ParKrylovFree,
1944 hypre_ParKrylovCommInfo,
1945 hypre_ParKrylovCreateVector,
1946 hypre_ParKrylovCreateVectorArray,
1947 hypre_ParKrylovDestroyVector, hypre_ParKrylovMatvecCreate,
1948 hypre_ParKrylovMatvec,
1949 hypre_ParKrylovMatvecDestroy,
1950 hypre_ParKrylovInnerProd, hypre_ParKrylovCopyVector,
1951 hypre_ParKrylovClearVector,
1952 hypre_ParKrylovScaleVector, hypre_ParKrylovAxpy,
1953 hypre_ParKrylovIdentitySetup, hypre_ParKrylovIdentity );
1954 pcg_solver = hypre_GMRESCreate( gmres_functions );
1955
1956 hypre_GMRESSetTol(pcg_solver, tol);
1957 hypre_GMRESSetAbsoluteTol(pcg_solver, a_tol);
1958 hypre_GMRESSetKDim(pcg_solver, k_dim);
1959 hypre_GMRESSetStopCrit(pcg_solver, stop_crit);
1960 hypre_GMRESSetRelChange(pcg_solver, rel_change);
1961 hypre_GMRESSetLogging(pcg_solver, logging);
1962 hypre_GMRESSetPrintLevel(pcg_solver, sol_print_level);
1963 hypre_GMRESSetHybrid(pcg_solver,-1);
1964
1965 pcg_precond = NULL;
1966 }
1967
1968 hypre_GMRESSetMaxIter(pcg_solver, dscg_max_its);
1969 hypre_GMRESSetConvergenceFactorTol(pcg_solver, cf_tol);
1970 hypre_GMRESSetPrecond((void*) pcg_solver,
1971 (HYPRE_Int (*)(void*, void*, void*, void*)) HYPRE_ParCSRDiagScale,
1972 (HYPRE_Int (*)(void*, void*, void*, void*)) HYPRE_ParCSRDiagScaleSetup,
1973 (void*) pcg_precond);
1974
1975 hypre_GMRESSetup(pcg_solver, (void*) A, (void*) b, (void*) x);
1976 (AMGhybrid_data -> pcg_solver) = pcg_solver;
1977
1978 tt2 = hypre_MPI_Wtime();
1979 AMGhybrid_data->setup_time1 = tt2 - tt1;
1980
1981 /*---------------------------------------------------------------------
1982 * Solve with diagonal scaled GMRES
1983 *---------------------------------------------------------------------*/
1984 tt1 = tt2;
1985
1986 hypre_GMRESSolve(pcg_solver, (void*) A, (void*) b, (void*) x);
1987
1988 /*---------------------------------------------------------------------
1989 * Get information for GMRES
1990 *---------------------------------------------------------------------*/
1991 hypre_GMRESGetNumIterations(pcg_solver, &dscg_num_its);
1992 (AMGhybrid_data -> dscg_num_its) = dscg_num_its;
1993 hypre_GMRESGetFinalRelativeResidualNorm(pcg_solver, &res_norm);
1994
1995 hypre_GMRESGetConverged(pcg_solver, &converged);
1996
1997 tt2 = hypre_MPI_Wtime();
1998 AMGhybrid_data->solve_time1 = tt2 - tt1;
1999 }
2000 else if (solver_type == 3)
2001 {
2002 tt1 = hypre_MPI_Wtime();
2003
2004 if (pcg_solver == NULL)
2005 {
2006 bicgstab_functions =
2007 hypre_BiCGSTABFunctionsCreate(
2008 hypre_ParKrylovCreateVector,
2009 hypre_ParKrylovDestroyVector, hypre_ParKrylovMatvecCreate,
2010 hypre_ParKrylovMatvec,
2011 hypre_ParKrylovMatvecDestroy,
2012 hypre_ParKrylovInnerProd, hypre_ParKrylovCopyVector,
2013 hypre_ParKrylovClearVector,
2014 hypre_ParKrylovScaleVector, hypre_ParKrylovAxpy,
2015 hypre_ParKrylovCommInfo,
2016 hypre_ParKrylovIdentitySetup, hypre_ParKrylovIdentity );
2017 pcg_solver = hypre_BiCGSTABCreate( bicgstab_functions );
2018
2019 hypre_BiCGSTABSetTol(pcg_solver, tol);
2020 hypre_BiCGSTABSetAbsoluteTol(pcg_solver, a_tol);
2021 hypre_BiCGSTABSetStopCrit(pcg_solver, stop_crit);
2022 hypre_BiCGSTABSetLogging(pcg_solver, logging);
2023 hypre_BiCGSTABSetPrintLevel(pcg_solver, sol_print_level);
2024 hypre_BiCGSTABSetHybrid(pcg_solver,-1);
2025
2026 pcg_precond = NULL;
2027 }
2028
2029 hypre_BiCGSTABSetMaxIter(pcg_solver, dscg_max_its);
2030 hypre_BiCGSTABSetConvergenceFactorTol(pcg_solver, cf_tol);
2031 hypre_BiCGSTABSetPrecond((void*) pcg_solver,
2032 (HYPRE_Int (*)(void*, void*, void*, void*)) HYPRE_ParCSRDiagScale,
2033 (HYPRE_Int (*)(void*, void*, void*, void*)) HYPRE_ParCSRDiagScaleSetup,
2034 (void*) pcg_precond);
2035
2036 hypre_BiCGSTABSetup(pcg_solver, (void*) A, (void*) b, (void*) x);
2037 (AMGhybrid_data -> pcg_solver) = pcg_solver;
2038
2039 tt2 = hypre_MPI_Wtime();
2040 AMGhybrid_data->setup_time1 = tt2 - tt1;
2041
2042 /*---------------------------------------------------------------------
2043 * Solve with diagonal scaled BiCGSTAB
2044 *---------------------------------------------------------------------*/
2045 tt1 = tt2;
2046
2047 hypre_BiCGSTABSolve(pcg_solver, (void*) A, (void*) b, (void*) x);
2048
2049 /*---------------------------------------------------------------------
2050 * Get information for BiCGSTAB
2051 *---------------------------------------------------------------------*/
2052 hypre_BiCGSTABGetNumIterations(pcg_solver, &dscg_num_its);
2053 (AMGhybrid_data -> dscg_num_its) = dscg_num_its;
2054 hypre_BiCGSTABGetFinalRelativeResidualNorm(pcg_solver, &res_norm);
2055
2056 hypre_BiCGSTABGetConverged(pcg_solver, &converged);
2057
2058 tt2 = hypre_MPI_Wtime();
2059 AMGhybrid_data->solve_time1 = tt2 - tt1;
2060 }
2061 }
2062
2063 /*---------------------------------------------------------------------
2064 * If converged, done...
2065 *---------------------------------------------------------------------*/
2066 if (converged)
2067 {
2068 if (logging)
2069 {
2070 (AMGhybrid_data -> final_rel_res_norm) = res_norm;
2071 }
2072 }
2073 /*-----------------------------------------------------------------------
2074 * ... otherwise, use AMG+solver
2075 *-----------------------------------------------------------------------*/
2076 else
2077 {
2078 tt1 = hypre_MPI_Wtime();
2079
2080 /*--------------------------------------------------------------------
2081 * Free up previous PCG solver structure and set up a new one.
2082 *--------------------------------------------------------------------*/
2083 if (solver_type == 1)
2084 {
2085 hypre_PCGSetMaxIter(pcg_solver, pcg_max_its);
2086 hypre_PCGSetConvergenceFactorTol(pcg_solver, 0.0);
2087 hypre_PCGSetHybrid(pcg_solver, 0);
2088 }
2089 else if (solver_type == 2)
2090 {
2091 hypre_GMRESSetMaxIter(pcg_solver, pcg_max_its);
2092 hypre_GMRESSetConvergenceFactorTol(pcg_solver, 0.0);
2093 hypre_GMRESSetHybrid(pcg_solver, 0);
2094 }
2095 else if (solver_type == 3)
2096 {
2097 hypre_BiCGSTABSetMaxIter(pcg_solver, pcg_max_its);
2098 hypre_BiCGSTABSetConvergenceFactorTol(pcg_solver, 0.0);
2099 hypre_BiCGSTABSetHybrid(pcg_solver, 0);
2100 }
2101
2102 /* Setup preconditioner */
2103 if (setup_type && pcg_default)
2104 {
2105 pcg_precond = hypre_BoomerAMGCreate();
2106 hypre_BoomerAMGSetMaxIter(pcg_precond, 1);
2107 hypre_BoomerAMGSetTol(pcg_precond, 0.0);
2108 hypre_BoomerAMGSetCoarsenType(pcg_precond, coarsen_type);
2109 hypre_BoomerAMGSetInterpType(pcg_precond, interp_type);
2110 hypre_BoomerAMGSetSetupType(pcg_precond, setup_type);
2111 hypre_BoomerAMGSetMeasureType(pcg_precond, measure_type);
2112 hypre_BoomerAMGSetStrongThreshold(pcg_precond, strong_threshold);
2113 hypre_BoomerAMGSetTruncFactor(pcg_precond, trunc_factor);
2114 hypre_BoomerAMGSetPMaxElmts(pcg_precond, pmax);
2115 hypre_BoomerAMGSetCycleType(pcg_precond, cycle_type);
2116 hypre_BoomerAMGSetPrintLevel(pcg_precond, pre_print_level);
2117 hypre_BoomerAMGSetMaxLevels(pcg_precond, max_levels);
2118 hypre_BoomerAMGSetMaxRowSum(pcg_precond, max_row_sum);
2119 hypre_BoomerAMGSetMaxCoarseSize(pcg_precond, max_coarse_size);
2120 hypre_BoomerAMGSetMinCoarseSize(pcg_precond, min_coarse_size);
2121 hypre_BoomerAMGSetSeqThreshold(pcg_precond, seq_threshold);
2122 hypre_BoomerAMGSetAggNumLevels(pcg_precond, agg_num_levels);
2123 hypre_BoomerAMGSetAggInterpType(pcg_precond, agg_interp_type);
2124 hypre_BoomerAMGSetNumPaths(pcg_precond, num_paths);
2125 hypre_BoomerAMGSetNumFunctions(pcg_precond, num_functions);
2126 hypre_BoomerAMGSetNodal(pcg_precond, nodal);
2127 hypre_BoomerAMGSetRelaxOrder(pcg_precond, relax_order);
2128 hypre_BoomerAMGSetKeepTranspose(pcg_precond, keepT);
2129 hypre_BoomerAMGSetNonGalerkTol(pcg_precond, nongalerk_num_tol, nongalerkin_tol);
2130 if (grid_relax_type)
2131 {
2132 boom_grt = hypre_CTAlloc(HYPRE_Int, 4, HYPRE_MEMORY_HOST);
2133 for (i=0; i < 4; i++)
2134 {
2135 boom_grt[i] = grid_relax_type[i];
2136 }
2137 hypre_BoomerAMGSetGridRelaxType(pcg_precond, boom_grt);
2138 }
2139 else
2140 {
2141 boom_grt = hypre_CTAlloc(HYPRE_Int, 4, HYPRE_MEMORY_HOST);
2142 boom_grt[0] = 3;
2143 boom_grt[1] = 13;
2144 boom_grt[2] = 14;
2145 boom_grt[3] = 9;
2146 hypre_BoomerAMGSetGridRelaxType(pcg_precond, boom_grt);
2147 }
2148
2149 hypre_ParAMGDataUserCoarseRelaxType((hypre_ParAMGData *) pcg_precond) = boom_grt[3];
2150 hypre_ParAMGDataUserRelaxType((hypre_ParAMGData *) pcg_precond) = boom_grt[0];
2151
2152 if (relax_weight)
2153 {
2154 boom_rlxw = hypre_CTAlloc(HYPRE_Real, max_levels, HYPRE_MEMORY_HOST);
2155 for (i=0; i < max_levels; i++)
2156 boom_rlxw[i] = relax_weight[i];
2157 hypre_BoomerAMGSetRelaxWeight(pcg_precond, boom_rlxw);
2158 }
2159 if (omega)
2160 {
2161 boom_omega = hypre_CTAlloc(HYPRE_Real, max_levels, HYPRE_MEMORY_HOST);
2162 for (i=0; i < max_levels; i++)
2163 boom_omega[i] = omega[i];
2164 hypre_BoomerAMGSetOmega(pcg_precond, boom_omega);
2165 }
2166 if (num_grid_sweeps)
2167 {
2168 boom_ngs = hypre_CTAlloc(HYPRE_Int, 4, HYPRE_MEMORY_HOST);
2169 for (i=0; i < 4; i++)
2170 boom_ngs[i] = num_grid_sweeps[i];
2171 hypre_BoomerAMGSetNumGridSweeps(pcg_precond, boom_ngs);
2172 if (grid_relax_points)
2173 {
2174 boom_grp = hypre_CTAlloc(HYPRE_Int*, 4, HYPRE_MEMORY_HOST);
2175 for (i=0; i < 4; i++)
2176 {
2177 boom_grp[i] = hypre_CTAlloc(HYPRE_Int, num_grid_sweeps[i], HYPRE_MEMORY_HOST);
2178 for (j=0; j < num_grid_sweeps[i]; j++)
2179 boom_grp[i][j] = grid_relax_points[i][j];
2180 }
2181 hypre_BoomerAMGSetGridRelaxPoints(pcg_precond, boom_grp);
2182 }
2183 }
2184 if (dof_func)
2185 {
2186 boom_dof_func = hypre_CTAlloc(HYPRE_Int, num_variables, HYPRE_MEMORY_HOST);
2187 for (i=0; i < num_variables; i++)
2188 boom_dof_func[i] = dof_func[i];
2189 hypre_BoomerAMGSetDofFunc(pcg_precond, boom_dof_func);
2190 }
2191 pcg_precond_solve = (HYPRE_Int (*)(void*, void*, void*, void*)) hypre_BoomerAMGSolve;
2192 pcg_precond_setup = (HYPRE_Int (*)(void*, void*, void*, void*)) hypre_BoomerAMGSetup;
2193 (AMGhybrid_data -> pcg_precond_setup) = pcg_precond_setup;
2194 (AMGhybrid_data -> pcg_precond_solve) = pcg_precond_solve;
2195 (AMGhybrid_data -> pcg_precond) = pcg_precond;
2196 /*(AMGhybrid_data -> pcg_default) = 0;*/
2197 /*(AMGhybrid_data -> setup_type) = 0;*/
2198 }
2199 else
2200 {
2201 pcg_precond = (AMGhybrid_data -> pcg_precond);
2202 pcg_precond_solve = (AMGhybrid_data -> pcg_precond_solve);
2203 pcg_precond_setup = (AMGhybrid_data -> pcg_precond_setup);
2204 hypre_BoomerAMGSetSetupType(pcg_precond, setup_type);
2205 }
2206
2207 /* Complete setup of solver+AMG */
2208 if (solver_type == 1)
2209 {
2210 hypre_PCGSetPrecond((void*) pcg_solver,
2211 (HYPRE_Int (*)(void*, void*, void*, void*)) pcg_precond_solve,
2212 (HYPRE_Int (*)(void*, void*, void*, void*)) pcg_precond_setup,
2213 (void*) pcg_precond);
2214
2215 hypre_PCGSetup(pcg_solver, (void*) A, (void*) b, (void*) x);
2216
2217 tt2 = hypre_MPI_Wtime();
2218 AMGhybrid_data->setup_time2 = tt2 - tt1;
2219
2220 /* Solve */
2221 tt1 = tt2;
2222
2223 hypre_PCGSolve(pcg_solver, (void*) A, (void*) b, (void*) x);
2224
2225 /* Get information from PCG that is always logged in AMGhybrid solver*/
2226 hypre_PCGGetNumIterations(pcg_solver, &pcg_num_its);
2227 (AMGhybrid_data -> pcg_num_its) = pcg_num_its;
2228 if (logging)
2229 {
2230 hypre_PCGGetFinalRelativeResidualNorm(pcg_solver, &res_norm);
2231 (AMGhybrid_data -> final_rel_res_norm) = res_norm;
2232 }
2233
2234 tt2 = hypre_MPI_Wtime();
2235 AMGhybrid_data->solve_time2 = tt2 - tt1;
2236 }
2237 else if (solver_type == 2)
2238 {
2239 hypre_GMRESSetPrecond((void*) pcg_solver,
2240 (HYPRE_Int (*)(void*, void*, void*, void*)) pcg_precond_solve,
2241 (HYPRE_Int (*)(void*, void*, void*, void*)) pcg_precond_setup,
2242 (void*) pcg_precond);
2243
2244 hypre_GMRESSetup(pcg_solver, (void*) A, (void*) b, (void*) x);
2245
2246 tt2 = hypre_MPI_Wtime();
2247 AMGhybrid_data->setup_time2 = tt2 - tt1;
2248
2249 /* Solve */
2250 tt1 = tt2;
2251
2252 hypre_GMRESSolve(pcg_solver, (void*) A, (void*) b, (void*) x);
2253
2254 /* Get information from GMRES that is always logged in AMGhybrid solver*/
2255 hypre_GMRESGetNumIterations(pcg_solver, &pcg_num_its);
2256 (AMGhybrid_data -> pcg_num_its) = pcg_num_its;
2257 if (logging)
2258 {
2259 hypre_GMRESGetFinalRelativeResidualNorm(pcg_solver, &res_norm);
2260 (AMGhybrid_data -> final_rel_res_norm) = res_norm;
2261 }
2262
2263 tt2 = hypre_MPI_Wtime();
2264 AMGhybrid_data->solve_time2 = tt2 - tt1;
2265 }
2266 else if (solver_type == 3)
2267 {
2268 hypre_BiCGSTABSetPrecond((void*) pcg_solver,
2269 (HYPRE_Int (*)(void*, void*, void*, void*)) pcg_precond_solve,
2270 (HYPRE_Int (*)(void*, void*, void*, void*)) pcg_precond_setup,
2271 (void*) pcg_precond);
2272
2273 hypre_BiCGSTABSetup(pcg_solver, (void*) A, (void*) b, (void*) x);
2274
2275 tt2 = hypre_MPI_Wtime();
2276 AMGhybrid_data->setup_time2 = tt2 - tt1;
2277
2278 /* Solve */
2279 tt1 = tt2;
2280
2281 hypre_BiCGSTABSolve(pcg_solver, (void*) A, (void*) b, (void*) x);
2282
2283 /* Get information from BiCGSTAB that is always logged in AMGhybrid solver*/
2284 hypre_BiCGSTABGetNumIterations(pcg_solver, &pcg_num_its);
2285 (AMGhybrid_data -> pcg_num_its) = pcg_num_its;
2286 if (logging)
2287 {
2288 hypre_BiCGSTABGetFinalRelativeResidualNorm(pcg_solver, &res_norm);
2289 (AMGhybrid_data -> final_rel_res_norm) = res_norm;
2290 }
2291
2292 tt2 = hypre_MPI_Wtime();
2293 AMGhybrid_data->solve_time2 = tt2 - tt1;
2294 }
2295 }
2296
2297 return hypre_error_flag;
2298 }
2299
2300