Hex  1.0
Hydrogen-electron collision solver
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
preconditioners.h
Go to the documentation of this file.
1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * *\
2  * *
3  * / / / / __ \ \ / / *
4  * / /__ / / / _ \ \ \/ / *
5  * / ___ / | |/_/ / /\ \ *
6  * / / / / \_\ / / \ \ *
7  * *
8  * Jakub Benda (c) 2014 *
9  * Charles University in Prague *
10  * *
11 \* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
12 
13 #ifndef HEX_PRECONDITIONERS
14 #define HEX_PRECONDITIONERS
15 
16 #include <tuple>
17 
18 #include "arrays.h"
19 #include "input.h"
20 #include "matrix.h"
21 #include "radial.h"
22 #include "parallel.h"
23 
24 #ifndef NO_OPENCL
25  #include "opencl.h"
26 #endif
27 
59 cArray iChol (cArrayView const & A, lArrayView const & I, lArrayView const & P);
60 
105 SymDiaMatrix DIC (SymDiaMatrix const & A);
106 
115 SymDiaMatrix SSOR (SymDiaMatrix const & A);
116 
126 CooMatrix SPAI (SymDiaMatrix const & A, const iArrayView diagonals);
127 
138 {
139  public:
140 
144  virtual RadialIntegrals const & rad() const = 0;
145 
153  virtual void setup () = 0;
154 
161  virtual void update (double E) = 0;
162 
166  virtual void rhs (cArrayView chi, int ienergy, int instate) const = 0;
167 
174  virtual void multiply (const cArrayView p, cArrayView q) const = 0;
175 
186  virtual void precondition (const cArrayView r, cArrayView z) const = 0;
187 };
188 
201 {
202  public:
203 
204  static const std::string name;
205  static const std::string description;
206 
208  Parallel const & par,
209  InputFile const & inp,
210  std::vector<std::pair<int,int>> const & ll,
211  Bspline const & bspline,
212  CommandLine const & cmd
213  ) : PreconditionerBase(), cmd_(cmd), par_(par), inp_(inp), l1_l2_(ll),
214  s_bspline_(bspline), s_rad_(s_bspline_) {}
215 
216  virtual RadialIntegrals const & rad () const { return s_rad_; }
217 
218  virtual void setup ();
219  virtual void update (double E);
220  virtual void rhs (cArrayView chi, int ienergy, int instate) const;
221  virtual void multiply (const cArrayView p, cArrayView q) const;
222  virtual void precondition (const cArrayView r, cArrayView z) const;
223 
224  protected:
225 
226  // command line switches
227  CommandLine const & cmd_;
228 
229  // parallel environment
230  Parallel const & par_;
231 
232  // input parameters
233  InputFile const & inp_;
234 
235  // coupled states
236  std::vector<std::pair<int,int>> const & l1_l2_;
237 
238  // diagonal blocks in DIA format (these will be used in matrix multiplication)
239  mutable std::vector<SymDiaMatrix> dia_blocks_;
240 
241  // B-spline environment for the solution
243 
244  // radial integrals for the solution
246 
247  // Kronecker products
251 };
252 
261 {
262  public:
263 
264  static const std::string name;
265  static const std::string description;
266 
268  Parallel const & par,
269  InputFile const & inp,
270  std::vector<std::pair<int,int>> const & ll,
271  Bspline const & bspline,
272  CommandLine const & cmd
273  ) : NoPreconditioner(par, inp, ll, bspline, cmd) {}
274 
275  // reuse parent definitions
276  virtual RadialIntegrals const & rad () const { return NoPreconditioner::rad(); }
277  virtual void setup () { return NoPreconditioner::setup(); }
278  virtual void update (double E) { return NoPreconditioner::update(E); }
279  virtual void rhs (cArrayView chi, int ienergy, int instate) const { NoPreconditioner::rhs(chi, ienergy, instate); }
280  virtual void multiply (const cArrayView p, cArrayView q) const { NoPreconditioner::multiply(p, q); }
281 
282  // declare own definitions
283  virtual void precondition (const cArrayView r, cArrayView z) const;
284 
285  // inner CG callbacks
286  virtual void CG_mmul (int iblock, const cArrayView p, cArrayView q) const;
287  virtual void CG_prec (int iblock, const cArrayView r, cArrayView z) const;
288 };
289 
297 #ifndef NO_OPENCL
299 {
300  public:
301 
302  static const std::string name;
303  static const std::string description;
304 
306  Parallel const & par,
307  InputFile const & inp,
308  std::vector<std::pair<int,int>> const & ll,
309  Bspline const & bspline,
310  CommandLine const & cmd
311  ) : NoPreconditioner(par, inp, ll, bspline, cmd) {}
312 
313  // reuse parent definitions
314  virtual RadialIntegrals const & rad () const { return NoPreconditioner::rad(); }
315  virtual void setup ();
316  virtual void update (double E);
317  virtual void rhs (cArrayView chi, int ienergy, int instate) const { NoPreconditioner::rhs(chi, ienergy, instate); }
318  virtual void multiply (const cArrayView p, cArrayView q) const { NoPreconditioner::multiply(p, q); }
319 
320  // declare own definitions
321  virtual void precondition (const cArrayView r, cArrayView z) const;
322 
323  private:
324 
325  // diagonal blocks
326  mutable std::vector<CsrMatrix> csr_blocks_;
327  mutable std::vector<cArray> block_;
328 
329  // OpenCL environment
330  cl_platform_id platform_;
331  cl_device_id device_;
332  cl_context context_;
333  cl_command_queue queue_;
334  cl_program program_;
335 
336  // size of a workgroup
337  size_t Nlocal_;
338 
339  // computational kernels
340  cl_kernel mmul_;
341  cl_kernel amul_;
342  cl_kernel axby_;
343  cl_kernel vnrm_;
344  cl_kernel norm_;
345  cl_kernel spro_;
346 };
347 #endif
348 
356 {
357  public:
358 
359  static const std::string name;
360  static const std::string description;
361 
363  Parallel const & par,
364  InputFile const & inp,
365  std::vector<std::pair<int,int>> const & ll,
366  Bspline const & bspline,
367  CommandLine const & cmd
368  ) : CGPreconditioner(par, inp, ll, bspline, cmd) {}
369 
370  // reuse parent definitions
371  virtual RadialIntegrals const & rad () const { return CGPreconditioner::rad(); }
372  virtual void rhs (cArrayView chi, int ienergy, int instate) const { CGPreconditioner::rhs(chi, ienergy, instate); }
373  virtual void multiply (const cArrayView p, cArrayView q) const { CGPreconditioner::multiply(p, q); }
374  virtual void precondition (const cArrayView r, cArrayView z) const { CGPreconditioner::precondition(r, z); }
375 
376  // declare own definitions
377  virtual void setup ();
378  virtual void update (double E);
379 
380  // inner CG callback (needed by parent)
381  virtual void CG_prec (int iblock, const cArrayView r, cArrayView z) const { z = invd_[iblock] * r; }
382 
383  protected:
384 
385  // inverse diagonals for every block
387 };
388 
396 {
397  public:
398 
399  static const std::string name;
400  static const std::string description;
401 
403  Parallel const & par,
404  InputFile const & inp,
405  std::vector<std::pair<int,int>> const & ll,
406  Bspline const & bspline,
407  CommandLine const & cmd
408  ) : CGPreconditioner(par, inp, ll, bspline, cmd) {}
409 
410  // reuse parent definitions
411  virtual RadialIntegrals const & rad () const { return CGPreconditioner::rad(); }
412  virtual void rhs (cArrayView chi, int ienergy, int instate) const { CGPreconditioner::rhs(chi, ienergy, instate); }
413  virtual void multiply (const cArrayView p, cArrayView q) const { CGPreconditioner::multiply(p, q); }
414  virtual void precondition (const cArrayView r, cArrayView z) const { CGPreconditioner::precondition(r, z); }
415 
416  // declare own definitions
417  virtual void setup ();
418  virtual void update (double E);
419 
420  // inner CG callback (needed by parent)
421  virtual void CG_prec (int iblock, const cArrayView r, cArrayView z) const;
422 
423  protected:
424 
425  // inverse diagonals for every block
426  mutable std::vector<SymDiaMatrix> SSOR_;
427 };
428 
437 {
438  public:
439 
440  static const std::string name;
441  static const std::string description;
442 
444  Parallel const & par,
445  InputFile const & inp,
446  std::vector<std::pair<int,int>> const & ll,
447  Bspline const & bspline,
448  CommandLine const & cmd
449  ) : CGPreconditioner(par, inp, ll, bspline, cmd), droptol_(cmd.droptol) {}
450 
451  // reuse parent definitions
452  virtual RadialIntegrals const & rad () const { return CGPreconditioner::rad(); }
453  virtual void multiply (const cArrayView p, cArrayView q) const { CGPreconditioner::multiply(p,q); }
454  virtual void rhs (cArrayView chi, int ienergy, int instate) const { CGPreconditioner::rhs(chi,ienergy,instate); }
455  virtual void precondition (const cArrayView r, cArrayView z) const { CGPreconditioner::precondition(r,z); }
456 
457  // declare own definitions
458  virtual void setup ();
459  virtual void update (double E);
460 
461  // inner CG callback (needed by parent)
462  virtual void CG_prec (int iblock, const cArrayView r, cArrayView z) const;
463 
464  protected:
465 
466  // drop tolarance for the factorizations
467  double droptol_;
468 
469  // diagonal CSR block for every coupled state
470  mutable std::vector<CsrMatrix> csr_blocks_;
471 
472  // LU decompositions of the CSR blocks
473  mutable std::vector<CsrMatrix::LUft> lu_;
474 };
475 
484 {
485  public:
486 
487  static const std::string name;
488  static const std::string description;
489 
491  Parallel const & par,
492  InputFile const & inp,
493  std::vector<std::pair<int,int>> const & ll,
494  Bspline const & bspline,
495  CommandLine const & cmd
496  ) : CGPreconditioner(par, inp, ll, bspline, cmd) {}
497 
498  // reuse parent definitions
499  virtual RadialIntegrals const & rad () const { return CGPreconditioner::rad(); }
500  virtual void multiply (const cArrayView p, cArrayView q) const { CGPreconditioner::multiply(p,q); }
501  virtual void rhs (cArrayView chi, int ienergy, int instate) const { CGPreconditioner::rhs(chi,ienergy,instate); }
502  virtual void precondition (const cArrayView r, cArrayView z) const { CGPreconditioner::precondition(r,z); }
503 
504  // declare own definitions
505  virtual void setup ();
506  virtual void update (double E);
507 
508  // inner CG callback (needed by parent)
509  virtual void CG_prec (int iblock, const cArrayView r, cArrayView z) const
510  {
511  z = DIC_[iblock].upperSolve( DIC_[iblock].dot( DIC_[iblock].lowerSolve(r), diagonal ) );
512  }
513 
514  private:
515 
516  std::vector<SymDiaMatrix> DIC_;
517 };
518 
527 #ifndef NO_LAPACK
529 {
530  public:
531 
532  static const std::string name;
533  static const std::string description;
534 
536  Parallel const & par,
537  InputFile const & inp,
538  std::vector<std::pair<int,int>> const & ll,
539  Bspline const & bspline,
540  CommandLine const & cmd
541  ) : CGPreconditioner(par, inp, ll, bspline, cmd) {}
542 
543  // reuse parent definitions
544  virtual RadialIntegrals const & rad () const { return CGPreconditioner::rad(); }
545  virtual void rhs (cArrayView chi, int ienergy, int instate) const { CGPreconditioner::rhs(chi,ienergy,instate); }
546  virtual void multiply (const cArrayView p, cArrayView q) const { CGPreconditioner::multiply(p,q); }
547  virtual void precondition (const cArrayView r, cArrayView z) const { CGPreconditioner::precondition(r,z); }
548 
549  // declare own definitions
550  virtual void setup ();
551  virtual void update (double E);
552 
553  // inner CG callback (needed by parent)
554  virtual void CG_prec (int iblock, const cArrayView r, cArrayView z) const
555  {
556  z = spai_[iblock].dot(r);
557  }
558 
559  private:
560 
561  // SPAIs for every diagonal block
562  std::vector<CsrMatrix> spai_;
563 };
564 #endif
565 
627 {
628  public:
629 
630  static const std::string name;
631  static const std::string description;
632 
633  // constructor
635  Parallel const & par,
636  InputFile const & inp,
637  std::vector<std::pair<int,int>> const & ll,
638  Bspline const & s_bspline,
639  CommandLine const & cmd
640  ) : SSORCGPreconditioner(par,inp,ll,s_bspline, cmd),
641  p_bspline_ (
642  1, // use first order for preconditioner basis
643  sorted_unique(s_bspline.rknots(), 1), // allow just one consecutive occurence
644  s_bspline.ECStheta(), // copy rotation point
645  sorted_unique(s_bspline.cknots(), 1) // allow just one consecutive occurence
646  ),
648  {}
649 
650  // reuse parent definitions
651  RadialIntegrals const & rad () const { return SSORCGPreconditioner::rad(); }
653 
654  // declare own definitions
655  void setup();
656  void update (double E);
657  void rhs (cArrayView chi, int ienergy, int instate) const;
658  void multiply (const cArrayView p, cArrayView q) const;
659 
660  // inner CG callback
661  virtual void CG_prec (int iblock, const cArrayView r, cArrayView z) const;
662 
663  protected:
664 
665  // compute transfer matrix
666  void computeSigma_();
667 
668  // compute contribution to an element of the transfer matrix
669  Complex computeSigma_iknot_ (int qord, int is, int iknots, int ip, int iknotp) const;
670 
671  // transfer overlap matrix (s|p) multiplied from left by S⁻¹ (in s-basis)
673 
674  // transfer overlap matrix (p|s) multiplied from left by S⁻¹ (in p-basis)
676 
680 
681  // preconditioner objects
682 
683  // - B-spline environment for the preconditioning
685 
686  // - radial integrals for the preconditioning
688 
689  // diagonal blocks in the preconditioner basis
690  std::vector<CsrMatrix> p_csr_;
691 
692  // factorization of p_csr
693  std::vector<CsrMatrix::LUft> p_lu_;
694 
695  // - Kronecker products
698 
699  // integrator
701 };
702 
709 {
710  public:
711 
712  static const std::string name;
713  static const std::string description;
714 
715  // constructor
717  Parallel const & par,
718  InputFile const & inp,
719  std::vector<std::pair<int,int>> const & ll,
720  Bspline const & bspline,
721  CommandLine const & cmd
722  );
723 
724  // destructor
726 
727  // use data from the highest order
728  RadialIntegrals const & rad() const
729  {
730  return p_.back()->rad();
731  }
732 
733  // setup all levels
734  void setup ()
735  {
736  for (PreconditionerBase* p : p_)
737  p->setup();
738  }
739 
740  // update all levels
741  void update (double E)
742  {
743  for (PreconditionerBase* p : p_)
744  p->update(E);
745  }
746 
747  // use RHS from the highest preconditioner
748  void rhs (cArrayView chi, int ienergy, int instate) const
749  {
750  p_.back()->rhs(chi, ienergy, instate);
751  }
752 
753  // multiply by the matrix of the highest order
754  void multiply (const cArrayView p, cArrayView q) const
755  {
756  p_.back()->multiply(p,q);
757  }
758 
759  // execute the preconditioner
760  void precondition (const cArrayView r, cArrayView z) const;
761 
762  private:
763 
764  // array of per-level preconditioners
765  std::vector<PreconditionerBase*> p_;
766 };
767 
791 {
792  public:
793 
801  typedef std::tuple <
802  ILUCGPreconditioner, // Solve diagonal blocks by drop-tolerance incomplete LU factorization.
803  NoPreconditioner, // No preconditioner.
804  CGPreconditioner, // Solve diagonal blocks by non-preconditioned CG iterations.
805  JacobiCGPreconditioner, // Solve diagonal blocks by Jacobi-preconditioned CG iterations.
806 #ifndef NO_OPENCL
807  GPUCGPreconditioner, // Solve diagonal blocks by Jacobi-preconditioned CG iterations (GPU variant).
808 #endif
809  SSORCGPreconditioner // Solve diagonal blocks by SSOR-preconditioned CG iterations.
810 // MultiresPreconditioner, // Multi-resolution preconditioner.
811 // SPAICGPreconditioner, // Solve diagonal blocks by SPAI-preconditioned CG iterations.
812 // TwoLevelPreconditioner, // Solve diagonal blocks by two-level precondtioned CG iterations.
813 // DICCGPreconditioner // Diagonal incomplete Cholesky factorization + CG iterations.
815 
819  static size_t size ()
820  {
821  return std::tuple_size<AvailableTypes>::value;
822  }
823 
831  template <int i = 0> static inline typename std::enable_if<(i < std::tuple_size<AvailableTypes>::value), PreconditionerBase*>::type choose (
832  Parallel const & par, InputFile const & inp, std::vector<std::pair<int,int>> const & ll, Bspline const & bspline, CommandLine const & cmd
833  ){
834  // check if i-th preconditioner is the right one
835  // - Yes : return new pointer of its type
836  // - No : try the next preconditioner
837  return (i == cmd.preconditioner) ? new typename std::tuple_element<i,AvailableTypes>::type (par, inp, ll, bspline, cmd) : choose<i+1>(par, inp, ll, bspline, cmd);
838  }
839  template <int i = 0> static inline typename std::enable_if<(i == std::tuple_size<AvailableTypes>::value), PreconditionerBase*>::type choose (
840  Parallel const & par, InputFile const & inp, std::vector<std::pair<int,int>> const & ll, Bspline const & bspline, CommandLine const & cmd
841  ){
842  // we visited all available preconditioners without success : return NULL pointer
843  return nullptr;
844  }
846 
859  template <int i = 0> static inline typename std::enable_if<(i < std::tuple_size<AvailableTypes>::value),int>::type findByName (std::string name)
860  {
861  // check if i-th preconditioner is the right one
862  // - Yes : return new pointer of its type
863  // - No : try the next preconditioner
864  return (name == std::tuple_element<i,AvailableTypes>::type::name) ? i : findByName<i+1>(name);
865  }
866  template <int i = 0> static inline typename std::enable_if<(i == std::tuple_size<AvailableTypes>::value),int>::type findByName (std::string name)
867  {
868  // we visited all available preconditioners without success : return invalid index (-1)
869  return -1;
870  }
872 
877  template <int i = 0> static inline typename std::enable_if<(i < std::tuple_size<AvailableTypes>::value),std::string>::type name (unsigned idx)
878  {
879  return (i == idx) ? std::tuple_element<i,AvailableTypes>::type::name : name<i+1>(idx);
880  }
881  template <int i = 0> static inline typename std::enable_if<(i == std::tuple_size<AvailableTypes>::value),std::string>::type name (unsigned idx)
882  {
883  return "";
884  }
886 
891  template <int i = 0> static inline typename std::enable_if<(i < std::tuple_size<AvailableTypes>::value),std::string>::type description (unsigned idx)
892  {
893  return (i == idx) ? std::tuple_element<i,AvailableTypes>::type::description : description<i+1>(idx);
894  }
895  template <int i = 0> static inline typename std::enable_if<(i == std::tuple_size<AvailableTypes>::value),std::string>::type description (unsigned idx)
896  {
897  return "";
898  }
900 };
901 #endif
ColMatrix< Complex > SpsSigma
Definition: preconditioners.h:679
virtual RadialIntegrals const & rad() const =0
Get radial integrals.
virtual void CG_prec(int iblock, const cArrayView r, cArrayView z) const
Definition: preconditioners.h:381
Array view.
Definition: arrays.h:186
virtual void precondition(const cArrayView r, cArrayView z) const
Precondition the equation.
Definition: preconditioners.h:374
RadialIntegrals const & rad() const
Get radial integrals.
Definition: preconditioners.h:651
static const std::string description
Definition: preconditioners.h:441
Definition: radial.h:76
static std::enable_if<(i< std::tuple_size< AvailableTypes >::value), std::string >::type name(unsigned idx)
Get name of the preconditioner.
Definition: preconditioners.h:877
static const std::string name
Definition: preconditioners.h:630
virtual void rhs(cArrayView chi, int ienergy, int instate) const
Return the right-hand side.
Definition: preconditioners.h:454
SymDiaMatrix S_kron_Mm2_
Definition: preconditioners.h:248
void multiply(const cArrayView p, cArrayView q) const
Multiply by the matrix equation.
std::tuple< ILUCGPreconditioner, NoPreconditioner, CGPreconditioner, JacobiCGPreconditioner, GPUCGPreconditioner, SSORCGPreconditioner > AvailableTypes
List of available preconditioners.
Definition: preconditioners.h:814
virtual void update(double E)
Update the preconditioner for the next energy.
Definition: preconditioners.cpp:1095
RadialIntegrals s_rad_
Definition: preconditioners.h:245
Two-resolution preconditioner.
Definition: preconditioners.h:626
static const std::string description
Definition: preconditioners.h:360
static std::enable_if<(i< std::tuple_size< AvailableTypes >::value), PreconditionerBase * >::type choose(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
Return pointer to new preconditioner object.
Definition: preconditioners.h:831
void update(double E)
Update the preconditioner for the next energy.
std::vector< SymDiaMatrix > dia_blocks_
Definition: preconditioners.h:239
RadialIntegrals p_rad_
Definition: preconditioners.h:687
virtual void precondition(const cArrayView r, cArrayView z) const =0
Precondition the equation.
SymDiaMatrix p_half_D_minus_Mm1_tr_kron_S_
Definition: preconditioners.h:696
virtual RadialIntegrals const & rad() const
Get radial integrals.
Definition: preconditioners.h:216
DICCGPreconditioner(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
Definition: preconditioners.h:490
static const std::string description
Definition: preconditioners.h:205
virtual void rhs(cArrayView chi, int ienergy, int instate) const
Return the right-hand side.
Definition: preconditioners.h:279
virtual RadialIntegrals const & rad() const
Get radial integrals.
Definition: preconditioners.h:276
CommandLine const & cmd_
Definition: preconditioners.h:227
SymDiaMatrix p_S_kron_Mm2_
Definition: preconditioners.h:696
DIC-preconditioned CG-based preconditioner.
Definition: preconditioners.h:483
virtual void precondition(const cArrayView r, cArrayView z) const
Precondition the equation.
Definition: preconditioners.h:547
SSOR-preconditioned CG-based preconditioner.
Definition: preconditioners.h:395
virtual void precondition(const cArrayView r, cArrayView z) const
Precondition the equation.
Definition: preconditioners.h:414
Bspline p_bspline_
Definition: preconditioners.h:684
static const std::string name
Definition: preconditioners.h:532
NumberArray< T > sorted_unique(const ArrayView< T > v, int n=1)
Drop all redundant repetitions from sorted array.
Definition: arrays.h:1935
ColMatrix< Complex > psSigma
Definition: preconditioners.h:679
virtual void setup()
Initialize the preconditioner.
virtual void update(double E)
Update the preconditioner for the next energy.
Definition: preconditioners.cpp:860
SymDiaMatrix DIC(SymDiaMatrix const &A)
DIC preconditioner.
Definition: preconditioners.cpp:94
Multi-resolution (V-cycle) preconditioner.
Definition: preconditioners.h:708
void rhs(cArrayView chi, int ienergy, int instate) const
Return the right-hand side.
SymDiaMatrix p_Mm2_kron_S_
Definition: preconditioners.h:696
virtual void rhs(cArrayView chi, int ienergy, int instate) const =0
Return the right-hand side.
std::vector< CsrMatrix > p_csr_
Definition: preconditioners.h:690
virtual void CG_prec(int iblock, const cArrayView r, cArrayView z) const
Definition: preconditioners.cpp:754
A comfortable number array class.
Definition: arrays.h:171
ILU-preconditioned CG-based preconditioner.
Definition: preconditioners.h:436
static const std::string name
Definition: preconditioners.h:302
virtual void update(double E)
Update the preconditioner for the next energy.
MultiresPreconditioner(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
virtual void update(double E)=0
Update the preconditioner for the next energy.
static const std::string name
Definition: preconditioners.h:359
Solution driver without actual preconditioner.
Definition: preconditioners.h:200
B-spline environment.
Definition: bspline.h:35
static const std::string description
Definition: preconditioners.h:400
static const std::string name
Definition: preconditioners.h:204
RowMatrix< Complex > spSigma_
Definition: preconditioners.h:672
ILUCGPreconditioner(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
Definition: preconditioners.h:443
SymDiaMatrix p_S_kron_half_D_minus_Mm1_tr_
Definition: preconditioners.h:696
void update(double E)
Update the preconditioner for the next energy.
Definition: preconditioners.h:741
std::vector< CsrMatrix > csr_blocks_
Definition: preconditioners.h:470
Gauss-Legendre quadrature.
Definition: gauss.h:33
virtual void multiply(const cArrayView p, cArrayView q) const
Multiply by the matrix equation.
Definition: preconditioners.h:373
cArrays invd_
Definition: preconditioners.h:386
static const std::string description
Definition: preconditioners.h:265
virtual void setup()
Initialize the preconditioner.
Definition: preconditioners.cpp:1047
virtual void CG_prec(int iblock, const cArrayView r, cArrayView z) const
Definition: preconditioners.cpp:1120
Parallel const & par_
Definition: preconditioners.h:230
virtual void precondition(const cArrayView r, cArrayView z) const
Precondition the equation.
Definition: preconditioners.h:502
GaussLegendre g_
Definition: preconditioners.h:700
virtual void CG_prec(int iblock, const cArrayView r, cArrayView z) const
void setup()
Initialize the preconditioner.
virtual void update(double E)
Update the preconditioner for the next energy.
Definition: preconditioners.cpp:405
static const std::string name
Definition: preconditioners.h:440
void precondition(const cArrayView r, cArrayView z) const
Precondition the equation.
Definition: preconditioners.h:652
ColMatrix< Complex > SspSigma
Definition: preconditioners.h:678
static std::enable_if<(i< std::tuple_size< AvailableTypes >::value), int >::type findByName(std::string name)
Find preconditioner by name.
Definition: preconditioners.h:859
virtual void precondition(const cArrayView r, cArrayView z) const
Precondition the equation.
Definition: preconditioners.h:455
virtual void precondition(const cArrayView r, cArrayView z) const
Precondition the equation.
Definition: preconditioners.cpp:693
virtual void multiply(const cArrayView p, cArrayView q) const
Multiply by the matrix equation.
Definition: preconditioners.h:500
static std::enable_if<(i< std::tuple_size< AvailableTypes >::value), std::string >::type description(unsigned idx)
Get description of the preconditioner.
Definition: preconditioners.h:891
Command line parameters.
Definition: input.h:34
static const std::string description
Definition: preconditioners.h:533
virtual RadialIntegrals const & rad() const
Get radial integrals.
Definition: preconditioners.h:371
SymDiaMatrix SSOR(SymDiaMatrix const &A)
SSOR preconditioner.
Definition: preconditioners.cpp:175
ColMatrix< Complex > spSigma
DEBUG.
Definition: preconditioners.h:678
virtual RadialIntegrals const & rad() const
Get radial integrals.
Definition: preconditioners.h:499
virtual void CG_prec(int iblock, const cArrayView r, cArrayView z) const
Definition: preconditioners.h:509
virtual void rhs(cArrayView chi, int ienergy, int instate) const
Return the right-hand side.
Definition: preconditioners.h:317
static size_t size()
Number of available preconditioners.
Definition: preconditioners.h:819
virtual void setup()
Initialize the preconditioner.
Definition: preconditioners.h:277
InputFile const & inp_
Definition: preconditioners.h:233
virtual void CG_mmul(int iblock, const cArrayView p, cArrayView q) const
Definition: preconditioners.cpp:735
int preconditioner
Preconditioner to use. See Preconditioners::AvailableTypes for available types.
Definition: input.h:99
std::vector< CsrMatrix::LUft > lu_
Definition: preconditioners.h:473
virtual void multiply(const cArrayView p, cArrayView q) const
Multiply by the matrix equation.
Definition: preconditioners.h:413
void setup()
Initialize the preconditioner.
Definition: preconditioners.h:734
SymDiaMatrix half_D_minus_Mm1_tr_kron_S_
Definition: preconditioners.h:248
static const std::string description
Definition: preconditioners.h:631
Definition: matrix.h:518
cArray iChol(cArrayView const &A, lArrayView const &I, lArrayView const &P)
Sparse incomplete Cholesky decomposition.
virtual RadialIntegrals const & rad() const
Get radial integrals.
Definition: preconditioners.h:452
Jacobi-preconditioned CG-based preconditioner.
Definition: preconditioners.h:355
virtual void setup()
Initialize the preconditioner.
virtual void rhs(cArrayView chi, int ienergy, int instate) const
Return the right-hand side.
Definition: preconditioners.h:545
CG iteration-based preconditioner (GPU variant).
Definition: preconditioners.h:298
void rhs(cArrayView chi, int ienergy, int instate) const
Return the right-hand side.
Definition: preconditioners.h:748
NoPreconditioner(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
Definition: preconditioners.h:207
CG iteration-based preconditioner.
Definition: preconditioners.h:260
Preconditioner traits.
Definition: preconditioners.h:790
Input parameters.
Definition: input.h:147
std::vector< std::pair< int, int > > const & l1_l2_
Definition: preconditioners.h:236
RowMatrix< Complex > psSigma_
Definition: preconditioners.h:675
CGPreconditioner(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
Definition: preconditioners.h:267
virtual void multiply(const cArrayView p, cArrayView q) const
Multiply by the matrix equation.
Definition: preconditioners.cpp:601
static const std::string name
Definition: preconditioners.h:399
Complex computeSigma_iknot_(int qord, int is, int iknots, int ip, int iknotp) const
virtual void rhs(cArrayView chi, int ienergy, int instate) const
Return the right-hand side.
Definition: preconditioners.h:372
void precondition(const cArrayView r, cArrayView z) const
Precondition the equation.
MPI info.
Definition: parallel.h:29
virtual void multiply(const cArrayView p, cArrayView q) const
Multiply by the matrix equation.
Definition: preconditioners.h:453
std::vector< CsrMatrix::LUft > p_lu_
Definition: preconditioners.h:693
virtual RadialIntegrals const & rad() const
Get radial integrals.
Definition: preconditioners.h:544
SymDiaMatrix Mm2_kron_S_
Definition: preconditioners.h:248
virtual void precondition(const cArrayView r, cArrayView z) const
Precondition the equation.
Definition: preconditioners.cpp:896
std::vector< SymDiaMatrix > SSOR_
Definition: preconditioners.h:426
static const std::string description
Definition: preconditioners.h:713
SymDiaMatrix S_kron_half_D_minus_Mm1_tr_
Definition: preconditioners.h:248
SPAI-preconditioned CG-based preconditioner.
Definition: preconditioners.h:528
SymDiaMatrix S_kron_Mm1_tr_
Definition: preconditioners.h:248
static const std::string name
Definition: preconditioners.h:264
virtual void update(double E)
Update the preconditioner for the next energy.
Definition: preconditioners.h:278
virtual void CG_prec(int iblock, const cArrayView r, cArrayView z) const
Definition: preconditioners.h:554
virtual void multiply(const cArrayView p, cArrayView q) const
Multiply by the matrix equation.
Definition: preconditioners.h:280
virtual void multiply(const cArrayView p, cArrayView q) const =0
Multiply by the matrix equation.
RadialIntegrals const & rad() const
Get radial integrals.
Definition: preconditioners.h:728
static const std::string name
Definition: preconditioners.h:712
TwoLevelPreconditioner(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &s_bspline, CommandLine const &cmd)
Definition: preconditioners.h:634
virtual void multiply(const cArrayView p, cArrayView q) const
Multiply by the matrix equation.
Definition: preconditioners.h:318
virtual RadialIntegrals const & rad() const
Get radial integrals.
Definition: preconditioners.h:411
virtual RadialIntegrals const & rad() const
Get radial integrals.
Definition: preconditioners.h:314
virtual void setup()=0
Initialize the preconditioner.
virtual void CG_prec(int iblock, const cArrayView r, cArrayView z) const
Definition: preconditioners.cpp:1214
CooMatrix SPAI(SymDiaMatrix const &A, iArrayView diagonals)
SPAI preconditioner.
Definition: preconditioners.cpp:273
Complex COO matrix.
Definition: matrix.h:1120
virtual void update(double E)
Update the preconditioner for the next energy.
Definition: preconditioners.cpp:1055
virtual void rhs(cArrayView chi, int ienergy, int instate) const
Return the right-hand side.
Definition: preconditioners.h:412
virtual void update(double E)
Update the preconditioner for the next energy.
virtual void setup()
Initialize the preconditioner.
Definition: preconditioners.cpp:1146
std::complex< double > Complex
Definition: complex.h:20
Symmetric diagonal matrix.
Definition: matrix.h:1547
GPUCGPreconditioner(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
Definition: preconditioners.h:305
double droptol_
Definition: preconditioners.h:467
static const std::string description
Definition: preconditioners.h:488
void multiply(const cArrayView p, cArrayView q) const
Multiply by the matrix equation.
Definition: preconditioners.h:754
SymDiaMatrix S_kron_S_
Definition: preconditioners.h:248
Preconditioner template.
Definition: preconditioners.h:137
virtual void precondition(const cArrayView r, cArrayView z) const
Precondition the equation.
Definition: preconditioners.cpp:701
virtual void update(double E)
Update the preconditioner for the next energy.
Definition: preconditioners.cpp:1155
SSORCGPreconditioner(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
Definition: preconditioners.h:402
virtual void setup()
Initialize the preconditioner.
Definition: preconditioners.cpp:362
virtual void rhs(cArrayView chi, int ienergy, int instate) const
Return the right-hand side.
Definition: preconditioners.h:501
SymDiaMatrix Mm1_tr_kron_S_
Definition: preconditioners.h:248
virtual void setup()
Initialize the preconditioner.
Definition: preconditioners.cpp:773
virtual void multiply(const cArrayView p, cArrayView q) const
Multiply by the matrix equation.
Definition: preconditioners.h:546
SPAICGPreconditioner(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
Definition: preconditioners.h:535
Bspline s_bspline_
Definition: preconditioners.h:242
virtual void rhs(cArrayView chi, int ienergy, int instate) const
Return the right-hand side.
Definition: preconditioners.cpp:478
static const std::string description
Definition: preconditioners.h:303
virtual void setup()
Initialize the preconditioner.
Definition: preconditioners.cpp:1086
SymDiaMatrix half_D_minus_Mm1_tr_
Definition: preconditioners.h:248
JacobiCGPreconditioner(Parallel const &par, InputFile const &inp, std::vector< std::pair< int, int >> const &ll, Bspline const &bspline, CommandLine const &cmd)
Definition: preconditioners.h:362
static const std::string name
Definition: preconditioners.h:487
SymDiaMatrix p_S_kron_S_
Definition: preconditioners.h:696