source: CIVL/examples/mpi-omp/AMG2013/krylov/pcg.c

main
Last change on this file was ea777aa, checked in by Alex Wilton <awilton@…>, 3 years ago

Moved examples, include, build_default.properties, common.xml, and README out from dev.civl.com into the root of the repo.

git-svn-id: ​svn://vsl.cis.udel.edu/civl/trunk@5704 fb995dde-84ed-4084-dfe6-e5aef3e2452c
  • Property mode set to 100644
File size: 30.0 KB
Line 
1/*BHEADER**********************************************************************
2 * Copyright (c) 2008, Lawrence Livermore National Security, LLC.
3 * Produced at the Lawrence Livermore National Laboratory.
4 * This file is part of HYPRE. See file COPYRIGHT for details.
5 *
6 * HYPRE is free software; you can redistribute it and/or modify it under the
7 * terms of the GNU Lesser General Public License (as published by the Free
8 * Software Foundation) version 2.1 dated February 1999.
9 *
10 * $Revision: 2.4 $
11 ***********************************************************************EHEADER*/
12
13
14
15
16/******************************************************************************
17 *
18 * Preconditioned conjugate gradient (Omin) functions
19 *
20 *****************************************************************************/
21
22/* This was based on the pcg.c formerly in struct_ls, with
23 changes (GetPrecond and stop_crit) for compatibility with the pcg.c
24 in parcsr_ls and elsewhere. Incompatibilities with the
25 parcsr_ls version:
26 - logging is different; no attempt has been made to be the same
27 - treatment of b=0 in Ax=b is different: this returns x=0; the parcsr
28 version iterates with a special stopping criterion
29*/
30
31#include "krylov.h"
32#include "utilities.h"
33
34/*--------------------------------------------------------------------------
35 * hypre_PCGFunctionsCreate
36 *--------------------------------------------------------------------------*/
37
38hypre_PCGFunctions *
39hypre_PCGFunctionsCreate(
40 char * (*CAlloc) ( int count, int elt_size ),
41 int (*Free) ( char *ptr ),
42 int (*CommInfo) ( void *A, int *my_id, int *num_procs ),
43 void * (*CreateVector) ( void *vector ),
44 int (*DestroyVector) ( void *vector ),
45 void * (*MatvecCreate) ( void *A, void *x ),
46 int (*Matvec) ( void *matvec_data, double alpha, void *A,
47 void *x, double beta, void *y ),
48 int (*MatvecDestroy) ( void *matvec_data ),
49 double (*InnerProd) ( void *x, void *y ),
50 int (*CopyVector) ( void *x, void *y ),
51 int (*ClearVector) ( void *x ),
52 int (*ScaleVector) ( double alpha, void *x ),
53 int (*Axpy) ( double alpha, void *x, void *y ),
54 int (*PrecondSetup) ( void *vdata, void *A, void *b, void *x ),
55 int (*Precond) ( void *vdata, void *A, void *b, void *x )
56 )
57{
58 hypre_PCGFunctions * pcg_functions;
59 pcg_functions = (hypre_PCGFunctions *)
60 CAlloc( 1, sizeof(hypre_PCGFunctions) );
61
62 pcg_functions->CAlloc = CAlloc;
63 pcg_functions->Free = Free;
64 pcg_functions->CommInfo = CommInfo;
65 pcg_functions->CreateVector = CreateVector;
66 pcg_functions->DestroyVector = DestroyVector;
67 pcg_functions->MatvecCreate = MatvecCreate;
68 pcg_functions->Matvec = Matvec;
69 pcg_functions->MatvecDestroy = MatvecDestroy;
70 pcg_functions->InnerProd = InnerProd;
71 pcg_functions->CopyVector = CopyVector;
72 pcg_functions->ClearVector = ClearVector;
73 pcg_functions->ScaleVector = ScaleVector;
74 pcg_functions->Axpy = Axpy;
75/* default preconditioner must be set here but can be changed later... */
76 pcg_functions->precond_setup = PrecondSetup;
77 pcg_functions->precond = Precond;
78
79 return pcg_functions;
80}
81
82/*--------------------------------------------------------------------------
83 * hypre_PCGCreate
84 *--------------------------------------------------------------------------*/
85
86void *
87hypre_PCGCreate( hypre_PCGFunctions *pcg_functions )
88{
89 hypre_PCGData *pcg_data;
90
91 pcg_data = hypre_CTAllocF(hypre_PCGData, 1, pcg_functions);
92
93 pcg_data -> functions = pcg_functions;
94
95 /* set defaults */
96 (pcg_data -> tol) = 1.0e-06;
97 (pcg_data -> atolf) = 0.0;
98 (pcg_data -> cf_tol) = 0.0;
99 (pcg_data -> max_iter) = 1000;
100 (pcg_data -> two_norm) = 0;
101 (pcg_data -> rel_change) = 0;
102 (pcg_data -> stop_crit) = 0;
103 (pcg_data -> converged) = 0;
104 (pcg_data -> owns_matvec_data ) = 1;
105 (pcg_data -> matvec_data) = NULL;
106 (pcg_data -> precond_data) = NULL;
107 (pcg_data -> print_level) = 0;
108 (pcg_data -> logging) = 0;
109 (pcg_data -> norms) = NULL;
110 (pcg_data -> rel_norms) = NULL;
111 (pcg_data -> p) = NULL;
112 (pcg_data -> s) = NULL;
113 (pcg_data -> r) = NULL;
114
115 return (void *) pcg_data;
116}
117
118/*--------------------------------------------------------------------------
119 * hypre_PCGDestroy
120 *--------------------------------------------------------------------------*/
121
122int
123hypre_PCGDestroy( void *pcg_vdata )
124{
125 hypre_PCGData *pcg_data = pcg_vdata;
126 hypre_PCGFunctions *pcg_functions = pcg_data->functions;
127 int ierr = 0;
128
129 if (pcg_data)
130 {
131 if ( (pcg_data -> norms) != NULL )
132 {
133 hypre_TFreeF( pcg_data -> norms, pcg_functions );
134 pcg_data -> norms = NULL;
135 }
136 if ( (pcg_data -> rel_norms) != NULL )
137 {
138 hypre_TFreeF( pcg_data -> rel_norms, pcg_functions );
139 pcg_data -> rel_norms = NULL;
140 }
141 if ( pcg_data -> matvec_data != NULL && pcg_data->owns_matvec_data )
142 {
143 (*(pcg_functions->MatvecDestroy))(pcg_data -> matvec_data);
144 pcg_data -> matvec_data = NULL;
145 }
146 if ( pcg_data -> p != NULL )
147 {
148 (*(pcg_functions->DestroyVector))(pcg_data -> p);
149 pcg_data -> p = NULL;
150 }
151 if ( pcg_data -> s != NULL )
152 {
153 (*(pcg_functions->DestroyVector))(pcg_data -> s);
154 pcg_data -> s = NULL;
155 }
156 if ( pcg_data -> r != NULL )
157 {
158 (*(pcg_functions->DestroyVector))(pcg_data -> r);
159 pcg_data -> r = NULL;
160 }
161 hypre_TFreeF( pcg_data, pcg_functions );
162 hypre_TFreeF( pcg_functions, pcg_functions );
163 }
164
165 return(ierr);
166}
167
168/*--------------------------------------------------------------------------
169 * hypre_PCGGetResidual
170 *--------------------------------------------------------------------------*/
171
172int hypre_PCGGetResidual( void *pcg_vdata, void **residual )
173{
174 /* returns a pointer to the residual vector */
175 int ierr = 0;
176 hypre_PCGData *pcg_data = pcg_vdata;
177 *residual = pcg_data->r;
178 return ierr;
179}
180
181/*--------------------------------------------------------------------------
182 * hypre_PCGSetup
183 *--------------------------------------------------------------------------*/
184
185int
186hypre_PCGSetup( void *pcg_vdata,
187 void *A,
188 void *b,
189 void *x )
190{
191 hypre_PCGData *pcg_data = pcg_vdata;
192 hypre_PCGFunctions *pcg_functions = pcg_data->functions;
193 int max_iter = (pcg_data -> max_iter);
194 int (*precond_setup)(void *vdata, void *A, void *b, void *x) = (pcg_functions -> precond_setup);
195 void *precond_data = (pcg_data -> precond_data);
196 int ierr = 0;
197
198 (pcg_data -> A) = A;
199
200 /*--------------------------------------------------
201 * The arguments for CreateVector are important to
202 * maintain consistency between the setup and
203 * compute phases of matvec and the preconditioner.
204 *--------------------------------------------------*/
205
206 if ( pcg_data -> p != NULL )
207 (*(pcg_functions->DestroyVector))(pcg_data -> p);
208 (pcg_data -> p) = (*(pcg_functions->CreateVector))(x);
209
210 if ( pcg_data -> s != NULL )
211 (*(pcg_functions->DestroyVector))(pcg_data -> s);
212 (pcg_data -> s) = (*(pcg_functions->CreateVector))(x);
213
214 if ( pcg_data -> r != NULL )
215 (*(pcg_functions->DestroyVector))(pcg_data -> r);
216 (pcg_data -> r) = (*(pcg_functions->CreateVector))(b);
217
218 if ( pcg_data -> matvec_data != NULL && pcg_data->owns_matvec_data )
219 (*(pcg_functions->MatvecDestroy))(pcg_data -> matvec_data);
220 (pcg_data -> matvec_data) = (*(pcg_functions->MatvecCreate))(A, x);
221
222 ierr = precond_setup(precond_data, A, b, x);
223
224 /*-----------------------------------------------------
225 * Allocate space for log info
226 *-----------------------------------------------------*/
227
228 if ( (pcg_data->logging)>0 || (pcg_data->print_level)>0 )
229 {
230 if ( (pcg_data -> norms) != NULL )
231 hypre_TFreeF( pcg_data -> norms, pcg_functions );
232 (pcg_data -> norms) = hypre_CTAllocF( double, max_iter + 1,
233 pcg_functions);
234
235 if ( (pcg_data -> rel_norms) != NULL )
236 hypre_TFreeF( pcg_data -> rel_norms, pcg_functions );
237 (pcg_data -> rel_norms) = hypre_CTAllocF( double, max_iter + 1,
238 pcg_functions );
239 }
240
241 return ierr;
242}
243
244/*--------------------------------------------------------------------------
245 * hypre_PCGSolve
246 *--------------------------------------------------------------------------
247 *
248 * We use the following convergence test as the default (see Ashby, Holst,
249 * Manteuffel, and Saylor):
250 *
251 * ||e||_A ||r||_C
252 * ------- <= [kappa_A(C*A)]^(1/2) ------- < tol
253 * ||x||_A ||b||_C
254 *
255 * where we let (for the time being) kappa_A(CA) = 1.
256 * We implement the test as:
257 *
258 * gamma = <C*r,r>/<C*b,b> < (tol^2) = eps
259 *
260 *--------------------------------------------------------------------------*/
261
262int
263hypre_PCGSolve( void *pcg_vdata,
264 void *A,
265 void *b,
266 void *x )
267{
268 hypre_PCGData *pcg_data = pcg_vdata;
269 hypre_PCGFunctions *pcg_functions = pcg_data->functions;
270
271 double tol = (pcg_data -> tol);
272 double atolf = (pcg_data -> atolf);
273 double cf_tol = (pcg_data -> cf_tol);
274 int max_iter = (pcg_data -> max_iter);
275 int two_norm = (pcg_data -> two_norm);
276 int rel_change = (pcg_data -> rel_change);
277 int stop_crit = (pcg_data -> stop_crit);
278/*
279 int converged = (pcg_data -> converged);
280*/
281 void *p = (pcg_data -> p);
282 void *s = (pcg_data -> s);
283 void *r = (pcg_data -> r);
284 void *matvec_data = (pcg_data -> matvec_data);
285 int (*precond)(void*, void*, void*, void*) = (pcg_functions -> precond);
286 void *precond_data = (pcg_data -> precond_data);
287 int print_level = (pcg_data -> print_level);
288 int logging = (pcg_data -> logging);
289 double *norms = (pcg_data -> norms);
290 double *rel_norms = (pcg_data -> rel_norms);
291
292 double alpha, beta;
293 double gamma, gamma_old;
294 double bi_prod, i_prod, eps;
295 double pi_prod, xi_prod;
296 double ieee_check = 0.;
297
298 double i_prod_0;
299 double cf_ave_0 = 0.0;
300 double cf_ave_1 = 0.0;
301 double weight;
302 double ratio;
303
304 double guard_zero_residual, sdotp;
305
306 int i = 0;
307 int ierr = 0;
308 int my_id, num_procs;
309
310 (pcg_data -> converged) = 0;
311
312 (*(pcg_functions->CommInfo))(A,&my_id,&num_procs);
313
314 /*-----------------------------------------------------------------------
315 * With relative change convergence test on, it is possible to attempt
316 * another iteration with a zero residual. This causes the parameter
317 * alpha to go NaN. The guard_zero_residual parameter is to circumvent
318 * this. Perhaps it should be set to something non-zero (but small).
319 *-----------------------------------------------------------------------*/
320
321 guard_zero_residual = 0.0;
322
323 /*-----------------------------------------------------------------------
324 * Start pcg solve
325 *-----------------------------------------------------------------------*/
326
327 /* compute eps */
328 if (two_norm)
329 {
330 /* bi_prod = <b,b> */
331 bi_prod = (*(pcg_functions->InnerProd))(b, b);
332 if (print_level > 1 && my_id == 0)
333 printf("<b,b>: %e\n",bi_prod);
334 }
335 else
336 {
337 /* bi_prod = <C*b,b> */
338 (*(pcg_functions->ClearVector))(p);
339 precond(precond_data, A, b, p);
340 bi_prod = (*(pcg_functions->InnerProd))(p, b);
341 if (print_level > 1 && my_id == 0)
342 printf("<C*b,b>: %e\n",bi_prod);
343 };
344
345 /* Since it is does not diminish performance, attempt to return an error flag
346 and notify users when they supply bad input. */
347 if (bi_prod != 0.) ieee_check = bi_prod/bi_prod; /* INF -> NaN conversion */
348 if (ieee_check != ieee_check)
349 {
350 /* ...INFs or NaNs in input can make ieee_check a NaN. This test
351 for ieee_check self-equality works on all IEEE-compliant compilers/
352 machines, c.f. page 8 of "Lecture Notes on the Status of IEEE 754"
353 by W. Kahan, May 31, 1996. Currently (July 2002) this paper may be
354 found at http://HTTP.CS.Berkeley.EDU/~wkahan/ieee754status/IEEE754.PDF */
355 if (print_level > 0 || logging > 0)
356 {
357 printf("\n\nERROR detected by Hypre ... BEGIN\n");
358 printf("ERROR -- hypre_PCGSolve: INFs and/or NaNs detected in input.\n");
359 printf("User probably placed non-numerics in supplied b.\n");
360 printf("Returning error flag += 101. Program not terminated.\n");
361 printf("ERROR detected by Hypre ... END\n\n\n");
362 }
363 ierr += 101;
364 return ierr;
365 }
366
367 eps = tol*tol;
368 if ( bi_prod > 0.0 ) {
369 if ( stop_crit && !rel_change && atolf<=0 ) { /* pure absolute tolerance */
370 eps = eps / bi_prod;
371 /* Note: this section is obsolete. Aside from backwards comatability
372 concerns, we could delete the stop_crit parameter and related code,
373 using tol & atolf instead. */
374 }
375 else if ( atolf>0 ) /* mixed relative and absolute tolerance */
376 bi_prod += atolf;
377 }
378 else /* bi_prod==0.0: the rhs vector b is zero */
379 {
380 /* Set x equal to zero and return */
381 (*(pcg_functions->CopyVector))(b, x);
382 if (logging>0 || print_level>0)
383 {
384 norms[0] = 0.0;
385 rel_norms[i] = 0.0;
386 }
387 ierr = 0;
388 return ierr;
389 /* In this case, for the original parcsr pcg, the code would take special
390 action to force iterations even though the exact value was known. */
391 };
392
393 /* r = b - Ax */
394 (*(pcg_functions->CopyVector))(b, r);
395 (*(pcg_functions->Matvec))(matvec_data, -1.0, A, x, 1.0, r);
396
397 /* p = C*r */
398 (*(pcg_functions->ClearVector))(p);
399 precond(precond_data, A, r, p);
400
401 /* gamma = <r,p> */
402 gamma = (*(pcg_functions->InnerProd))(r,p);
403
404 /* Since it is does not diminish performance, attempt to return an error flag
405 and notify users when they supply bad input. */
406 if (gamma != 0.) ieee_check = gamma/gamma; /* INF -> NaN conversion */
407 if (ieee_check != ieee_check)
408 {
409 /* ...INFs or NaNs in input can make ieee_check a NaN. This test
410 for ieee_check self-equality works on all IEEE-compliant compilers/
411 machines, c.f. page 8 of "Lecture Notes on the Status of IEEE 754"
412 by W. Kahan, May 31, 1996. Currently (July 2002) this paper may be
413 found at http://HTTP.CS.Berkeley.EDU/~wkahan/ieee754status/IEEE754.PDF */
414 if (print_level > 0 || logging > 0)
415 {
416 printf("\n\nERROR detected by Hypre ... BEGIN\n");
417 printf("ERROR -- hypre_PCGSolve: INFs and/or NaNs detected in input.\n");
418 printf("User probably placed non-numerics in supplied A or x_0.\n");
419 printf("Returning error flag += 101. Program not terminated.\n");
420 printf("ERROR detected by Hypre ... END\n\n\n");
421 }
422 ierr += 101;
423 return ierr;
424 }
425
426 /* Set initial residual norm */
427 if ( logging>0 || print_level > 0 || cf_tol > 0.0 )
428 {
429 if (two_norm)
430 i_prod_0 = (*(pcg_functions->InnerProd))(r,r);
431 else
432 i_prod_0 = gamma;
433
434 if ( logging>0 || print_level>0 ) norms[0] = sqrt(i_prod_0);
435 }
436 if ( print_level > 1 && my_id==0 ) /* formerly for par_csr only */
437 {
438 printf("\n\n");
439 if (two_norm)
440 {
441 if ( stop_crit && !rel_change && atolf==0 ) { /* pure absolute tolerance */
442 printf("Iters ||r||_2 conv.rate\n");
443 printf("----- ------------ ---------\n");
444 }
445 else {
446 printf("Iters ||r||_2 conv.rate ||r||_2/||b||_2\n");
447 printf("----- ------------ --------- ------------ \n");
448 }
449 }
450 else /* !two_norm */
451 {
452 printf("Iters ||r||_C ||r||_C/||b||_C\n");
453 printf("----- ------------ ------------ \n");
454 }
455 }
456
457 while ((i+1) <= max_iter)
458 {
459 i++;
460
461 /* s = A*p */
462 (*(pcg_functions->Matvec))(matvec_data, 1.0, A, p, 0.0, s);
463
464 /* alpha = gamma / <s,p> */
465 sdotp = (*(pcg_functions->InnerProd))(s, p);
466 if ( sdotp==0.0 )
467 {
468 ++ierr;
469 if (i==1) i_prod=i_prod_0;
470 break;
471 }
472 alpha = gamma / sdotp;
473
474 gamma_old = gamma;
475
476 /* x = x + alpha*p */
477 (*(pcg_functions->Axpy))(alpha, p, x);
478
479 /* r = r - alpha*s */
480 (*(pcg_functions->Axpy))(-alpha, s, r);
481
482 /* s = C*r */
483 (*(pcg_functions->ClearVector))(s);
484 precond(precond_data, A, r, s);
485
486 /* gamma = <r,s> */
487 gamma = (*(pcg_functions->InnerProd))(r, s);
488
489 /* set i_prod for convergence test */
490 if (two_norm)
491 i_prod = (*(pcg_functions->InnerProd))(r,r);
492 else
493 i_prod = gamma;
494
495#if 0
496 if (two_norm)
497 printf("Iter (%d): ||r||_2 = %e, ||r||_2/||b||_2 = %e\n",
498 i, sqrt(i_prod), (bi_prod ? sqrt(i_prod/bi_prod) : 0));
499 else
500 printf("Iter (%d): ||r||_C = %e, ||r||_C/||b||_C = %e\n",
501 i, sqrt(i_prod), (bi_prod ? sqrt(i_prod/bi_prod) : 0));
502#endif
503
504 /* print norm info */
505 if ( logging>0 || print_level>0 )
506 {
507 norms[i] = sqrt(i_prod);
508 rel_norms[i] = bi_prod ? sqrt(i_prod/bi_prod) : 0;
509 }
510 if ( print_level > 1 && my_id==0 )
511 {
512 if (two_norm)
513 {
514 if ( stop_crit && !rel_change && atolf==0 ) { /* pure absolute tolerance */
515 printf("% 5d %e %f\n", i, norms[i],
516 norms[i]/norms[i-1] );
517 }
518 else
519 {
520 printf("% 5d %e %f %e\n", i, norms[i],
521 norms[i]/norms[i-1], rel_norms[i] );
522 }
523 }
524 else
525 {
526 printf("% 5d %e %f %e\n", i, norms[i],
527 norms[i]/norms[i-1], rel_norms[i] );
528 }
529 }
530
531
532 /* check for convergence */
533 if (i_prod / bi_prod < eps)
534 {
535 if (rel_change && i_prod > guard_zero_residual)
536 {
537 pi_prod = (*(pcg_functions->InnerProd))(p,p);
538 xi_prod = (*(pcg_functions->InnerProd))(x,x);
539 ratio = alpha*alpha*pi_prod/xi_prod;
540 if (ratio < eps)
541 {
542 (pcg_data -> converged) = 1;
543 break;
544 }
545 }
546 else
547 {
548 (pcg_data -> converged) = 1;
549 break;
550 }
551 }
552
553 if ( (gamma<1.0e-292) && ((-gamma)<1.0e-292) ) {
554 ierr = 1;
555 break;
556 }
557 /* ... gamma should be >=0. IEEE subnormal numbers are < 2**(-1022)=2.2e-308
558 (and >= 2**(-1074)=4.9e-324). So a gamma this small means we're getting
559 dangerously close to subnormal or zero numbers (usually if gamma is small,
560 so will be other variables). Thus further calculations risk a crash.
561 Such small gamma generally means no hope of progress anyway. */
562
563 /*--------------------------------------------------------------------
564 * Optional test to see if adequate progress is being made.
565 * The average convergence factor is recorded and compared
566 * against the tolerance 'cf_tol'. The weighting factor is
567 * intended to pay more attention to the test when an accurate
568 * estimate for average convergence factor is available.
569 *--------------------------------------------------------------------*/
570
571 if (cf_tol > 0.0)
572 {
573 cf_ave_0 = cf_ave_1;
574 if ( i_prod_0<1.0e-292 ) {
575 /* i_prod_0 is zero, or (almost) subnormal, yet i_prod wasn't small
576 enough to pass the convergence test. Therefore initial guess was good,
577 and we're just calculating garbage - time to bail out before the
578 next step, which will be a divide by zero (or close to it). */
579 ierr = 1;
580 break;
581 }
582 cf_ave_1 = pow( i_prod / i_prod_0, 1.0/(2.0*i));
583
584 weight = fabs(cf_ave_1 - cf_ave_0);
585 weight = weight / hypre_max(cf_ave_1, cf_ave_0);
586 weight = 1.0 - weight;
587#if 0
588 printf("I = %d: cf_new = %e, cf_old = %e, weight = %e\n",
589 i, cf_ave_1, cf_ave_0, weight );
590#endif
591 if (weight * cf_ave_1 > cf_tol) break;
592 }
593
594 /* beta = gamma / gamma_old */
595 beta = gamma / gamma_old;
596
597 /* p = s + beta p */
598 (*(pcg_functions->ScaleVector))(beta, p);
599 (*(pcg_functions->Axpy))(1.0, s, p);
600 }
601
602 if ( print_level > 1 && my_id==0 ) /* formerly for par_csr only */
603 printf("\n\n");
604
605 (pcg_data -> num_iterations) = i;
606 if (bi_prod > 0.0)
607 (pcg_data -> rel_residual_norm) = sqrt(i_prod/bi_prod);
608 else /* actually, we'll never get here... */
609 (pcg_data -> rel_residual_norm) = 0.0;
610
611 return ierr;
612}
613
614/*--------------------------------------------------------------------------
615 * hypre_PCGSetTol, hypre_PCGGetTol
616 *--------------------------------------------------------------------------*/
617
618int
619hypre_PCGSetTol( void *pcg_vdata,
620 double tol )
621{
622 hypre_PCGData *pcg_data = pcg_vdata;
623 int ierr = 0;
624
625 (pcg_data -> tol) = tol;
626
627 return ierr;
628}
629
630int
631hypre_PCGGetTol( void *pcg_vdata,
632 double * tol )
633{
634 hypre_PCGData *pcg_data = pcg_vdata;
635 int ierr = 0;
636
637 *tol = (pcg_data -> tol);
638
639 return ierr;
640}
641
642/*--------------------------------------------------------------------------
643 * hypre_PCGSetAbsoluteTolFactor, hypre_PCGGetAbsoluteTolFactor
644 *--------------------------------------------------------------------------*/
645
646int
647hypre_PCGSetAbsoluteTolFactor( void *pcg_vdata,
648 double atolf )
649{
650 hypre_PCGData *pcg_data = pcg_vdata;
651 int ierr = 0;
652
653 (pcg_data -> atolf) = atolf;
654
655 return ierr;
656}
657
658int
659hypre_PCGGetAbsoluteTolFactor( void *pcg_vdata,
660 double * atolf )
661{
662 hypre_PCGData *pcg_data = pcg_vdata;
663 int ierr = 0;
664
665 *atolf = (pcg_data -> atolf);
666
667 return ierr;
668}
669
670/*--------------------------------------------------------------------------
671 * hypre_PCGSetConvergenceFactorTol, hypre_PCGGetConvergenceFactorTol
672 *--------------------------------------------------------------------------*/
673
674int
675hypre_PCGSetConvergenceFactorTol( void *pcg_vdata,
676 double cf_tol )
677{
678 hypre_PCGData *pcg_data = pcg_vdata;
679 int ierr = 0;
680
681 (pcg_data -> cf_tol) = cf_tol;
682
683 return ierr;
684}
685
686int
687hypre_PCGGetConvergenceFactorTol( void *pcg_vdata,
688 double * cf_tol )
689{
690 hypre_PCGData *pcg_data = pcg_vdata;
691 int ierr = 0;
692
693 *cf_tol = (pcg_data -> cf_tol);
694
695 return ierr;
696}
697
698/*--------------------------------------------------------------------------
699 * hypre_PCGSetMaxIter, hypre_PCGGetMaxIter
700 *--------------------------------------------------------------------------*/
701
702int
703hypre_PCGSetMaxIter( void *pcg_vdata,
704 int max_iter )
705{
706 hypre_PCGData *pcg_data = pcg_vdata;
707 int ierr = 0;
708
709 (pcg_data -> max_iter) = max_iter;
710
711 return ierr;
712}
713
714int
715hypre_PCGGetMaxIter( void *pcg_vdata,
716 int * max_iter )
717{
718 hypre_PCGData *pcg_data = pcg_vdata;
719 int ierr = 0;
720
721 *max_iter = (pcg_data -> max_iter);
722
723 return ierr;
724}
725
726/*--------------------------------------------------------------------------
727 * hypre_PCGSetTwoNorm, hypre_PCGGetTwoNorm
728 *--------------------------------------------------------------------------*/
729
730int
731hypre_PCGSetTwoNorm( void *pcg_vdata,
732 int two_norm )
733{
734 hypre_PCGData *pcg_data = pcg_vdata;
735 int ierr = 0;
736
737 (pcg_data -> two_norm) = two_norm;
738
739 return ierr;
740}
741
742int
743hypre_PCGGetTwoNorm( void *pcg_vdata,
744 int * two_norm )
745{
746 hypre_PCGData *pcg_data = pcg_vdata;
747 int ierr = 0;
748
749 *two_norm = (pcg_data -> two_norm);
750
751 return ierr;
752}
753
754/*--------------------------------------------------------------------------
755 * hypre_PCGSetRelChange, hypre_PCGGetRelChange
756 *--------------------------------------------------------------------------*/
757
758int
759hypre_PCGSetRelChange( void *pcg_vdata,
760 int rel_change )
761{
762 hypre_PCGData *pcg_data = pcg_vdata;
763 int ierr = 0;
764
765 (pcg_data -> rel_change) = rel_change;
766
767 return ierr;
768}
769
770int
771hypre_PCGGetRelChange( void *pcg_vdata,
772 int * rel_change )
773{
774 hypre_PCGData *pcg_data = pcg_vdata;
775 int ierr = 0;
776
777 *rel_change = (pcg_data -> rel_change);
778
779 return ierr;
780}
781
782/*--------------------------------------------------------------------------
783 * hypre_PCGSetStopCrit, hypre_PCGGetStopCrit
784 *--------------------------------------------------------------------------*/
785
786int
787hypre_PCGSetStopCrit( void *pcg_vdata,
788 int stop_crit )
789{
790 hypre_PCGData *pcg_data = pcg_vdata;
791 int ierr = 0;
792
793 (pcg_data -> stop_crit) = stop_crit;
794
795 return ierr;
796}
797
798int
799hypre_PCGGetStopCrit( void *pcg_vdata,
800 int * stop_crit )
801{
802 hypre_PCGData *pcg_data = pcg_vdata;
803 int ierr = 0;
804
805 *stop_crit = (pcg_data -> stop_crit);
806
807 return ierr;
808}
809
810/*--------------------------------------------------------------------------
811 * hypre_PCGGetPrecond
812 *--------------------------------------------------------------------------*/
813
814int
815hypre_PCGGetPrecond( void *pcg_vdata,
816 HYPRE_Solver *precond_data_ptr )
817{
818 hypre_PCGData *pcg_data = pcg_vdata;
819 int ierr = 0;
820
821 *precond_data_ptr = (HYPRE_Solver)(pcg_data -> precond_data);
822
823 return ierr;
824}
825
826/*--------------------------------------------------------------------------
827 * hypre_PCGSetPrecond
828 *--------------------------------------------------------------------------*/
829
830int
831hypre_PCGSetPrecond( void *pcg_vdata,
832 int (*precond)(void *vdata, void *A, void *b, void *x),
833 int (*precond_setup)(void *vdata, void *A, void *b, void *x),
834 void *precond_data )
835{
836 hypre_PCGData *pcg_data = pcg_vdata;
837 hypre_PCGFunctions *pcg_functions = pcg_data->functions;
838 int ierr = 0;
839
840 (pcg_functions -> precond) = precond;
841 (pcg_functions -> precond_setup) = precond_setup;
842 (pcg_data -> precond_data) = precond_data;
843
844 return ierr;
845}
846
847/*--------------------------------------------------------------------------
848 * hypre_PCGSetPrintLevel, hypre_PCGGetPrintLevel
849 *--------------------------------------------------------------------------*/
850
851int
852hypre_PCGSetPrintLevel( void *pcg_vdata,
853 int level)
854{
855 hypre_PCGData *pcg_data = pcg_vdata;
856 int ierr = 0;
857
858 (pcg_data -> print_level) = level;
859
860 return ierr;
861}
862
863int
864hypre_PCGGetPrintLevel( void *pcg_vdata,
865 int * level)
866{
867 hypre_PCGData *pcg_data = pcg_vdata;
868 int ierr = 0;
869
870 *level = (pcg_data -> print_level);
871
872 return ierr;
873}
874
875/*--------------------------------------------------------------------------
876 * hypre_PCGSetLogging, hypre_PCGGetLogging
877 *--------------------------------------------------------------------------*/
878
879int
880hypre_PCGSetLogging( void *pcg_vdata,
881 int level)
882{
883 hypre_PCGData *pcg_data = pcg_vdata;
884 int ierr = 0;
885
886 (pcg_data -> logging) = level;
887
888 return ierr;
889}
890
891int
892hypre_PCGGetLogging( void *pcg_vdata,
893 int * level)
894{
895 hypre_PCGData *pcg_data = pcg_vdata;
896 int ierr = 0;
897
898 *level = (pcg_data -> logging);
899
900 return ierr;
901}
902
903/*--------------------------------------------------------------------------
904 * hypre_PCGGetNumIterations
905 *--------------------------------------------------------------------------*/
906
907int
908hypre_PCGGetNumIterations( void *pcg_vdata,
909 int *num_iterations )
910{
911 hypre_PCGData *pcg_data = pcg_vdata;
912 int ierr = 0;
913
914 *num_iterations = (pcg_data -> num_iterations);
915
916 return ierr;
917}
918
919/*--------------------------------------------------------------------------
920 * hypre_PCGGetConverged
921 *--------------------------------------------------------------------------*/
922
923int
924hypre_PCGGetConverged( void *pcg_vdata,
925 int *converged)
926{
927 hypre_PCGData *pcg_data = pcg_vdata;
928 int ierr = 0;
929
930 *converged = (pcg_data -> converged);
931
932 return ierr;
933}
934
935/*--------------------------------------------------------------------------
936 * hypre_PCGPrintLogging
937 *--------------------------------------------------------------------------*/
938
939int
940hypre_PCGPrintLogging( void *pcg_vdata,
941 int myid)
942{
943 hypre_PCGData *pcg_data = pcg_vdata;
944
945 int num_iterations = (pcg_data -> num_iterations);
946 int print_level = (pcg_data -> print_level);
947 double *norms = (pcg_data -> norms);
948 double *rel_norms = (pcg_data -> rel_norms);
949
950 int i;
951 int ierr = 0;
952
953 if (myid == 0)
954 {
955 if (print_level > 0)
956 {
957 for (i = 0; i < num_iterations; i++)
958 {
959 printf("Residual norm[%d] = %e ", i, norms[i]);
960 printf("Relative residual norm[%d] = %e\n", i, rel_norms[i]);
961 }
962 }
963 }
964
965 return ierr;
966}
967
968/*--------------------------------------------------------------------------
969 * hypre_PCGGetFinalRelativeResidualNorm
970 *--------------------------------------------------------------------------*/
971
972int
973hypre_PCGGetFinalRelativeResidualNorm( void *pcg_vdata,
974 double *relative_residual_norm )
975{
976 hypre_PCGData *pcg_data = pcg_vdata;
977
978 int ierr = 0;
979
980 double rel_residual_norm = (pcg_data -> rel_residual_norm);
981
982 *relative_residual_norm = rel_residual_norm;
983
984 return ierr;
985}
986
Note: See TracBrowser for help on using the repository browser.