/[escript]/trunk/paso/src/Solver_jacobi.c
ViewVC logotype

Contents of /trunk/paso/src/Solver_jacobi.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3094 - (show annotations)
Fri Aug 13 08:38:06 2010 UTC (9 years, 6 months ago) by gross
File MIME type: text/plain
File size: 3044 byte(s)
The MPI and sequational GAUSS_SEIDEL have been merged.
The couring and main diagonal pointer is now manged by the patternm which means that they are calculated once only even if the preconditioner is deleted.



1
2 /*******************************************************
3 *
4 * Copyright (c) 2003-2010 by University of Queensland
5 * Earth Systems Science Computational Center (ESSCC)
6 * http://www.uq.edu.au/esscc
7 *
8 * Primary Business: Queensland, Australia
9 * Licensed under the Open Software License version 3.0
10 * http://www.opensource.org/licenses/osl-3.0.php
11 *
12 *******************************************************/
13
14
15 /**************************************************************/
16
17 /* Paso: jacobi preconditioner: */
18
19 /**************************************************************/
20
21 /* Copyrights by ACcESS Australia 2003, 2004, 2005 */
22 /* Author: Lutz Gross, l.gross@uq.edu.au */
23
24 /**************************************************************/
25
26 #include "Paso.h"
27 #include "Solver.h"
28
29 /**************************************************************/
30
31 /* frees the Jacobi preconditioner */
32
33 void Paso_Solver_Jacobi_free(Paso_Solver_Jacobi * in) {
34 if (in!=NULL) {
35 MEMFREE(in->values);
36 MEMFREE(in->pivot);
37 MEMFREE(in);
38 }
39 }
40 Paso_Solver_Jacobi* Paso_Solver_getJacobi(Paso_SystemMatrix * A_p)
41 {
42 Paso_Solver_Jacobi* jac=Paso_Solver_getLocalJacobi(A_p->mainBlock);
43 if (Paso_MPIInfo_noError(A_p->mpi_info)) {
44 return jac;
45 } else {
46 return NULL;
47 }
48 }
49
50 /**************************************************************/
51
52 /* Jacobi precondioner set up */
53
54 Paso_Solver_Jacobi* Paso_Solver_getLocalJacobi(Paso_SparseMatrix * A_p) {
55 Paso_Solver_Jacobi* out=NULL;
56 dim_t n = A_p->numCols;
57 dim_t n_block=A_p->row_block_size;
58 dim_t block_size=A_p->block_size;
59
60 /* check matrix is square */
61 if (A_p->col_block_size !=A_p->row_block_size) {
62 Paso_setError(TYPE_ERROR, "Paso_Solver_getJacobi: Jacobi preconditioner square block size.");
63 return NULL;
64 }
65
66 out=MEMALLOC(1,Paso_Solver_Jacobi);
67 if (! Paso_checkPtr(out)) {
68 out->values = MEMALLOC( ((size_t) n) * ((size_t) block_size),double);
69 out->pivot = MEMALLOC( ((size_t) n) * ((size_t) n_block), index_t);
70
71
72 if (! ( Paso_checkPtr(out->values) || Paso_checkPtr(out->pivot) )) {
73 Paso_SparseMatrix_invMain(A_p, out->values, out->pivot );
74 }
75
76 }
77 if (Paso_noError()) {
78 return out;
79 } else {
80 Paso_Solver_Jacobi_free(out);
81 return NULL;
82 }
83 }
84 /**************************************************************/
85
86 /* applies Jacobi preconditioner */
87
88 /* should be called within a parallel region */
89 /* barrier synconization should be performed to make sure that the input vector available */
90
91 void Paso_Solver_solveJacobi(Paso_SystemMatrix * A_p, Paso_Solver_Jacobi * prec, double * x, double * b) {
92 Paso_Solver_solveLocalJacobi( A_p->mainBlock, prec, x, b) ;
93 return;
94 }
95
96 void Paso_Solver_solveLocalJacobi(Paso_SparseMatrix * A_p, Paso_Solver_Jacobi * prec, double * x, double * b) {
97 dim_t n = A_p->numCols;
98 dim_t n_block=A_p->row_block_size;
99 Paso_Solver_applyBlockDiagonalMatrix(n_block,n,prec->values,prec->pivot,x,b);
100 return;
101 }
102

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26