MBDyn-1.7.3
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups
linsol.cc
Go to the documentation of this file.
1 /* $Header: /var/cvs/mbdyn/mbdyn/mbdyn-1.0/libraries/libmbwrap/linsol.cc,v 1.82 2017/01/12 14:44:25 masarati Exp $ */
2 /*
3  * MBDyn (C) is a multibody analysis code.
4  * http://www.mbdyn.org
5  *
6  * Copyright (C) 1996-2017
7  *
8  * Pierangelo Masarati <masarati@aero.polimi.it>
9  * Paolo Mantegazza <mantegazza@aero.polimi.it>
10  *
11  * Dipartimento di Ingegneria Aerospaziale - Politecnico di Milano
12  * via La Masa, 34 - 20156 Milano, Italy
13  * http://www.aero.polimi.it
14  *
15  * Changing this copyright notice is forbidden.
16  *
17  * This program is free software; you can redistribute it and/or modify
18  * it under the terms of the GNU General Public License as published by
19  * the Free Software Foundation (version 2 of the License).
20  *
21  *
22  * This program is distributed in the hope that it will be useful,
23  * but WITHOUT ANY WARRANTY; without even the implied warranty of
24  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
25  * GNU General Public License for more details.
26  *
27  * You should have received a copy of the GNU General Public License
28  * along with this program; if not, write to the Free Software
29  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
30  */
31 
32 #include "mbconfig.h" /* This goes first in every *.c,*.cc file */
33 
34 #include "myassert.h"
35 #include "ac/sys_sysinfo.h"
36 
37 #include "parser.h"
38 
39 #include "spmapmh.h"
40 #include "ccmh.h"
41 #include "dirccmh.h"
42 #include "harwrap.h"
43 #include "mschwrap.h"
44 #include "y12wrap.h"
45 #include "umfpackwrap.h"
46 #include "kluwrap.h"
47 #include "parsuperluwrap.h"
48 #include "superluwrap.h"
49 #include "lapackwrap.h"
50 #include "taucswrap.h"
51 #include "naivewrap.h"
52 #include "parnaivewrap.h"
53 
54 
55 #include "linsol.h"
56 
57 /* solver data */
58 const LinSol::solver_t solver[] = {
59  { "Empty", NULL,
60  LinSol::EMPTY_SOLVER,
61  LinSol::SOLVER_FLAGS_NONE,
62  LinSol::SOLVER_FLAGS_NONE,
63  -1., -1. },
64  { "Harwell", NULL,
65  LinSol::HARWELL_SOLVER,
66  LinSol::SOLVER_FLAGS_NONE,
67  LinSol::SOLVER_FLAGS_NONE,
68  -1., -1. },
69  { "Lapack", NULL,
70  LinSol::LAPACK_SOLVER,
71  LinSol::SOLVER_FLAGS_NONE,
72  LinSol::SOLVER_FLAGS_NONE,
73  -1., -1. },
74  { "Meschach", NULL,
75  LinSol::MESCHACH_SOLVER,
76  LinSol::SOLVER_FLAGS_NONE,
77  LinSol::SOLVER_FLAGS_NONE,
78  -1., -1. },
79  { "Naive", NULL,
80  LinSol::NAIVE_SOLVER,
81  LinSol::SOLVER_FLAGS_ALLOWS_REVERSE_CUTHILL_MC_KEE |
82  LinSol::SOLVER_FLAGS_ALLOWS_COLAMD |
83  LinSol::SOLVER_FLAGS_ALLOWS_MMDATA |
84  LinSol::SOLVER_FLAGS_ALLOWS_MDAPLUSAT |
85  LinSol::SOLVER_FLAGS_ALLOWS_REVERSE_CUTHILL_MC_KEE |
86  LinSol::SOLVER_FLAGS_ALLOWS_KING |
87  LinSol::SOLVER_FLAGS_ALLOWS_SLOAN |
88  LinSol::SOLVER_FLAGS_ALLOWS_NESTED_DISSECTION |
89  LinSol::SOLVER_FLAGS_ALLOWS_MT_ASS |
90 #ifdef USE_NAIVE_MULTITHREAD
91  LinSol::SOLVER_FLAGS_ALLOWS_MT_FCT |
92 #endif /* USE_NAIVE_MULTITHREAD */
93  0,
94  LinSol::SOLVER_FLAGS_NONE,
95  1.e-8, -1. },
96  { "SuperLU", NULL,
97  LinSol::SUPERLU_SOLVER,
98  LinSol::SOLVER_FLAGS_ALLOWS_COLAMD |
99  LinSol::SOLVER_FLAGS_ALLOWS_MMDATA |
100  LinSol::SOLVER_FLAGS_ALLOWS_MAP |
101  LinSol::SOLVER_FLAGS_ALLOWS_CC |
102  LinSol::SOLVER_FLAGS_ALLOWS_DIR |
103 #ifdef USE_SUPERLU_MT
104  LinSol::SOLVER_FLAGS_ALLOWS_MT_FCT |
105 #endif /* USE_SUPERLU_MT */
106  0,
107  LinSol::SOLVER_FLAGS_ALLOWS_MAP|LinSol::SOLVER_FLAGS_ALLOWS_COLAMD,
108  1., -1. },
109  { "Taucs", NULL,
110  LinSol::TAUCS_SOLVER,
111  LinSol::SOLVER_FLAGS_ALLOWS_MAP|LinSol::SOLVER_FLAGS_ALLOWS_CC|LinSol::SOLVER_FLAGS_ALLOWS_DIR,
112  LinSol::SOLVER_FLAGS_ALLOWS_MAP,
113  -1., -1. },
114  { "Umfpack", "umfpack3",
115  LinSol::UMFPACK_SOLVER,
116  LinSol::SOLVER_FLAGS_ALLOWS_MAP|LinSol::SOLVER_FLAGS_ALLOWS_CC|LinSol::SOLVER_FLAGS_ALLOWS_DIR|LinSol::SOLVER_FLAGS_ALLOWS_MT_ASS,
117  LinSol::SOLVER_FLAGS_ALLOWS_MAP|LinSol::SOLVER_FLAGS_ALLOWS_MT_ASS,
118  .1,
119 #ifdef UMFPACK_DROPTOL
120  0.
121 #else // ! UMFPACK_DROPTOL
122  -1.
123 #endif // ! UMFPACK_DROPTOL
124  },
125  { "KLU", NULL,
126  LinSol::KLU_SOLVER,
127  LinSol::SOLVER_FLAGS_ALLOWS_MAP|LinSol::SOLVER_FLAGS_ALLOWS_CC|LinSol::SOLVER_FLAGS_ALLOWS_DIR|LinSol::SOLVER_FLAGS_ALLOWS_MT_ASS,
128  LinSol::SOLVER_FLAGS_ALLOWS_MAP|LinSol::SOLVER_FLAGS_ALLOWS_MT_ASS,
129  .1 },
130  { "Y12", NULL,
131  LinSol::Y12_SOLVER,
132  LinSol::SOLVER_FLAGS_ALLOWS_MAP|LinSol::SOLVER_FLAGS_ALLOWS_CC|LinSol::SOLVER_FLAGS_ALLOWS_DIR|LinSol::SOLVER_FLAGS_ALLOWS_MT_ASS,
133  LinSol::SOLVER_FLAGS_ALLOWS_MAP|LinSol::SOLVER_FLAGS_ALLOWS_MT_ASS,
134  -1., -1. },
135  { NULL, NULL,
136  LinSol::EMPTY_SOLVER,
137  LinSol::SOLVER_FLAGS_NONE,
138  LinSol::SOLVER_FLAGS_NONE,
139  -1., -1. }
140 };
141 
142 /*
143  * Default solver
144  */
145 LinSol::SolverType LinSol::defaultSolver =
146 #if defined(USE_UMFPACK)
147  LinSol::UMFPACK_SOLVER
148 #elif /* !USE_UMFPACK */ defined(USE_KLU)
149  LinSol::KLU_SOLVER
150 #else /* !USE_KLU */
151  // Naive always present
152  LinSol::NAIVE_SOLVER
153 #endif
154  ;
155 
156 LinSol::LinSol(void)
157 : currSolver(LinSol::defaultSolver),
158 solverFlags(0),
159 nThreads(1),
160 iWorkSpaceSize(0),
161 blockSize(0),
162 dPivotFactor(-1.),
163 dDropTolerance(0.),
164 iMaxIter(0) // Restore the original behavior by default
165 {
166  NO_OP;
167 }
168 
169 LinSol::~LinSol(void)
170 {
171  NO_OP;
172 }
173 
174 LinSol::SolverType
175 LinSol::GetSolver(void) const
176 {
177  return currSolver;
178 }
179 
180 bool
181 LinSol::SetSolver(LinSol::SolverType t, unsigned f)
182 {
183  /* check flags */
184  if ((::solver[t].s_flags & f) != f) {
185  return false;
186  }
187 
188  solverFlags = f;
189 
190  switch (t) {
191  case LinSol::UMFPACK_SOLVER:
192 #ifdef USE_UMFPACK
193  currSolver = t;
194  return true;
195 #endif /* USE_UMFPACK */
196 
197  case LinSol::KLU_SOLVER:
198 #ifdef USE_KLU
199  currSolver = t;
200  return true;
201 #endif /* USE_UMFPACK */
202 
203  case LinSol::SUPERLU_SOLVER:
204 #ifdef USE_SUPERLU
205  currSolver = t;
206  return true;
207 #endif /* USE_SUPERLU */
208 
209  case LinSol::LAPACK_SOLVER:
210 #ifdef USE_LAPACK
211  currSolver = t;
212  return true;
213 #endif /* USE_LAPACK */
214 
215  case LinSol::TAUCS_SOLVER:
216 #ifdef USE_TAUCS
217  currSolver = t;
218  return true;
219 #endif /* USE_TAUCS */
220 
221  case LinSol::Y12_SOLVER:
222 #ifdef USE_Y12
223  currSolver = t;
224  return true;
225 #endif /* USE_Y12 */
226 
227  case LinSol::HARWELL_SOLVER:
228 #ifdef USE_HARWELL
229  currSolver = t;
230  return true;
231 #endif /* USE_HARWELL */
232 
233  case LinSol::MESCHACH_SOLVER:
234 #ifdef USE_MESCHACH
235  currSolver = t;
236  return true;
237 #endif /* USE_MESCHACH */
238 
239  /* else */
240  silent_cerr(::solver[t].s_name << " unavailable" << std::endl);
241  return false;
242 
243  case LinSol::NAIVE_SOLVER:
244  currSolver = t;
245  return true;
246 
247  case LinSol::EMPTY_SOLVER:
248  currSolver = t;
249  return true;
250 
251  default:
252  return false;
253  }
254 }
255 
256 unsigned
257 LinSol::GetSolverFlags(void) const
258 {
259  return solverFlags;
260 }
261 
262 unsigned
263 LinSol::GetSolverFlags(SolverType t) const
264 {
265  return ::solver[t].s_flags;
266 }
267 
268 const char *const
269 LinSol::GetSolverName(void) const
270 {
271  return ::solver[currSolver].s_name;
272 }
273 
274 const char *const
275 LinSol::GetSolverName(SolverType t) const
276 {
277  return ::solver[t].s_name;
278 }
279 
280 bool
281 LinSol::SetSolverFlags(unsigned f)
282 {
283  if ((::solver[currSolver].s_flags & f) == f) {
284  solverFlags = f;
285  return true;
286  }
287 
288  return false;
289 }
290 
291 bool
292 LinSol::AddSolverFlags(unsigned f)
293 {
294  if ((::solver[currSolver].s_flags & f) == f) {
295  solverFlags |= f;
296  return true;
297  }
298 
299  return false;
300 }
301 
302 bool
303 LinSol::MaskSolverFlags(unsigned f)
304 {
305  if ((::solver[currSolver].s_flags & f) == f) {
306  solverFlags &= ~f;
307  return true;
308  }
309 
310  return false;
311 }
312 
313 bool
314 LinSol::SetNumThreads(unsigned nt)
315 {
316  if (GetSolverFlags(currSolver) & LinSol::SOLVER_FLAGS_ALLOWS_MT_FCT) {
317  if (nt == 0) {
318  solverFlags &= ~LinSol::SOLVER_FLAGS_ALLOWS_MT_FCT;
319 
320  } else {
321  solverFlags |= LinSol::SOLVER_FLAGS_ALLOWS_MT_FCT;
322  }
323 
324  nThreads = nt;
325 
326  return true;
327  }
328 
329  return false;
330 }
331 
332 unsigned
333 LinSol::GetNumThreads(void) const
334 {
335  return nThreads;
336 }
337 
338 integer
339 LinSol::iGetWorkSpaceSize(void) const
340 {
341  return iWorkSpaceSize;
342 }
343 
344 const doublereal&
345 LinSol::dGetPivotFactor(void) const
346 {
347  return dPivotFactor;
348 }
349 
350 const doublereal&
351 LinSol::dGetDropTolerance(void) const
352 {
353  return dDropTolerance;
354 }
355 
356 bool
357 LinSol::SetWorkSpaceSize(integer i)
358 {
359  switch (currSolver) {
360  case LinSol::Y12_SOLVER:
361  iWorkSpaceSize = i;
362  break;
363 
364  default:
365  return false;
366  }
367 
368  return true;
369 }
370 
371 bool
372 LinSol::SetPivotFactor(const doublereal& d)
373 {
374  if (::solver[currSolver].s_pivot_factor == -1.) {
375  return false;
376  }
377 
378  dPivotFactor = d;
379 
380  return true;
381 }
382 
383 bool
384 LinSol::SetDropTolerance(const doublereal& d)
385 {
386  if (::solver[currSolver].s_drop_tolerance == -1.) {
387  return false;
388  }
389 
390  dDropTolerance = d;
391 
392  return true;
393 }
394 
395 unsigned
396 LinSol::GetBlockSize(void) const
397 {
398  return blockSize;
399 }
400 
401 bool
402 LinSol::SetBlockSize(unsigned bs)
403 {
404  switch (currSolver) {
405  case LinSol::UMFPACK_SOLVER:
406  blockSize = bs;
407  break;
408 
409  default:
410  return false;
411  }
412 
413  return true;
414 }
415 
416 bool
417 LinSol::SetScale(const SolutionManager::ScaleOpt& s)
418 {
419  switch (currSolver) {
420  case LinSol::NAIVE_SOLVER:
421  case LinSol::UMFPACK_SOLVER:
422  case LinSol::KLU_SOLVER:
423  scale = s;
424  break;
425 
426  default:
427  return false;
428  }
429 
430  return true;
431 }
432 
433 integer
434 LinSol::GetMaxIterations(void) const
435 {
436  return iMaxIter;
437 }
438 
439 bool
440 LinSol::SetMaxIterations(integer iMaxIterations)
441 {
442  switch (currSolver) {
443  case LinSol::UMFPACK_SOLVER:
444  iMaxIter = iMaxIterations;
445  break;
446 
447  default:
448  return false;
449  }
450 
451  return true;
452 }
453 
454 SolutionManager *const
455 LinSol::GetSolutionManager(integer iNLD, integer iLWS) const
456 {
457  SolutionManager *pCurrSM = NULL;
458  const unsigned type = (solverFlags & LinSol::SOLVER_FLAGS_TYPE_MASK);
459  const unsigned perm = (solverFlags & LinSol::SOLVER_FLAGS_PERM_MASK);
460  const bool mt = (solverFlags & LinSol::SOLVER_FLAGS_ALLOWS_MT_FCT);
461 
462  /* silence warning */
463  if (mt) {
464  NO_OP;
465  }
466 
467  ASSERT((::solver[currSolver].s_flags & solverFlags) == solverFlags);
468 
469  if (iLWS == 0) {
470  iLWS = iWorkSpaceSize;
471  }
472 
473  switch (currSolver) {
474  case LinSol::Y12_SOLVER:
475 #ifdef USE_Y12
476  switch (type) {
477  case LinSol::SOLVER_FLAGS_ALLOWS_DIR: {
478  typedef Y12SparseCCSolutionManager<DirCColMatrixHandler<1> > CCSM;
479  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM,
480  CCSM(iNLD, iLWS, dPivotFactor));
481  break;
482  }
483 
484  case LinSol::SOLVER_FLAGS_ALLOWS_CC: {
485  typedef Y12SparseCCSolutionManager<CColMatrixHandler<1> > CCSM;
486  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM,
487  CCSM(iNLD, iLWS, dPivotFactor));
488  break;
489  }
490 
491  default:
492  SAFENEWWITHCONSTRUCTOR(pCurrSM,
493  Y12SparseSolutionManager,
494  Y12SparseSolutionManager(iNLD, iLWS,
495  dPivotFactor));
496  break;
497  }
498  break;
499 #else /* !USE_Y12 */
500  silent_cerr("Configure with --with-y12 "
501  "to enable Y12 solver" << std::endl);
502  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
503 #endif /* !USE_Y12 */
504 
505  case LinSol::SUPERLU_SOLVER:
506 #ifdef USE_SUPERLU
507 #ifdef USE_SUPERLU_MT
508  switch (type) {
509  case LinSol::SOLVER_FLAGS_ALLOWS_DIR: {
510  typedef ParSuperLUSparseCCSolutionManager<DirCColMatrixHandler<0> > CCSM;
511  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM,
512  CCSM(nThreads, iNLD, dPivotFactor));
513  break;
514  }
515 
516  case LinSol::SOLVER_FLAGS_ALLOWS_CC: {
517  typedef ParSuperLUSparseCCSolutionManager<CColMatrixHandler<0> > CCSM;
518  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM,
519  CCSM(nThreads, iNLD, dPivotFactor));
520  break;
521  }
522 
523  default:
524  SAFENEWWITHCONSTRUCTOR(pCurrSM,
525  ParSuperLUSparseSolutionManager,
526  ParSuperLUSparseSolutionManager(nThreads, iNLD,
527  dPivotFactor));
528  break;
529  }
530  break;
531 
532  if (nThreads == 1) {
533  silent_cerr("warning, using multithread SuperLU with only one thread; "
534  << std::endl);
535  }
536 #else /* !USE_SUPERLU_MT */
537  if (nThreads > 1 && mt) {
538  silent_cerr("multithread SuperLU solver support not compiled; "
539  << std::endl);
540  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
541  }
542  switch (type) {
543  case LinSol::SOLVER_FLAGS_ALLOWS_DIR: {
544  typedef SuperLUSparseCCSolutionManager<DirCColMatrixHandler<0> > CCSM;
545  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM,
546  CCSM(iNLD, dPivotFactor));
547  break;
548  }
549 
550  case LinSol::SOLVER_FLAGS_ALLOWS_CC: {
551  typedef SuperLUSparseCCSolutionManager<CColMatrixHandler<0> > CCSM;
552  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM,
553  CCSM(iNLD, dPivotFactor));
554  break;
555  }
556 
557  default:
558  SAFENEWWITHCONSTRUCTOR(pCurrSM,
559  SuperLUSparseSolutionManager,
560  SuperLUSparseSolutionManager(iNLD,
561  dPivotFactor));
562  break;
563  }
564 #endif /* !USE_SUPERLU_MT */
565  break;
566 #else /* !USE_SUPERLU */
567  silent_cerr("Configure with --with-superlu "
568  "to enable superlu solver" << std::endl);
569  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
570 #endif /* !USE_SUPERLU */
571 
572  case LinSol::MESCHACH_SOLVER:
573 #ifdef USE_MESCHACH
574  SAFENEWWITHCONSTRUCTOR(pCurrSM,
575  MeschachSparseSolutionManager,
576  MeschachSparseSolutionManager(iNLD, iLWS,
577  dPivotFactor));
578  break;
579 #else /* !USE_MESCHACH */
580  silent_cerr("Configure with --with-meschach "
581  "to enable Meschach solver" << std::endl);
582  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
583 #endif /* !USE_MESCHACH */
584 
585  case LinSol::LAPACK_SOLVER:
586 #ifdef USE_LAPACK
587  SAFENEWWITHCONSTRUCTOR(pCurrSM,
588  LapackSolutionManager,
589  LapackSolutionManager(iNLD, dPivotFactor));
590  break;
591 #else /* !USE_LAPACK */
592  silent_cerr("Configure with --with-lapack "
593  "to enable Lapack dense solver" << std::endl);
594  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
595 #endif /* !USE_LAPACK */
596 
597  case LinSol::TAUCS_SOLVER:
598 #ifdef USE_TAUCS
599  switch (type) {
600  case LinSol::SOLVER_FLAGS_ALLOWS_DIR: {
601  typedef TaucsSparseCCSolutionManager<DirCColMatrixHandler<0> > CCSM;
602  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM, CCSM(iNLD));
603  break;
604  }
605 
606  case LinSol::SOLVER_FLAGS_ALLOWS_CC: {
607  typedef TaucsSparseCCSolutionManager<CColMatrixHandler<0> > CCSM;
608  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM, CCSM(iNLD));
609  break;
610  }
611 
612  default:
613  SAFENEWWITHCONSTRUCTOR(pCurrSM,
614  TaucsSparseSolutionManager,
615  TaucsSparseSolutionManager(iNLD));
616  break;
617  }
618  break;
619 #else /* !USE_TAUCS */
620  silent_cerr("Configure with --with-taucs "
621  "to enable Taucs sparse solver" << std::endl);
622  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
623 #endif /* !USE_TAUCS */
624 
625  case LinSol::HARWELL_SOLVER:
626 #ifdef USE_HARWELL
627  SAFENEWWITHCONSTRUCTOR(pCurrSM,
628  HarwellSparseSolutionManager,
629  HarwellSparseSolutionManager(iNLD, iLWS,
630  dPivotFactor));
631  break;
632 #else /* !USE_HARWELL */
633  silent_cerr("Configure with --with-harwell "
634  "to enable Harwell solver" << std::endl);
635  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
636 #endif /* !USE_HARWELL */
637 
638  case LinSol::UMFPACK_SOLVER:
639 #ifdef USE_UMFPACK
640  {
641  switch (type) {
642  case LinSol::SOLVER_FLAGS_ALLOWS_DIR: {
643  typedef UmfpackSparseCCSolutionManager<DirCColMatrixHandler<0> > CCSM;
644  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM,
645  CCSM(iNLD, dPivotFactor, dDropTolerance, blockSize, scale, iMaxIter));
646  break;
647  }
648 
649  case LinSol::SOLVER_FLAGS_ALLOWS_CC: {
650  typedef UmfpackSparseCCSolutionManager<CColMatrixHandler<0> > CCSM;
651  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM,
652  CCSM(iNLD, dPivotFactor, dDropTolerance, blockSize, scale, iMaxIter));
653  break;
654  }
655 
656  default:
657  SAFENEWWITHCONSTRUCTOR(pCurrSM,
658  UmfpackSparseSolutionManager,
659  UmfpackSparseSolutionManager(iNLD,
660  dPivotFactor, dDropTolerance, blockSize, scale, iMaxIter));
661  break;
662  }
663  } break;
664 #else /* !USE_UMFPACK */
665  silent_cerr("Configure with --with-umfpack "
666  "to enable Umfpack solver" << std::endl);
667  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
668 #endif /* !USE_UMFPACK */
669 
670  case LinSol::KLU_SOLVER:
671 #ifdef USE_KLU
672  {
673  switch (type) {
674  case LinSol::SOLVER_FLAGS_ALLOWS_DIR: {
675  typedef KLUSparseCCSolutionManager<DirCColMatrixHandler<0> > CCSM;
676  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM,
677  CCSM(iNLD, dPivotFactor, scale));
678  break;
679  }
680 
681  case LinSol::SOLVER_FLAGS_ALLOWS_CC: {
682  typedef KLUSparseCCSolutionManager<CColMatrixHandler<0> > CCSM;
683  SAFENEWWITHCONSTRUCTOR(pCurrSM, CCSM,
684  CCSM(iNLD, dPivotFactor, scale));
685  break;
686  }
687 
688  default:
689  SAFENEWWITHCONSTRUCTOR(pCurrSM,
690  KLUSparseSolutionManager,
691  KLUSparseSolutionManager(iNLD,
692  dPivotFactor,
693  scale));
694  break;
695  }
696 
697  } break;
698 #else /* !USE_KLU */
699  silent_cerr("Configure with --with-klu "
700  "to enable KLU solver" << std::endl);
701  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
702 #endif /* !USE_KLU */
703 
704  case LinSol::NAIVE_SOLVER:
705  if (perm == LinSol::SOLVER_FLAGS_ALLOWS_COLAMD) {
706  if (nThreads == 1) {
707  SAFENEWWITHCONSTRUCTOR(pCurrSM,
708  NaiveSparsePermSolutionManager<Colamd_ordering>,
709  NaiveSparsePermSolutionManager<Colamd_ordering>(iNLD, dPivotFactor, scale));
710  } else {
711 #ifdef USE_NAIVE_MULTITHREAD
712  SAFENEWWITHCONSTRUCTOR(pCurrSM,
713  ParNaiveSparsePermSolutionManager,
714  ParNaiveSparsePermSolutionManager(nThreads, iNLD, dPivotFactor));
715 #else /* ! USE_NAIVE_MULTITHREAD */
716  silent_cerr("multithread naive solver support not compiled; "
717  "you can configure --enable-multithread-naive "
718  "on a linux ix86 to get it"
719  << std::endl);
720  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
721 #endif /* ! USE_NAIVE_MULTITHREAD */
722  }
723 #ifdef USE_BOOST
724 #ifdef HAVE_BOOST_GRAPH_CUTHILL_MCKEE_ORDERING_HPP
725  } else if (perm == LinSol::SOLVER_FLAGS_ALLOWS_REVERSE_CUTHILL_MC_KEE) {
726  if (nThreads == 1) {
727  SAFENEWWITHCONSTRUCTOR(pCurrSM,
728  NaiveSparsePermSolutionManager<rcmk_ordering>,
729  NaiveSparsePermSolutionManager<rcmk_ordering>(iNLD, dPivotFactor, scale));
730  } else {
731 #ifdef USE_NAIVE_MULTITHREAD
732 #if 0
733  SAFENEWWITHCONSTRUCTOR(pCurrSM,
734  ParNaiveSparsePermSolutionManager,
735  ParNaiveSparsePermSolutionManager(nThreads, iNLD, dPivotFactor));
736 #endif
737  silent_cerr("multithread naive solver with"
738  "reverse Cuthill-McKee permutation not"
739  "available yet. Patches welcome"
740  << std::endl);
741  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
742 #else /* ! USE_NAIVE_MULTITHREAD */
743  silent_cerr("multithread naive solver support not compiled; "
744  "you can configure --enable-multithread-naive "
745  "on a linux ix86 to get it"
746  << std::endl);
747  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
748 #endif /* ! USE_NAIVE_MULTITHREAD */
749  }
750 #endif /* HAVE_BOOST_GRAPH_CUTHILL_MCKEE_ORDERING_HPP */
751 
752 #if 0 /* ?!? */
753  } else if (perm == LinSol::SOLVER_FLAGS_ALLOWS_MMDATA) {
754  if (nThreads == 1) {
755  SAFENEWWITHCONSTRUCTOR(pCurrSM,
756  NaiveSparsePermSolutionManager<amd_ordering>,
757  NaiveSparsePermSolutionManager<amd_ordering>(iNLD, dPivotFactor, scale));
758  } else {
759 #ifdef USE_NAIVE_MULTITHREAD
760 #if 0
761  SAFENEWWITHCONSTRUCTOR(pCurrSM,
762  ParNaiveSparsePermSolutionManager,
763  ParNaiveSparsePermSolutionManager(nThreads, iNLD, dPivotFactor));
764 #endif
765  silent_cerr("multithread naive solver with"
766  "approximate minimum degree permutation not"
767  "available yet. Patches welcome"
768  << std::endl);
769  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
770 #else
771  silent_cerr("multithread naive solver support not compiled; "
772  "you can configure --enable-multithread-naive "
773  "on a linux ix86 to get it"
774  << std::endl);
775  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
776 #endif /* USE_NAIVE_MULTITHREAD */
777  }
778 #endif
779 
780 #ifdef HAVE_BOOST_GRAPH_KING_ORDERING_HPP
781  } else if (perm == LinSol::SOLVER_FLAGS_ALLOWS_KING) {
782  if (nThreads == 1) {
783  SAFENEWWITHCONSTRUCTOR(pCurrSM,
784  NaiveSparsePermSolutionManager<king_ordering>,
785  NaiveSparsePermSolutionManager<king_ordering>(iNLD, dPivotFactor, scale));
786  } else {
787 #ifdef USE_NAIVE_MULTITHREAD
788 #if 0
789  SAFENEWWITHCONSTRUCTOR(pCurrSM,
790  ParNaiveSparsePermSolutionManager,
791  ParNaiveSparsePermSolutionManager(nThreads, iNLD, dPivotFactor));
792 #endif
793  silent_cerr("multithread naive solver with"
794  "king permutation not"
795  "available yet. Patches welcome"
796  << std::endl);
797  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
798 #else /* ! USE_NAIVE_MULTITHREAD */
799  silent_cerr("multithread naive solver support not compiled; "
800  "you can configure --enable-multithread-naive "
801  "on a linux ix86 to get it"
802  << std::endl);
803  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
804 #endif /* ! USE_NAIVE_MULTITHREAD */
805  }
806 #endif /* HAVE_BOOST_GRAPH_KING_ORDERING_HPP */
807 
808 #ifdef HAVE_BOOST_GRAPH_SLOAN_ORDERING_HPP
809  } else if (perm == LinSol::SOLVER_FLAGS_ALLOWS_SLOAN) {
810  if (nThreads == 1) {
811  SAFENEWWITHCONSTRUCTOR(pCurrSM,
812  NaiveSparsePermSolutionManager<sloan_ordering>,
813  NaiveSparsePermSolutionManager<sloan_ordering>(iNLD, dPivotFactor, scale));
814  } else {
815 #ifdef USE_NAIVE_MULTITHREAD
816 #if 0
817  SAFENEWWITHCONSTRUCTOR(pCurrSM,
818  ParNaiveSparsePermSolutionManager,
819  ParNaiveSparsePermSolutionManager(nThreads, iNLD, dPivotFactor));
820 #endif
821  silent_cerr("multithread naive solver with"
822  "sloan permutation not"
823  "available yet. Patches welcome"
824  << std::endl);
825  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
826 #else /* ! USE_NAIVE_MULTITHREAD */
827  silent_cerr("multithread naive solver support not compiled; "
828  "you can configure --enable-multithread-naive "
829  "on a linux ix86 to get it"
830  << std::endl);
831  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
832 #endif /* ! USE_NAIVE_MULTITHREAD */
833  }
834 #endif /* HAVE_BOOST_GRAPH_SLOAN_ORDERING_HPP */
835 
836 #ifdef HAVE_BOOST_GRAPH_MINIMUM_DEGREE_ORDERING_HPP
837  } else if (perm == LinSol::SOLVER_FLAGS_ALLOWS_MDAPLUSAT) {
838  if (nThreads == 1) {
839  SAFENEWWITHCONSTRUCTOR(pCurrSM,
840  NaiveSparsePermSolutionManager<md_ordering>,
841  NaiveSparsePermSolutionManager<md_ordering>(iNLD, dPivotFactor, scale));
842  } else {
843 #ifdef USE_NAIVE_MULTITHREAD
844 #if 0
845  SAFENEWWITHCONSTRUCTOR(pCurrSM,
846  ParNaiveSparsePermSolutionManager,
847  ParNaiveSparsePermSolutionManager(nThreads, iNLD, dPivotFactor));
848 #endif
849  silent_cerr("multithread naive solver with"
850  "md permutation not"
851  "available yet. Patches welcome"
852  << std::endl);
853  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
854 #else /* ! USE_NAIVE_MULTITHREAD */
855  silent_cerr("multithread naive solver support not compiled; "
856  "you can configure --enable-multithread-naive "
857  "on a linux ix86 to get it"
858  << std::endl);
859  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
860 #endif /* ! USE_NAIVE_MULTITHREAD */
861  }
862 #endif /* HAVE_BOOST_GRAPH_MINIMUM_DEGREE_ORDERING_HPP */
863 #endif /* USE_BOOST */
864 
865  } else if (perm == LinSol::SOLVER_FLAGS_ALLOWS_NESTED_DISSECTION) {
866 #ifdef USE_METIS
867  if (nThreads == 1) {
868  SAFENEWWITHCONSTRUCTOR(pCurrSM,
869  NaiveSparsePermSolutionManager<metis_ordering>,
870  NaiveSparsePermSolutionManager<metis_ordering>(iNLD, dPivotFactor, scale));
871  } else {
872 #ifdef USE_NAIVE_MULTITHREAD
873 #if 0
874  SAFENEWWITHCONSTRUCTOR(pCurrSM,
875  ParNaiveSparsePermSolutionManager,
876  ParNaiveSparsePermSolutionManager(nThreads, iNLD, dPivotFactor));
877 #endif
878  silent_cerr("multithread naive solver with"
879  "nested dissection permutation not"
880  "available yet. Patches welcome"
881  << std::endl);
882  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
883 #else /* ! USE_NAIVE_MULTITHREAD */
884  silent_cerr("multithread naive solver support not compiled; "
885  "you can configure --enable-multithread-naive "
886  "on a linux ix86 to get it"
887  << std::endl);
888  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
889 #endif /* ! USE_NAIVE_MULTITHREAD */
890  }
891 #else /* ! USE_METIS */
892  silent_cerr("you should not get here("<< __FILE__ << ":" <<
893  __LINE__ << ")" << std::endl);
894  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
895 #endif /* ! USE_METIS */
896 
897  } else {
898  if (nThreads == 1) {
899  SAFENEWWITHCONSTRUCTOR(pCurrSM,
900  NaiveSparseSolutionManager,
901  NaiveSparseSolutionManager(iNLD, dPivotFactor, scale));
902  } else {
903 #ifdef USE_NAIVE_MULTITHREAD
904  SAFENEWWITHCONSTRUCTOR(pCurrSM,
905  ParNaiveSparseSolutionManager,
906  ParNaiveSparseSolutionManager(nThreads, iNLD, dPivotFactor));
907 #else
908  silent_cerr("multithread naive solver support not compiled; "
909  "you can configure --enable-multithread-naive "
910  "on a linux ix86 to get it"
911  << std::endl);
912  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
913 #endif /* USE_NAIVE_MULTITHREAD */
914  }
915  }
916  break;
917 
918  case LinSol::EMPTY_SOLVER:
919  break;
920 
921  default:
922  ASSERT(0);
923  throw ErrGeneric(MBDYN_EXCEPT_ARGS);
924 
925  }
926 
927  return pCurrSM;
928 }
929 
LinSol::SetScale
bool SetScale(const SolutionManager::ScaleOpt &scale)
Definition: linsol.cc:417
LinSol::defaultSolver
static SolverType defaultSolver
Definition: linsol.h:151
LinSol::SetPivotFactor
bool SetPivotFactor(const doublereal &d)
Definition: linsol.cc:372
LinSol::GetSolverFlags
unsigned GetSolverFlags(void) const
Definition: linsol.cc:257
LinSol::MaskSolverFlags
bool MaskSolverFlags(unsigned f)
Definition: linsol.cc:303
MBDYN_EXCEPT_ARGS
#define MBDYN_EXCEPT_ARGS
Definition: except.h:63
LinSol::dPivotFactor
doublereal dPivotFactor
Definition: linsol.h:127
LinSol::GetNumThreads
unsigned GetNumThreads(void) const
Definition: linsol.cc:333
LinSol::GetSolutionManager
SolutionManager *const GetSolutionManager(integer iNLD, integer iLWS=0) const
Definition: linsol.cc:455
LinSol::SetNumThreads
bool SetNumThreads(unsigned nt)
Definition: linsol.cc:314
NO_OP
#define NO_OP
Definition: myassert.h:74
LinSol::solver_t
Definition: linsol.h:84
LinSol::blockSize
unsigned blockSize
Definition: linsol.h:117
LinSol::LinSol
LinSol(void)
Definition: linsol.cc:156
LinSol
Definition: linsol.h:39
LinSol::dGetDropTolerance
const doublereal & dGetDropTolerance(void) const
Definition: linsol.cc:351
LinSol::nThreads
unsigned nThreads
Definition: linsol.h:103
LinSol::solver_t::s_name
const char *const s_name
Definition: linsol.h:85
LinSol::SetSolver
bool SetSolver(SolverType t, unsigned f=SOLVER_FLAGS_NONE)
Definition: linsol.cc:181
solver
const LinSol::solver_t solver[]
Definition: linsol.cc:58
LinSol::GetMaxIterations
integer GetMaxIterations(void) const
Definition: linsol.cc:434
LinSol::iMaxIter
integer iMaxIter
Definition: linsol.h:148
LinSol::SetSolverFlags
bool SetSolverFlags(unsigned f)
Definition: linsol.cc:281
LinSol::SetMaxIterations
bool SetMaxIterations(integer iMaxIter)
Definition: linsol.cc:440
LinSol::SetWorkSpaceSize
bool SetWorkSpaceSize(integer)
Definition: linsol.cc:357
LinSol::AddSolverFlags
bool AddSolverFlags(unsigned f)
Definition: linsol.cc:292
LinSol::~LinSol
virtual ~LinSol(void)
Definition: linsol.cc:169
LinSol::dGetPivotFactor
const doublereal & dGetPivotFactor(void) const
Definition: linsol.cc:345
LinSol::scale
SolutionManager::ScaleOpt scale
Definition: linsol.h:141
LinSol::iWorkSpaceSize
integer iWorkSpaceSize
Definition: linsol.h:110
LinSol::solverFlags
unsigned solverFlags
Definition: linsol.h:96
ASSERT
#define ASSERT(expression)
Definition: colamd.c:977
SAFENEWWITHCONSTRUCTOR
#define SAFENEWWITHCONSTRUCTOR(pnt, item, constructor)
Definition: mynewmem.h:698
LinSol::currSolver
SolverType currSolver
Definition: linsol.h:95
LinSol::SetBlockSize
bool SetBlockSize(unsigned bs)
Definition: linsol.cc:402
LinSol::GetSolver
SolverType GetSolver(void) const
Definition: linsol.cc:175
LinSol::dDropTolerance
doublereal dDropTolerance
Definition: linsol.h:134
LinSol::iGetWorkSpaceSize
integer iGetWorkSpaceSize(void) const
Definition: linsol.cc:339
LinSol::SolverType
SolverType
Definition: linsol.h:41
LinSol::SetDropTolerance
bool SetDropTolerance(const doublereal &d)
Definition: linsol.cc:384
LinSol::GetBlockSize
unsigned GetBlockSize(void) const
Definition: linsol.cc:396
LinSol::GetSolverName
const char *const GetSolverName(void) const
Definition: linsol.cc:269
doublereal
double doublereal
Definition: colamd.c:52
integer
long int integer
Definition: colamd.c:51
LinSol::solver_t::s_flags
unsigned s_flags
Definition: linsol.h:88