1 #include "BSprivate.h"
2
3 /*@ BSpar_solve - General solver of a system of equations preconditioned
4 by one of several preconditioners and using one of several
5 possible methods. The rhs can be a block of vectors.
6
7 Input Parameters:
8 . A - a sparse matrix
9 . fact_A - the incomplete factored version of A, if any
10 . comm_A - the communication structure for A
11 . rhs - the contiguous block of vectors forming the rhs
12 . residual - the final computed residual
13
14 Output Parameters:
15 . x - the contiguous block of vectors containing the solution (can
16 contain an initial guess if BSctx_set_guess() is set.
17
18 Now specified in procinfo context:
19 . pre_option - the preconditioner to use:
20 $ PRE_ICC: incomplete Cholesky factorization
21 $ PRE_ILU: incomplete LU factorization
22 $ PRE_SSOR: Successive over relaxation
23 $ PRE_BJACOBI: Block Jacobi
24 . err_tol - the tolerance to which to solve the problem
25 stop if the estimated norm of the residual divided by
26 the norm of the rhs is less than err_tol
27 . max_iter - the maximum number of iterations to take
28 . guess - if TRUE, then initialize out_x to 0, otherwise the program
29 assumes that out_x contains an initial guess
30 . procinfo - the usual processor stuff
31
32 Returns:
33 The number of iterations or a negative number indicating the number
34 of iterations prior to finding that the matrix (or preconditioner)
35 is not positive definite.
36
37 Notes:
38 The preconditioners must be computed prior to calling BSpar_solve.
39 For more information on the preconditioners, see the manual or
40 BSctx_set_pre().
41
42 @*/
BSpar_solve(BSpar_mat * A,BSpar_mat * fact_A,BScomm * comm_A,FLOAT * rhs,FLOAT * x,FLOAT * residual,BSprocinfo * procinfo)43 int BSpar_solve(BSpar_mat *A, BSpar_mat *fact_A, BScomm *comm_A,
44 FLOAT *rhs, FLOAT *x, FLOAT *residual, BSprocinfo *procinfo)
45 {
46 int num_iter;
47
48 /* Do some error checking */
49 if((A->icc_storage)&&(!procinfo->single)) {
50 if(comm_A->num_rhs!=procinfo->num_rhs) {
51 MY_SETERRCN(PAR_SOLVE_ERROR,"num_rhs in context does not agree num_rhs in comm\n");
52 }
53 }
54 if((procinfo->scaling)&&(A->scale_diag==NULL)) {
55 MY_SETERRCN(SCALING_ERROR,"Scaling set in context however matrix not scaled\n");
56 }
57
58 if ((procinfo->print) && (PSISROOT(procinfo))) {
59 BSctx_print(procinfo);
60 }
61
62 if(procinfo->method==GMRES) {
63 num_iter = BSpar_gmres(procinfo->num_rhs,A,fact_A,comm_A,rhs,x,
64 procinfo->preconditioner,procinfo->restart,
65 procinfo->residual_tol,procinfo->max_iterations,
66 residual,procinfo->guess,procinfo);
67 } else {
68 num_iter = BSpar_sym_solve(procinfo->num_rhs,A,fact_A,
69 comm_A,rhs,x,procinfo->preconditioner,
70 procinfo->residual_tol,procinfo->max_iterations,
71 residual,procinfo->guess,procinfo);
72 }
73 return(num_iter);
74 }
75