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 * GMRES gmres
11 *
12 *****************************************************************************/
13
14 #include "krylov.h"
15 #include "_hypre_utilities.h"
16
17 /*--------------------------------------------------------------------------
18 * hypre_GMRESFunctionsCreate
19 *--------------------------------------------------------------------------*/
20
21 hypre_GMRESFunctions *
hypre_GMRESFunctionsCreate(void * (* CAlloc)(size_t count,size_t elt_size,HYPRE_MemoryLocation location),HYPRE_Int (* Free)(void * ptr),HYPRE_Int (* CommInfo)(void * A,HYPRE_Int * my_id,HYPRE_Int * num_procs),void * (* CreateVector)(void * vector),void * (* CreateVectorArray)(HYPRE_Int size,void * vectors),HYPRE_Int (* DestroyVector)(void * vector),void * (* MatvecCreate)(void * A,void * x),HYPRE_Int (* Matvec)(void * matvec_data,HYPRE_Complex alpha,void * A,void * x,HYPRE_Complex beta,void * y),HYPRE_Int (* MatvecDestroy)(void * matvec_data),HYPRE_Real (* InnerProd)(void * x,void * y),HYPRE_Int (* CopyVector)(void * x,void * y),HYPRE_Int (* ClearVector)(void * x),HYPRE_Int (* ScaleVector)(HYPRE_Complex alpha,void * x),HYPRE_Int (* Axpy)(HYPRE_Complex alpha,void * x,void * y),HYPRE_Int (* PrecondSetup)(void * vdata,void * A,void * b,void * x),HYPRE_Int (* Precond)(void * vdata,void * A,void * b,void * x))22 hypre_GMRESFunctionsCreate(
23 void * (*CAlloc) ( size_t count, size_t elt_size, HYPRE_MemoryLocation location ),
24 HYPRE_Int (*Free) ( void *ptr ),
25 HYPRE_Int (*CommInfo) ( void *A, HYPRE_Int *my_id,
26 HYPRE_Int *num_procs ),
27 void * (*CreateVector) ( void *vector ),
28 void * (*CreateVectorArray) ( HYPRE_Int size, void *vectors ),
29 HYPRE_Int (*DestroyVector) ( void *vector ),
30 void * (*MatvecCreate) ( void *A, void *x ),
31 HYPRE_Int (*Matvec) ( void *matvec_data, HYPRE_Complex alpha, void *A,
32 void *x, HYPRE_Complex beta, void *y ),
33 HYPRE_Int (*MatvecDestroy) ( void *matvec_data ),
34 HYPRE_Real (*InnerProd) ( void *x, void *y ),
35 HYPRE_Int (*CopyVector) ( void *x, void *y ),
36 HYPRE_Int (*ClearVector) ( void *x ),
37 HYPRE_Int (*ScaleVector) ( HYPRE_Complex alpha, void *x ),
38 HYPRE_Int (*Axpy) ( HYPRE_Complex alpha, void *x, void *y ),
39 HYPRE_Int (*PrecondSetup) ( void *vdata, void *A, void *b, void *x ),
40 HYPRE_Int (*Precond) ( void *vdata, void *A, void *b, void *x )
41 )
42 {
43 hypre_GMRESFunctions * gmres_functions;
44 gmres_functions = (hypre_GMRESFunctions *)
45 CAlloc( 1, sizeof(hypre_GMRESFunctions), HYPRE_MEMORY_HOST );
46
47 gmres_functions->CAlloc = CAlloc;
48 gmres_functions->Free = Free;
49 gmres_functions->CommInfo = CommInfo; /* not in PCGFunctionsCreate */
50 gmres_functions->CreateVector = CreateVector;
51 gmres_functions->CreateVectorArray = CreateVectorArray; /* not in PCGFunctionsCreate */
52 gmres_functions->DestroyVector = DestroyVector;
53 gmres_functions->MatvecCreate = MatvecCreate;
54 gmres_functions->Matvec = Matvec;
55 gmres_functions->MatvecDestroy = MatvecDestroy;
56 gmres_functions->InnerProd = InnerProd;
57 gmres_functions->CopyVector = CopyVector;
58 gmres_functions->ClearVector = ClearVector;
59 gmres_functions->ScaleVector = ScaleVector;
60 gmres_functions->Axpy = Axpy;
61 /* default preconditioner must be set here but can be changed later... */
62 gmres_functions->precond_setup = PrecondSetup;
63 gmres_functions->precond = Precond;
64
65 return gmres_functions;
66 }
67
68 /*--------------------------------------------------------------------------
69 * hypre_GMRESCreate
70 *--------------------------------------------------------------------------*/
71
72 void *
hypre_GMRESCreate(hypre_GMRESFunctions * gmres_functions)73 hypre_GMRESCreate( hypre_GMRESFunctions *gmres_functions )
74 {
75 hypre_GMRESData *gmres_data;
76
77 HYPRE_ANNOTATE_FUNC_BEGIN;
78
79 gmres_data = hypre_CTAllocF(hypre_GMRESData, 1, gmres_functions, HYPRE_MEMORY_HOST);
80 gmres_data->functions = gmres_functions;
81
82 /* set defaults */
83 (gmres_data -> k_dim) = 5;
84 (gmres_data -> tol) = 1.0e-06; /* relative residual tol */
85 (gmres_data -> cf_tol) = 0.0;
86 (gmres_data -> a_tol) = 0.0; /* abs. residual tol */
87 (gmres_data -> min_iter) = 0;
88 (gmres_data -> max_iter) = 1000;
89 (gmres_data -> rel_change) = 0;
90 (gmres_data -> skip_real_r_check) = 0;
91 (gmres_data -> stop_crit) = 0; /* rel. residual norm - this is obsolete!*/
92 (gmres_data -> converged) = 0;
93 (gmres_data -> hybrid) = 0;
94 (gmres_data -> precond_data) = NULL;
95 (gmres_data -> print_level) = 0;
96 (gmres_data -> logging) = 0;
97 (gmres_data -> p) = NULL;
98 (gmres_data -> r) = NULL;
99 (gmres_data -> w) = NULL;
100 (gmres_data -> w_2) = NULL;
101 (gmres_data -> matvec_data) = NULL;
102 (gmres_data -> norms) = NULL;
103 (gmres_data -> log_file_name) = NULL;
104
105 HYPRE_ANNOTATE_FUNC_END;
106
107 return (void *) gmres_data;
108 }
109
110 /*--------------------------------------------------------------------------
111 * hypre_GMRESDestroy
112 *--------------------------------------------------------------------------*/
113
114 HYPRE_Int
hypre_GMRESDestroy(void * gmres_vdata)115 hypre_GMRESDestroy( void *gmres_vdata )
116 {
117 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
118 HYPRE_Int i;
119
120 HYPRE_ANNOTATE_FUNC_BEGIN;
121 if (gmres_data)
122 {
123 hypre_GMRESFunctions *gmres_functions = gmres_data->functions;
124 if ( (gmres_data->logging>0) || (gmres_data->print_level) > 0 )
125 {
126 if ( (gmres_data -> norms) != NULL )
127 hypre_TFreeF( gmres_data -> norms, gmres_functions );
128 }
129
130 if ( (gmres_data -> matvec_data) != NULL )
131 (*(gmres_functions->MatvecDestroy))(gmres_data -> matvec_data);
132
133 if ( (gmres_data -> r) != NULL )
134 (*(gmres_functions->DestroyVector))(gmres_data -> r);
135 if ( (gmres_data -> w) != NULL )
136 (*(gmres_functions->DestroyVector))(gmres_data -> w);
137 if ( (gmres_data -> w_2) != NULL )
138 (*(gmres_functions->DestroyVector))(gmres_data -> w_2);
139
140
141 if ( (gmres_data -> p) != NULL )
142 {
143 for (i = 0; i < (gmres_data -> k_dim+1); i++)
144 {
145 if ( (gmres_data -> p)[i] != NULL )
146 {
147 (*(gmres_functions->DestroyVector))( (gmres_data -> p) [i]);
148 }
149 }
150 hypre_TFreeF( gmres_data->p, gmres_functions );
151 }
152 hypre_TFreeF( gmres_data, gmres_functions );
153 hypre_TFreeF( gmres_functions, gmres_functions );
154 }
155 HYPRE_ANNOTATE_FUNC_END;
156
157 return hypre_error_flag;
158 }
159
160 /*--------------------------------------------------------------------------
161 * hypre_GMRESGetResidual
162 *--------------------------------------------------------------------------*/
163
hypre_GMRESGetResidual(void * gmres_vdata,void ** residual)164 HYPRE_Int hypre_GMRESGetResidual( void *gmres_vdata, void **residual )
165 {
166
167 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
168 *residual = gmres_data->r;
169
170 return hypre_error_flag;
171 }
172
173 /*--------------------------------------------------------------------------
174 * hypre_GMRESSetup
175 *--------------------------------------------------------------------------*/
176
177 HYPRE_Int
hypre_GMRESSetup(void * gmres_vdata,void * A,void * b,void * x)178 hypre_GMRESSetup( void *gmres_vdata,
179 void *A,
180 void *b,
181 void *x )
182 {
183 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
184 hypre_GMRESFunctions *gmres_functions = gmres_data->functions;
185
186 HYPRE_Int k_dim = (gmres_data -> k_dim);
187 HYPRE_Int max_iter = (gmres_data -> max_iter);
188 HYPRE_Int (*precond_setup)(void*,void*,void*,void*) = (gmres_functions->precond_setup);
189 void *precond_data = (gmres_data -> precond_data);
190
191 HYPRE_Int rel_change = (gmres_data -> rel_change);
192
193 HYPRE_ANNOTATE_FUNC_BEGIN;
194
195 (gmres_data -> A) = A;
196
197 /*--------------------------------------------------
198 * The arguments for NewVector are important to
199 * maintain consistency between the setup and
200 * compute phases of matvec and the preconditioner.
201 *--------------------------------------------------*/
202
203 if ((gmres_data -> p) == NULL)
204 (gmres_data -> p) = (void**)(*(gmres_functions->CreateVectorArray))(k_dim+1,x);
205 if ((gmres_data -> r) == NULL)
206 (gmres_data -> r) = (*(gmres_functions->CreateVector))(b);
207 if ((gmres_data -> w) == NULL)
208 (gmres_data -> w) = (*(gmres_functions->CreateVector))(b);
209
210 if (rel_change)
211 {
212 if ((gmres_data -> w_2) == NULL)
213 (gmres_data -> w_2) = (*(gmres_functions->CreateVector))(b);
214 }
215
216
217 if ((gmres_data -> matvec_data) == NULL)
218 (gmres_data -> matvec_data) = (*(gmres_functions->MatvecCreate))(A, x);
219
220 precond_setup(precond_data, A, b, x);
221
222 /*-----------------------------------------------------
223 * Allocate space for log info
224 *-----------------------------------------------------*/
225
226 if ( (gmres_data->logging)>0 || (gmres_data->print_level) > 0 )
227 {
228 if ((gmres_data -> norms) != NULL)
229 hypre_TFreeF(gmres_data -> norms,gmres_functions);
230 (gmres_data -> norms) = hypre_CTAllocF(HYPRE_Real, max_iter + 1,gmres_functions, HYPRE_MEMORY_HOST);
231 }
232 if ( (gmres_data->print_level) > 0 )
233 {
234 if ((gmres_data -> log_file_name) == NULL)
235 {
236 (gmres_data -> log_file_name) = (char*)"gmres.out.log";
237 }
238 }
239
240 HYPRE_ANNOTATE_FUNC_END;
241
242 return hypre_error_flag;
243 }
244
245 /*--------------------------------------------------------------------------
246 * hypre_GMRESSolve
247 *-------------------------------------------------------------------------*/
248
249 HYPRE_Int
hypre_GMRESSolve(void * gmres_vdata,void * A,void * b,void * x)250 hypre_GMRESSolve(void *gmres_vdata,
251 void *A,
252 void *b,
253 void *x)
254 {
255 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
256 hypre_GMRESFunctions *gmres_functions = gmres_data->functions;
257 HYPRE_Int k_dim = (gmres_data -> k_dim);
258 HYPRE_Int min_iter = (gmres_data -> min_iter);
259 HYPRE_Int max_iter = (gmres_data -> max_iter);
260 HYPRE_Int rel_change = (gmres_data -> rel_change);
261 HYPRE_Int skip_real_r_check = (gmres_data -> skip_real_r_check);
262 HYPRE_Int hybrid = (gmres_data -> hybrid);
263 HYPRE_Real r_tol = (gmres_data -> tol);
264 HYPRE_Real cf_tol = (gmres_data -> cf_tol);
265 HYPRE_Real a_tol = (gmres_data -> a_tol);
266 void *matvec_data = (gmres_data -> matvec_data);
267 void *r = (gmres_data -> r);
268 void *w = (gmres_data -> w);
269 /* note: w_2 is only allocated if rel_change = 1 */
270 void *w_2 = (gmres_data -> w_2);
271
272 void **p = (gmres_data -> p);
273
274
275 HYPRE_Int (*precond)(void*,void*,void*,void*) = (gmres_functions -> precond);
276 HYPRE_Int *precond_data = (HYPRE_Int*) (gmres_data -> precond_data);
277
278 HYPRE_Int print_level = (gmres_data -> print_level);
279 HYPRE_Int logging = (gmres_data -> logging);
280
281 HYPRE_Real *norms = (gmres_data -> norms);
282 /* not used yet char *log_file_name = (gmres_data -> log_file_name);*/
283 /* FILE *fp; */
284
285 HYPRE_Int break_value = 0;
286 HYPRE_Int i, j, k;
287 HYPRE_Real *rs, **hh, *c, *s, *rs_2;
288 HYPRE_Int iter;
289 HYPRE_Int my_id, num_procs;
290 HYPRE_Real epsilon, gamma, t, r_norm, b_norm, den_norm, x_norm;
291 HYPRE_Real w_norm;
292
293 HYPRE_Real epsmac = 1.e-16;
294 HYPRE_Real ieee_check = 0.;
295
296 HYPRE_Real guard_zero_residual;
297 HYPRE_Real cf_ave_0 = 0.0;
298 HYPRE_Real cf_ave_1 = 0.0;
299 HYPRE_Real weight;
300 HYPRE_Real r_norm_0;
301 HYPRE_Real relative_error = 1.0;
302
303 HYPRE_Int rel_change_passed = 0, num_rel_change_check = 0;
304
305 HYPRE_Real real_r_norm_old, real_r_norm_new;
306
307 HYPRE_ANNOTATE_FUNC_BEGIN;
308
309 (gmres_data -> converged) = 0;
310 /*-----------------------------------------------------------------------
311 * With relative change convergence test on, it is possible to attempt
312 * another iteration with a zero residual. This causes the parameter
313 * alpha to go NaN. The guard_zero_residual parameter is to circumvent
314 * this. Perhaps it should be set to something non-zero (but small).
315 *-----------------------------------------------------------------------*/
316 guard_zero_residual = 0.0;
317
318 (*(gmres_functions->CommInfo))(A,&my_id,&num_procs);
319 if ( logging>0 || print_level>0 )
320 {
321 norms = (gmres_data -> norms);
322 }
323
324 /* initialize work arrays */
325 rs = hypre_CTAllocF(HYPRE_Real,k_dim+1,gmres_functions, HYPRE_MEMORY_HOST);
326 c = hypre_CTAllocF(HYPRE_Real,k_dim,gmres_functions, HYPRE_MEMORY_HOST);
327 s = hypre_CTAllocF(HYPRE_Real,k_dim,gmres_functions, HYPRE_MEMORY_HOST);
328 if (rel_change)
329 {
330 rs_2 = hypre_CTAllocF(HYPRE_Real,k_dim+1,gmres_functions, HYPRE_MEMORY_HOST);
331 }
332 hh = hypre_CTAllocF(HYPRE_Real*,k_dim+1,gmres_functions, HYPRE_MEMORY_HOST);
333 for (i=0; i < k_dim+1; i++)
334 {
335 hh[i] = hypre_CTAllocF(HYPRE_Real,k_dim,gmres_functions, HYPRE_MEMORY_HOST);
336 }
337
338 (*(gmres_functions->CopyVector))(b,p[0]);
339
340 /* compute initial residual */
341 (*(gmres_functions->Matvec))(matvec_data,-1.0, A, x, 1.0, p[0]);
342
343 b_norm = sqrt((*(gmres_functions->InnerProd))(b,b));
344 real_r_norm_old = b_norm;
345
346 /* Since it is does not diminish performance, attempt to return an error flag
347 and notify users when they supply bad input. */
348 if (b_norm != 0.)
349 {
350 ieee_check = b_norm/b_norm; /* INF -> NaN conversion */
351 }
352 if (ieee_check != ieee_check)
353 {
354 /* ...INFs or NaNs in input can make ieee_check a NaN. This test
355 for ieee_check self-equality works on all IEEE-compliant compilers/
356 machines, c.f. page 8 of "Lecture Notes on the Status of IEEE 754"
357 by W. Kahan, May 31, 1996. Currently (July 2002) this paper may be
358 found at http://HTTP.CS.Berkeley.EDU/~wkahan/ieee754status/IEEE754.PDF */
359 if (logging > 0 || print_level > 0)
360 {
361 hypre_printf("\n\nERROR detected by Hypre ... BEGIN\n");
362 hypre_printf("ERROR -- hypre_GMRESSolve: INFs and/or NaNs detected in input.\n");
363 hypre_printf("User probably placed non-numerics in supplied b.\n");
364 hypre_printf("Returning error flag += 101. Program not terminated.\n");
365 hypre_printf("ERROR detected by Hypre ... END\n\n\n");
366 }
367 hypre_error(HYPRE_ERROR_GENERIC);
368 HYPRE_ANNOTATE_FUNC_END;
369
370 return hypre_error_flag;
371 }
372
373 r_norm = sqrt((*(gmres_functions->InnerProd))(p[0],p[0]));
374 r_norm_0 = r_norm;
375
376 /* Since it is does not diminish performance, attempt to return an error flag
377 and notify users when they supply bad input. */
378 if (r_norm != 0.)
379 {
380 ieee_check = r_norm/r_norm; /* INF -> NaN conversion */
381 }
382 if (ieee_check != ieee_check)
383 {
384 /* ...INFs or NaNs in input can make ieee_check a NaN. This test
385 for ieee_check self-equality works on all IEEE-compliant compilers/
386 machines, c.f. page 8 of "Lecture Notes on the Status of IEEE 754"
387 by W. Kahan, May 31, 1996. Currently (July 2002) this paper may be
388 found at http://HTTP.CS.Berkeley.EDU/~wkahan/ieee754status/IEEE754.PDF */
389 if (logging > 0 || print_level > 0)
390 {
391 hypre_printf("\n\nERROR detected by Hypre ... BEGIN\n");
392 hypre_printf("ERROR -- hypre_GMRESSolve: INFs and/or NaNs detected in input.\n");
393 hypre_printf("User probably placed non-numerics in supplied A or x_0.\n");
394 hypre_printf("Returning error flag += 101. Program not terminated.\n");
395 hypre_printf("ERROR detected by Hypre ... END\n\n\n");
396 }
397 hypre_error(HYPRE_ERROR_GENERIC);
398 HYPRE_ANNOTATE_FUNC_END;
399
400 return hypre_error_flag;
401 }
402
403 if ( logging>0 || print_level > 0)
404 {
405 norms[0] = r_norm;
406 if ( print_level>1 && my_id == 0 )
407 {
408 hypre_printf("L2 norm of b: %e\n", b_norm);
409 if (b_norm == 0.0)
410 {
411 hypre_printf("Rel_resid_norm actually contains the residual norm\n");
412 }
413 hypre_printf("Initial L2 norm of residual: %e\n", r_norm);
414 }
415 }
416 iter = 0;
417
418 if (b_norm > 0.0)
419 {
420 /* convergence criterion |r_i|/|b| <= accuracy if |b| > 0 */
421 den_norm= b_norm;
422 }
423 else
424 {
425 /* convergence criterion |r_i|/|r0| <= accuracy if |b| = 0 */
426 den_norm= r_norm;
427 }
428
429
430 /* convergence criteria: |r_i| <= max( a_tol, r_tol * den_norm)
431 den_norm = |r_0| or |b|
432 note: default for a_tol is 0.0, so relative residual criteria is used unless
433 user specifies a_tol, or sets r_tol = 0.0, which means absolute
434 tol only is checked */
435
436 epsilon = hypre_max(a_tol,r_tol*den_norm);
437
438 /* so now our stop criteria is |r_i| <= epsilon */
439
440 if ( print_level>1 && my_id == 0 )
441 {
442 if (b_norm > 0.0)
443 {
444 hypre_printf("=============================================\n\n");
445 hypre_printf("Iters resid.norm conv.rate rel.res.norm\n");
446 hypre_printf("----- ------------ ---------- ------------\n");
447 }
448 else
449 {
450 hypre_printf("=============================================\n\n");
451 hypre_printf("Iters resid.norm conv.rate\n");
452 hypre_printf("----- ------------ ----------\n");
453 }
454 }
455
456 /* once the rel. change check has passed, we do not want to check it again */
457 rel_change_passed = 0;
458
459 /* outer iteration cycle */
460 while (iter < max_iter)
461 {
462 /* initialize first term of hessenberg system */
463
464 rs[0] = r_norm;
465 if (r_norm == 0.0)
466 {
467 hypre_TFreeF(c,gmres_functions);
468 hypre_TFreeF(s,gmres_functions);
469 hypre_TFreeF(rs,gmres_functions);
470 if (rel_change) hypre_TFreeF(rs_2,gmres_functions);
471 for (i=0; i < k_dim+1; i++) hypre_TFreeF(hh[i],gmres_functions);
472 hypre_TFreeF(hh,gmres_functions);
473 HYPRE_ANNOTATE_FUNC_END;
474
475 return hypre_error_flag;
476 }
477
478 /* see if we are already converged and
479 should print the final norm and exit */
480 if (r_norm <= epsilon && iter >= min_iter)
481 {
482 if (!rel_change) /* shouldn't exit after no iterations if
483 * relative change is on*/
484 {
485 (*(gmres_functions->CopyVector))(b,r);
486 (*(gmres_functions->Matvec))(matvec_data,-1.0,A,x,1.0,r);
487 r_norm = sqrt((*(gmres_functions->InnerProd))(r,r));
488 if (r_norm <= epsilon)
489 {
490 if ( print_level>1 && my_id == 0)
491 {
492 hypre_printf("\n\n");
493 hypre_printf("Final L2 norm of residual: %e\n\n", r_norm);
494 }
495 break;
496 }
497 else
498 {
499 if ( print_level>0 && my_id == 0)
500 {
501 hypre_printf("false convergence 1\n");
502 }
503 }
504 }
505 }
506
507 t = 1.0 / r_norm;
508 (*(gmres_functions->ScaleVector))(t,p[0]);
509 i = 0;
510
511 /***RESTART CYCLE (right-preconditioning) ***/
512 while (i < k_dim && iter < max_iter)
513 {
514 i++;
515 iter++;
516 (*(gmres_functions->ClearVector))(r);
517 precond(precond_data, A, p[i-1], r);
518 (*(gmres_functions->Matvec))(matvec_data, 1.0, A, r, 0.0, p[i]);
519 /* modified Gram_Schmidt */
520 for (j=0; j < i; j++)
521 {
522 hh[j][i-1] = (*(gmres_functions->InnerProd))(p[j],p[i]);
523 (*(gmres_functions->Axpy))(-hh[j][i-1],p[j],p[i]);
524 }
525 t = sqrt((*(gmres_functions->InnerProd))(p[i],p[i]));
526 hh[i][i-1] = t;
527 if (t != 0.0)
528 {
529 t = 1.0/t;
530 (*(gmres_functions->ScaleVector))(t,p[i]);
531 }
532 /* done with modified Gram_schmidt and Arnoldi step.
533 update factorization of hh */
534 for (j = 1; j < i; j++)
535 {
536 t = hh[j-1][i-1];
537 hh[j-1][i-1] = s[j-1]*hh[j][i-1] + c[j-1]*t;
538 hh[j][i-1] = -s[j-1]*t + c[j-1]*hh[j][i-1];
539 }
540 t= hh[i][i-1]*hh[i][i-1];
541 t+= hh[i-1][i-1]*hh[i-1][i-1];
542 gamma = sqrt(t);
543 if (gamma == 0.0)
544 {
545 gamma = epsmac;
546 }
547 c[i-1] = hh[i-1][i-1]/gamma;
548 s[i-1] = hh[i][i-1]/gamma;
549 rs[i] = -hh[i][i-1]*rs[i-1];
550 rs[i]/= gamma;
551 rs[i-1] = c[i-1]*rs[i-1];
552 /* determine residual norm */
553 hh[i-1][i-1] = s[i-1]*hh[i][i-1] + c[i-1]*hh[i-1][i-1];
554 r_norm = fabs(rs[i]);
555
556 /* print ? */
557 if ( print_level>0 )
558 {
559 norms[iter] = r_norm;
560 if ( print_level>1 && my_id == 0 )
561 {
562 if (b_norm > 0.0)
563 {
564 hypre_printf("% 5d %e %f %e\n", iter,
565 norms[iter],norms[iter]/norms[iter-1],
566 norms[iter]/b_norm);
567 }
568 else
569 {
570 hypre_printf("% 5d %e %f\n", iter, norms[iter],
571 norms[iter]/norms[iter-1]);
572 }
573 }
574 }
575 /*convergence factor tolerance */
576 if (cf_tol > 0.0)
577 {
578 cf_ave_0 = cf_ave_1;
579 cf_ave_1 = pow( r_norm / r_norm_0, 1.0/(2.0*iter));
580
581 weight = fabs(cf_ave_1 - cf_ave_0);
582 weight = weight / hypre_max(cf_ave_1, cf_ave_0);
583 weight = 1.0 - weight;
584 #if 0
585 hypre_printf("I = %d: cf_new = %e, cf_old = %e, weight = %e\n",
586 i, cf_ave_1, cf_ave_0, weight );
587 #endif
588 if (weight * cf_ave_1 > cf_tol)
589 {
590 break_value = 1;
591 break;
592 }
593 }
594 /* should we exit the restart cycle? (conv. check) */
595 if (r_norm <= epsilon && iter >= min_iter)
596 {
597 if (rel_change && !rel_change_passed)
598 {
599
600 /* To decide whether to break here: to actually
601 determine the relative change requires the approx
602 solution (so a triangular solve) and a
603 precond. solve - so if we have to do this many
604 times, it will be expensive...(unlike cg where is
605 is relatively straightforward)
606
607 previously, the intent (there was a bug), was to
608 exit the restart cycle based on the residual norm
609 and check the relative change outside the cycle.
610 Here we will check the relative here as we don't
611 want to exit the restart cycle prematurely */
612
613 for (k=0; k<i; k++)
614 {
615 /* extra copy of rs so we don't need to change the later solve */
616 rs_2[k] = rs[k];
617 }
618
619 /* solve tri. system*/
620 rs_2[i-1] = rs_2[i-1]/hh[i-1][i-1];
621 for (k = i-2; k >= 0; k--)
622 {
623 t = 0.0;
624 for (j = k+1; j < i; j++)
625 {
626 t -= hh[k][j]*rs_2[j];
627 }
628 t+= rs_2[k];
629 rs_2[k] = t/hh[k][k];
630 }
631
632 (*(gmres_functions->CopyVector))(p[i-1],w);
633 (*(gmres_functions->ScaleVector))(rs_2[i-1],w);
634 for (j = i-2; j >=0; j--)
635 {
636 (*(gmres_functions->Axpy))(rs_2[j], p[j], w);
637 }
638 (*(gmres_functions->ClearVector))(r);
639 /* find correction (in r) */
640 precond(precond_data, A, w, r);
641 /* copy current solution (x) to w (don't want to over-write x)*/
642 (*(gmres_functions->CopyVector))(x,w);
643
644 /* add the correction */
645 (*(gmres_functions->Axpy))(1.0,r,w);
646
647 /* now w is the approx solution - get the norm*/
648 x_norm = sqrt( (*(gmres_functions->InnerProd))(w,w) );
649
650 if ( !(x_norm <= guard_zero_residual ))
651 /* don't divide by zero */
652 { /* now get x_i - x_i-1 */
653
654 if (num_rel_change_check)
655 {
656 /* have already checked once so we can avoid another precond.
657 solve */
658 (*(gmres_functions->CopyVector))(w, r);
659 (*(gmres_functions->Axpy))(-1.0, w_2, r);
660 /* now r contains x_i - x_i-1*/
661
662 /* save current soln w in w_2 for next time */
663 (*(gmres_functions->CopyVector))(w, w_2);
664 }
665 else
666 {
667 /* first time to check rel change*/
668
669 /* first save current soln w in w_2 for next time */
670 (*(gmres_functions->CopyVector))(w, w_2);
671
672 /* for relative change take x_(i-1) to be
673 x + M^{-1}[sum{j=0..i-2} rs_j p_j ].
674 Now
675 x_i - x_{i-1}= {x + M^{-1}[sum{j=0..i-1} rs_j p_j ]}
676 - {x + M^{-1}[sum{j=0..i-2} rs_j p_j ]}
677 = M^{-1} rs_{i-1}{p_{i-1}} */
678
679 (*(gmres_functions->ClearVector))(w);
680 (*(gmres_functions->Axpy))(rs_2[i-1], p[i-1], w);
681 (*(gmres_functions->ClearVector))(r);
682 /* apply the preconditioner */
683 precond(precond_data, A, w, r);
684 /* now r contains x_i - x_i-1 */
685 }
686 /* find the norm of x_i - x_i-1 */
687 w_norm = sqrt( (*(gmres_functions->InnerProd))(r,r) );
688 relative_error = w_norm/x_norm;
689 if (relative_error <= r_tol)
690 {
691 rel_change_passed = 1;
692 break;
693 }
694 }
695 else
696 {
697 rel_change_passed = 1;
698 break;
699
700 }
701 num_rel_change_check++;
702 }
703 else /* no relative change */
704 {
705 break;
706 }
707 }
708 } /*** end of restart cycle ***/
709
710 /* now compute solution, first solve upper triangular system */
711
712 if (break_value)
713 {
714 break;
715 }
716
717 rs[i-1] = rs[i-1]/hh[i-1][i-1];
718 for (k = i-2; k >= 0; k--)
719 {
720 t = 0.0;
721 for (j = k+1; j < i; j++)
722 {
723 t -= hh[k][j]*rs[j];
724 }
725 t += rs[k];
726 rs[k] = t/hh[k][k];
727 }
728
729 (*(gmres_functions->CopyVector))(p[i-1],w);
730 (*(gmres_functions->ScaleVector))(rs[i-1],w);
731 for (j = i-2; j >=0; j--)
732 (*(gmres_functions->Axpy))(rs[j], p[j], w);
733
734 (*(gmres_functions->ClearVector))(r);
735 /* find correction (in r) */
736 precond(precond_data, A, w, r);
737
738 /* update current solution x (in x) */
739 (*(gmres_functions->Axpy))(1.0,r,x);
740
741 /* check for convergence by evaluating the actual residual */
742 if (r_norm <= epsilon && iter >= min_iter)
743 {
744 if (skip_real_r_check)
745 {
746 (gmres_data -> converged) = 1;
747 break;
748 }
749
750 /* calculate actual residual norm*/
751 (*(gmres_functions->CopyVector))(b,r);
752 (*(gmres_functions->Matvec))(matvec_data,-1.0,A,x,1.0,r);
753 real_r_norm_new = r_norm = sqrt( (*(gmres_functions->InnerProd))(r,r) );
754
755 if (r_norm <= epsilon)
756 {
757 if (rel_change && !rel_change_passed) /* calculate the relative change */
758 {
759 /* calculate the norm of the solution */
760 x_norm = sqrt( (*(gmres_functions->InnerProd))(x,x) );
761
762 if ( !(x_norm <= guard_zero_residual ))
763 /* don't divide by zero */
764 {
765 /* for relative change take x_(i-1) to be
766 x + M^{-1}[sum{j=0..i-2} rs_j p_j ].
767 Now
768 x_i - x_{i-1}= {x + M^{-1}[sum{j=0..i-1} rs_j p_j ]}
769 - {x + M^{-1}[sum{j=0..i-2} rs_j p_j ]}
770 = M^{-1} rs_{i-1}{p_{i-1}} */
771 (*(gmres_functions->ClearVector))(w);
772 (*(gmres_functions->Axpy))(rs[i-1], p[i-1], w);
773 (*(gmres_functions->ClearVector))(r);
774 /* apply the preconditioner */
775 precond(precond_data, A, w, r);
776 /* find the norm of x_i - x_i-1 */
777 w_norm = sqrt( (*(gmres_functions->InnerProd))(r,r) );
778 relative_error= w_norm/x_norm;
779 if ( relative_error < r_tol )
780 {
781 (gmres_data -> converged) = 1;
782 if ( print_level>1 && my_id == 0 )
783 {
784 hypre_printf("\n\n");
785 hypre_printf("Final L2 norm of residual: %e\n\n", r_norm);
786 }
787 break;
788 }
789 }
790 else
791 {
792 (gmres_data -> converged) = 1;
793 if ( print_level>1 && my_id == 0 )
794 {
795 hypre_printf("\n\n");
796 hypre_printf("Final L2 norm of residual: %e\n\n", r_norm);
797 }
798 break;
799 }
800
801 }
802 else /* don't need to check rel. change */
803 {
804 if ( print_level>1 && my_id == 0 )
805 {
806 hypre_printf("\n\n");
807 hypre_printf("Final L2 norm of residual: %e\n\n", r_norm);
808 }
809 (gmres_data -> converged) = 1;
810 break;
811 }
812 }
813 else /* conv. has not occurred, according to true residual */
814 {
815 /* exit if the real residual norm has not decreased */
816 if (real_r_norm_new >= real_r_norm_old)
817 {
818 if (print_level > 1 && my_id == 0)
819 {
820 hypre_printf("\n\n");
821 hypre_printf("Final L2 norm of residual: %e\n\n", r_norm);
822 }
823 (gmres_data -> converged) = 1;
824 break;
825 }
826
827 /* report discrepancy between real/GMRES residuals and restart */
828 if ( print_level>0 && my_id == 0)
829 {
830 hypre_printf("false convergence 2, L2 norm of residual: %e\n", r_norm);
831 }
832 (*(gmres_functions->CopyVector))(r,p[0]);
833 i = 0;
834 real_r_norm_old = real_r_norm_new;
835 }
836 } /* end of convergence check */
837
838 /* compute residual vector and continue loop */
839 for (j=i ; j > 0; j--)
840 {
841 rs[j-1] = -s[j-1]*rs[j];
842 rs[j] = c[j-1]*rs[j];
843 }
844
845 if (i) (*(gmres_functions->Axpy))(rs[i]-1.0,p[i],p[i]);
846 for (j=i-1 ; j > 0; j--)
847 (*(gmres_functions->Axpy))(rs[j],p[j],p[i]);
848
849 if (i)
850 {
851 (*(gmres_functions->Axpy))(rs[0]-1.0,p[0],p[0]);
852 (*(gmres_functions->Axpy))(1.0,p[i],p[0]);
853 }
854 } /* END of iteration while loop */
855
856
857 if ( print_level>1 && my_id == 0 )
858 {
859 hypre_printf("\n\n");
860 }
861
862 (gmres_data -> num_iterations) = iter;
863
864 if (b_norm > 0.0)
865 {
866 (gmres_data -> rel_residual_norm) = r_norm/b_norm;
867 }
868
869 if (b_norm == 0.0)
870 {
871 (gmres_data -> rel_residual_norm) = r_norm;
872 }
873
874 if (iter >= max_iter && r_norm > epsilon && epsilon > 0 && hybrid != -1)
875 {
876 hypre_error(HYPRE_ERROR_CONV);
877 }
878
879 hypre_TFreeF(c, gmres_functions);
880 hypre_TFreeF(s, gmres_functions);
881 hypre_TFreeF(rs, gmres_functions);
882
883 if (rel_change)
884 {
885 hypre_TFreeF(rs_2,gmres_functions);
886 }
887
888 for (i=0; i < k_dim+1; i++)
889 {
890 hypre_TFreeF(hh[i],gmres_functions);
891 }
892
893 hypre_TFreeF(hh, gmres_functions);
894
895 HYPRE_ANNOTATE_FUNC_END;
896
897 return hypre_error_flag;
898 }
899
900 /*--------------------------------------------------------------------------
901 * hypre_GMRESSetKDim, hypre_GMRESGetKDim
902 *--------------------------------------------------------------------------*/
903
904 HYPRE_Int
hypre_GMRESSetKDim(void * gmres_vdata,HYPRE_Int k_dim)905 hypre_GMRESSetKDim( void *gmres_vdata,
906 HYPRE_Int k_dim )
907 {
908 hypre_GMRESData *gmres_data =(hypre_GMRESData *) gmres_vdata;
909
910
911 (gmres_data -> k_dim) = k_dim;
912
913 return hypre_error_flag;
914
915 }
916
917 HYPRE_Int
hypre_GMRESGetKDim(void * gmres_vdata,HYPRE_Int * k_dim)918 hypre_GMRESGetKDim( void *gmres_vdata,
919 HYPRE_Int *k_dim )
920 {
921 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
922
923
924 *k_dim = (gmres_data -> k_dim);
925
926 return hypre_error_flag;
927 }
928
929 /*--------------------------------------------------------------------------
930 * hypre_GMRESSetTol, hypre_GMRESGetTol
931 *--------------------------------------------------------------------------*/
932
933 HYPRE_Int
hypre_GMRESSetTol(void * gmres_vdata,HYPRE_Real tol)934 hypre_GMRESSetTol( void *gmres_vdata,
935 HYPRE_Real tol )
936 {
937 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
938
939
940 (gmres_data -> tol) = tol;
941
942 return hypre_error_flag;
943 }
944
945 HYPRE_Int
hypre_GMRESGetTol(void * gmres_vdata,HYPRE_Real * tol)946 hypre_GMRESGetTol( void *gmres_vdata,
947 HYPRE_Real *tol )
948 {
949 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
950
951
952 *tol = (gmres_data -> tol);
953
954 return hypre_error_flag;
955 }
956 /*--------------------------------------------------------------------------
957 * hypre_GMRESSetAbsoluteTol, hypre_GMRESGetAbsoluteTol
958 *--------------------------------------------------------------------------*/
959
960 HYPRE_Int
hypre_GMRESSetAbsoluteTol(void * gmres_vdata,HYPRE_Real a_tol)961 hypre_GMRESSetAbsoluteTol( void *gmres_vdata,
962 HYPRE_Real a_tol )
963 {
964 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
965
966
967 (gmres_data -> a_tol) = a_tol;
968
969 return hypre_error_flag;
970 }
971
972 HYPRE_Int
hypre_GMRESGetAbsoluteTol(void * gmres_vdata,HYPRE_Real * a_tol)973 hypre_GMRESGetAbsoluteTol( void *gmres_vdata,
974 HYPRE_Real *a_tol )
975 {
976 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
977
978
979 *a_tol = (gmres_data -> a_tol);
980
981 return hypre_error_flag;
982 }
983 /*--------------------------------------------------------------------------
984 * hypre_GMRESSetConvergenceFactorTol, hypre_GMRESGetConvergenceFactorTol
985 *--------------------------------------------------------------------------*/
986
987 HYPRE_Int
hypre_GMRESSetConvergenceFactorTol(void * gmres_vdata,HYPRE_Real cf_tol)988 hypre_GMRESSetConvergenceFactorTol( void *gmres_vdata,
989 HYPRE_Real cf_tol )
990 {
991 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
992
993
994 (gmres_data -> cf_tol) = cf_tol;
995
996 return hypre_error_flag;
997 }
998
999 HYPRE_Int
hypre_GMRESGetConvergenceFactorTol(void * gmres_vdata,HYPRE_Real * cf_tol)1000 hypre_GMRESGetConvergenceFactorTol( void *gmres_vdata,
1001 HYPRE_Real *cf_tol )
1002 {
1003 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1004
1005
1006 *cf_tol = (gmres_data -> cf_tol);
1007
1008 return hypre_error_flag;
1009 }
1010
1011 /*--------------------------------------------------------------------------
1012 * hypre_GMRESSetMinIter, hypre_GMRESGetMinIter
1013 *--------------------------------------------------------------------------*/
1014
1015 HYPRE_Int
hypre_GMRESSetMinIter(void * gmres_vdata,HYPRE_Int min_iter)1016 hypre_GMRESSetMinIter( void *gmres_vdata,
1017 HYPRE_Int min_iter )
1018 {
1019 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1020
1021
1022 (gmres_data -> min_iter) = min_iter;
1023
1024 return hypre_error_flag;
1025 }
1026
1027 HYPRE_Int
hypre_GMRESGetMinIter(void * gmres_vdata,HYPRE_Int * min_iter)1028 hypre_GMRESGetMinIter( void *gmres_vdata,
1029 HYPRE_Int *min_iter )
1030 {
1031 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1032
1033
1034 *min_iter = (gmres_data -> min_iter);
1035
1036 return hypre_error_flag;
1037 }
1038
1039 /*--------------------------------------------------------------------------
1040 * hypre_GMRESSetMaxIter, hypre_GMRESGetMaxIter
1041 *--------------------------------------------------------------------------*/
1042
1043 HYPRE_Int
hypre_GMRESSetMaxIter(void * gmres_vdata,HYPRE_Int max_iter)1044 hypre_GMRESSetMaxIter( void *gmres_vdata,
1045 HYPRE_Int max_iter )
1046 {
1047 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1048
1049
1050 (gmres_data -> max_iter) = max_iter;
1051
1052 return hypre_error_flag;
1053 }
1054
1055 HYPRE_Int
hypre_GMRESGetMaxIter(void * gmres_vdata,HYPRE_Int * max_iter)1056 hypre_GMRESGetMaxIter( void *gmres_vdata,
1057 HYPRE_Int *max_iter )
1058 {
1059 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1060
1061
1062 *max_iter = (gmres_data -> max_iter);
1063
1064 return hypre_error_flag;
1065 }
1066
1067 /*--------------------------------------------------------------------------
1068 * hypre_GMRESSetRelChange, hypre_GMRESGetRelChange
1069 *--------------------------------------------------------------------------*/
1070
1071 HYPRE_Int
hypre_GMRESSetRelChange(void * gmres_vdata,HYPRE_Int rel_change)1072 hypre_GMRESSetRelChange( void *gmres_vdata,
1073 HYPRE_Int rel_change )
1074 {
1075 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1076
1077
1078 (gmres_data -> rel_change) = rel_change;
1079
1080 return hypre_error_flag;
1081 }
1082
1083 HYPRE_Int
hypre_GMRESGetRelChange(void * gmres_vdata,HYPRE_Int * rel_change)1084 hypre_GMRESGetRelChange( void *gmres_vdata,
1085 HYPRE_Int *rel_change )
1086 {
1087 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1088
1089
1090 *rel_change = (gmres_data -> rel_change);
1091
1092 return hypre_error_flag;
1093 }
1094
1095 /*--------------------------------------------------------------------------
1096 * hypre_GMRESSetSkipRealResidualCheck, hypre_GMRESGetSkipRealResidualCheck
1097 *--------------------------------------------------------------------------*/
1098
1099 HYPRE_Int
hypre_GMRESSetSkipRealResidualCheck(void * gmres_vdata,HYPRE_Int skip_real_r_check)1100 hypre_GMRESSetSkipRealResidualCheck( void *gmres_vdata,
1101 HYPRE_Int skip_real_r_check )
1102 {
1103 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1104
1105 (gmres_data -> skip_real_r_check) = skip_real_r_check;
1106
1107 return hypre_error_flag;
1108 }
1109
1110 HYPRE_Int
hypre_GMRESGetSkipRealResidualCheck(void * gmres_vdata,HYPRE_Int * skip_real_r_check)1111 hypre_GMRESGetSkipRealResidualCheck( void *gmres_vdata,
1112 HYPRE_Int *skip_real_r_check)
1113 {
1114 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1115
1116 *skip_real_r_check = (gmres_data -> skip_real_r_check);
1117
1118 return hypre_error_flag;
1119 }
1120
1121 /*--------------------------------------------------------------------------
1122 * hypre_GMRESSetStopCrit, hypre_GMRESGetStopCrit
1123 *
1124 * OBSOLETE
1125 *--------------------------------------------------------------------------*/
1126
1127 HYPRE_Int
hypre_GMRESSetStopCrit(void * gmres_vdata,HYPRE_Int stop_crit)1128 hypre_GMRESSetStopCrit( void *gmres_vdata,
1129 HYPRE_Int stop_crit )
1130 {
1131 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1132
1133
1134 (gmres_data -> stop_crit) = stop_crit;
1135
1136 return hypre_error_flag;
1137 }
1138
1139 HYPRE_Int
hypre_GMRESGetStopCrit(void * gmres_vdata,HYPRE_Int * stop_crit)1140 hypre_GMRESGetStopCrit( void *gmres_vdata,
1141 HYPRE_Int *stop_crit )
1142 {
1143 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1144
1145
1146 *stop_crit = (gmres_data -> stop_crit);
1147
1148 return hypre_error_flag;
1149 }
1150
1151 /*--------------------------------------------------------------------------
1152 * hypre_GMRESSetPrecond
1153 *--------------------------------------------------------------------------*/
1154
1155 HYPRE_Int
hypre_GMRESSetPrecond(void * gmres_vdata,HYPRE_Int (* precond)(void *,void *,void *,void *),HYPRE_Int (* precond_setup)(void *,void *,void *,void *),void * precond_data)1156 hypre_GMRESSetPrecond( void *gmres_vdata,
1157 HYPRE_Int (*precond)(void*,void*,void*,void*),
1158 HYPRE_Int (*precond_setup)(void*,void*,void*,void*),
1159 void *precond_data )
1160 {
1161 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1162 hypre_GMRESFunctions *gmres_functions = gmres_data->functions;
1163
1164
1165 (gmres_functions -> precond) = precond;
1166 (gmres_functions -> precond_setup) = precond_setup;
1167 (gmres_data -> precond_data) = precond_data;
1168
1169 return hypre_error_flag;
1170 }
1171
1172 /*--------------------------------------------------------------------------
1173 * hypre_GMRESGetPrecond
1174 *--------------------------------------------------------------------------*/
1175
1176 HYPRE_Int
hypre_GMRESGetPrecond(void * gmres_vdata,HYPRE_Solver * precond_data_ptr)1177 hypre_GMRESGetPrecond( void *gmres_vdata,
1178 HYPRE_Solver *precond_data_ptr )
1179 {
1180 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1181
1182
1183 *precond_data_ptr = (HYPRE_Solver)(gmres_data -> precond_data);
1184
1185 return hypre_error_flag;
1186 }
1187
1188 /*--------------------------------------------------------------------------
1189 * hypre_GMRESSetPrintLevel, hypre_GMRESGetPrintLevel
1190 *--------------------------------------------------------------------------*/
1191
1192 HYPRE_Int
hypre_GMRESSetPrintLevel(void * gmres_vdata,HYPRE_Int level)1193 hypre_GMRESSetPrintLevel( void *gmres_vdata,
1194 HYPRE_Int level)
1195 {
1196 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1197
1198
1199 (gmres_data -> print_level) = level;
1200
1201 return hypre_error_flag;
1202 }
1203
1204 HYPRE_Int
hypre_GMRESGetPrintLevel(void * gmres_vdata,HYPRE_Int * level)1205 hypre_GMRESGetPrintLevel( void *gmres_vdata,
1206 HYPRE_Int *level)
1207 {
1208 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1209
1210
1211 *level = (gmres_data -> print_level);
1212
1213 return hypre_error_flag;
1214 }
1215
1216 /*--------------------------------------------------------------------------
1217 * hypre_GMRESSetLogging, hypre_GMRESGetLogging
1218 *--------------------------------------------------------------------------*/
1219
1220 HYPRE_Int
hypre_GMRESSetLogging(void * gmres_vdata,HYPRE_Int level)1221 hypre_GMRESSetLogging( void *gmres_vdata,
1222 HYPRE_Int level)
1223 {
1224 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1225
1226 (gmres_data -> logging) = level;
1227
1228 return hypre_error_flag;
1229 }
1230
1231 HYPRE_Int
hypre_GMRESGetLogging(void * gmres_vdata,HYPRE_Int * level)1232 hypre_GMRESGetLogging( void *gmres_vdata,
1233 HYPRE_Int *level)
1234 {
1235 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1236
1237 *level = (gmres_data -> logging);
1238
1239 return hypre_error_flag;
1240 }
1241
1242 HYPRE_Int
hypre_GMRESSetHybrid(void * gmres_vdata,HYPRE_Int level)1243 hypre_GMRESSetHybrid( void *gmres_vdata,
1244 HYPRE_Int level)
1245 {
1246 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1247
1248 (gmres_data -> hybrid) = level;
1249
1250 return hypre_error_flag;
1251 }
1252
1253 /*--------------------------------------------------------------------------
1254 * hypre_GMRESGetNumIterations
1255 *--------------------------------------------------------------------------*/
1256
1257 HYPRE_Int
hypre_GMRESGetNumIterations(void * gmres_vdata,HYPRE_Int * num_iterations)1258 hypre_GMRESGetNumIterations( void *gmres_vdata,
1259 HYPRE_Int *num_iterations )
1260 {
1261 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1262
1263
1264 *num_iterations = (gmres_data -> num_iterations);
1265
1266 return hypre_error_flag;
1267 }
1268
1269 /*--------------------------------------------------------------------------
1270 * hypre_GMRESGetConverged
1271 *--------------------------------------------------------------------------*/
1272
1273 HYPRE_Int
hypre_GMRESGetConverged(void * gmres_vdata,HYPRE_Int * converged)1274 hypre_GMRESGetConverged( void *gmres_vdata,
1275 HYPRE_Int *converged )
1276 {
1277 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1278
1279
1280 *converged = (gmres_data -> converged);
1281
1282 return hypre_error_flag;
1283 }
1284
1285 /*--------------------------------------------------------------------------
1286 * hypre_GMRESGetFinalRelativeResidualNorm
1287 *--------------------------------------------------------------------------*/
1288
1289 HYPRE_Int
hypre_GMRESGetFinalRelativeResidualNorm(void * gmres_vdata,HYPRE_Real * relative_residual_norm)1290 hypre_GMRESGetFinalRelativeResidualNorm( void *gmres_vdata,
1291 HYPRE_Real *relative_residual_norm )
1292 {
1293 hypre_GMRESData *gmres_data = (hypre_GMRESData *)gmres_vdata;
1294
1295
1296 *relative_residual_norm = (gmres_data -> rel_residual_norm);
1297
1298 return hypre_error_flag;
1299 }
1300