escript
Revision_
|
00001 00002 /***************************************************************************** 00003 * 00004 * Copyright (c) 2003-2014 by University of Queensland 00005 * http://www.uq.edu.au 00006 * 00007 * Primary Business: Queensland, Australia 00008 * Licensed under the Open Software License version 3.0 00009 * http://www.opensource.org/licenses/osl-3.0.php 00010 * 00011 * Development until 2012 by Earth Systems Science Computational Center (ESSCC) 00012 * Development 2012-2013 by School of Earth Sciences 00013 * Development from 2014 by Centre for Geoscience Computing (GeoComp) 00014 * 00015 *****************************************************************************/ 00016 00017 #ifndef __PASO_PRECONDITIONER_H__ 00018 #define __PASO_PRECONDITIONER_H__ 00019 00020 #include "SystemMatrix.h" 00021 #include "BOOMERAMG.h" 00022 00023 #define PRECONDITIONER_NO_ERROR 0 00024 #define PRECONDITIONER_MAXITER_REACHED 1 00025 #define PRECONDITIONER_INPUT_ERROR -1 00026 #define PRECONDITIONER_MEMORY_ERROR -9 00027 #define PRECONDITIONER_BREAKDOWN -10 00028 #define PRECONDITIONER_NEGATIVE_NORM_ERROR -11 00029 #define PRECONDITIONER_DIVERGENCE -12 00030 00031 namespace paso { 00032 00033 struct MergedSolver; 00034 struct Preconditioner; 00035 typedef boost::shared_ptr<Preconditioner> Preconditioner_ptr; 00036 typedef boost::shared_ptr<const Preconditioner> const_Preconditioner_ptr; 00037 00038 struct Preconditioner_Smoother; 00039 struct Preconditioner_AMG_Root; 00040 struct Solver_ILU; 00041 struct Solver_RILU; 00042 00043 // general preconditioner interface 00044 struct Preconditioner 00045 { 00046 dim_t type; 00047 dim_t sweeps; 00049 Preconditioner_Smoother* jacobi; 00051 Preconditioner_Smoother* gs; 00053 Preconditioner_AMG_Root* amg; 00055 Solver_ILU* ilu; 00057 Solver_RILU* rilu; 00058 }; 00059 00060 void Preconditioner_free(Preconditioner*); 00061 Preconditioner* Preconditioner_alloc(SystemMatrix_ptr A, Options* options); 00062 void Preconditioner_solve(Preconditioner* prec, SystemMatrix_ptr A, double*, double*); 00063 00064 00065 // GAUSS SEIDEL & Jacobi 00066 struct Preconditioner_LocalSmoother 00067 { 00068 bool Jacobi; 00069 double* diag; 00070 double* buffer; 00071 index_t* pivot; 00072 }; 00073 00074 struct Preconditioner_Smoother 00075 { 00076 Preconditioner_LocalSmoother* localSmoother; 00077 bool is_local; 00078 }; 00079 00080 void Preconditioner_Smoother_free(Preconditioner_Smoother * in); 00081 void Preconditioner_LocalSmoother_free(Preconditioner_LocalSmoother * in); 00082 00083 Preconditioner_Smoother* Preconditioner_Smoother_alloc( 00084 SystemMatrix_ptr A, bool jacobi, bool is_local, bool verbose); 00085 00086 Preconditioner_LocalSmoother* Preconditioner_LocalSmoother_alloc( 00087 SparseMatrix_ptr A, bool jacobi, bool verbose); 00088 00089 void Preconditioner_Smoother_solve(SystemMatrix_ptr A, 00090 Preconditioner_Smoother* gs, double* x, const double* b, 00091 dim_t sweeps, bool x_is_initial); 00092 00093 void Preconditioner_LocalSmoother_solve(SparseMatrix_ptr A, 00094 Preconditioner_LocalSmoother* gs, double* x, const double* b, 00095 dim_t sweeps, bool x_is_initial); 00096 00097 err_t Preconditioner_Smoother_solve_byTolerance(SystemMatrix_ptr A, 00098 Preconditioner_Smoother* gs, double* x, const double* b, 00099 double atol, dim_t* sweeps, bool x_is_initial); 00100 00101 void Preconditioner_LocalSmoother_Sweep(SparseMatrix_ptr A, 00102 Preconditioner_LocalSmoother* gs, double* x); 00103 00104 void Preconditioner_LocalSmoother_Sweep_sequential( 00105 SparseMatrix_ptr A, Preconditioner_LocalSmoother* gs, 00106 double* x); 00107 00108 void Preconditioner_LocalSmoother_Sweep_tiled(SparseMatrix_ptr A, 00109 Preconditioner_LocalSmoother* gs, double* x); 00110 00111 void Preconditioner_LocalSmoother_Sweep_colored(SparseMatrix_ptr A, 00112 Preconditioner_LocalSmoother* gs, double* x); 00113 00114 00115 typedef enum 00116 { 00117 PASO_AMG_UNDECIDED=-1, 00118 PASO_AMG_IN_F=0, 00119 PASO_AMG_IN_C=1 00120 } AMGBlockSelect; 00121 00123 struct Preconditioner_AMG 00124 { 00125 dim_t level; 00127 SystemMatrix_ptr A_C; 00129 SystemMatrix_ptr P; 00131 SystemMatrix_ptr R; 00132 00133 Preconditioner_Smoother* Smoother; 00134 dim_t post_sweeps; 00135 dim_t pre_sweeps; 00137 dim_t options_smoother; 00139 bool verbose; 00141 index_t reordering; 00143 dim_t refinements; 00145 double* r; 00147 double* x_C; 00149 double* b_C; 00151 MergedSolver* merged_solver; 00152 Preconditioner_AMG* AMG_C; 00153 }; 00154 00155 void Preconditioner_AMG_free(Preconditioner_AMG* in); 00156 Preconditioner_AMG* Preconditioner_AMG_alloc(SystemMatrix_ptr A, dim_t level, Options* options); 00157 void Preconditioner_AMG_solve(SystemMatrix_ptr A, Preconditioner_AMG* amg, double* x, double* b); 00158 void Preconditioner_AMG_setStrongConnections(SystemMatrix_ptr A, dim_t* degree_S, index_t* offset_S, index_t* S, double theta, double tau); 00159 void Preconditioner_AMG_setStrongConnections_Block(SystemMatrix_ptr A, dim_t* degree_S, index_t* offset_S, index_t* S, double theta, double tau); 00160 SystemMatrix_ptr Preconditioner_AMG_getProlongation(SystemMatrix_ptr A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, dim_t n_C, index_t* counter_C, index_t interpolation_method); 00161 void Preconditioner_AMG_setClassicProlongation(SystemMatrix_ptr P, SystemMatrix_ptr A, const index_t* offset_S, const dim_t* degree_S, const index_t* S,const index_t *counter_C); 00162 void Preconditioner_AMG_setClassicProlongation_Block(SystemMatrix_ptr P, SystemMatrix_ptr A, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t* counter_C); 00163 void Preconditioner_AMG_setDirectProlongation(SystemMatrix_ptr P, SystemMatrix_ptr A, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t* counter_C); 00164 void Preconditioner_AMG_setDirectProlongation_Block(SystemMatrix_ptr P, SystemMatrix_ptr A, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t* counter_C); 00165 double Preconditioner_AMG_getCoarseLevelSparsity(const Preconditioner_AMG* in); 00166 dim_t Preconditioner_AMG_getNumCoarseUnknowns(const Preconditioner_AMG* in); 00167 index_t Preconditioner_AMG_getMaxLevel(const Preconditioner_AMG* in); 00168 void Preconditioner_AMG_transposeStrongConnections(dim_t n, const dim_t* degree_S, const index_t* offset_S, const index_t* S, dim_t nT, dim_t* degree_ST, index_t* offset_ST, index_t* ST); 00169 void Preconditioner_AMG_CIJPCoarsening(dim_t n, dim_t my_n, 00170 AMGBlockSelect* split_marker, const dim_t* degree_S, 00171 const index_t* offset_S, const index_t* S, const dim_t* degree_ST, 00172 const index_t* offset_ST, const index_t* ST, 00173 Connector_ptr col_connector, const_Distribution_ptr col_dist); 00174 00175 SystemMatrix_ptr Preconditioner_AMG_getRestriction(SystemMatrix_ptr P); 00176 00177 SystemMatrix_ptr Preconditioner_AMG_buildInterpolationOperator( 00178 SystemMatrix_ptr A, SystemMatrix_ptr P, SystemMatrix_ptr R); 00179 00180 SystemMatrix_ptr Preconditioner_AMG_buildInterpolationOperatorBlock( 00181 SystemMatrix_ptr A, SystemMatrix_ptr P, SystemMatrix_ptr R); 00182 00183 SparseMatrix_ptr Preconditioner_AMG_mergeSystemMatrix(SystemMatrix_ptr A); 00184 00185 void Preconditioner_AMG_mergeSolve(Preconditioner_AMG* amg); 00186 00188 struct Preconditioner_LocalAMG 00189 { 00190 dim_t level; 00191 SparseMatrix_ptr A_C; /* coarse level matrix */ 00192 SparseMatrix_ptr P; /* prolongation n x n_C*/ 00193 SparseMatrix_ptr R; /* restriction n_C x n */ 00194 00195 Preconditioner_LocalSmoother* Smoother; 00196 dim_t post_sweeps; 00197 dim_t pre_sweeps; 00198 index_t reordering; /* applied reordering in direct solver */ 00199 dim_t refinements; /* number of refinements in direct solver (typically =0) */ 00200 double* r; /* buffer for residual */ 00201 double* x_C; /* solution of coarse level system */ 00202 double* b_C; /* right hand side of coarse level system */ 00203 struct Preconditioner_LocalAMG * AMG_C; 00204 }; 00205 00206 void Preconditioner_LocalAMG_free(Preconditioner_LocalAMG * in); 00207 Preconditioner_LocalAMG* Preconditioner_LocalAMG_alloc(SparseMatrix_ptr A, dim_t level, Options* options); 00208 void Preconditioner_LocalAMG_solve(SparseMatrix_ptr A, Preconditioner_LocalAMG * amg, double * x, double * b); 00209 00210 void Preconditioner_LocalAMG_RungeStuebenSearch(const dim_t n, const index_t* offset, const dim_t* degree, const index_t* S, AMGBlockSelect*split_marker, const bool usePanel); 00211 void Preconditioner_LocalAMG_setStrongConnections_Block(SparseMatrix_ptr A, dim_t *degree, index_t *S, const double theta, const double tau); 00212 void Preconditioner_LocalAMG_setStrongConnections(SparseMatrix_ptr A, dim_t *degree, index_t *S, const double theta, const double tau); 00213 00214 SparseMatrix_ptr Preconditioner_LocalAMG_getProlongation( 00215 SparseMatrix_ptr A_p, const index_t* offset_S, 00216 const dim_t* degree_S, const index_t* S, dim_t n_C, 00217 const index_t* counter_C, index_t interpolation_method); 00218 00219 void Preconditioner_LocalAMG_setDirectProlongation_Block(SparseMatrix_ptr P_p, const_SparseMatrix_ptr A_p, const index_t *counter_C); 00220 00221 void Preconditioner_LocalAMG_setDirectProlongation(SparseMatrix_ptr P_p, const_SparseMatrix_ptr A_p, const index_t *counter_C); 00222 void Preconditioner_LocalAMG_setClassicProlongation(SparseMatrix_ptr P_p, SparseMatrix_ptr A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t *counter_C); 00223 void Preconditioner_LocalAMG_setClassicProlongation_Block(SparseMatrix_ptr P_p, SparseMatrix_ptr A_p, const index_t* offset_S, const dim_t* degree_S, const index_t* S, const index_t *counter_C); 00224 index_t Preconditioner_LocalAMG_getMaxLevel(const Preconditioner_LocalAMG * in); 00225 double Preconditioner_LocalAMG_getCoarseLevelSparsity(const Preconditioner_LocalAMG * in); 00226 dim_t Preconditioner_LocalAMG_getNumCoarseUnknowns(const Preconditioner_LocalAMG * in); 00227 void Preconditioner_LocalAMG_enforceFFConnectivity(const dim_t n, const index_t* offset_S, const dim_t* degree_S, const index_t* S, AMGBlockSelect*split_marker); 00228 00229 00230 struct Preconditioner_AMG_Root 00231 { 00232 bool is_local; 00233 Preconditioner_AMG* amg; 00234 Preconditioner_LocalAMG* localamg; 00235 Preconditioner_BoomerAMG* boomeramg; 00236 dim_t sweeps; 00237 Preconditioner_Smoother* amgsubstitute; 00238 }; 00239 00240 Preconditioner_AMG_Root* Preconditioner_AMG_Root_alloc(SystemMatrix_ptr A, Options* options); 00241 void Preconditioner_AMG_Root_free(Preconditioner_AMG_Root * in); 00242 void Preconditioner_AMG_Root_solve(SystemMatrix_ptr A, Preconditioner_AMG_Root * amg, double * x, double * b); 00243 00245 struct Solver_ILU 00246 { 00247 double* factors; 00248 }; 00249 00251 struct Solver_RILU 00252 { 00253 dim_t n; 00254 dim_t n_block; 00255 dim_t n_F; 00256 dim_t n_C; 00257 double* inv_A_FF; 00258 index_t* A_FF_pivot; 00259 SparseMatrix_ptr A_FC; 00260 SparseMatrix_ptr A_CF; 00261 index_t* rows_in_F; 00262 index_t* rows_in_C; 00263 index_t* mask_F; 00264 index_t* mask_C; 00265 double* x_F; 00266 double* b_F; 00267 double* x_C; 00268 double* b_C; 00269 Solver_RILU* RILU_of_Schur; 00270 }; 00271 00272 void Solver_ILU_free(Solver_ILU * in); 00273 Solver_ILU* Solver_getILU(SparseMatrix_ptr A, bool verbose); 00274 void Solver_solveILU(SparseMatrix_ptr A, Solver_ILU* ilu, double* x, const double* b); 00275 00276 void Solver_RILU_free(Solver_RILU* in); 00277 Solver_RILU* Solver_getRILU(SparseMatrix_ptr A, bool verbose); 00278 void Solver_solveRILU(Solver_RILU* rilu, double* x, double* b); 00279 00280 void Solver_updateIncompleteSchurComplement(SparseMatrix_ptr A_CC, 00281 SparseMatrix_ptr A_CF, double* invA_FF, index_t* A_FF_pivot, 00282 SparseMatrix_ptr A_FC); 00283 00284 } // namespace paso 00285 00286 #endif // __PASO_PRECONDITIONER_H__ 00287