/[escript]/trunk/escriptcore/src/SolverOptions.cpp
ViewVC logotype

Contents of /trunk/escriptcore/src/SolverOptions.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6144 - (show annotations)
Wed Apr 6 05:25:13 2016 UTC (2 years, 7 months ago) by caltinay
File size: 28488 byte(s)
last round of namespacing defines.

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 Apache License, version 2.0
9 * http://www.apache.org/licenses/LICENSE-2.0
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 "SolverOptions.h"
18 #include "EsysException.h"
19
20 #include <boost/python.hpp>
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 ValueError("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 ValueError("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 ValueError("setting num_level to non-int value");
294 num_level = i;
295 } else if (name == "num_inner_iter") {
296 if (!ib)
297 throw ValueError("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 ValueError("setting time to non-double value");
302 cum_time += time = d;
303 } else if (name == "set_up_time") {
304 if (!db)
305 throw ValueError("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 ValueError("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 ValueError("setting residual_norm to non-double value");
314 residual_norm = d;
315 } else if (name == "converged") {
316 if (!bb)
317 throw ValueError("setting converged to non-bool value");
318 converged = b;
319 } else if (name == "time_step_backtracking_used") {
320 if (!bb)
321 throw ValueError("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 ValueError("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 ValueError("setting num_coarse_unknowns to non-int value");
330 num_coarse_unknowns = i;
331 } else {
332 throw ValueError(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 ValueError(std::string("unknown diagnostic item: ") + 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 ValueError("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 ValueError("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 /*
411 #ifdef ESYS_MPI
412 throw ValueError("AMG preconditioner is not supported in MPI builds");
413 break;
414 #endif
415 */
416 case SO_PRECONDITIONER_AMLI:
417 case SO_PRECONDITIONER_BOOMERAMG:
418 case SO_PRECONDITIONER_GAUSS_SEIDEL:
419 case SO_PRECONDITIONER_ILU0:
420 case SO_PRECONDITIONER_ILUT:
421 case SO_PRECONDITIONER_JACOBI:
422 case SO_PRECONDITIONER_NONE:
423 case SO_PRECONDITIONER_REC_ILU:
424 case SO_PRECONDITIONER_RILU:
425 this->preconditioner = preconditioner;
426 break;
427 default:
428 throw ValueError("unknown preconditioner");
429 }
430 }
431
432 SolverOptions SolverBuddy::getPreconditioner() const
433 {
434 return preconditioner;
435 }
436
437 void SolverBuddy::setSmoother(int s)
438 {
439 SolverOptions smoother = static_cast<SolverOptions>(s);
440 if (smoother != SO_PRECONDITIONER_JACOBI &&
441 smoother != SO_PRECONDITIONER_GAUSS_SEIDEL) {
442 throw ValueError("unknown smoother");
443 }
444 this->smoother = smoother;
445 }
446
447 SolverOptions SolverBuddy::getSmoother() const
448 {
449 return smoother;
450 }
451
452 void SolverBuddy::setSolverMethod(int method)
453 {
454 SolverOptions meth = static_cast<SolverOptions>(method);
455 switch(meth) {
456 case SO_DEFAULT:
457 case SO_METHOD_BICGSTAB:
458 case SO_METHOD_CGLS:
459 case SO_METHOD_CGS:
460 case SO_METHOD_CHOLEVSKY:
461 case SO_METHOD_CR:
462 case SO_METHOD_GMRES:
463 case SO_METHOD_HRZ_LUMPING:
464 case SO_METHOD_ITERATIVE:
465 case SO_METHOD_LSQR:
466 case SO_METHOD_MINRES:
467 case SO_METHOD_NONLINEAR_GMRES:
468 case SO_METHOD_PCG:
469 case SO_METHOD_PRES20:
470 case SO_METHOD_ROWSUM_LUMPING:
471 case SO_METHOD_TFQMR:
472 this->method = meth;
473 break;
474 case SO_METHOD_DIRECT:
475 #ifdef ESYS_HAVE_UMFPACK
476 this->method = meth;
477 break;
478 #elif defined ESYS_HAVE_TRILINOS
479 this->method = meth;
480 break;
481 #elif defined ESYS_HAVE_MKL
482 this->method = meth;
483 break;
484 #elif defined PASTIX
485 this->method = meth;
486 break;
487 #else
488 throw ValueError("Cannot use DIRECT solver method, the running "
489 "escript was not compiled with a direct solver enabled");
490 #endif
491 default:
492 throw ValueError("unknown solver method");
493 }
494 }
495
496 SolverOptions SolverBuddy::getSolverMethod() const
497 {
498 return method;
499 }
500
501 void SolverBuddy::setSolverTarget(int target)
502 {
503 SolverOptions targ = static_cast<SolverOptions>(target);
504 switch (targ) {
505 case SO_TARGET_CPU:
506 case SO_TARGET_GPU:
507 this->target = targ;
508 break;
509 default:
510 throw ValueError("unknown solver target");
511 }
512 }
513
514 SolverOptions SolverBuddy::getSolverTarget() const
515 {
516 return target;
517 }
518
519 void SolverBuddy::setPackage(int package)
520 {
521 SolverOptions pack = static_cast<SolverOptions>(package);
522 switch (pack) {
523 case SO_DEFAULT:
524 case SO_PACKAGE_CUSP:
525 case SO_PACKAGE_MKL:
526 case SO_PACKAGE_PASO:
527 case SO_PACKAGE_PASTIX:
528 case SO_PACKAGE_SUPER_LU:
529 case SO_PACKAGE_TRILINOS:
530 case SO_PACKAGE_UMFPACK:
531 this->package = pack;
532 break;
533 default:
534 throw ValueError("unknown solver package");
535 }
536 }
537
538 SolverOptions SolverBuddy::getPackage() const
539 {
540 return package;
541 }
542
543 void SolverBuddy::setReordering(int ordering)
544 {
545 SolverOptions ord = static_cast<SolverOptions>(ordering);
546 switch (ordering) {
547 case SO_REORDERING_DEFAULT:
548 case SO_REORDERING_MINIMUM_FILL_IN:
549 case SO_REORDERING_NESTED_DISSECTION:
550 case SO_REORDERING_NONE:
551 reordering = ord;
552 break;
553 default:
554 throw ValueError("unknown reordering strategy");
555 }
556 }
557
558 SolverOptions SolverBuddy::getReordering() const
559 {
560 return reordering;
561 }
562
563 void SolverBuddy::setRestart(int restart)
564 {
565 if (restart < 0)
566 throw ValueError("restart must be non-negative.");
567
568 this->restart = restart;
569 }
570
571 int SolverBuddy::getRestart() const
572 {
573 return restart;
574 }
575
576 int SolverBuddy::_getRestartForC() const
577 {
578 int r = getRestart();
579 if (r == 0)
580 return -1;
581 return r;
582 }
583
584 void SolverBuddy::setDiagonalDominanceThreshold(double value)
585 {
586 if (value < 0. || value > 1.)
587 throw ValueError("Diagonal dominance threshold must be between 0 and 1.");
588 diagonal_dominance_threshold = value;
589 }
590
591 double SolverBuddy::getDiagonalDominanceThreshold() const
592 {
593 return diagonal_dominance_threshold;
594 }
595
596 void SolverBuddy::setTruncation(int truncation)
597 {
598 if (truncation < 1)
599 throw ValueError("truncation must be positive.");
600 this->truncation = truncation;
601 }
602
603 int SolverBuddy::getTruncation() const
604 {
605 return truncation;
606 }
607
608 void SolverBuddy::setInnerIterMax(int iter_max)
609 {
610 if (iter_max < 1)
611 throw ValueError("maximum number of inner iteration must be positive.");
612 inner_iter_max = iter_max;
613 }
614
615 int SolverBuddy::getInnerIterMax() const
616 {
617 return inner_iter_max;
618 }
619
620 void SolverBuddy::setIterMax(int iter_max)
621 {
622 if (iter_max < 1)
623 throw ValueError("maximum number of iteration steps must be positive.");
624 this->iter_max = iter_max;
625 }
626
627 int SolverBuddy::getIterMax() const
628 {
629 return iter_max;
630 }
631
632 void SolverBuddy::setLevelMax(int level_max)
633 {
634 if (level_max < 0)
635 throw ValueError("maximum number of coarsening levels must be non-negative.");
636 this->level_max = level_max;
637 }
638
639 int SolverBuddy::getLevelMax() const
640 {
641 return level_max;
642 }
643
644 void SolverBuddy::setCycleType(int cycle_type)
645 {
646 this->cycle_type = cycle_type;
647 }
648
649 int SolverBuddy::getCycleType() const
650 {
651 return cycle_type;
652 }
653
654 void SolverBuddy::setCoarseningThreshold(double theta)
655 {
656 if (theta < 0. || theta > 1.)
657 throw ValueError("threshold must be between 0 and 1.");
658 coarsening_threshold = theta;
659 }
660
661 double SolverBuddy::getCoarseningThreshold() const
662 {
663 return coarsening_threshold;
664 }
665
666 void SolverBuddy::setNumSweeps(int sweeps)
667 {
668 if (sweeps < 1)
669 throw ValueError("number of sweeps must be positive.");
670 this->sweeps = sweeps;
671 }
672
673 int SolverBuddy::getNumSweeps() const
674 {
675 return sweeps;
676 }
677
678 void SolverBuddy::setNumPreSweeps(int sweeps)
679 {
680 if (sweeps < 1)
681 throw ValueError("number of pre-sweeps must be positive.");
682 pre_sweeps = sweeps;
683 }
684
685 int SolverBuddy::getNumPreSweeps() const
686 {
687 return pre_sweeps;
688 }
689
690 void SolverBuddy::setNumPostSweeps(int sweeps)
691 {
692 if (sweeps < 1)
693 throw ValueError("number of post-sweeps must be positive.");
694 post_sweeps = sweeps;
695 }
696
697 int SolverBuddy::getNumPostSweeps() const
698 {
699 return post_sweeps;
700 }
701
702 void SolverBuddy::setTolerance(double rtol)
703 {
704 if (rtol < 0. || rtol > 1.)
705 throw ValueError("tolerance must be between 0 and 1.");
706 tolerance = rtol;
707 }
708
709 double SolverBuddy::getTolerance() const
710 {
711 return tolerance;
712 }
713
714 void SolverBuddy::setAbsoluteTolerance(double atol)
715 {
716 if (atol < 0.)
717 throw ValueError("absolute tolerance must be non-negative.");
718 absolute_tolerance = atol;
719 }
720
721 double SolverBuddy::getAbsoluteTolerance() const
722 {
723 return absolute_tolerance;
724 }
725
726 void SolverBuddy::setInnerTolerance(double rtol)
727 {
728 if (rtol <= 0. || rtol > 1.)
729 throw ValueError("tolerance must be positive and less than or equal to 1.");
730 inner_tolerance = rtol;
731 }
732
733 double SolverBuddy::getInnerTolerance() const
734 {
735 return inner_tolerance;
736 }
737
738 void SolverBuddy::setDropTolerance(double drop_tol)
739 {
740 if (drop_tol < 0. || drop_tol > 1.)
741 throw ValueError("drop tolerance must be between 0 and 1.");
742 drop_tolerance = drop_tol;
743 }
744
745 double SolverBuddy::getDropTolerance() const
746 {
747 return drop_tolerance;
748 }
749
750 void SolverBuddy::setDropStorage(double storage)
751 {
752 if (storage < 1.)
753 throw ValueError("allowed storage increase must be greater than or equal to 1.");
754 drop_storage = storage;
755 }
756
757 double SolverBuddy::getDropStorage() const
758 {
759 return drop_storage;
760 }
761
762 void SolverBuddy::setRelaxationFactor(double factor)
763 {
764 if (factor < 0.)
765 throw ValueError("relaxation factor must be non-negative.");
766 relaxation = factor;
767 }
768
769 double SolverBuddy::getRelaxationFactor() const
770 {
771 return relaxation;
772 }
773
774 bool SolverBuddy::isSymmetric() const
775 {
776 return symmetric;
777 }
778
779 void SolverBuddy::setSymmetryOn()
780 {
781 symmetric = true;
782 }
783
784 void SolverBuddy::setSymmetryOff()
785 {
786 symmetric = false;
787 }
788
789 void SolverBuddy::setSymmetry(bool flag)
790 {
791 if (flag)
792 setSymmetryOn();
793 else
794 setSymmetryOff();
795 }
796
797 bool SolverBuddy::isVerbose() const
798 {
799 return verbose;
800 }
801
802 void SolverBuddy::setVerbosityOn()
803 {
804 verbose = true;
805 }
806
807 void SolverBuddy::setVerbosityOff()
808 {
809 verbose = false;
810 }
811
812 void SolverBuddy::setVerbosity(bool verbose)
813 {
814 if (verbose)
815 setVerbosityOn();
816 else
817 setVerbosityOff();
818 }
819
820 bool SolverBuddy::adaptInnerTolerance() const
821 {
822 return adapt_inner_tolerance;
823 }
824
825 void SolverBuddy::setInnerToleranceAdaptionOn()
826 {
827 adapt_inner_tolerance = true;
828 }
829
830 void SolverBuddy::setInnerToleranceAdaptionOff()
831 {
832 adapt_inner_tolerance = false;
833 }
834
835 void SolverBuddy::setInnerToleranceAdaption(bool adapt)
836 {
837 if (adapt)
838 setInnerToleranceAdaptionOn();
839 else
840 setInnerToleranceAdaptionOff();
841 }
842
843 bool SolverBuddy::acceptConvergenceFailure() const
844 {
845 return accept_convergence_failure;
846 }
847
848 void SolverBuddy::setAcceptanceConvergenceFailureOn()
849 {
850 accept_convergence_failure = true;
851 }
852
853 void SolverBuddy::setAcceptanceConvergenceFailureOff()
854 {
855 accept_convergence_failure = false;
856 }
857
858 void SolverBuddy::setAcceptanceConvergenceFailure(bool accept)
859 {
860 if (accept)
861 setAcceptanceConvergenceFailureOn();
862 else
863 setAcceptanceConvergenceFailureOff();
864 }
865
866 bool SolverBuddy::useLocalPreconditioner() const
867 {
868 return use_local_preconditioner;
869 }
870
871 void SolverBuddy::setLocalPreconditionerOn()
872 {
873 use_local_preconditioner = true;
874 }
875
876 void SolverBuddy::setLocalPreconditionerOff()
877 {
878 use_local_preconditioner=false;
879 }
880
881 void SolverBuddy::setLocalPreconditioner(bool use)
882 {
883 if (use)
884 setLocalPreconditionerOn();
885 else
886 setLocalPreconditionerOff();
887 }
888
889 void SolverBuddy::setMinCoarseMatrixSparsity(double sparsity)
890 {
891 if (sparsity < 0. || sparsity > 1.)
892 throw ValueError("sparsity must be between 0 and 1.");
893 min_sparsity = sparsity;
894 }
895
896 double SolverBuddy::getMinCoarseMatrixSparsity() const
897 {
898 return min_sparsity;
899 }
900
901 void SolverBuddy::setNumRefinements(int refinements)
902 {
903 if (refinements < 0)
904 throw ValueError("number of refinements must be non-negative.");
905 this->refinements = refinements;
906 }
907
908 int SolverBuddy::getNumRefinements() const
909 {
910 return refinements;
911 }
912
913 void SolverBuddy::setNumCoarseMatrixRefinements(int refinements)
914 {
915 if (refinements < 0)
916 throw ValueError("number of coarse matrix refinements must be non-negative.");
917 coarse_refinements = refinements;
918 }
919
920 int SolverBuddy::getNumCoarseMatrixRefinements() const
921 {
922 return coarse_refinements;
923 }
924
925 bool SolverBuddy::usePanel() const
926 {
927 return use_panel;
928 }
929
930 void SolverBuddy::setUsePanelOn()
931 {
932 use_panel = true;
933 }
934
935 void SolverBuddy::setUsePanelOff()
936 {
937 use_panel = false;
938 }
939
940 void SolverBuddy::setUsePanel(bool use)
941 {
942 if (use)
943 setUsePanelOn();
944 else
945 setUsePanelOff();
946 }
947
948 void SolverBuddy::setAMGInterpolation(int method)
949 {
950 SolverOptions meth = static_cast<SolverOptions>(method);
951 switch (meth) {
952 case SO_INTERPOLATION_CLASSIC:
953 case SO_INTERPOLATION_CLASSIC_WITH_FF_COUPLING:
954 case SO_INTERPOLATION_DIRECT:
955 amg_interpolation_method = meth;
956 break;
957 default:
958 throw ValueError("unknown AMG interpolation method");
959 }
960 }
961
962 SolverOptions SolverBuddy::getAMGInterpolation() const
963 {
964 return amg_interpolation_method;
965 }
966
967 void SolverBuddy::setODESolver(int method)
968 {
969 SolverOptions ode = static_cast<SolverOptions>(method);
970 switch (ode) {
971 case SO_ODESOLVER_BACKWARD_EULER:
972 case SO_ODESOLVER_CRANK_NICOLSON:
973 case SO_ODESOLVER_LINEAR_CRANK_NICOLSON:
974 ode_solver = ode;
975 break;
976 default:
977 throw ValueError("unknown ODE solver method");
978 }
979 }
980
981 SolverOptions SolverBuddy::getODESolver() const
982 {
983 return ode_solver;
984 }
985
986 } // namespace escript
987

  ViewVC Help
Powered by ViewVC 1.1.26