/[escript]/branches/trilinos_from_5897/escriptcore/src/SolverOptions.cpp
ViewVC logotype

Contents of /branches/trilinos_from_5897/escriptcore/src/SolverOptions.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 5914 - (show annotations)
Wed Feb 10 06:15:12 2016 UTC (2 years, 4 months ago) by caltinay
File size: 28956 byte(s)
Added support for more Trilinos preconditioners and solvers.
Need to think about how to deal with the wealth of options!

1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2016 by The University of Queensland
5 * http://www.uq.edu.au
6 *
7 * Primary Business: Queensland, Australia
8 * Licensed under the Open Software License version 3.0
9 * http://www.opensource.org/licenses/osl-3.0.php
10 *
11 * Development until 2012 by Earth Systems Science Computational Center (ESSCC)
12 * Development 2012-2013 by School of Earth Sciences
13 * Development from 2014 by Centre for Geoscience Computing (GeoComp)
14 *
15 *****************************************************************************/
16
17 #include <boost/python.hpp>
18
19 #include "SolverOptions.h"
20 #include "SolverOptionsException.h"
21
22 namespace bp = boost::python;
23
24 template <class R>
25 bool convert(bp::object bpo, R& result)
26 {
27 bool b = bp::extract<R>(bpo).check();
28 if (b)
29 result = bp::extract<R>(bpo);
30 return b;
31 }
32
33 namespace escript {
34
35 SolverBuddy::SolverBuddy() :
36 target(SO_TARGET_CPU),
37 package(SO_DEFAULT),
38 method(SO_DEFAULT),
39 preconditioner(SO_PRECONDITIONER_JACOBI),
40 ode_solver(SO_ODESOLVER_LINEAR_CRANK_NICOLSON),
41 smoother(SO_PRECONDITIONER_GAUSS_SEIDEL),
42 reordering(SO_REORDERING_DEFAULT),
43 coarsening(SO_DEFAULT),
44 amg_interpolation_method(SO_INTERPOLATION_DIRECT),
45 level_max(100),
46 coarsening_threshold(0.25),
47 sweeps(1),
48 pre_sweeps(1),
49 post_sweeps(1),
50 tolerance(1e-8),
51 absolute_tolerance(0.),
52 inner_tolerance(0.9),
53 drop_tolerance(0.01),
54 drop_storage(2.),
55 iter_max(100000),
56 inner_iter_max(10),
57 truncation(20),
58 restart(0),
59 symmetric(false),
60 verbose(false),
61 adapt_inner_tolerance(true),
62 accept_convergence_failure(false),
63 min_coarse_matrix_size(500),
64 relaxation(0.3),
65 use_local_preconditioner(false),
66 min_sparsity(0.05),
67 refinements(2),
68 coarse_refinements(2),
69 use_panel(true),
70 diagonal_dominance_threshold(0.5),
71 cycle_type(1)
72 {
73 resetDiagnostics(true);
74 }
75
76 std::string SolverBuddy::getSummary() const
77 {
78 std::stringstream out;
79 out << "Solver Package = " << getName(getPackage()) << std::endl
80 << "Solver target = " << getName(getSolverTarget()) << std::endl
81 << "Verbosity = " << isVerbose() << std::endl
82 << "Accept failed convergence = " << acceptConvergenceFailure()
83 << std::endl
84 << "Relative tolerance = " << getTolerance() << std::endl
85 << "Absolute tolerance = " << getAbsoluteTolerance() << std::endl
86 << "Symmetric problem = " << isSymmetric() << std::endl
87 << "Maximum number of iteration steps = " << getIterMax() << std::endl
88 << "Inner tolerance = " << getInnerTolerance() << std::endl
89 << "Adapt innner tolerance = " << adaptInnerTolerance() << std::endl;
90
91 if (getPackage() == SO_DEFAULT || getPackage() == SO_PACKAGE_PASO ||
92 getPackage() == SO_PACKAGE_CUSP ||
93 getPackage() == SO_PACKAGE_TRILINOS) {
94 out << "Solver method = " << getName(getSolverMethod()) << std::endl;
95 if (getSolverMethod() == SO_METHOD_GMRES) {
96 out << "Truncation = " << getTruncation() << std::endl
97 << "Restart = " << getRestart() << std::endl;
98 } else if (getSolverMethod() == SO_PRECONDITIONER_AMG) {
99 out << "Number of pre / post sweeps = " << getNumPreSweeps()
100 << " / " << getNumPostSweeps() << ", " << getNumSweeps()
101 << std::endl
102 << "Maximum number of levels = " << getLevelMax() << std::endl
103 << "Coarsening threshold = " << getCoarseningThreshold()
104 << std::endl
105 << "Coarsening method = " << getName(getCoarsening())
106 << std::endl;
107 }
108 out << "Preconditioner = " << getName(getPreconditioner()) << std::endl
109 << "Apply preconditioner locally = " << useLocalPreconditioner()
110 << std::endl;
111 switch (getPreconditioner()) {
112 case SO_PRECONDITIONER_AMG:
113 out << "Maximum number of levels = " << getLevelMax()
114 << std::endl
115 << "Coarsening threshold = " << getCoarseningThreshold()
116 << std::endl
117 << "Minimal sparsity on coarsest level = "
118 << getMinCoarseMatrixSparsity() << std::endl
119 << "Smoother = " << getName(getSmoother()) << std::endl
120 << "Minimum size of the coarsest level matrix = "
121 << getMinCoarseMatrixSize() << std::endl
122 << "Number of pre / post sweeps = " << getNumPreSweeps()
123 << " / " << getNumPostSweeps() << ", " << getNumSweeps()
124 << std::endl
125 << "Number of refinement steps in coarsest level solver = "
126 << getNumCoarseMatrixRefinements() << std::endl
127 << "Use node panel = " << usePanel() << std::endl
128 << "Interpolation = " << getName(getAMGInterpolation())
129 << std::endl
130 << "Threshold for diagonal dominant rows = "
131 << getDiagonalDominanceThreshold() << std::endl;
132 break;
133 case SO_PRECONDITIONER_AMLI:
134 out << "Maximum number of levels = " << getLevelMax()
135 << std::endl
136 << "Coarsening method = " << getName(getCoarsening())
137 << std::endl
138 << "Coarsening threshold = " << getMinCoarseMatrixSize()
139 << std::endl
140 << "Minimum size of the coarsest level matrix = "
141 << getCoarseningThreshold() << std::endl
142 << "Number of pre / post sweeps = " << getNumPreSweeps()
143 << " / " << getNumPostSweeps() << ", " << getNumSweeps()
144 << std::endl;
145 break;
146 case SO_PRECONDITIONER_BOOMERAMG:
147 out << "Maximum number of levels = " << getLevelMax()
148 << std::endl
149 << "Coarsening threshold = " << getCoarseningThreshold()
150 << std::endl
151 << "Threshold for diagonal dominant rows = "
152 << getDiagonalDominanceThreshold() << std::endl
153 << "Coarsening method = " << getName(getCoarsening())
154 << std::endl
155 << "V-cycle (1) or W-cyle (2) = " << getCycleType()
156 << std::endl
157 << "Number of pre / post sweeps = " << getNumPreSweeps()
158 << " / " << getNumPostSweeps() << ", " << getNumSweeps()
159 << std::endl
160 << "Smoother = " << getName(getSmoother()) << std::endl;
161 break;
162 case SO_PRECONDITIONER_GAUSS_SEIDEL:
163 out << "Number of sweeps = " << getNumSweeps() << std::endl;
164 break;
165 case SO_PRECONDITIONER_ILUT:
166 out << "Drop tolerance = " << getDropTolerance() << std::endl;
167 out << "Storage increase = " << getDropStorage() << std::endl;
168 break;
169 case SO_PRECONDITIONER_RILU:
170 out << "Relaxation factor = " << getRelaxationFactor()
171 << std::endl;
172 break;
173 default:
174 break;
175 } // preconditioner switch
176 out << "ODE solver = " << getName(getODESolver()) << std::endl;
177 }
178 return out.str();
179 }
180
181 const char* SolverBuddy::getName(int key) const
182 {
183 switch (static_cast<SolverOptions>(key)) {
184 case SO_DEFAULT: return "DEFAULT";
185 case SO_TARGET_CPU: return "CPU";
186 case SO_TARGET_GPU: return "GPU";
187
188 case SO_PACKAGE_CUSP: return "CUSP";
189 case SO_PACKAGE_MKL: return "MKL";
190 case SO_PACKAGE_PASO: return "PASO";
191 case SO_PACKAGE_PASTIX: return "PASTIX";
192 case SO_PACKAGE_SUPER_LU: return "SUPER_LU";
193 case SO_PACKAGE_TRILINOS: return "TRILINOS";
194 case SO_PACKAGE_UMFPACK: return "UMFPACK";
195
196 case SO_METHOD_BICGSTAB: return "BICGSTAB";
197 case SO_METHOD_CGLS: return "CGLS";
198 case SO_METHOD_CGS: return "CGS";
199 case SO_METHOD_CHOLEVSKY: return "CHOLEVSKY";
200 case SO_METHOD_CR: return "CR";
201 case SO_METHOD_DIRECT: return "DIRECT";
202 case SO_METHOD_GMRES: return "GMRES";
203 case SO_METHOD_HRZ_LUMPING: return "HRZ_LUMPING";
204 case SO_METHOD_ITERATIVE: return "ITERATIVE";
205 case SO_METHOD_LSQR: return "LSQR";
206 case SO_METHOD_MINRES: return "MINRES";
207 case SO_METHOD_NONLINEAR_GMRES: return "NONLINEAR_GMRES";
208 case SO_METHOD_PCG: return "PCG";
209 case SO_METHOD_PRES20: return "PRES20";
210 case SO_METHOD_ROWSUM_LUMPING: return "ROWSUM_LUMPING";
211 case SO_METHOD_TFQMR: return "TFQMR";
212
213 case SO_PRECONDITIONER_AMG: return "AMG";
214 case SO_PRECONDITIONER_AMLI: return "AMLI";
215 case SO_PRECONDITIONER_BOOMERAMG: return "BOOMERAMG";
216 case SO_PRECONDITIONER_GAUSS_SEIDEL: return "GAUSS_SEIDEL";
217 case SO_PRECONDITIONER_ILU0: return "ILU0";
218 case SO_PRECONDITIONER_ILUT: return "ILUT";
219 case SO_PRECONDITIONER_JACOBI: return "JACOBI";
220 case SO_PRECONDITIONER_NONE: return "NO_PRECONDITIONER";
221 case SO_PRECONDITIONER_REC_ILU: return "REC_ILU";
222 case SO_PRECONDITIONER_RILU: return "RILU";
223
224 case SO_ODESOLVER_BACKWARD_EULER: return "BACKWARD_EULER";
225 case SO_ODESOLVER_CRANK_NICOLSON: return "CRANK_NICOLSON";
226 case SO_ODESOLVER_LINEAR_CRANK_NICOLSON: return "LINEAR_CRANK_NICOLSON";
227
228 case SO_INTERPOLATION_CLASSIC: return "CLASSIC_INTERPOLATION";
229 case SO_INTERPOLATION_CLASSIC_WITH_FF_COUPLING:
230 return "CLASSIC_INTERPOLATION_WITH_FF";
231 case SO_INTERPOLATION_DIRECT: return "DIRECT_INTERPOLATION";
232
233 case SO_COARSENING_AGGREGATION: return "AGGREGATION_COARSENING";
234 case SO_COARSENING_CIJP: return "CIJP_COARSENING";
235 case SO_COARSENING_CIJP_FIXED_RANDOM:
236 return "CIJP_FIXED_RANDOM_COARSENING";
237 case SO_COARSENING_FALGOUT: return "FALGOUT_COARSENING";
238 case SO_COARSENING_HMIS: return "HMIS_COARSENING";
239 case SO_COARSENING_PMIS: return "PMIS_COARSENING";
240 case SO_COARSENING_RUGE_STUEBEN: return "RUGE_STUEBEN_COARSENING";
241 case SO_COARSENING_STANDARD: return "STANDARD_COARSENING";
242 case SO_COARSENING_YAIR_SHAPIRA: return "YAIR_SHAPIRA_COARSENING";
243
244 case SO_REORDERING_DEFAULT: return "DEFAULT_REORDERING";
245 case SO_REORDERING_MINIMUM_FILL_IN: return "MINIMUM_FILL_IN";
246 case SO_REORDERING_NESTED_DISSECTION: return "NESTED_DISSECTION";
247 case SO_REORDERING_NONE: return "NO_REORDERING";
248 default:
249 throw SolverOptionsException("getName() invalid option given");
250 }
251 return "invalid option";
252 }
253
254 void SolverBuddy::resetDiagnostics(bool all)
255 {
256 num_iter = 0;
257 num_level = 0;
258 num_inner_iter = 0;
259 time = 0.;
260 set_up_time = 0.;
261 net_time = 0.;
262 residual_norm = 0.;
263 converged = false;
264 preconditioner_size = -1;
265 time_step_backtracking_used = false;
266 coarse_level_sparsity = 0;
267 num_coarse_unknowns = 0;
268 if (all) {
269 cum_num_inner_iter = 0;
270 cum_num_iter = 0;
271 cum_time = 0.;
272 cum_set_up_time = 0.;
273 cum_net_time = 0.;
274 }
275 }
276
277 void SolverBuddy::updateDiagnostics(std::string name, const bp::object& value)
278 {
279
280 int i=0;
281 double d=0; // to keep older compilers happy
282 bool b=false;
283 bool ib = convert<int>(value, i);
284 bool db = convert<double>(value, d);
285 bool bb = convert<bool>(value, b);
286
287 if (name == "num_iter") {
288 if (!ib)
289 throw SolverOptionsException("setting num_iter to non-int value");
290 cum_num_iter += num_iter = i;
291 } else if (name == "num_level") {
292 if (!ib)
293 throw SolverOptionsException("setting num_level to non-int value");
294 num_level = i;
295 } else if (name == "num_inner_iter") {
296 if (!ib)
297 throw SolverOptionsException("setting num_inner_iter to non-int value");
298 cum_num_inner_iter += num_inner_iter = i;
299 } else if (name == "time") {
300 if (!db)
301 throw SolverOptionsException("setting time to non-double value");
302 cum_time += time = d;
303 } else if (name == "set_up_time") {
304 if (!db)
305 throw SolverOptionsException("setting set_up_time to non-double value");
306 cum_set_up_time += set_up_time = d;
307 } else if (name == "net_time") {
308 if (!db)
309 throw SolverOptionsException("setting net_time to non-double value");
310 cum_net_time += net_time = d;
311 } else if (name == "residual_norm") {
312 if (!db)
313 throw SolverOptionsException("setting residual_norm to non-double value");
314 residual_norm = d;
315 } else if (name == "converged") {
316 if (!bb)
317 throw SolverOptionsException("setting converged to non-bool value");
318 converged = b;
319 } else if (name == "time_step_backtracking_used") {
320 if (!bb)
321 throw SolverOptionsException("setting time_step_backtracking_used to non-bool value");
322 time_step_backtracking_used = b;
323 } else if (name == "coarse_level_sparsity") {
324 if (!db)
325 throw SolverOptionsException("setting coarse_level_sparsity to non-double value");
326 coarse_level_sparsity = d;
327 } else if (name == "num_coarse_unknowns") {
328 if (!ib)
329 throw SolverOptionsException("setting num_coarse_unknowns to non-int value");
330 num_coarse_unknowns = i;
331 } else {
332 throw SolverOptionsException(std::string("Unknown diagnostic: ") + name);
333 }
334 }
335
336 double SolverBuddy::getDiagnostics(const std::string name) const
337 {
338 if (name == "num_iter") return num_iter;
339 else if (name == "cum_num_iter") return cum_num_iter;
340 else if (name == "num_level") return num_level;
341 else if (name == "num_inner_iter") return num_inner_iter;
342 else if (name == "cum_num_inner_iter") return cum_num_inner_iter;
343 else if (name == "time") return time;
344 else if (name == "cum_time") return cum_time;
345 else if (name == "set_up_time") return set_up_time;
346 else if (name == "cum_set_up_time") return cum_set_up_time;
347 else if (name == "net_time") return net_time;
348 else if (name == "cum_net_time") return cum_net_time;
349 else if (name == "residual_norm") return residual_norm;
350 else if (name == "converged") return converged;
351 else if (name == "preconditioner_size") return preconditioner_size;
352 else if (name == "time_step_backtracking_used")
353 return time_step_backtracking_used;
354 else if (name == "coarse_level_sparsity") return coarse_level_sparsity;
355 else if (name == "num_coarse_unknowns") return num_coarse_unknowns;
356
357 throw SolverOptionsException(std::string("unknown diagnostic item: ")
358 + name);
359 }
360
361 bool SolverBuddy::hasConverged() const
362 {
363 return converged;
364 }
365
366 void SolverBuddy::setCoarsening(int method)
367 {
368 SolverOptions meth = static_cast<SolverOptions>(method);
369 switch (meth) {
370 case SO_DEFAULT:
371 case SO_COARSENING_AGGREGATION:
372 case SO_COARSENING_CIJP:
373 case SO_COARSENING_CIJP_FIXED_RANDOM:
374 case SO_COARSENING_FALGOUT:
375 case SO_COARSENING_HMIS:
376 case SO_COARSENING_PMIS:
377 case SO_COARSENING_RUGE_STUEBEN:
378 case SO_COARSENING_STANDARD:
379 case SO_COARSENING_YAIR_SHAPIRA:
380 coarsening = meth;
381 break;
382 default:
383 throw SolverOptionsException("unknown coarsening method");
384 }
385 }
386
387 SolverOptions SolverBuddy::getCoarsening() const
388 {
389 return coarsening;
390 }
391
392 void SolverBuddy::setMinCoarseMatrixSize(int size)
393 {
394 if (size < 0) {
395 throw SolverOptionsException("minimum size of the coarsest level "
396 "matrix must be non-negative.");
397 }
398 min_coarse_matrix_size = size;
399 }
400
401 int SolverBuddy::getMinCoarseMatrixSize() const
402 {
403 return min_coarse_matrix_size;
404 }
405
406 void SolverBuddy::setPreconditioner(int precon)
407 {
408 SolverOptions preconditioner = static_cast<SolverOptions>(precon);
409 switch(preconditioner) {
410 case SO_PRECONDITIONER_AMG:
411 /*
412 #ifdef ESYS_MPI
413 throw SolverOptionsException("AMG preconditioner is not supported in MPI builds");
414 break;
415 #endif
416 */
417 case SO_PRECONDITIONER_AMLI:
418 case SO_PRECONDITIONER_BOOMERAMG:
419 case SO_PRECONDITIONER_GAUSS_SEIDEL:
420 case SO_PRECONDITIONER_ILU0:
421 case SO_PRECONDITIONER_ILUT:
422 case SO_PRECONDITIONER_JACOBI:
423 case SO_PRECONDITIONER_NONE:
424 case SO_PRECONDITIONER_REC_ILU:
425 case SO_PRECONDITIONER_RILU:
426 this->preconditioner = preconditioner;
427 break;
428 default:
429 throw SolverOptionsException("unknown preconditioner");
430 }
431 }
432
433 SolverOptions SolverBuddy::getPreconditioner() const
434 {
435 return preconditioner;
436 }
437
438 void SolverBuddy::setSmoother(int s)
439 {
440 SolverOptions smoother = static_cast<SolverOptions>(s);
441 if (smoother != SO_PRECONDITIONER_JACOBI &&
442 smoother != SO_PRECONDITIONER_GAUSS_SEIDEL) {
443 throw SolverOptionsException("unknown smoother");
444 }
445 this->smoother = smoother;
446 }
447
448 SolverOptions SolverBuddy::getSmoother() const
449 {
450 return smoother;
451 }
452
453 void SolverBuddy::setSolverMethod(int method)
454 {
455 SolverOptions meth = static_cast<SolverOptions>(method);
456 switch(meth) {
457 case SO_DEFAULT:
458 case SO_METHOD_BICGSTAB:
459 case SO_METHOD_CGLS:
460 case SO_METHOD_CGS:
461 case SO_METHOD_CHOLEVSKY:
462 case SO_METHOD_CR:
463 case SO_METHOD_GMRES:
464 case SO_METHOD_HRZ_LUMPING:
465 case SO_METHOD_ITERATIVE:
466 case SO_METHOD_LSQR:
467 case SO_METHOD_MINRES:
468 case SO_METHOD_NONLINEAR_GMRES:
469 case SO_METHOD_PCG:
470 case SO_METHOD_PRES20:
471 case SO_METHOD_ROWSUM_LUMPING:
472 case SO_METHOD_TFQMR:
473 this->method = meth;
474 break;
475 case SO_METHOD_DIRECT:
476 #ifdef USE_UMFPACK
477 this->method = meth;
478 break;
479 #elif defined MKL
480 this->method = meth;
481 break;
482 #elif defined PASTIX
483 this->method = meth;
484 break;
485 #else
486 throw SolverOptionsException("Cannot use DIRECT solver method, the running escript was not compiled with a direct solver enabled");
487 #endif
488 default:
489 throw SolverOptionsException("unknown solver method");
490 }
491 }
492
493 SolverOptions SolverBuddy::getSolverMethod() const
494 {
495 return method;
496 }
497
498 void SolverBuddy::setSolverTarget(int target)
499 {
500 SolverOptions targ = static_cast<SolverOptions>(target);
501 switch (targ) {
502 case SO_TARGET_CPU:
503 case SO_TARGET_GPU:
504 this->target = targ;
505 break;
506 default:
507 throw SolverOptionsException("unknown solver target");
508 }
509 }
510
511 SolverOptions SolverBuddy::getSolverTarget() const
512 {
513 return target;
514 }
515
516 void SolverBuddy::setPackage(int package)
517 {
518 SolverOptions pack = static_cast<SolverOptions>(package);
519 switch (pack) {
520 case SO_DEFAULT:
521 case SO_PACKAGE_CUSP:
522 case SO_PACKAGE_MKL:
523 case SO_PACKAGE_PASO:
524 case SO_PACKAGE_PASTIX:
525 case SO_PACKAGE_SUPER_LU:
526 case SO_PACKAGE_TRILINOS:
527 case SO_PACKAGE_UMFPACK:
528 this->package = pack;
529 break;
530 default:
531 throw SolverOptionsException("unknown solver package");
532 }
533 }
534
535 SolverOptions SolverBuddy::getPackage() const
536 {
537 return package;
538 }
539
540 void SolverBuddy::setReordering(int ordering)
541 {
542 SolverOptions ord = static_cast<SolverOptions>(ordering);
543 switch (ordering) {
544 case SO_REORDERING_DEFAULT:
545 case SO_REORDERING_MINIMUM_FILL_IN:
546 case SO_REORDERING_NESTED_DISSECTION:
547 case SO_REORDERING_NONE:
548 reordering = ord;
549 break;
550 default:
551 throw SolverOptionsException("unknown reordering strategy");
552 }
553 }
554
555 SolverOptions SolverBuddy::getReordering() const
556 {
557 return reordering;
558 }
559
560 void SolverBuddy::setRestart(int restart)
561 {
562 if (restart < 0)
563 throw SolverOptionsException("restart must be non-negative.");
564
565 this->restart = restart;
566 }
567
568 int SolverBuddy::getRestart() const
569 {
570 return restart;
571 }
572
573 int SolverBuddy::_getRestartForC() const
574 {
575 int r = getRestart();
576 if (r == 0)
577 return -1;
578 return r;
579 }
580
581 void SolverBuddy::setDiagonalDominanceThreshold(double value)
582 {
583 if (value < 0. || value > 1.)
584 throw SolverOptionsException("Diagonal dominance threshold must be between 0 and 1.");
585 diagonal_dominance_threshold = value;
586 }
587
588 double SolverBuddy::getDiagonalDominanceThreshold() const
589 {
590 return diagonal_dominance_threshold;
591 }
592
593 void SolverBuddy::setTruncation(int truncation)
594 {
595 if (truncation < 1)
596 throw SolverOptionsException("truncation must be positive.");
597 this->truncation = truncation;
598 }
599
600 int SolverBuddy::getTruncation() const
601 {
602 return truncation;
603 }
604
605 void SolverBuddy::setInnerIterMax(int iter_max)
606 {
607 if (iter_max < 1)
608 throw SolverOptionsException("maximum number of inner iteration must be positive.");
609 inner_iter_max = iter_max;
610 }
611
612 int SolverBuddy::getInnerIterMax() const
613 {
614 return inner_iter_max;
615 }
616
617 void SolverBuddy::setIterMax(int iter_max)
618 {
619 if (iter_max < 1)
620 throw SolverOptionsException("maximum number of iteration steps must be positive.");
621 this->iter_max = iter_max;
622 }
623
624 int SolverBuddy::getIterMax() const
625 {
626 return iter_max;
627 }
628
629 void SolverBuddy::setLevelMax(int level_max)
630 {
631 if (level_max < 0)
632 throw SolverOptionsException("maximum number of coarsening levels must be non-negative.");
633 this->level_max = level_max;
634 }
635
636 int SolverBuddy::getLevelMax() const
637 {
638 return level_max;
639 }
640
641 void SolverBuddy::setCycleType(int cycle_type)
642 {
643 this->cycle_type = cycle_type;
644 }
645
646 int SolverBuddy::getCycleType() const
647 {
648 return cycle_type;
649 }
650
651 void SolverBuddy::setCoarseningThreshold(double theta)
652 {
653 if (theta < 0. || theta > 1.)
654 throw SolverOptionsException("threshold must be between 0 and 1.");
655 coarsening_threshold = theta;
656 }
657
658 double SolverBuddy::getCoarseningThreshold() const
659 {
660 return coarsening_threshold;
661 }
662
663 void SolverBuddy::setNumSweeps(int sweeps)
664 {
665 if (sweeps < 1)
666 throw SolverOptionsException("number of sweeps must be positive.");
667 this->sweeps = sweeps;
668 }
669
670 int SolverBuddy::getNumSweeps() const
671 {
672 return sweeps;
673 }
674
675 void SolverBuddy::setNumPreSweeps(int sweeps)
676 {
677 if (sweeps < 1)
678 throw SolverOptionsException("number of pre-sweeps must be positive.");
679 pre_sweeps = sweeps;
680 }
681
682 int SolverBuddy::getNumPreSweeps() const
683 {
684 return pre_sweeps;
685 }
686
687 void SolverBuddy::setNumPostSweeps(int sweeps)
688 {
689 if (sweeps < 1)
690 throw SolverOptionsException("number of post-sweeps must be positive.");
691 post_sweeps = sweeps;
692 }
693
694 int SolverBuddy::getNumPostSweeps() const
695 {
696 return post_sweeps;
697 }
698
699 void SolverBuddy::setTolerance(double rtol)
700 {
701 if (rtol < 0. || rtol > 1.)
702 throw SolverOptionsException("tolerance must be between 0 and 1.");
703 tolerance = rtol;
704 }
705
706 double SolverBuddy::getTolerance() const
707 {
708 return tolerance;
709 }
710
711 void SolverBuddy::setAbsoluteTolerance(double atol)
712 {
713 if (atol < 0.)
714 throw SolverOptionsException("absolute tolerance must be non-negative.");
715 absolute_tolerance = atol;
716 }
717
718 double SolverBuddy::getAbsoluteTolerance() const
719 {
720 return absolute_tolerance;
721 }
722
723 void SolverBuddy::setInnerTolerance(double rtol)
724 {
725 if (rtol <= 0. || rtol > 1.)
726 throw SolverOptionsException("tolerance must be positive and less than or equal to 1.");
727 inner_tolerance = rtol;
728 }
729
730 double SolverBuddy::getInnerTolerance() const
731 {
732 return inner_tolerance;
733 }
734
735 void SolverBuddy::setDropTolerance(double drop_tol)
736 {
737 if (drop_tol < 0. || drop_tol > 1.)
738 throw SolverOptionsException("drop tolerance must be between 0 and 1.");
739 drop_tolerance = drop_tol;
740 }
741
742 double SolverBuddy::getDropTolerance() const
743 {
744 return drop_tolerance;
745 }
746
747 void SolverBuddy::setDropStorage(double storage)
748 {
749 if (storage < 1.)
750 throw SolverOptionsException("allowed storage increase must be greater than or equal to 1.");
751 drop_storage = storage;
752 }
753
754 double SolverBuddy::getDropStorage() const
755 {
756 return drop_storage;
757 }
758
759 void SolverBuddy::setRelaxationFactor(double factor)
760 {
761 if (factor < 0.)
762 throw SolverOptionsException("relaxation factor must be non-negative.");
763 relaxation = factor;
764 }
765
766 double SolverBuddy::getRelaxationFactor() const
767 {
768 return relaxation;
769 }
770
771 bool SolverBuddy::isSymmetric() const
772 {
773 return symmetric;
774 }
775
776 void SolverBuddy::setSymmetryOn()
777 {
778 symmetric = true;
779 }
780
781 void SolverBuddy::setSymmetryOff()
782 {
783 symmetric = false;
784 }
785
786 void SolverBuddy::setSymmetry(bool flag)
787 {
788 if (flag)
789 setSymmetryOn();
790 else
791 setSymmetryOff();
792 }
793
794 bool SolverBuddy::isVerbose() const
795 {
796 return verbose;
797 }
798
799 void SolverBuddy::setVerbosityOn()
800 {
801 verbose = true;
802 }
803
804 void SolverBuddy::setVerbosityOff()
805 {
806 verbose = false;
807 }
808
809 void SolverBuddy::setVerbosity(bool verbose)
810 {
811 if (verbose)
812 setVerbosityOn();
813 else
814 setVerbosityOff();
815 }
816
817 bool SolverBuddy::adaptInnerTolerance() const
818 {
819 return adapt_inner_tolerance;
820 }
821
822 void SolverBuddy::setInnerToleranceAdaptionOn()
823 {
824 adapt_inner_tolerance = true;
825 }
826
827 void SolverBuddy::setInnerToleranceAdaptionOff()
828 {
829 adapt_inner_tolerance = false;
830 }
831
832 void SolverBuddy::setInnerToleranceAdaption(bool adapt)
833 {
834 if (adapt)
835 setInnerToleranceAdaptionOn();
836 else
837 setInnerToleranceAdaptionOff();
838 }
839
840 bool SolverBuddy::acceptConvergenceFailure() const
841 {
842 return accept_convergence_failure;
843 }
844
845 void SolverBuddy::setAcceptanceConvergenceFailureOn()
846 {
847 accept_convergence_failure = true;
848 }
849
850 void SolverBuddy::setAcceptanceConvergenceFailureOff()
851 {
852 accept_convergence_failure = false;
853 }
854
855 void SolverBuddy::setAcceptanceConvergenceFailure(bool accept)
856 {
857 if (accept)
858 setAcceptanceConvergenceFailureOn();
859 else
860 setAcceptanceConvergenceFailureOff();
861 }
862
863 bool SolverBuddy::useLocalPreconditioner() const
864 {
865 return use_local_preconditioner;
866 }
867
868 void SolverBuddy::setLocalPreconditionerOn()
869 {
870 use_local_preconditioner = true;
871 }
872
873 void SolverBuddy::setLocalPreconditionerOff()
874 {
875 use_local_preconditioner=false;
876 }
877
878 void SolverBuddy::setLocalPreconditioner(bool use)
879 {
880 if (use)
881 setLocalPreconditionerOn();
882 else
883 setLocalPreconditionerOff();
884 }
885
886 void SolverBuddy::setMinCoarseMatrixSparsity(double sparsity)
887 {
888 if (sparsity < 0. || sparsity > 1.)
889 throw SolverOptionsException("sparsity must be between 0 and 1.");
890 min_sparsity = sparsity;
891 }
892
893 double SolverBuddy::getMinCoarseMatrixSparsity() const
894 {
895 return min_sparsity;
896 }
897
898 void SolverBuddy::setNumRefinements(int refinements)
899 {
900 if (refinements < 0)
901 throw SolverOptionsException("number of refinements must be non-negative.");
902 this->refinements = refinements;
903 }
904
905 int SolverBuddy::getNumRefinements() const
906 {
907 return refinements;
908 }
909
910 void SolverBuddy::setNumCoarseMatrixRefinements(int refinements)
911 {
912 if (refinements < 0)
913 throw SolverOptionsException("number of coarse matrix refinements must be non-negative.");
914 coarse_refinements = refinements;
915 }
916
917 int SolverBuddy::getNumCoarseMatrixRefinements() const
918 {
919 return coarse_refinements;
920 }
921
922 bool SolverBuddy::usePanel() const
923 {
924 return use_panel;
925 }
926
927 void SolverBuddy::setUsePanelOn()
928 {
929 use_panel = true;
930 }
931
932 void SolverBuddy::setUsePanelOff()
933 {
934 use_panel = false;
935 }
936
937 void SolverBuddy::setUsePanel(bool use)
938 {
939 if (use)
940 setUsePanelOn();
941 else
942 setUsePanelOff();
943 }
944
945 void SolverBuddy::setAMGInterpolation(int method)
946 {
947 SolverOptions meth = static_cast<SolverOptions>(method);
948 switch (meth) {
949 case SO_INTERPOLATION_CLASSIC:
950 case SO_INTERPOLATION_CLASSIC_WITH_FF_COUPLING:
951 case SO_INTERPOLATION_DIRECT:
952 amg_interpolation_method = meth;
953 break;
954 default:
955 throw SolverOptionsException("unknown AMG interpolation method");
956 }
957 }
958
959 SolverOptions SolverBuddy::getAMGInterpolation() const
960 {
961 return amg_interpolation_method;
962 }
963
964 void SolverBuddy::setODESolver(int method)
965 {
966 SolverOptions ode = static_cast<SolverOptions>(method);
967 switch (ode) {
968 case SO_ODESOLVER_BACKWARD_EULER:
969 case SO_ODESOLVER_CRANK_NICOLSON:
970 case SO_ODESOLVER_LINEAR_CRANK_NICOLSON:
971 ode_solver = ode;
972 break;
973 default:
974 throw SolverOptionsException("unknown ODE solver method");
975 }
976 }
977
978 SolverOptions SolverBuddy::getODESolver() const
979 {
980 return ode_solver;
981 }
982
983 } // namespace escript
984

  ViewVC Help
Powered by ViewVC 1.1.26