Reference documentation for deal.II version GIT 4efb66ecd0 2023-02-07 13:45:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
cuda_precondition.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
18 
20 
21 namespace
22 {
29  template <typename Number>
30  cusparseStatus_t
31  cusparseXcsrilu02(cusparseHandle_t handle,
32  int m,
33  int nnz,
34  const cusparseMatDescr_t descrA,
35  Number * csrValA_valM,
36  const int * csrRowPtrA,
37  const int * csrColIndA,
38  csrilu02Info_t info,
39  cusparseSolvePolicy_t policy,
40  void * pBuffer)
41  {
43  return CUSPARSE_STATUS_INVALID_VALUE;
44  }
45 
46  template <>
47  cusparseStatus_t
48  cusparseXcsrilu02<float>(cusparseHandle_t handle,
49  int m,
50  int nnz,
51  const cusparseMatDescr_t descrA,
52  float * csrValA_valM,
53  const int * csrRowPtrA,
54  const int * csrColIndA,
55  csrilu02Info_t info,
56  cusparseSolvePolicy_t policy,
57  void * pBuffer)
58  {
59  return cusparseScsrilu02(handle,
60  m,
61  nnz,
62  descrA,
63  csrValA_valM,
64  csrRowPtrA,
65  csrColIndA,
66  info,
67  policy,
68  pBuffer);
69  }
70 
71  template <>
72  cusparseStatus_t
73  cusparseXcsrilu02<double>(cusparseHandle_t handle,
74  int m,
75  int nnz,
76  const cusparseMatDescr_t descrA,
77  double * csrValA_valM,
78  const int * csrRowPtrA,
79  const int * csrColIndA,
80  csrilu02Info_t info,
81  cusparseSolvePolicy_t policy,
82  void * pBuffer)
83  {
84  return cusparseDcsrilu02(handle,
85  m,
86  nnz,
87  descrA,
88  csrValA_valM,
89  csrRowPtrA,
90  csrColIndA,
91  info,
92  policy,
93  pBuffer);
94  }
95 
96  /*
97  template <>
98  cusparseStatus_t
99  cusparseXcsrilu02<cuComplex>(cusparseHandle_t handle,
100  int m,
101  int nnz,
102  const cusparseMatDescr_t descrA,
103  cuComplex * csrValA_valM,
104  const int * csrRowPtrA,
105  const int * csrColIndA,
106  csrilu02Info_t info,
107  cusparseSolvePolicy_t policy,
108  void * pBuffer)
109  {
110  return cusparseCcsrilu02(handle,
111  m,
112  nnz,
113  descrA,
114  csrValA_valM,
115  csrRowPtrA,
116  csrColIndA,
117  info,
118  policy,
119  pBuffer);
120  }
121 
122  template <>
123  cusparseStatus_t
124  cusparseXcsrilu02<cuDoubleComplex>(cusparseHandle_t handle,
125  int m,
126  int nnz,
127  const cusparseMatDescr_t descrA,
128  cuDoubleComplex * csrValA_valM,
129  const int * csrRowPtrA,
130  const int * csrColIndA,
131  csrilu02Info_t info,
132  cusparseSolvePolicy_t policy,
133  void * pBuffer)
134  {
135  return cusparseZcsrilu02(handle,
136  m,
137  nnz,
138  descrA,
139  csrValA_valM,
140  csrRowPtrA,
141  csrColIndA,
142  info,
143  policy,
144  pBuffer);
145  }
146  */
147 
148 
149 
156  template <typename Number>
157  cusparseStatus_t
158  cusparseXcsrilu02_analysis(cusparseHandle_t /*handle*/,
159  int /*m*/,
160  int /*nnz*/,
161  const cusparseMatDescr_t /*descrA*/,
162  const Number * /*csrValA*/,
163  const int * /*csrRowPtrA*/,
164  const int * /*csrColIndA*/,
165  csrilu02Info_t /*info*/,
166  cusparseSolvePolicy_t /*policy*/,
167  void * /*pBuffer*/)
168  {
169  AssertThrow(false, ExcNotImplemented());
170  return CUSPARSE_STATUS_INVALID_VALUE;
171  }
172 
173  template <>
174  cusparseStatus_t
175  cusparseXcsrilu02_analysis<float>(cusparseHandle_t handle,
176  int m,
177  int nnz,
178  const cusparseMatDescr_t descrA,
179  const float * csrValA,
180  const int * csrRowPtrA,
181  const int * csrColIndA,
182  csrilu02Info_t info,
183  cusparseSolvePolicy_t policy,
184  void * pBuffer)
185  {
186  return cusparseScsrilu02_analysis(handle,
187  m,
188  nnz,
189  descrA,
190  csrValA,
191  csrRowPtrA,
192  csrColIndA,
193  info,
194  policy,
195  pBuffer);
196  }
197 
198  template <>
199  cusparseStatus_t
200  cusparseXcsrilu02_analysis<double>(cusparseHandle_t handle,
201  int m,
202  int nnz,
203  const cusparseMatDescr_t descrA,
204  const double * csrValA,
205  const int * csrRowPtrA,
206  const int * csrColIndA,
207  csrilu02Info_t info,
208  cusparseSolvePolicy_t policy,
209  void * pBuffer)
210  {
211  return cusparseDcsrilu02_analysis(handle,
212  m,
213  nnz,
214  descrA,
215  csrValA,
216  csrRowPtrA,
217  csrColIndA,
218  info,
219  policy,
220  pBuffer);
221  }
222 
223  /*
224  template <>
225  cusparseStatus_t
226  cusparseXcsrilu02_analysis<cuComplex>(cusparseHandle_t handle,
227  int m,
228  int nnz,
229  const cusparseMatDescr_t descrA,
230  const cuComplex * csrValA,
231  const int * csrRowPtrA,
232  const int * csrColIndA,
233  csrilu02Info_t info,
234  cusparseSolvePolicy_t policy,
235  void * pBuffer)
236  {
237  return cusparseCcsrilu02_analysis(handle,
238  m,
239  nnz,
240  descrA,
241  csrValA,
242  csrRowPtrA,
243  csrColIndA,
244  info,
245  policy,
246  pBuffer);
247  }
248 
249  template <>
250  cusparseStatus_t
251  cusparseXcsrilu02_analysis<cuDoubleComplex>(cusparseHandle_t handle,
252  int m,
253  int nnz,
254  const cusparseMatDescr_t descrA,
255  const cuDoubleComplex * csrValA,
256  const int * csrRowPtrA,
257  const int * csrColIndA,
258  csrilu02Info_t info,
259  cusparseSolvePolicy_t policy,
260  void * pBuffer)
261  {
262  return cusparseZcsrilu02_analysis(handle,
263  m,
264  nnz,
265  descrA,
266  csrValA,
267  csrRowPtrA,
268  csrColIndA,
269  info,
270  policy,
271  pBuffer);
272  }
273 */
274 
275 
276 
283  template <typename Number>
284  cusparseStatus_t
285  cusparseXcsrilu02_bufferSize(cusparseHandle_t /*handle*/,
286  int /*m*/,
287  int /*nnz*/,
288  const cusparseMatDescr_t /*descrA*/,
289  Number * /*csrValA*/,
290  const int * /*csrRowPtrA*/,
291  const int * /*csrColIndA*/,
292  csrilu02Info_t /*info*/,
293  int * /*pBufferSizeInBytes*/)
294  {
295  AssertThrow(false, ExcNotImplemented());
296  return CUSPARSE_STATUS_INVALID_VALUE;
297  }
298 
299  template <>
300  cusparseStatus_t
301  cusparseXcsrilu02_bufferSize<float>(cusparseHandle_t handle,
302  int m,
303  int nnz,
304  const cusparseMatDescr_t descrA,
305  float * csrValA,
306  const int * csrRowPtrA,
307  const int * csrColIndA,
308  csrilu02Info_t info,
309  int *pBufferSizeInBytes)
310  {
311  return cusparseScsrilu02_bufferSize(handle,
312  m,
313  nnz,
314  descrA,
315  csrValA,
316  csrRowPtrA,
317  csrColIndA,
318  info,
319  pBufferSizeInBytes);
320  }
321 
322  template <>
323  cusparseStatus_t
324  cusparseXcsrilu02_bufferSize<double>(cusparseHandle_t handle,
325  int m,
326  int nnz,
327  const cusparseMatDescr_t descrA,
328  double * csrValA,
329  const int * csrRowPtrA,
330  const int * csrColIndA,
331  csrilu02Info_t info,
332  int *pBufferSizeInBytes)
333  {
334  return cusparseDcsrilu02_bufferSize(handle,
335  m,
336  nnz,
337  descrA,
338  csrValA,
339  csrRowPtrA,
340  csrColIndA,
341  info,
342  pBufferSizeInBytes);
343  }
344 
345  /*
346  template <>
347  cusparseStatus_t
348  cusparseXcsrilu02_bufferSize<cuComplex>(cusparseHandle_t handle,
349  int m,
350  int nnz,
351  const cusparseMatDescr_t descrA,
352  cuComplex * csrValA,
353  const int * csrRowPtrA,
354  const int * csrColIndA,
355  csrilu02Info_t info,
356  int *pBufferSizeInBytes)
357  {
358  return cusparseCcsrilu02_bufferSize(handle,
359  m,
360  nnz,
361  descrA,
362  csrValA,
363  csrRowPtrA,
364  csrColIndA,
365  info,
366  pBufferSizeInBytes);
367  }
368 
369  template <>
370  cusparseStatus_t
371  cusparseXcsrilu02_bufferSize<cuDoubleComplex>(cusparseHandle_t handle, int
372  m, int nnz, const cusparseMatDescr_t descrA,
373  cuDoubleComplex *csrValA,
374  const int * csrRowPtrA,
375  const int * csrColIndA,
376  csrilu02Info_t info,
377  int *pBufferSizeInBytes)
378  {
379  return cusparseZcsrilu02_bufferSize(handle,
380  m,
381  nnz,
382  descrA,
383  csrValA,
384  csrRowPtrA,
385  csrColIndA,
386  info,
387  pBufferSizeInBytes);
388  }
389 */
390 
391 
392 
399  template <typename Number>
400  cusparseStatus_t
401  cusparseXcsric02(cusparseHandle_t /*handle*/,
402  int /*m*/,
403  int /*nnz*/,
404  const cusparseMatDescr_t /*descrA*/,
405  Number * /*csrValA_valM*/,
406  const int * /*csrRowPtrA*/,
407  const int * /*csrColIndA*/,
408  csric02Info_t /*info*/,
409  cusparseSolvePolicy_t /*policy*/,
410  void * /*pBuffer*/)
411  {
412  AssertThrow(false, ExcNotImplemented());
413  return CUSPARSE_STATUS_INVALID_VALUE;
414  }
415 
416  template <>
417  cusparseStatus_t
418  cusparseXcsric02<float>(cusparseHandle_t handle,
419  int m,
420  int nnz,
421  const cusparseMatDescr_t descrA,
422  float * csrValA_valM,
423  const int * csrRowPtrA,
424  const int * csrColIndA,
425  csric02Info_t info,
426  cusparseSolvePolicy_t policy,
427  void * pBuffer)
428  {
429  return cusparseScsric02(handle,
430  m,
431  nnz,
432  descrA,
433  csrValA_valM,
434  csrRowPtrA,
435  csrColIndA,
436  info,
437  policy,
438  pBuffer);
439  }
440 
441  template <>
442  cusparseStatus_t
443  cusparseXcsric02<double>(cusparseHandle_t handle,
444  int m,
445  int nnz,
446  const cusparseMatDescr_t descrA,
447  double * csrValA_valM,
448  const int * csrRowPtrA,
449  const int * csrColIndA,
450  csric02Info_t info,
451  cusparseSolvePolicy_t policy,
452  void * pBuffer)
453  {
454  return cusparseDcsric02(handle,
455  m,
456  nnz,
457  descrA,
458  csrValA_valM,
459  csrRowPtrA,
460  csrColIndA,
461  info,
462  policy,
463  pBuffer);
464  }
465 
466  /*
467  template <>
468  cusparseStatus_t
469  cusparseXcsric02<cuComplex>(cusparseHandle_t handle,
470  int m,
471  int nnz,
472  const cusparseMatDescr_t descrA,
473  cuComplex * csrValA_valM,
474  const int * csrRowPtrA,
475  const int * csrColIndA,
476  csric02Info_t info,
477  cusparseSolvePolicy_t policy,
478  void * pBuffer)
479  {
480  return cusparseCcsric02(handle,
481  m,
482  nnz,
483  descrA,
484  csrValA_valM,
485  csrRowPtrA,
486  csrColIndA,
487  info,
488  policy,
489  pBuffer);
490  }
491 
492  template <>
493  cusparseStatus_t
494  cusparseXcsric02<cuDoubleComplex>(cusparseHandle_t handle,
495  int m,
496  int nnz,
497  const cusparseMatDescr_t descrA,
498  cuDoubleComplex * csrValA_valM,
499  const int * csrRowPtrA,
500  const int * csrColIndA,
501  csric02Info_t info,
502  cusparseSolvePolicy_t policy,
503  void * pBuffer)
504  {
505  return cusparseZcsric02(handle,
506  m,
507  nnz,
508  descrA,
509  csrValA_valM,
510  csrRowPtrA,
511  csrColIndA,
512  info,
513  policy,
514  pBuffer);
515  }
516  */
517 
518 
519 
526  template <typename Number>
527  cusparseStatus_t
528  cusparseXcsrsv2_solve(cusparseHandle_t /*handle*/,
529  cusparseOperation_t /*transA*/,
530  int /*m*/,
531  int /*nnz*/,
532  const Number * /*alpha*/,
533  const cusparseMatDescr_t /*descra*/,
534  const Number * /*csrValA*/,
535  const int * /*csrRowPtrA*/,
536  const int * /*csrColIndA*/,
537  csrsv2Info_t /*info*/,
538  const Number * /*x*/,
539  Number * /*y*/,
540  cusparseSolvePolicy_t /*policy*/,
541  void * /*pBuffer*/)
542  {
543  AssertThrow(false, ExcNotImplemented());
544  return CUSPARSE_STATUS_INVALID_VALUE;
545  }
546 
547  template <>
548  cusparseStatus_t
549  cusparseXcsrsv2_solve<float>(cusparseHandle_t handle,
550  cusparseOperation_t transA,
551  int m,
552  int nnz,
553  const float * alpha,
554  const cusparseMatDescr_t descra,
555  const float * csrValA,
556  const int * csrRowPtrA,
557  const int * csrColIndA,
558  csrsv2Info_t info,
559  const float * x,
560  float * y,
561  cusparseSolvePolicy_t policy,
562  void * pBuffer)
563  {
564  return cusparseScsrsv2_solve(handle,
565  transA,
566  m,
567  nnz,
568  alpha,
569  descra,
570  csrValA,
571  csrRowPtrA,
572  csrColIndA,
573  info,
574  x,
575  y,
576  policy,
577  pBuffer);
578  }
579 
580  template <>
581  cusparseStatus_t
582  cusparseXcsrsv2_solve<double>(cusparseHandle_t handle,
583  cusparseOperation_t transA,
584  int m,
585  int nnz,
586  const double * alpha,
587  const cusparseMatDescr_t descra,
588  const double * csrValA,
589  const int * csrRowPtrA,
590  const int * csrColIndA,
591  csrsv2Info_t info,
592  const double * x,
593  double * y,
594  cusparseSolvePolicy_t policy,
595  void * pBuffer)
596  {
597  return cusparseDcsrsv2_solve(handle,
598  transA,
599  m,
600  nnz,
601  alpha,
602  descra,
603  csrValA,
604  csrRowPtrA,
605  csrColIndA,
606  info,
607  x,
608  y,
609  policy,
610  pBuffer);
611  }
612 
613  /*
614  template <>
615  cusparseStatus_t
616  cusparseXcsrsv2_solve<cuComplex>(cusparseHandle_t handle,
617  cusparseOperation_t transA,
618  int m,
619  int nnz,
620  const cuComplex * alpha,
621  const cusparseMatDescr_t descra,
622  const cuComplex * csrValA,
623  const int * csrRowPtrA,
624  const int * csrColIndA,
625  csrsv2Info_t info,
626  const cuComplex * x,
627  cuComplex * y,
628  cusparseSolvePolicy_t policy,
629  void * pBuffer)
630  {
631  return cusparseCcsrsv2_solve(handle,
632  transA,
633  m,
634  nnz,
635  alpha,
636  descra,
637  csrValA,
638  csrRowPtrA,
639  csrColIndA,
640  info,
641  x,
642  y,
643  policy,
644  pBuffer);
645  }
646 
647  template <>
648  cusparseStatus_t
649  cusparseXcsrsv2_solve<cuDoubleComplex>(cusparseHandle_t handle,
650  cusparseOperation_t transA,
651  int m,
652  int nnz,
653  const cuDoubleComplex * alpha,
654  const cusparseMatDescr_t descra,
655  const cuDoubleComplex * csrValA,
656  const int * csrRowPtrA,
657  const int * csrColIndA,
658  csrsv2Info_t info,
659  const cuDoubleComplex * x,
660  cuDoubleComplex * y,
661  cusparseSolvePolicy_t policy,
662  void * pBuffer)
663  {
664  return cusparseZcsrsv2_solve(handle,
665  transA,
666  m,
667  nnz,
668  alpha,
669  descra,
670  csrValA,
671  csrRowPtrA,
672  csrColIndA,
673  info,
674  x,
675  y,
676  policy,
677  pBuffer);
678  }
679 */
680 
681 
682 
689  template <typename Number>
690  cusparseStatus_t
691  cusparseXcsrsv2_analysis(cusparseHandle_t /*handle*/,
692  cusparseOperation_t /*transA*/,
693  int /*m*/,
694  int /*nnz*/,
695  const cusparseMatDescr_t /*descrA*/,
696  const Number * /*csrValA*/,
697  const int * /*csrRowPtrA*/,
698  const int * /*csrColIndA*/,
699  csrsv2Info_t /*info*/,
700  cusparseSolvePolicy_t /*policy*/,
701  void * /*pBuffer*/)
702  {
703  AssertThrow(false, ExcNotImplemented());
704  return CUSPARSE_STATUS_INVALID_VALUE;
705  }
706 
707  template <>
708  cusparseStatus_t
709  cusparseXcsrsv2_analysis<float>(cusparseHandle_t handle,
710  cusparseOperation_t transA,
711  int m,
712  int nnz,
713  const cusparseMatDescr_t descrA,
714  const float * csrValA,
715  const int * csrRowPtrA,
716  const int * csrColIndA,
717  csrsv2Info_t info,
718  cusparseSolvePolicy_t policy,
719  void * pBuffer)
720  {
721  return cusparseScsrsv2_analysis(handle,
722  transA,
723  m,
724  nnz,
725  descrA,
726  csrValA,
727  csrRowPtrA,
728  csrColIndA,
729  info,
730  policy,
731  pBuffer);
732  }
733 
734  template <>
735  cusparseStatus_t
736  cusparseXcsrsv2_analysis<double>(cusparseHandle_t handle,
737  cusparseOperation_t transA,
738  int m,
739  int nnz,
740  const cusparseMatDescr_t descrA,
741  const double * csrValA,
742  const int * csrRowPtrA,
743  const int * csrColIndA,
744  csrsv2Info_t info,
745  cusparseSolvePolicy_t policy,
746  void * pBuffer)
747  {
748  return cusparseDcsrsv2_analysis(handle,
749  transA,
750  m,
751  nnz,
752  descrA,
753  csrValA,
754  csrRowPtrA,
755  csrColIndA,
756  info,
757  policy,
758  pBuffer);
759  }
760 
761  /*
762  template <>
763  cusparseStatus_t
764  cusparseXcsrsv2_analysis<cuComplex>(cusparseHandle_t handle,
765  cusparseOperation_t transA,
766  int m,
767  int nnz,
768  const cusparseMatDescr_t descrA,
769  const cuComplex * csrValA,
770  const int * csrRowPtrA,
771  const int * csrColIndA,
772  csrsv2Info_t info,
773  cusparseSolvePolicy_t policy,
774  void * pBuffer)
775  {
776  return cusparseCcsrsv2_analysis(handle,
777  transA,
778  m,
779  nnz,
780  descrA,
781  csrValA,
782  csrRowPtrA,
783  csrColIndA,
784  info,
785  policy,
786  pBuffer);
787  }
788 
789  template <>
790  cusparseStatus_t
791  cusparseXcsrsv2_analysis<cuDoubleComplex>(cusparseHandle_t handle,
792  cusparseOperation_t transA,
793  int m,
794  int nnz,
795  const cusparseMatDescr_t descrA,
796  const cuDoubleComplex * csrValA,
797  const int * csrRowPtrA,
798  const int * csrColIndA,
799  csrsv2Info_t info,
800  cusparseSolvePolicy_t policy,
801  void * pBuffer)
802  {
803  return cusparseZcsrsv2_analysis(handle,
804  transA,
805  m,
806  nnz,
807  descrA,
808  csrValA,
809  csrRowPtrA,
810  csrColIndA,
811  info,
812  policy,
813  pBuffer);
814  }
815 */
816 
817 
824  template <typename Number>
825  cusparseStatus_t
826  cusparseXcsric02_analysis(cusparseHandle_t /*handle*/,
827  int /*m*/,
828  int /*nnz*/,
829  const cusparseMatDescr_t /*descrA*/,
830  const Number * /*csrValA*/,
831  const int * /*csrRowPtrA*/,
832  const int * /*csrColIndA*/,
833  csric02Info_t /*info*/,
834  cusparseSolvePolicy_t /*policy*/,
835  void * /*pBuffer*/)
836  {
837  AssertThrow(false, ExcNotImplemented());
838  return CUSPARSE_STATUS_INVALID_VALUE;
839  }
840 
841  template <>
842  cusparseStatus_t
843  cusparseXcsric02_analysis<float>(cusparseHandle_t handle,
844  int m,
845  int nnz,
846  const cusparseMatDescr_t descrA,
847  const float * csrValA,
848  const int * csrRowPtrA,
849  const int * csrColIndA,
850  csric02Info_t info,
851  cusparseSolvePolicy_t policy,
852  void * pBuffer)
853  {
854  return cusparseScsric02_analysis(handle,
855  m,
856  nnz,
857  descrA,
858  csrValA,
859  csrRowPtrA,
860  csrColIndA,
861  info,
862  policy,
863  pBuffer);
864  }
865 
866  template <>
867  cusparseStatus_t
868  cusparseXcsric02_analysis<double>(cusparseHandle_t handle,
869  int m,
870  int nnz,
871  const cusparseMatDescr_t descrA,
872  const double * csrValA,
873  const int * csrRowPtrA,
874  const int * csrColIndA,
875  csric02Info_t info,
876  cusparseSolvePolicy_t policy,
877  void * pBuffer)
878  {
879  return cusparseDcsric02_analysis(handle,
880  m,
881  nnz,
882  descrA,
883  csrValA,
884  csrRowPtrA,
885  csrColIndA,
886  info,
887  policy,
888  pBuffer);
889  }
890 
891  /*
892  template <>
893  cusparseStatus_t
894  cusparseXcsric02_analysis<cuComplex>(cusparseHandle_t handle,
895  int m,
896  int nnz,
897  const cusparseMatDescr_t descrA,
898  const cuComplex * csrValA,
899  const int * csrRowPtrA,
900  const int * csrColIndA,
901  csric02Info_t info,
902  cusparseSolvePolicy_t policy,
903  void * pBuffer)
904  {
905  return cusparseCcsric02_analysis(handle,
906  m,
907  nnz,
908  descrA,
909  csrValA,
910  csrRowPtrA,
911  csrColIndA,
912  info,
913  policy,
914  pBuffer);
915  }
916 
917  template <>
918  cusparseStatus_t
919  cusparseXcsric02_analysis<cuDoubleComplex>(cusparseHandle_t handle,
920  int m,
921  int nnz,
922  const cusparseMatDescr_t descrA,
923  const cuDoubleComplex * csrValA,
924  const int * csrRowPtrA,
925  const int * csrColIndA,
926  csric02Info_t info,
927  cusparseSolvePolicy_t policy,
928  void * pBuffer)
929  {
930  return cusparseZcsric02_analysis(handle,
931  m,
932  nnz,
933  descrA,
934  csrValA,
935  csrRowPtrA,
936  csrColIndA,
937  info,
938  policy,
939  pBuffer);
940  }
941 */
942 
943 
944 
951  template <typename Number>
952  cusparseStatus_t
953  cusparseXcsrsv2_bufferSize(cusparseHandle_t /*handle*/,
954  cusparseOperation_t /*transA*/,
955  int /*m*/,
956  int /*nnz*/,
957  const cusparseMatDescr_t /*descrA*/,
958  Number * /*csrValA*/,
959  const int * /*csrRowPtrA*/,
960  const int * /*csrColIndA*/,
961  csrsv2Info_t /*info*/,
962  int * /*pBufferSizeInBytes*/)
963  {
964  AssertThrow(false, ExcNotImplemented());
965  return CUSPARSE_STATUS_INVALID_VALUE;
966  }
967 
968  template <>
969  cusparseStatus_t
970  cusparseXcsrsv2_bufferSize<float>(cusparseHandle_t handle,
971  cusparseOperation_t transA,
972  int m,
973  int nnz,
974  const cusparseMatDescr_t descrA,
975  float * csrValA,
976  const int * csrRowPtrA,
977  const int * csrColIndA,
978  csrsv2Info_t info,
979  int * pBufferSizeInBytes)
980  {
981  return cusparseScsrsv2_bufferSize(handle,
982  transA,
983  m,
984  nnz,
985  descrA,
986  csrValA,
987  csrRowPtrA,
988  csrColIndA,
989  info,
990  pBufferSizeInBytes);
991  }
992 
993  template <>
994  cusparseStatus_t
995  cusparseXcsrsv2_bufferSize<double>(cusparseHandle_t handle,
996  cusparseOperation_t transA,
997  int m,
998  int nnz,
999  const cusparseMatDescr_t descrA,
1000  double * csrValA,
1001  const int * csrRowPtrA,
1002  const int * csrColIndA,
1003  csrsv2Info_t info,
1004  int *pBufferSizeInBytes)
1005  {
1006  return cusparseDcsrsv2_bufferSize(handle,
1007  transA,
1008  m,
1009  nnz,
1010  descrA,
1011  csrValA,
1012  csrRowPtrA,
1013  csrColIndA,
1014  info,
1015  pBufferSizeInBytes);
1016  }
1017 
1018  /*
1019  template <>
1020  cusparseStatus_t
1021  cusparseXcsrsv2_bufferSize<cuComplex>(cusparseHandle_t handle,
1022  cusparseOperation_t transA,
1023  int m,
1024  int nnz,
1025  const cusparseMatDescr_t descrA,
1026  cuComplex * csrValA,
1027  const int * csrRowPtrA,
1028  const int * csrColIndA,
1029  csrsv2Info_t info,
1030  int *pBufferSizeInBytes)
1031  {
1032  return cusparseCcsrsv2_bufferSize(handle,
1033  transA,
1034  m,
1035  nnz,
1036  descrA,
1037  csrValA,
1038  csrRowPtrA,
1039  csrColIndA,
1040  info,
1041  pBufferSizeInBytes);
1042  }
1043 
1044  template <>
1045  cusparseStatus_t
1046  cusparseXcsrsv2_bufferSize<cuDoubleComplex>(cusparseHandle_t handle,
1047  cusparseOperation_t transA,
1048  int m,
1049  int nnz,
1050  const cusparseMatDescr_t descrA,
1051  cuDoubleComplex * csrValA,
1052  const int * csrRowPtrA,
1053  const int * csrColIndA,
1054  csrsv2Info_t info,
1055  int * pBufferSizeInBytes)
1056  {
1057  return cusparseZcsrsv2_bufferSize(handle,
1058  transA,
1059  m,
1060  nnz,
1061  descrA,
1062  csrValA,
1063  csrRowPtrA,
1064  csrColIndA,
1065  info,
1066  pBufferSizeInBytes);
1067  }
1068 */
1069 
1070 
1071 
1078  template <typename Number>
1079  cusparseStatus_t
1080  cusparseXcsric02_bufferSize(cusparseHandle_t /*handle*/,
1081  int /*m*/,
1082  int /*nnz*/,
1083  const cusparseMatDescr_t /*descrA*/,
1084  Number * /*csrValA*/,
1085  const int * /*csrRowPtrA*/,
1086  const int * /*csrColIndA*/,
1087  csric02Info_t /*info*/,
1088  int * /*pBufferSizeInBytes*/)
1089  {
1090  AssertThrow(false, ExcNotImplemented());
1091  return CUSPARSE_STATUS_INVALID_VALUE;
1092  }
1093 
1094  template <>
1095  cusparseStatus_t
1096  cusparseXcsric02_bufferSize<float>(cusparseHandle_t handle,
1097  int m,
1098  int nnz,
1099  const cusparseMatDescr_t descrA,
1100  float * csrValA,
1101  const int * csrRowPtrA,
1102  const int * csrColIndA,
1103  csric02Info_t info,
1104  int *pBufferSizeInBytes)
1105  {
1106  return cusparseScsric02_bufferSize(handle,
1107  m,
1108  nnz,
1109  descrA,
1110  csrValA,
1111  csrRowPtrA,
1112  csrColIndA,
1113  info,
1114  pBufferSizeInBytes);
1115  }
1116 
1117  template <>
1118  cusparseStatus_t
1119  cusparseXcsric02_bufferSize<double>(cusparseHandle_t handle,
1120  int m,
1121  int nnz,
1122  const cusparseMatDescr_t descrA,
1123  double * csrValA,
1124  const int * csrRowPtrA,
1125  const int * csrColIndA,
1126  csric02Info_t info,
1127  int *pBufferSizeInBytes)
1128  {
1129  return cusparseDcsric02_bufferSize(handle,
1130  m,
1131  nnz,
1132  descrA,
1133  csrValA,
1134  csrRowPtrA,
1135  csrColIndA,
1136  info,
1137  pBufferSizeInBytes);
1138  }
1139 
1140  /*
1141  template <>
1142  cusparseStatus_t
1143  cusparseXcsric02_bufferSize<cuComplex>(cusparseHandle_t handle,
1144  int m,
1145  int nnz,
1146  const cusparseMatDescr_t descrA,
1147  cuComplex * csrValA,
1148  const int * csrRowPtrA,
1149  const int * csrColIndA,
1150  csric02Info_t info,
1151  int *pBufferSizeInBytes)
1152  {
1153  return cusparseCcsric02_bufferSize(handle,
1154  m,
1155  nnz,
1156  descrA,
1157  csrValA,
1158  csrRowPtrA,
1159  csrColIndA,
1160  info,
1161  pBufferSizeInBytes);
1162  }
1163 
1164  template <>
1165  cusparseStatus_t
1166  cusparseXcsric02_bufferSize<cuDoubleComplex>(cusparseHandle_t handle,
1167  int m,
1168  int nnz,
1169  const cusparseMatDescr_t descrA,
1170  cuDoubleComplex * csrValA,
1171  const int * csrRowPtrA,
1172  const int * csrColIndA,
1173  csric02Info_t info,
1174  int * pBufferSizeInBytes)
1175  {
1176  return cusparseZcsric02_bufferSize(handle,
1177  m,
1178  nnz,
1179  descrA,
1180  csrValA,
1181  csrRowPtrA,
1182  csrColIndA,
1183  info,
1184  pBufferSizeInBytes);
1185  }
1186  */
1187 } // namespace
1188 
1189 namespace CUDAWrappers
1190 {
1191  template <typename Number>
1193  bool use_level_analysis_)
1194  : use_level_analysis(use_level_analysis_)
1195  {}
1196 
1197 
1198 
1199  template <typename Number>
1202  , P_val_dev(nullptr, Utilities::CUDA::delete_device_data<Number>)
1203  , P_row_ptr_dev(nullptr)
1204  , P_column_index_dev(nullptr)
1205  , tmp_dev(nullptr, Utilities::CUDA::delete_device_data<Number>)
1206  , buffer_dev(nullptr, Utilities::CUDA::delete_device_data<void>)
1207  , policy_L(CUSPARSE_SOLVE_POLICY_USE_LEVEL)
1208  , policy_Lt(CUSPARSE_SOLVE_POLICY_USE_LEVEL)
1209  , policy_M(CUSPARSE_SOLVE_POLICY_USE_LEVEL)
1210  , n_rows(0)
1211  , n_nonzero_elements(0)
1212  {
1213  // step 1: create a descriptor which contains
1214  // - matrix M is base-0
1215  // - matrix L is base-0
1216  // - matrix L is lower triangular
1217  // - matrix L has non-unit diagonal
1218  cusparseStatus_t status = cusparseCreateMatDescr(&descr_M);
1219  AssertCusparse(status);
1220  status = cusparseSetMatIndexBase(descr_M, CUSPARSE_INDEX_BASE_ZERO);
1221  AssertCusparse(status);
1222  status = cusparseSetMatType(descr_M, CUSPARSE_MATRIX_TYPE_GENERAL);
1223  AssertCusparse(status);
1224 
1225  status = cusparseCreateMatDescr(&descr_L);
1226  AssertCusparse(status);
1227  status = cusparseSetMatIndexBase(descr_L, CUSPARSE_INDEX_BASE_ZERO);
1228  AssertCusparse(status);
1229  status = cusparseSetMatType(descr_L, CUSPARSE_MATRIX_TYPE_GENERAL);
1230  AssertCusparse(status);
1231  status = cusparseSetMatFillMode(descr_L, CUSPARSE_FILL_MODE_LOWER);
1232  AssertCusparse(status);
1233  status = cusparseSetMatDiagType(descr_L, CUSPARSE_DIAG_TYPE_NON_UNIT);
1234  AssertCusparse(status);
1235 
1236  // step 2: create a empty info structure
1237  // we need one info for csric02 and two info's for csrsv2
1238  status = cusparseCreateCsric02Info(&info_M);
1239  AssertCusparse(status);
1240  status = cusparseCreateCsrsv2Info(&info_L);
1241  AssertCusparse(status);
1242  status = cusparseCreateCsrsv2Info(&info_Lt);
1243  AssertCusparse(status);
1244  }
1245 
1246 
1247 
1248  template <typename Number>
1250  {
1251  // step 8: free resources
1252  cusparseStatus_t status = cusparseDestroyMatDescr(descr_M);
1253  AssertNothrowCusparse(status);
1254 
1255  status = cusparseDestroyMatDescr(descr_L);
1256  AssertNothrowCusparse(status);
1257 
1258  status = cusparseDestroyCsric02Info(info_M);
1259  AssertNothrowCusparse(status);
1260 
1261  status = cusparseDestroyCsrsv2Info(info_L);
1262  AssertNothrowCusparse(status);
1263 
1264  status = cusparseDestroyCsrsv2Info(info_Lt);
1265  AssertNothrowCusparse(status);
1266  }
1267 
1268 
1269 
1270  template <typename Number>
1271  void
1273  const AdditionalData &additional_data)
1274  {
1275  if (additional_data.use_level_analysis)
1276  {
1277  policy_L = CUSPARSE_SOLVE_POLICY_USE_LEVEL;
1278  policy_Lt = CUSPARSE_SOLVE_POLICY_USE_LEVEL;
1279  policy_M = CUSPARSE_SOLVE_POLICY_USE_LEVEL;
1280  }
1281  else
1282  {
1283  policy_L = CUSPARSE_SOLVE_POLICY_NO_LEVEL;
1284  policy_Lt = CUSPARSE_SOLVE_POLICY_NO_LEVEL;
1285  policy_M = CUSPARSE_SOLVE_POLICY_NO_LEVEL;
1286  }
1287 
1288  n_rows = A.m();
1289  n_nonzero_elements = A.n_nonzero_elements();
1290  AssertDimension(A.m(), A.n());
1291 
1292  matrix_pointer = &A;
1293  const auto cusparse_matrix = A.get_cusparse_matrix();
1294  const Number *const A_val_dev = std::get<0>(cusparse_matrix);
1295 
1296  // create a copy of the matrix entries since the algorithm works in-place.
1297  P_val_dev.reset(
1298  Utilities::CUDA::allocate_device_data<Number>(n_nonzero_elements));
1299  cudaError_t cuda_status = cudaMemcpy(P_val_dev.get(),
1300  A_val_dev,
1301  n_nonzero_elements * sizeof(Number),
1302  cudaMemcpyDeviceToDevice);
1303 
1304  P_column_index_dev = std::get<1>(cusparse_matrix);
1305  P_row_ptr_dev = std::get<2>(cusparse_matrix);
1306  const cusparseMatDescr_t mat_descr = std::get<3>(cusparse_matrix);
1307 
1308  // initialize an internal buffer we need later on
1309  tmp_dev.reset(Utilities::CUDA::allocate_device_data<Number>(n_rows));
1310 
1311  // step 3: query how much memory used in csric02 and csrsv2, and allocate
1312  // the buffer
1313  int BufferSize_M;
1314  cusparseStatus_t status = cusparseXcsric02_bufferSize(cusparse_handle,
1315  n_rows,
1316  n_nonzero_elements,
1317  descr_M,
1318  P_val_dev.get(),
1319  P_row_ptr_dev,
1320  P_column_index_dev,
1321  info_M,
1322  &BufferSize_M);
1323  AssertCusparse(status);
1324 
1325  int BufferSize_L;
1326  status = cusparseXcsrsv2_bufferSize(cusparse_handle,
1327  CUSPARSE_OPERATION_NON_TRANSPOSE,
1328  n_rows,
1329  n_nonzero_elements,
1330  descr_L,
1331  P_val_dev.get(),
1332  P_row_ptr_dev,
1333  P_column_index_dev,
1334  info_L,
1335  &BufferSize_L);
1336  AssertCusparse(status);
1337 
1338  int BufferSize_Lt;
1339  status = cusparseXcsrsv2_bufferSize(cusparse_handle,
1340  CUSPARSE_OPERATION_TRANSPOSE,
1341  n_rows,
1342  n_nonzero_elements,
1343  descr_L,
1344  P_val_dev.get(),
1345  P_row_ptr_dev,
1346  P_column_index_dev,
1347  info_Lt,
1348  &BufferSize_Lt);
1349  AssertCusparse(status);
1350 
1351  const int BufferSize =
1352  std::max(BufferSize_M, std::max(BufferSize_L, BufferSize_Lt));
1353  // workaround: since allocate_device_data needs a type, we pass char
1354  // which is required to have size 1.
1355  buffer_dev.reset(static_cast<void *>(
1356  Utilities::CUDA::allocate_device_data<char>(BufferSize / sizeof(char))));
1357 
1358  // step 4: perform analysis of incomplete Cholesky on M
1359  // perform analysis of triangular solve on L
1360  // perform analysis of triangular solve on L'
1361  // The lower triangular part of M has the same sparsity pattern as L, so
1362  // we can do analysis of csric02 and csrsv2 simultaneously.
1363 
1364  status = cusparseXcsric02_analysis(cusparse_handle,
1365  n_rows,
1366  n_nonzero_elements,
1367  descr_M,
1368  P_val_dev.get(),
1369  P_row_ptr_dev,
1370  P_column_index_dev,
1371  info_M,
1372  policy_M,
1373  buffer_dev.get());
1374  AssertCusparse(status);
1375 
1376  int structural_zero;
1377  status =
1378  cusparseXcsric02_zeroPivot(cusparse_handle, info_M, &structural_zero);
1379  AssertCusparse(status);
1380 
1381  status = cusparseXcsrsv2_analysis(cusparse_handle,
1382  CUSPARSE_OPERATION_TRANSPOSE,
1383  n_rows,
1384  n_nonzero_elements,
1385  descr_L,
1386  P_val_dev.get(),
1387  P_row_ptr_dev,
1388  P_column_index_dev,
1389  info_Lt,
1390  policy_Lt,
1391  buffer_dev.get());
1392  AssertCusparse(status);
1393 
1394  status = cusparseXcsrsv2_analysis(cusparse_handle,
1395  CUSPARSE_OPERATION_NON_TRANSPOSE,
1396  n_rows,
1397  n_nonzero_elements,
1398  descr_L,
1399  P_val_dev.get(),
1400  P_row_ptr_dev,
1401  P_column_index_dev,
1402  info_L,
1403  policy_L,
1404  buffer_dev.get());
1405  AssertCusparse(status);
1406 
1407  // step 5: M = L * L'
1408  status = cusparseXcsric02(cusparse_handle,
1409  n_rows,
1410  n_nonzero_elements,
1411  descr_M,
1412  P_val_dev.get(),
1413  P_row_ptr_dev,
1414  P_column_index_dev,
1415  info_M,
1416  policy_M,
1417  buffer_dev.get());
1418  AssertCusparse(status);
1419 
1420  int numerical_zero;
1421  status =
1422  cusparseXcsric02_zeroPivot(cusparse_handle, info_M, &numerical_zero);
1423  AssertCusparse(status);
1424  }
1425 
1426 
1427 
1428  template <typename Number>
1429  void
1433  {
1434  Assert(P_val_dev != nullptr, ExcNotInitialized());
1435  Assert(P_row_ptr_dev != nullptr, ExcNotInitialized());
1436  Assert(P_column_index_dev != nullptr, ExcNotInitialized());
1437  AssertDimension(dst.size(), static_cast<unsigned int>(n_rows));
1438  AssertDimension(src.size(), static_cast<unsigned int>(n_rows));
1439  Assert(tmp_dev != nullptr, ExcInternalError());
1440 
1441  const Number *const src_dev = src.get_values();
1442  Number *const dst_dev = dst.get_values();
1443  // step 6: solve L*z = alpha*x
1444  const Number alpha = internal::NumberType<Number>::value(1.);
1445  cusparseStatus_t status =
1446  cusparseXcsrsv2_solve(cusparse_handle,
1447  CUSPARSE_OPERATION_NON_TRANSPOSE,
1448  n_rows,
1449  n_nonzero_elements,
1450  &alpha,
1451  descr_L,
1452  P_val_dev.get(),
1453  P_row_ptr_dev,
1454  P_column_index_dev,
1455  info_L,
1456  src_dev,
1457  tmp_dev.get(),
1458  policy_L,
1459  buffer_dev.get());
1460  AssertCusparse(status);
1461 
1462  // step 7: solve L'*y = alpha*z
1463  status = cusparseXcsrsv2_solve(cusparse_handle,
1464  CUSPARSE_OPERATION_TRANSPOSE,
1465  n_rows,
1466  n_nonzero_elements,
1467  &alpha,
1468  descr_L,
1469  P_val_dev.get(),
1470  P_row_ptr_dev,
1471  P_column_index_dev,
1472  info_Lt,
1473  tmp_dev.get(),
1474  dst_dev,
1475  policy_Lt,
1476  buffer_dev.get());
1477  AssertCusparse(status);
1478  }
1479 
1480 
1481 
1482  template <typename Number>
1483  void
1487  {
1488  // the constructed preconditioner is symmetric
1489  vmult(dst, src);
1490  }
1491 
1492 
1493 
1494  template <typename Number>
1496  bool use_level_analysis_)
1497  : use_level_analysis(use_level_analysis_)
1498  {}
1499 
1500 
1501 
1502  template <typename Number>
1504  const Utilities::CUDA::Handle &handle)
1506  , P_val_dev(nullptr, Utilities::CUDA::delete_device_data<Number>)
1507  , P_row_ptr_dev(nullptr)
1508  , P_column_index_dev(nullptr)
1509  , tmp_dev(nullptr, Utilities::CUDA::delete_device_data<Number>)
1510  , buffer_dev(nullptr, Utilities::CUDA::delete_device_data<void>)
1511  , policy_L(CUSPARSE_SOLVE_POLICY_USE_LEVEL)
1512  , policy_U(CUSPARSE_SOLVE_POLICY_USE_LEVEL)
1513  , policy_M(CUSPARSE_SOLVE_POLICY_USE_LEVEL)
1514  , n_rows(0)
1515  , n_nonzero_elements(0)
1516  {
1517  // step 1: create a descriptor which contains
1518  // - matrix M is base-0
1519  // - matrix L is base-0
1520  // - matrix L is lower triangular
1521  // - matrix L has unit diagonal
1522  // - matrix U is base-0
1523  // - matrix U is upper triangular
1524  // - matrix U has non-unit diagonal
1525  cusparseStatus_t status = cusparseCreateMatDescr(&descr_M);
1526  AssertCusparse(status);
1527  status = cusparseSetMatIndexBase(descr_M, CUSPARSE_INDEX_BASE_ZERO);
1528  AssertCusparse(status);
1529  status = cusparseSetMatType(descr_M, CUSPARSE_MATRIX_TYPE_GENERAL);
1530  AssertCusparse(status);
1531 
1532  status = cusparseCreateMatDescr(&descr_L);
1533  AssertCusparse(status);
1534  status = cusparseSetMatIndexBase(descr_L, CUSPARSE_INDEX_BASE_ZERO);
1535  AssertCusparse(status);
1536  status = cusparseSetMatType(descr_L, CUSPARSE_MATRIX_TYPE_GENERAL);
1537  AssertCusparse(status);
1538  status = cusparseSetMatFillMode(descr_L, CUSPARSE_FILL_MODE_LOWER);
1539  AssertCusparse(status);
1540  status = cusparseSetMatDiagType(descr_L, CUSPARSE_DIAG_TYPE_UNIT);
1541  AssertCusparse(status);
1542 
1543  status = cusparseCreateMatDescr(&descr_U);
1544  AssertCusparse(status);
1545  status = cusparseSetMatIndexBase(descr_U, CUSPARSE_INDEX_BASE_ZERO);
1546  AssertCusparse(status);
1547  status = cusparseSetMatType(descr_U, CUSPARSE_MATRIX_TYPE_GENERAL);
1548  AssertCusparse(status);
1549  status = cusparseSetMatFillMode(descr_U, CUSPARSE_FILL_MODE_UPPER);
1550  AssertCusparse(status);
1551  status = cusparseSetMatDiagType(descr_U, CUSPARSE_DIAG_TYPE_NON_UNIT);
1552  AssertCusparse(status);
1553 
1554  // step 2: create a empty info structure
1555  // we need one info for csrilu02 and two info's for csrsv2
1556  status = cusparseCreateCsrilu02Info(&info_M);
1557  AssertCusparse(status);
1558  status = cusparseCreateCsrsv2Info(&info_L);
1559  AssertCusparse(status);
1560  status = cusparseCreateCsrsv2Info(&info_U);
1561  AssertCusparse(status);
1562  }
1563 
1564 
1565 
1566  template <typename Number>
1568  {
1569  // step 8: free resources
1570  cusparseStatus_t status = cusparseDestroyMatDescr(descr_M);
1571  AssertNothrowCusparse(status);
1572 
1573  status = cusparseDestroyMatDescr(descr_L);
1574  AssertNothrowCusparse(status);
1575 
1576  status = cusparseDestroyMatDescr(descr_U);
1577  AssertNothrowCusparse(status);
1578 
1579  status = cusparseDestroyCsrilu02Info(info_M);
1580  AssertNothrowCusparse(status);
1581 
1582  status = cusparseDestroyCsrsv2Info(info_L);
1583  AssertNothrowCusparse(status);
1584 
1585  status = cusparseDestroyCsrsv2Info(info_U);
1586  AssertNothrowCusparse(status);
1587  }
1588 
1589 
1590 
1591  template <typename Number>
1592  void
1594  const AdditionalData &additional_data)
1595  {
1596  if (additional_data.use_level_analysis)
1597  {
1598  policy_L = CUSPARSE_SOLVE_POLICY_USE_LEVEL;
1599  policy_U = CUSPARSE_SOLVE_POLICY_USE_LEVEL;
1600  policy_M = CUSPARSE_SOLVE_POLICY_USE_LEVEL;
1601  }
1602  else
1603  {
1604  policy_L = CUSPARSE_SOLVE_POLICY_NO_LEVEL;
1605  policy_U = CUSPARSE_SOLVE_POLICY_NO_LEVEL;
1606  policy_M = CUSPARSE_SOLVE_POLICY_NO_LEVEL;
1607  }
1608 
1609  matrix_pointer = &A;
1610  n_rows = A.m();
1611  n_nonzero_elements = A.n_nonzero_elements();
1612  AssertDimension(A.m(), A.n());
1613 
1614  const auto cusparse_matrix = A.get_cusparse_matrix();
1615  const Number *const A_val_dev = std::get<0>(cusparse_matrix);
1616 
1617  // create a copy of the matrix entries since the algorithm works in-place.
1618  P_val_dev.reset(
1619  Utilities::CUDA::allocate_device_data<Number>(n_nonzero_elements));
1620  cudaError_t cuda_status = cudaMemcpy(P_val_dev.get(),
1621  A_val_dev,
1622  n_nonzero_elements * sizeof(Number),
1623  cudaMemcpyDeviceToDevice);
1624 
1625  P_column_index_dev = std::get<1>(cusparse_matrix);
1626  P_row_ptr_dev = std::get<2>(cusparse_matrix);
1627  const cusparseMatDescr_t mat_descr = std::get<3>(cusparse_matrix);
1628 
1629  // initialize an internal buffer we need later on
1630  tmp_dev.reset(Utilities::CUDA::allocate_device_data<Number>(n_rows));
1631 
1632  // step 3: query how much memory used in csrilu02 and csrsv2, and allocate
1633  // the buffer
1634  int BufferSize_M;
1635  cusparseStatus_t status = cusparseXcsrilu02_bufferSize(cusparse_handle,
1636  n_rows,
1637  n_nonzero_elements,
1638  descr_M,
1639  P_val_dev.get(),
1640  P_row_ptr_dev,
1641  P_column_index_dev,
1642  info_M,
1643  &BufferSize_M);
1644  AssertCusparse(status);
1645 
1646  int BufferSize_L;
1647  status = cusparseXcsrsv2_bufferSize(cusparse_handle,
1648  CUSPARSE_OPERATION_NON_TRANSPOSE,
1649  n_rows,
1650  n_nonzero_elements,
1651  descr_L,
1652  P_val_dev.get(),
1653  P_row_ptr_dev,
1654  P_column_index_dev,
1655  info_L,
1656  &BufferSize_L);
1657  AssertCusparse(status);
1658 
1659  int BufferSize_U;
1660  status = cusparseXcsrsv2_bufferSize(cusparse_handle,
1661  CUSPARSE_OPERATION_NON_TRANSPOSE,
1662  n_rows,
1663  n_nonzero_elements,
1664  descr_U,
1665  P_val_dev.get(),
1666  P_row_ptr_dev,
1667  P_column_index_dev,
1668  info_U,
1669  &BufferSize_U);
1670  AssertCusparse(status);
1671 
1672  const int BufferSize =
1673  std::max(BufferSize_M, std::max(BufferSize_L, BufferSize_U));
1674  // workaround: since allocate_device_data needs a type, we pass char
1675  // which is required to have size 1.
1676  buffer_dev.reset(static_cast<void *>(
1677  Utilities::CUDA::allocate_device_data<char>(BufferSize / sizeof(char))));
1678 
1679  // step 4: perform analysis of incomplete Cholesky on M
1680  // perform analysis of triangular solve on L
1681  // perform analysis of triangular solve on U
1682  // The lower(upper) triangular part of M has the same sparsity pattern as
1683  // L(U), we can do analysis of csrilu0 and csrsv2 simultaneously.
1684 
1685  status = cusparseXcsrilu02_analysis(cusparse_handle,
1686  n_rows,
1687  n_nonzero_elements,
1688  descr_M,
1689  P_val_dev.get(),
1690  P_row_ptr_dev,
1691  P_column_index_dev,
1692  info_M,
1693  policy_M,
1694  buffer_dev.get());
1695  AssertCusparse(status);
1696 
1697  int structural_zero;
1698  status =
1699  cusparseXcsrilu02_zeroPivot(cusparse_handle, info_M, &structural_zero);
1700  AssertCusparse(status);
1701 
1702  status = cusparseXcsrsv2_analysis(cusparse_handle,
1703  CUSPARSE_OPERATION_NON_TRANSPOSE,
1704  n_rows,
1705  n_nonzero_elements,
1706  descr_L,
1707  P_val_dev.get(),
1708  P_row_ptr_dev,
1709  P_column_index_dev,
1710  info_L,
1711  policy_L,
1712  buffer_dev.get());
1713  AssertCusparse(status);
1714 
1715  status = cusparseXcsrsv2_analysis(cusparse_handle,
1716  CUSPARSE_OPERATION_NON_TRANSPOSE,
1717  n_rows,
1718  n_nonzero_elements,
1719  descr_U,
1720  P_val_dev.get(),
1721  P_row_ptr_dev,
1722  P_column_index_dev,
1723  info_U,
1724  policy_U,
1725  buffer_dev.get());
1726 
1727  // step 5: M = L * U
1728  status = cusparseXcsrilu02(cusparse_handle,
1729  n_rows,
1730  n_nonzero_elements,
1731  descr_M,
1732  P_val_dev.get(),
1733  P_row_ptr_dev,
1734  P_column_index_dev,
1735  info_M,
1736  policy_M,
1737  buffer_dev.get());
1738  AssertCusparse(status);
1739 
1740  int numerical_zero;
1741  status =
1742  cusparseXcsrilu02_zeroPivot(cusparse_handle, info_M, &numerical_zero);
1743  AssertCusparse(status);
1744  }
1745 
1746 
1747 
1748  template <typename Number>
1749  void
1753  {
1754  Assert(P_val_dev != nullptr, ExcNotInitialized());
1755  Assert(P_row_ptr_dev != nullptr, ExcNotInitialized());
1756  Assert(P_column_index_dev != nullptr, ExcNotInitialized());
1757  AssertDimension(dst.size(), static_cast<unsigned int>(n_rows));
1758  AssertDimension(src.size(), static_cast<unsigned int>(n_rows));
1759  Assert(tmp_dev != nullptr, ExcInternalError());
1760 
1761  const Number *const src_dev = src.get_values();
1762  Number *const dst_dev = dst.get_values();
1763 
1764  // step 6: solve L*z = alpha*x
1765  const Number alpha = internal::NumberType<Number>::value(1.);
1766  cusparseStatus_t status =
1767  cusparseXcsrsv2_solve(cusparse_handle,
1768  CUSPARSE_OPERATION_NON_TRANSPOSE,
1769  n_rows,
1770  n_nonzero_elements,
1771  &alpha,
1772  descr_L,
1773  P_val_dev.get(),
1774  P_row_ptr_dev,
1775  P_column_index_dev,
1776  info_L,
1777  src_dev,
1778  tmp_dev.get(),
1779  policy_L,
1780  buffer_dev.get());
1781  AssertCusparse(status);
1782 
1783  // step 7: solve U*y = alpha*z
1784  status = cusparseXcsrsv2_solve(cusparse_handle,
1785  CUSPARSE_OPERATION_NON_TRANSPOSE,
1786  n_rows,
1787  n_nonzero_elements,
1788  &alpha,
1789  descr_U,
1790  P_val_dev.get(),
1791  P_row_ptr_dev,
1792  P_column_index_dev,
1793  info_U,
1794  tmp_dev.get(),
1795  dst_dev,
1796  policy_U,
1797  buffer_dev.get());
1798  AssertCusparse(status);
1799  }
1800 
1801 
1802 
1803  template <typename Number>
1804  void
1807  const LinearAlgebra::CUDAWrappers::Vector<Number> & /*src*/) const
1808  {
1809  Assert(false, ExcNotImplemented());
1810  }
1811 
1812 
1813 
1814  // explicit instantiations
1815  template class PreconditionIC<float>;
1816  template class PreconditionIC<double>;
1817  // template class PreconditionIC<cuComplex>;
1818  // template class PreconditionIC<cuDoubleComplex>;
1819  template class PreconditionILU<float>;
1820  template class PreconditionILU<double>;
1821  // template class PreconditionILU<cuComplex>;
1822  // template class PreconditionILU<cuDoubleComplex>;
1823 } // namespace CUDAWrappers
1824 
cusparseSolvePolicy_t policy_M
std::unique_ptr< Number[], void(*)(Number *)> P_val_dev
void Tvmult(LinearAlgebra::CUDAWrappers::Vector< Number > &dst, const LinearAlgebra::CUDAWrappers::Vector< Number > &src) const
std::unique_ptr< Number[], void(*)(Number *)> tmp_dev
void initialize(const SparseMatrix< Number > &matrix, const AdditionalData &additional_data=AdditionalData())
cusparseSolvePolicy_t policy_Lt
cusparseSolvePolicy_t policy_L
void vmult(LinearAlgebra::CUDAWrappers::Vector< Number > &dst, const LinearAlgebra::CUDAWrappers::Vector< Number > &src) const
PreconditionIC(const Utilities::CUDA::Handle &handle)
std::unique_ptr< void, void(*)(void *)> buffer_dev
cusparseSolvePolicy_t policy_M
void vmult(LinearAlgebra::CUDAWrappers::Vector< Number > &dst, const LinearAlgebra::CUDAWrappers::Vector< Number > &src) const
std::unique_ptr< void, void(*)(void *)> buffer_dev
void initialize(const SparseMatrix< Number > &matrix, const AdditionalData &additional_data=AdditionalData())
cusparseSolvePolicy_t policy_U
PreconditionILU(const Utilities::CUDA::Handle &handle)
std::unique_ptr< Number[], void(*)(Number *)> P_val_dev
std::unique_ptr< Number[], void(*)(Number *)> tmp_dev
void Tvmult(LinearAlgebra::CUDAWrappers::Vector< Number > &dst, const LinearAlgebra::CUDAWrappers::Vector< Number > &src) const
cusparseSolvePolicy_t policy_L
virtual size_type size() const override
Definition: cuda_vector.h:381
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:461
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:462
#define AssertCusparse(error_code)
Definition: exceptions.h:1990
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcNotInitialized()
#define Assert(cond, exc)
Definition: exceptions.h:1583
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1756
#define AssertNothrowCusparse(error_code)
Definition: exceptions.h:2020
#define AssertThrow(cond, exc)
Definition: exceptions.h:1672
static const char A
Default CUDA
Definition: memory_space.h:55
void delete_device_data(Number *device_ptr) noexcept
Definition: cuda.h:121
static constexpr DEAL_II_HOST_DEVICE_ALWAYS_INLINE const T & value(const T &t)
Definition: numbers.h:733