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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 5504 - (show annotations)
Wed Mar 4 22:58:13 2015 UTC (4 years, 1 month ago) by jfenwick
File size: 28524 byte(s)
Again with a more up to date copy


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

  ViewVC Help
Powered by ViewVC 1.1.26