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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 6799 - (show annotations)
Mon Mar 25 05:53:58 2019 UTC (3 weeks, 4 days ago) by aellery
File size: 24892 byte(s)
I have rewritten the solverbuddy. Briefly:

1. The remaining AMG code has been removed from PASO.
2. If Trilinos is available, escript will now use it by default.
3. eScript will use a direct solver by default, (if one is available,) when solving 2 dimensional problems and an iterative solver, by default, when solving 3 dimensional problems. This can be changed by a user by manually specifying which solver to use.
4. There is a new option available, setHermitian(), that allows a user to specify when a coefficient matrix is Hermitian.
5. Symmetry information is always passed onto the Trilinos solver when this information is relevant.
6. All tests have been updated, when relevant, to reflect these changes.
7. I fixed a couple of undocumented bugs.


1
2 /*****************************************************************************
3 *
4 * Copyright (c) 2003-2018 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 package(SO_DEFAULT),
37 method(SO_DEFAULT),
38 preconditioner(SO_PRECONDITIONER_JACOBI),
39 ode_solver(SO_ODESOLVER_LINEAR_CRANK_NICOLSON),
40 reordering(SO_REORDERING_DEFAULT),
41 sweeps(1),
42 tolerance(1e-8),
43 absolute_tolerance(0.),
44 inner_tolerance(0.9),
45 drop_tolerance(0.0005),
46 drop_storage(2.),
47 iter_max(100000),
48 inner_iter_max(10),
49 truncation(20),
50 restart(0),
51 is_complex(false),
52 symmetric(false),
53 hermitian(false),
54 verbose(false),
55 adapt_inner_tolerance(true),
56 accept_convergence_failure(false),
57 relaxation(0.3),
58 use_local_preconditioner(false),
59 refinements(2),
60 dim(2)
61 {
62 resetDiagnostics(true);
63 setup();
64 }
65
66 SolverBuddy::~SolverBuddy()
67 {
68 }
69
70
71 std::string SolverBuddy::getSummary() const
72 {
73 std::stringstream out;
74 out << "Solver Package = " << getName(getPackage()) << std::endl
75 << "Verbosity = " << isVerbose() << std::endl
76 << "Accept failed convergence = " << acceptConvergenceFailure()
77 << std::endl
78 << "Relative tolerance = " << getTolerance() << std::endl
79 << "Absolute tolerance = " << getAbsoluteTolerance() << std::endl
80 << "Symmetric problem = " << isSymmetric() << std::endl
81 << "Hermitian problem = " << isHermitian() << std::endl
82 << "Maximum number of iteration steps = " << getIterMax() << std::endl
83 << "Inner tolerance = " << getInnerTolerance() << std::endl
84 << "Adapt innner tolerance = " << adaptInnerTolerance() << std::endl;
85
86 if (getPackage() == SO_DEFAULT || getPackage() == SO_PACKAGE_PASO ||
87 getPackage() == SO_PACKAGE_TRILINOS) {
88 out << "Solver method = " << getName(getSolverMethod()) << std::endl;
89 if (getSolverMethod() == SO_METHOD_GMRES) {
90 out << "Truncation = " << getTruncation() << std::endl
91 << "Restart = " << getRestart() << std::endl;
92 }
93 out << "Preconditioner = " << getName(getPreconditioner()) << std::endl
94 << "Apply preconditioner locally = " << useLocalPreconditioner()
95 << std::endl;
96 switch (getPreconditioner()) {
97 case SO_PRECONDITIONER_GAUSS_SEIDEL:
98 out << "Number of sweeps = " << getNumSweeps() << std::endl;
99 break;
100 case SO_PRECONDITIONER_ILUT:
101 out << "Drop tolerance = " << getDropTolerance() << std::endl;
102 out << "Storage increase = " << getDropStorage() << std::endl;
103 break;
104 case SO_PRECONDITIONER_RILU:
105 out << "Relaxation factor = " << getRelaxationFactor()
106 << std::endl;
107 break;
108 default:
109 break;
110 } // preconditioner switch
111 out << "ODE solver = " << getName(getODESolver()) << std::endl;
112 }
113 return out.str();
114 }
115
116 const char* SolverBuddy::getName(int key) const
117 {
118 switch (static_cast<SolverOptions>(key)) {
119 case SO_DEFAULT: return "DEFAULT";
120
121 case SO_PACKAGE_MKL: return "MKL";
122 case SO_PACKAGE_PASO: return "PASO";
123 case SO_PACKAGE_TRILINOS: return "TRILINOS";
124 case SO_PACKAGE_UMFPACK: return "UMFPACK";
125
126 case SO_METHOD_BICGSTAB: return "BICGSTAB";
127 case SO_METHOD_CGLS: return "CGLS";
128 case SO_METHOD_CGS: return "CGS";
129 case SO_METHOD_CHOLEVSKY: return "CHOLEVSKY";
130 case SO_METHOD_CR: return "CR";
131 case SO_METHOD_DIRECT: return "DIRECT";
132 case SO_METHOD_DIRECT_MUMPS: return "DIRECT_MUMPS";
133 case SO_METHOD_DIRECT_PARDISO: return "DIRECT_PARDISO";
134 case SO_METHOD_DIRECT_SUPERLU: return "DIRECT_SUPERLU";
135 case SO_METHOD_DIRECT_TRILINOS: return "DIRECT_TRILINOS";
136 case SO_METHOD_GMRES: return "GMRES";
137 case SO_METHOD_HRZ_LUMPING: return "HRZ_LUMPING";
138 case SO_METHOD_ITERATIVE: return "ITERATIVE";
139 case SO_METHOD_LSQR: return "LSQR";
140 case SO_METHOD_MINRES: return "MINRES";
141 case SO_METHOD_NONLINEAR_GMRES: return "NONLINEAR_GMRES";
142 case SO_METHOD_PCG: return "PCG";
143 case SO_METHOD_PRES20: return "PRES20";
144 case SO_METHOD_ROWSUM_LUMPING: return "ROWSUM_LUMPING";
145 case SO_METHOD_TFQMR: return "TFQMR";
146
147 case SO_PRECONDITIONER_AMG: return "AMG";
148 case SO_PRECONDITIONER_GAUSS_SEIDEL: return "GAUSS_SEIDEL";
149 case SO_PRECONDITIONER_ILU0: return "ILU0";
150 case SO_PRECONDITIONER_ILUT: return "ILUT";
151 case SO_PRECONDITIONER_JACOBI: return "JACOBI";
152 case SO_PRECONDITIONER_NONE: return "NO_PRECONDITIONER";
153 case SO_PRECONDITIONER_REC_ILU: return "REC_ILU";
154 case SO_PRECONDITIONER_RILU: return "RILU";
155
156 case SO_ODESOLVER_BACKWARD_EULER: return "BACKWARD_EULER";
157 case SO_ODESOLVER_CRANK_NICOLSON: return "CRANK_NICOLSON";
158 case SO_ODESOLVER_LINEAR_CRANK_NICOLSON: return "LINEAR_CRANK_NICOLSON";
159
160 case SO_INTERPOLATION_CLASSIC: return "CLASSIC_INTERPOLATION";
161 case SO_INTERPOLATION_CLASSIC_WITH_FF_COUPLING:
162 return "CLASSIC_INTERPOLATION_WITH_FF";
163 case SO_INTERPOLATION_DIRECT: return "DIRECT_INTERPOLATION";
164
165
166 case SO_REORDERING_DEFAULT: return "DEFAULT_REORDERING";
167 case SO_REORDERING_MINIMUM_FILL_IN: return "MINIMUM_FILL_IN";
168 case SO_REORDERING_NESTED_DISSECTION: return "NESTED_DISSECTION";
169 case SO_REORDERING_NONE: return "NO_REORDERING";
170 default:
171 throw ValueError("getName() invalid option given");
172 }
173 return "invalid option";
174 }
175
176 void SolverBuddy::resetDiagnostics(bool all)
177 {
178 num_iter = 0;
179 num_level = 0;
180 num_inner_iter = 0;
181 time = 0.;
182 set_up_time = 0.;
183 net_time = 0.;
184 residual_norm = 0.;
185 converged = false;
186 preconditioner_size = -1;
187 time_step_backtracking_used = false;
188 coarse_level_sparsity = 0;
189 num_coarse_unknowns = 0;
190 if (all) {
191 cum_num_inner_iter = 0;
192 cum_num_iter = 0;
193 cum_time = 0.;
194 cum_set_up_time = 0.;
195 cum_net_time = 0.;
196 }
197 }
198
199 void SolverBuddy::setup(){
200 setPackage(SO_DEFAULT);
201 setSolverMethod(SO_DEFAULT);
202 setPreconditioner(SO_PRECONDITIONER_JACOBI);
203 setODESolver(SO_ODESOLVER_LINEAR_CRANK_NICOLSON);
204 setReordering(SO_REORDERING_DEFAULT);
205 }
206
207 void SolverBuddy::updateDiagnostics(const std::string& name, bool value)
208 {
209 if (name == "converged") {
210 converged = value;
211 } else if (name == "time_step_backtracking_used") {
212 time_step_backtracking_used = value;
213 } else {
214 throw ValueError(std::string("Unknown diagnostic: ") + name);
215 }
216 }
217
218 void SolverBuddy::updateDiagnostics(const std::string& name, int value)
219 {
220 if (name == "num_iter") {
221 cum_num_iter += num_iter = value;
222 } else if (name == "num_level") {
223 num_level = value;
224 } else if (name == "num_inner_iter") {
225 cum_num_inner_iter += num_inner_iter = value;
226 } else if (name == "num_coarse_unknowns") {
227 num_coarse_unknowns = value;
228 } else {
229 throw ValueError(std::string("Unknown diagnostic: ") + name);
230 }
231 }
232
233 void SolverBuddy::updateDiagnostics(const std::string& name, double value)
234 {
235 if (name == "time") {
236 cum_time += time = value;
237 } else if (name == "set_up_time") {
238 cum_set_up_time += set_up_time = value;
239 } else if (name == "net_time") {
240 cum_net_time += net_time = value;
241 } else if (name == "residual_norm") {
242 residual_norm = value;
243 } else if (name == "coarse_level_sparsity") {
244 coarse_level_sparsity = value;
245 } else {
246 throw ValueError(std::string("Unknown diagnostic: ") + name);
247 }
248 }
249
250 void SolverBuddy::updateDiagnosticsPy(const std::string& name,
251 const bp::object& value)
252 {
253 int i=0;
254 double d=0; // to keep older compilers happy
255 bool b=false;
256 bool ib = convert<int>(value, i);
257 bool db = convert<double>(value, d);
258 bool bb = convert<bool>(value, b);
259
260 if (name == "num_iter") {
261 if (!ib)
262 throw ValueError("setting num_iter to non-int value");
263 cum_num_iter += num_iter = i;
264 } else if (name == "num_level") {
265 if (!ib)
266 throw ValueError("setting num_level to non-int value");
267 num_level = i;
268 } else if (name == "num_inner_iter") {
269 if (!ib)
270 throw ValueError("setting num_inner_iter to non-int value");
271 cum_num_inner_iter += num_inner_iter = i;
272 } else if (name == "time") {
273 if (!db)
274 throw ValueError("setting time to non-double value");
275 cum_time += time = d;
276 } else if (name == "set_up_time") {
277 if (!db)
278 throw ValueError("setting set_up_time to non-double value");
279 cum_set_up_time += set_up_time = d;
280 } else if (name == "net_time") {
281 if (!db)
282 throw ValueError("setting net_time to non-double value");
283 cum_net_time += net_time = d;
284 } else if (name == "residual_norm") {
285 if (!db)
286 throw ValueError("setting residual_norm to non-double value");
287 residual_norm = d;
288 } else if (name == "converged") {
289 if (!bb)
290 throw ValueError("setting converged to non-bool value");
291 converged = b;
292 } else if (name == "time_step_backtracking_used") {
293 if (!bb)
294 throw ValueError("setting time_step_backtracking_used to non-bool value");
295 time_step_backtracking_used = b;
296 } else if (name == "coarse_level_sparsity") {
297 if (!db)
298 throw ValueError("setting coarse_level_sparsity to non-double value");
299 coarse_level_sparsity = d;
300 } else if (name == "num_coarse_unknowns") {
301 if (!ib)
302 throw ValueError("setting num_coarse_unknowns to non-int value");
303 num_coarse_unknowns = i;
304 } else {
305 throw ValueError(std::string("Unknown diagnostic: ") + name);
306 }
307 }
308
309 double SolverBuddy::getDiagnostics(const std::string name) const
310 {
311 if (name == "num_iter") return num_iter;
312 else if (name == "cum_num_iter") return cum_num_iter;
313 else if (name == "num_level") return num_level;
314 else if (name == "num_inner_iter") return num_inner_iter;
315 else if (name == "cum_num_inner_iter") return cum_num_inner_iter;
316 else if (name == "time") return time;
317 else if (name == "cum_time") return cum_time;
318 else if (name == "set_up_time") return set_up_time;
319 else if (name == "cum_set_up_time") return cum_set_up_time;
320 else if (name == "net_time") return net_time;
321 else if (name == "cum_net_time") return cum_net_time;
322 else if (name == "residual_norm") return residual_norm;
323 else if (name == "converged") return converged;
324 else if (name == "preconditioner_size") return preconditioner_size;
325 else if (name == "time_step_backtracking_used")
326 return time_step_backtracking_used;
327 else if (name == "coarse_level_sparsity") return coarse_level_sparsity;
328 else if (name == "num_coarse_unknowns") return num_coarse_unknowns;
329
330 throw ValueError(std::string("unknown diagnostic item: ") + name);
331 }
332
333 bool SolverBuddy::hasConverged() const
334 {
335 return converged;
336 }
337
338 void SolverBuddy::setPreconditioner(int precon)
339 {
340 SolverOptions preconditioner = static_cast<SolverOptions>(precon);
341 switch(preconditioner) {
342 case SO_PRECONDITIONER_AMG:
343 #ifndef ESYS_HAVE_TRILINOS
344 throw ValueError("escript was not compiled with Trilinos enabled");
345 #endif
346 case SO_PRECONDITIONER_GAUSS_SEIDEL:
347 case SO_PRECONDITIONER_JACOBI: // This is the default preconditioner in ifpack2
348 case SO_PRECONDITIONER_ILU0:
349 case SO_PRECONDITIONER_ILUT:
350 case SO_PRECONDITIONER_NONE:
351 case SO_PRECONDITIONER_REC_ILU:
352 case SO_PRECONDITIONER_RILU:
353 this->preconditioner = preconditioner;
354 break;
355 default:
356 throw ValueError("unknown preconditioner");
357 }
358 }
359
360 SolverOptions SolverBuddy::getPreconditioner() const
361 {
362 return preconditioner;
363 }
364
365 void SolverBuddy::setSolverMethod(int method)
366 {
367 SolverOptions meth = static_cast<SolverOptions>(method);
368
369 // Possible error: User has told escript to use a direct solver but no direct solver is available
370 #if !defined(ESYS_HAVE_TRILINOS) && !defined(ESYS_HAVE_UMFPACK) && !defined(ESYS_HAVE_MKL)
371 if(meth == SO_METHOD_DIRECT || meth == SO_METHOD_DIRECT_MUMPS
372 || meth == SO_METHOD_DIRECT_PARDISO || meth == SO_METHOD_DIRECT_SUPERLU
373 || meth == SO_METHOD_DIRECT_TRILINOS) {
374 throw ValueError("Cannot use DIRECT solver method. eScript was not compiled with a direct solver.");
375 }
376 #endif
377
378 switch(meth) {
379 // Default settings:
380 // If we're solving a 2d problem, default to the trilinos direct solver
381 // If trilinos is not available, default to the paso direct solver
382 // If trilinos is not available and there is no paso direct solver, fall back on GMRES
383 // In 3D, always use an iterative method
384
385 case SO_DEFAULT:
386 if(dim == 2){
387 #ifdef ESYS_HAVE_TRILINOS
388 this->method = SO_METHOD_DIRECT_TRILINOS;
389 #elif defined(ESYS_HAVE_UMFPACK) || defined(ESYS_HAVE_MKL)
390 this->method = SO_METHOD_DIRECT;
391 #else
392 this->method = SO_METHOD_GMRES;
393 #endif
394 break;
395 } else {
396 this->method = SO_METHOD_GMRES;
397 break;
398 }
399
400 case SO_METHOD_DIRECT:
401 #ifdef ESYS_HAVE_TRILINOS
402 // Default to trilinos
403 this->method = SO_METHOD_DIRECT_TRILINOS;
404 break;
405 #else
406 this->method = SO_METHOD_DIRECT;
407 break;
408 #endif
409 case SO_METHOD_DIRECT_MUMPS:
410 case SO_METHOD_DIRECT_PARDISO:
411 case SO_METHOD_DIRECT_SUPERLU:
412 case SO_METHOD_DIRECT_TRILINOS:
413 #ifdef ESYS_HAVE_TRILINOS
414 // We're assuming that trilinos was compiled with MUMPS, PARDISO, SUPERLU etc.
415 this->method = meth;
416 break;
417 #else
418 throw ValueError("escript was not compiled with Trilinos enabled");
419 #endif
420 //These methods are in paso
421 case SO_METHOD_ITERATIVE:
422 case SO_METHOD_BICGSTAB:
423 case SO_METHOD_CHOLEVSKY:
424 case SO_METHOD_CGS:
425 case SO_METHOD_GMRES:
426 case SO_METHOD_HRZ_LUMPING:
427 case SO_METHOD_TFQMR:
428 case SO_METHOD_MINRES:
429 case SO_METHOD_NONLINEAR_GMRES:
430 case SO_METHOD_PCG:
431 case SO_METHOD_PRES20:
432 case SO_METHOD_ROWSUM_LUMPING:
433 this->method = meth;
434 break;
435 //These methods are in Trilinos but not in paso
436 case SO_METHOD_CGLS:
437 case SO_METHOD_CR:
438 case SO_METHOD_LSQR:
439 #ifdef ESYS_HAVE_TRILINOS
440 this->package = SO_PACKAGE_TRILINOS;
441 this->method = meth;
442 break;
443 #else
444 throw ValueError("escript was not compiled with Trilinos enabled");
445 #endif
446 default:
447 throw ValueError("unknown solver method");
448 }
449 }
450
451 SolverOptions SolverBuddy::getSolverMethod() const
452 {
453 return method;
454 }
455
456 void SolverBuddy::setPackage(int package)
457 {
458 SolverOptions pack = static_cast<SolverOptions>(package);
459 switch (pack) {
460 case SO_DEFAULT:
461 // Default to trilinos if it is available, else use PASO
462 // This should always work because escript cannot be compiled
463 // unless at least one of these is available
464 #ifdef ESYS_HAVE_TRILINOS
465 this->package = SO_PACKAGE_TRILINOS;
466 break;
467 #else
468 this->package = SO_PACKAGE_PASO;
469 break;
470 #endif
471 // Set to PASO iff escript was compiled with PASO
472 case SO_PACKAGE_PASO:
473 #ifdef ESYS_HAVE_PASO
474 this->package = SO_PACKAGE_PASO;
475 break;
476 #else
477 throw ValueError("escript was not compiled with PASO enabled");
478 #endif
479 // Set to Trilinos iff escript was compiled with Trilinos
480 case SO_PACKAGE_TRILINOS:
481 #ifdef ESYS_HAVE_TRILINOS
482 this->package = SO_PACKAGE_TRILINOS;
483 break;
484 #else
485 throw ValueError("escript was not compiled with Trilinos enabled");
486 #endif
487 // Set to MKL iff escript was compiled with MKL
488 case SO_PACKAGE_MKL:
489 #ifdef ESYS_HAVE_MKL
490 this->package = SO_PACKAGE_MKL;
491 break;
492 #else
493 throw ValueError("escript was not compiled with MKL enabled");
494 #endif
495 // Set to Umfpack iff escript was compiled with Umfpack
496 case SO_PACKAGE_UMFPACK:
497 #ifdef ESYS_HAVE_UMFPACK
498 this->package = SO_PACKAGE_UMFPACK;
499 break;
500 #else
501 throw ValueError("escript was not compiled with UMFPACK enabled");
502 #endif
503 default:
504 throw ValueError("unknown solver package");
505 }
506 }
507
508 SolverOptions SolverBuddy::getPackage() const
509 {
510 return package;
511 }
512
513 void SolverBuddy::setReordering(int ordering)
514 {
515 SolverOptions ord = static_cast<SolverOptions>(ordering);
516 switch (ordering) {
517 case SO_REORDERING_DEFAULT:
518 case SO_REORDERING_MINIMUM_FILL_IN:
519 case SO_REORDERING_NESTED_DISSECTION:
520 case SO_REORDERING_NONE:
521 reordering = ord;
522 break;
523 default:
524 throw ValueError("unknown reordering strategy");
525 }
526 }
527
528 SolverOptions SolverBuddy::getReordering() const
529 {
530 return reordering;
531 }
532
533 void SolverBuddy::setRestart(int restart)
534 {
535 if (restart < 0)
536 throw ValueError("restart must be non-negative.");
537
538 this->restart = restart;
539 }
540
541 int SolverBuddy::getRestart() const
542 {
543 return restart;
544 }
545
546 int SolverBuddy::_getRestartForC() const
547 {
548 int r = getRestart();
549 if (r == 0)
550 return -1;
551 return r;
552 }
553
554 void SolverBuddy::setTruncation(int truncation)
555 {
556 if (truncation < 1)
557 throw ValueError("truncation must be positive.");
558 this->truncation = truncation;
559 }
560
561 int SolverBuddy::getTruncation() const
562 {
563 return truncation;
564 }
565
566 void SolverBuddy::setInnerIterMax(int iter_max)
567 {
568 if (iter_max < 1)
569 throw ValueError("maximum number of inner iteration must be positive.");
570 inner_iter_max = iter_max;
571 }
572
573 int SolverBuddy::getInnerIterMax() const
574 {
575 return inner_iter_max;
576 }
577
578 void SolverBuddy::setIterMax(int iter_max)
579 {
580 if (iter_max < 1)
581 throw ValueError("maximum number of iteration steps must be positive.");
582 this->iter_max = iter_max;
583 }
584
585 int SolverBuddy::getIterMax() const
586 {
587 return iter_max;
588 }
589
590 void SolverBuddy::setNumSweeps(int sweeps)
591 {
592 if (sweeps < 1)
593 throw ValueError("number of sweeps must be positive.");
594 this->sweeps = sweeps;
595 }
596
597 int SolverBuddy::getNumSweeps() const
598 {
599 return sweeps;
600 }
601
602 void SolverBuddy::setTolerance(double rtol)
603 {
604 if (rtol < 0. || rtol > 1.)
605 throw ValueError("tolerance must be between 0 and 1.");
606 tolerance = rtol;
607 }
608
609 double SolverBuddy::getTolerance() const
610 {
611 return tolerance;
612 }
613
614 void SolverBuddy::setAbsoluteTolerance(double atol)
615 {
616 if (atol < 0.)
617 throw ValueError("absolute tolerance must be non-negative.");
618 absolute_tolerance = atol;
619 }
620
621 double SolverBuddy::getAbsoluteTolerance() const
622 {
623 return absolute_tolerance;
624 }
625
626 void SolverBuddy::setInnerTolerance(double rtol)
627 {
628 if (rtol <= 0. || rtol > 1.)
629 throw ValueError("tolerance must be positive and less than or equal to 1.");
630 inner_tolerance = rtol;
631 }
632
633 double SolverBuddy::getInnerTolerance() const
634 {
635 return inner_tolerance;
636 }
637
638 void SolverBuddy::setDropTolerance(double drop_tol)
639 {
640 if (drop_tol < 0. || drop_tol > 1.)
641 throw ValueError("drop tolerance must be between 0 and 1.");
642 drop_tolerance = drop_tol;
643 }
644
645 double SolverBuddy::getDropTolerance() const
646 {
647 return drop_tolerance;
648 }
649
650 void SolverBuddy::setDropStorage(double storage)
651 {
652 if (storage < 1.)
653 throw ValueError("allowed storage increase must be greater than or equal to 1.");
654 drop_storage = storage;
655 }
656
657 double SolverBuddy::getDropStorage() const
658 {
659 return drop_storage;
660 }
661
662 void SolverBuddy::setRelaxationFactor(double factor)
663 {
664 if (factor < 0.)
665 throw ValueError("relaxation factor must be non-negative.");
666 relaxation = factor;
667 }
668
669 double SolverBuddy::getRelaxationFactor() const
670 {
671 return relaxation;
672 }
673
674 bool SolverBuddy::isComplex() const
675 {
676 return is_complex;
677 }
678
679 void SolverBuddy::setComplex(bool flag)
680 {
681 is_complex = flag;
682 }
683
684 bool SolverBuddy::isSymmetric() const
685 {
686 return symmetric;
687 }
688
689 void SolverBuddy::setSymmetryOn()
690 {
691 symmetric = true;
692 }
693
694 void SolverBuddy::setSymmetryOff()
695 {
696 symmetric = false;
697 }
698
699 void SolverBuddy::setSymmetry(bool flag)
700 {
701 if (flag)
702 setSymmetryOn();
703 else
704 setSymmetryOff();
705 }
706
707 bool SolverBuddy::isHermitian() const
708 {
709 return hermitian;
710 }
711
712 void SolverBuddy::setHermitianOn()
713 {
714 hermitian = true;
715 }
716
717 void SolverBuddy::setHermitianOff()
718 {
719 hermitian = false;
720 }
721
722 void SolverBuddy::setHermitian(bool flag)
723 {
724 if (flag)
725 setHermitianOn();
726 else
727 setHermitianOff();
728 }
729
730 bool SolverBuddy::isVerbose() const
731 {
732 return verbose;
733 }
734
735 void SolverBuddy::setVerbosityOn()
736 {
737 verbose = true;
738 }
739
740 void SolverBuddy::setVerbosityOff()
741 {
742 verbose = false;
743 }
744
745 void SolverBuddy::setVerbosity(bool verbose)
746 {
747 if (verbose)
748 setVerbosityOn();
749 else
750 setVerbosityOff();
751 }
752
753 bool SolverBuddy::adaptInnerTolerance() const
754 {
755 return adapt_inner_tolerance;
756 }
757
758 void SolverBuddy::setInnerToleranceAdaptionOn()
759 {
760 adapt_inner_tolerance = true;
761 }
762
763 void SolverBuddy::setInnerToleranceAdaptionOff()
764 {
765 adapt_inner_tolerance = false;
766 }
767
768 void SolverBuddy::setInnerToleranceAdaption(bool adapt)
769 {
770 if (adapt)
771 setInnerToleranceAdaptionOn();
772 else
773 setInnerToleranceAdaptionOff();
774 }
775
776 bool SolverBuddy::acceptConvergenceFailure() const
777 {
778 return accept_convergence_failure;
779 }
780
781 void SolverBuddy::setAcceptanceConvergenceFailureOn()
782 {
783 accept_convergence_failure = true;
784 }
785
786 void SolverBuddy::setAcceptanceConvergenceFailureOff()
787 {
788 accept_convergence_failure = false;
789 }
790
791 void SolverBuddy::setAcceptanceConvergenceFailure(bool accept)
792 {
793 if (accept)
794 setAcceptanceConvergenceFailureOn();
795 else
796 setAcceptanceConvergenceFailureOff();
797 }
798
799 bool SolverBuddy::useLocalPreconditioner() const
800 {
801 return use_local_preconditioner;
802 }
803
804 void SolverBuddy::setLocalPreconditionerOn()
805 {
806 use_local_preconditioner = true;
807 }
808
809 void SolverBuddy::setLocalPreconditionerOff()
810 {
811 use_local_preconditioner=false;
812 }
813
814 void SolverBuddy::setLocalPreconditioner(bool use)
815 {
816 if (use)
817 setLocalPreconditionerOn();
818 else
819 setLocalPreconditionerOff();
820 }
821
822 void SolverBuddy::setNumRefinements(int refinements)
823 {
824 if (refinements < 0)
825 throw ValueError("number of refinements must be non-negative.");
826 this->refinements = refinements;
827 }
828
829 int SolverBuddy::getNumRefinements() const
830 {
831 return refinements;
832 }
833
834 void SolverBuddy::setODESolver(int method)
835 {
836 SolverOptions ode = static_cast<SolverOptions>(method);
837 switch (ode) {
838 case SO_ODESOLVER_BACKWARD_EULER:
839 case SO_ODESOLVER_CRANK_NICOLSON:
840 case SO_ODESOLVER_LINEAR_CRANK_NICOLSON:
841 ode_solver = ode;
842 break;
843 default:
844 throw ValueError("unknown ODE solver method");
845 }
846 }
847
848 SolverOptions SolverBuddy::getODESolver() const
849 {
850 return ode_solver;
851 }
852
853 void SolverBuddy::setTrilinosParameter(const std::string& name,
854 const bp::object& value)
855 {
856 #ifdef ESYS_HAVE_TRILINOS
857 if (value == "true" || value == "false"){
858 trilinosParams[name] = (value == "true" ? true : false);
859 } else {
860 trilinosParams[name] = value;
861 }
862 #endif
863 }
864
865 bp::dict SolverBuddy::getTrilinosParameters() const
866 {
867 return trilinosParams;
868 }
869
870 void SolverBuddy::setDim(int dim)
871 {
872 if (dim != 2 && dim != 3)
873 throw ValueError("Dimension must be either 2 or 3.");
874 this->dim = dim;
875 }
876
877 int SolverBuddy::getDim()
878 {
879 return dim;
880 }
881
882 } // namespace escript
883

  ViewVC Help
Powered by ViewVC 1.1.26