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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3482 - (show annotations)
Wed Mar 23 04:06:52 2011 UTC (8 years, 5 months ago) by gross
File MIME type: text/plain
File size: 4944 byte(s)
some work on AMG MPI
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: AMG sets-up */
18
19 /**************************************************************/
20
21 /* Author: Lutz Gross, l.gross@uq.edu.au */
22
23 /**************************************************************/
24
25 #include "Paso.h"
26 #include "Preconditioner.h"
27
28 /***********************************************************************************/
29
30 /* free space */
31
32 void Paso_Preconditioner_AMG_Root_free(Paso_Preconditioner_AMG_Root* in) {
33 if (in!=NULL) {
34 Paso_Preconditioner_AMG_free(in->amg);
35 Paso_Preconditioner_LocalAMG_free(in->localamg);
36 Paso_Preconditioner_Smoother_free(in->amgsubstitute);
37 MEMFREE(in);
38 }
39 }
40
41 Paso_Preconditioner_AMG_Root* Paso_Preconditioner_AMG_Root_alloc(Paso_SystemMatrix* A,Paso_Options* options) {
42
43 Paso_Preconditioner_AMG_Root* prec=MEMALLOC(1,Paso_Preconditioner_AMG_Root);
44
45 if (! Esys_checkPtr(prec)) {
46
47 prec->amg=NULL;
48 prec->localamg=NULL;
49 prec->amgsubstitute=NULL;
50 prec->is_local=( A->mpi_info->size == 1 ) | options->use_local_preconditioner;
51 prec->is_local=FALSE;
52 /* prec->is_local=TRUE; */
53
54 if (prec->is_local) {
55 prec->localamg=Paso_Preconditioner_LocalAMG_alloc(A->mainBlock,1,options);
56 Esys_MPIInfo_noError(A->mpi_info);
57 } else {
58 prec->amg=Paso_Preconditioner_AMG_alloc(A,1,options);
59 }
60 if ( Esys_noError() ) {
61 if (options->verbose) {
62 if ( (prec->localamg != NULL) || ( prec->amg != NULL) ) {
63 printf("Paso_Preconditioner_AMG_Root: Smoother is ");
64 if (options->smoother == PASO_JACOBI) {
65 printf("Jacobi");
66 } else {
67 printf("Gauss-Seidel");
68 }
69 printf(" with %d/%d pre/post sweeps",options->pre_sweeps, options->post_sweeps);
70 if (options->interpolation_method == PASO_CLASSIC_INTERPOLATION) {
71 printf( " and classical interpolation.\n");
72 } else if (options->interpolation_method == PASO_CLASSIC_INTERPOLATION_WITH_FF_COUPLING) {
73 printf( " and classical interpolation with enforced FF coupling.\n");
74 } else {
75 printf( " and direct interpolation.\n");
76 }
77 } else {
78 printf("Paso_Preconditioner_AMG_Root: no coarsening constructed.\n");
79 }
80 }
81
82 if (prec->localamg != NULL) {
83 options->num_level=Paso_Preconditioner_LocalAMG_getMaxLevel(prec->localamg);
84 options->coarse_level_sparsity=Paso_Preconditioner_LocalAMG_getCoarseLevelSparsity(prec->localamg);
85 options->num_coarse_unknowns=Paso_Preconditioner_LocalAMG_getNumCoarseUnknwons(prec->localamg);
86 } else if ( prec->amg != NULL) {
87 options->num_level=Paso_Preconditioner_AMG_getMaxLevel(prec->amg);
88 options->coarse_level_sparsity=Paso_Preconditioner_AMG_getCoarseLevelSparsity(prec->amg);
89 options->num_coarse_unknowns=Paso_Preconditioner_AMG_getNumCoarseUnknwons(prec->amg);
90 } else {
91 prec->sweeps=options->sweeps;
92 prec->amgsubstitute=Paso_Preconditioner_Smoother_alloc(A, (options->smoother == PASO_JACOBI), prec->is_local, options->verbose);
93 options->num_level=0;
94 if (options->verbose) {
95 if (options->smoother == PASO_JACOBI) {
96 printf("Paso_Preconditioner: Jacobi(%d) preconditioner is used.\n",prec->sweeps);
97 } else {
98 printf("Paso_Preconditioner: Gauss-Seidel(%d) preconditioner is used.\n",prec->sweeps);
99 }
100 }
101 }
102 }
103 }
104 if (! Esys_noError() ){
105 Paso_Preconditioner_AMG_Root_free(prec);
106 return NULL;
107 } else {
108 return prec;
109 }
110 }
111
112 /* applies the preconditioner */
113 /* has to be called within a parallel reqion */
114 /* barrier synchronization is performed before the evaluation to make sure that the input vector is available */
115 void Paso_Preconditioner_AMG_Root_solve(Paso_SystemMatrix* A, Paso_Preconditioner_AMG_Root* prec, double* x,double* b)
116 {
117 if (prec->localamg != NULL) {
118 Paso_Preconditioner_LocalAMG_solve(A->mainBlock, prec->localamg,x,b);
119 } else if ( prec->amg != NULL) {
120 Paso_Preconditioner_AMG_solve(A, prec->amg,x,b);
121 } else {
122 Paso_Preconditioner_Smoother_solve(A, prec->amgsubstitute,x,b,prec->sweeps,FALSE);
123 }
124 }

  ViewVC Help
Powered by ViewVC 1.1.26