-
Notifications
You must be signed in to change notification settings - Fork 182
/
HYPRE_parcsr_ls.h
5005 lines (4375 loc) · 185 KB
/
HYPRE_parcsr_ls.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
/******************************************************************************
* Copyright (c) 1998 Lawrence Livermore National Security, LLC and other
* HYPRE Project Developers. See the top-level COPYRIGHT file for details.
*
* SPDX-License-Identifier: (Apache-2.0 OR MIT)
******************************************************************************/
#ifndef HYPRE_PARCSR_LS_HEADER
#define HYPRE_PARCSR_LS_HEADER
#include "HYPRE_utilities.h"
#include "HYPRE_seq_mv.h"
#include "HYPRE_parcsr_mv.h"
#include "HYPRE_IJ_mv.h"
#include "HYPRE_lobpcg.h"
#ifdef __cplusplus
extern "C" {
#endif
/*--------------------------------------------------------------------------
*--------------------------------------------------------------------------*/
/**
* @defgroup ParCSRSolvers ParCSR Solvers
*
* Linear solvers for sparse matrix systems. These solvers use matrix/vector
* storage schemes that are taylored for general sparse matrix systems.
*
* @{
**/
/*--------------------------------------------------------------------------
*--------------------------------------------------------------------------*/
/**
* @name ParCSR Solvers
*
* @{
**/
/**
* The solver object.
**/
typedef HYPRE_Int (*HYPRE_PtrToParSolverFcn)(HYPRE_Solver,
HYPRE_ParCSRMatrix,
HYPRE_ParVector,
HYPRE_ParVector);
#ifndef HYPRE_MODIFYPC
#define HYPRE_MODIFYPC
typedef HYPRE_Int (*HYPRE_PtrToModifyPCFcn)(HYPRE_Solver,
HYPRE_Int,
HYPRE_Real);
#endif
/**@}*/
/*--------------------------------------------------------------------------
*--------------------------------------------------------------------------*/
/**
* @name ParCSR BoomerAMG Solver and Preconditioner
*
* Parallel unstructured algebraic multigrid solver and preconditioner
*
* @{
**/
/**
* Create a solver object.
**/
HYPRE_Int HYPRE_BoomerAMGCreate(HYPRE_Solver *solver);
/**
* Destroy a solver object.
**/
HYPRE_Int HYPRE_BoomerAMGDestroy(HYPRE_Solver solver);
/**
* Set up the BoomerAMG solver or preconditioner.
* If used as a preconditioner, this function should be passed
* to the iterative solver \e SetPrecond function.
*
* @param solver [IN] object to be set up.
* @param A [IN] ParCSR matrix used to construct the solver/preconditioner.
* @param b Ignored by this function.
* @param x Ignored by this function.
**/
HYPRE_Int HYPRE_BoomerAMGSetup(HYPRE_Solver solver,
HYPRE_ParCSRMatrix A,
HYPRE_ParVector b,
HYPRE_ParVector x);
/**
* Solve the system or apply AMG as a preconditioner.
* If used as a preconditioner, this function should be passed
* to the iterative solver \e SetPrecond function.
*
* @param solver [IN] solver or preconditioner object to be applied.
* @param A [IN] ParCSR matrix, matrix of the linear system to be solved
* @param b [IN] right hand side of the linear system to be solved
* @param x [OUT] approximated solution of the linear system to be solved
**/
HYPRE_Int HYPRE_BoomerAMGSolve(HYPRE_Solver solver,
HYPRE_ParCSRMatrix A,
HYPRE_ParVector b,
HYPRE_ParVector x);
/**
* Solve the transpose system \f$A^T x = b\f$ or apply AMG as a preconditioner
* to the transpose system . Note that this function should only be used
* when preconditioning CGNR with BoomerAMG. It can only be used with
* Jacobi smoothing (relax_type 0 or 7) and without CF smoothing,
* i.e relax_order needs to be set to 0.
* If used as a preconditioner, this function should be passed
* to the iterative solver \e SetPrecond function.
*
* @param solver [IN] solver or preconditioner object to be applied.
* @param A [IN] ParCSR matrix
* @param b [IN] right hand side of the linear system to be solved
* @param x [OUT] approximated solution of the linear system to be solved
**/
HYPRE_Int HYPRE_BoomerAMGSolveT(HYPRE_Solver solver,
HYPRE_ParCSRMatrix A,
HYPRE_ParVector b,
HYPRE_ParVector x);
/**
* Recovers old default for coarsening and interpolation, i.e Falgout
* coarsening and untruncated modified classical interpolation.
* This option might be preferred for 2 dimensional problems.
**/
HYPRE_Int HYPRE_BoomerAMGSetOldDefault(HYPRE_Solver solver);
/**
* Returns the residual.
**/
HYPRE_Int HYPRE_BoomerAMGGetResidual(HYPRE_Solver solver,
HYPRE_ParVector *residual);
/**
* Returns the number of iterations taken.
**/
HYPRE_Int HYPRE_BoomerAMGGetNumIterations(HYPRE_Solver solver,
HYPRE_Int *num_iterations);
/**
* Returns cumulative num of nonzeros for A and P operators
**/
HYPRE_Int HYPRE_BoomerAMGGetCumNnzAP(HYPRE_Solver solver,
HYPRE_Real *cum_nnz_AP);
/**
* Activates cumulative num of nonzeros for A and P operators.
* Needs to be set to a positive number for activation.
**/
HYPRE_Int HYPRE_BoomerAMGSetCumNnzAP(HYPRE_Solver solver,
HYPRE_Real cum_nnz_AP);
/**
* Returns the norm of the final relative residual.
**/
HYPRE_Int HYPRE_BoomerAMGGetFinalRelativeResidualNorm(HYPRE_Solver solver,
HYPRE_Real *rel_resid_norm);
/**
* (Optional) Sets the size of the system of PDEs, if using the systems version.
* The default is 1, i.e. a scalar system.
**/
HYPRE_Int HYPRE_BoomerAMGSetNumFunctions(HYPRE_Solver solver,
HYPRE_Int num_functions);
/**
* (Optional) Sets the mapping that assigns the function to each variable,
* if using the systems version. If no assignment is made and the number of
* functions is k > 1, the mapping generated is (0,1,...,k-1,0,1,...,k-1,...).
**/
HYPRE_Int HYPRE_BoomerAMGSetDofFunc(HYPRE_Solver solver,
HYPRE_Int *dof_func);
/**
* (Optional) Set the type convergence checking
* 0: (default) norm(r)/norm(b), or norm(r) when b == 0
* 1: nomr(r) / norm(r_0)
**/
HYPRE_Int HYPRE_BoomerAMGSetConvergeType(HYPRE_Solver solver,
HYPRE_Int type);
/**
* (Optional) Set the convergence tolerance, if BoomerAMG is used
* as a solver. If it is used as a preconditioner, it should be set to 0.
* The default is 1.e-6.
**/
HYPRE_Int HYPRE_BoomerAMGSetTol(HYPRE_Solver solver,
HYPRE_Real tol);
/**
* (Optional) Sets maximum number of iterations, if BoomerAMG is used
* as a solver. If it is used as a preconditioner, it should be set to 1.
* The default is 20.
**/
HYPRE_Int HYPRE_BoomerAMGSetMaxIter(HYPRE_Solver solver,
HYPRE_Int max_iter);
/**
* (Optional)
**/
HYPRE_Int HYPRE_BoomerAMGSetMinIter(HYPRE_Solver solver,
HYPRE_Int min_iter);
/**
* (Optional) Sets maximum size of coarsest grid.
* The default is 9.
**/
HYPRE_Int HYPRE_BoomerAMGSetMaxCoarseSize(HYPRE_Solver solver,
HYPRE_Int max_coarse_size);
/**
* (Optional) Sets minimum size of coarsest grid.
* The default is 1.
**/
HYPRE_Int HYPRE_BoomerAMGSetMinCoarseSize(HYPRE_Solver solver,
HYPRE_Int min_coarse_size);
/**
* (Optional) Sets maximum number of multigrid levels.
* The default is 25.
**/
HYPRE_Int HYPRE_BoomerAMGSetMaxLevels(HYPRE_Solver solver,
HYPRE_Int max_levels);
/**
* (Optional) Sets cut factor for choosing isolated points
* during coarsening according to the rows' density. The default is 0.
* If nnzrow > coarsen_cut_factor*avg_nnzrow, where avg_nnzrow is the
* average number of nonzeros per row of the global matrix, holds for
* a given row, it is set as fine, and interpolation weights are not computed.
**/
HYPRE_Int HYPRE_BoomerAMGSetCoarsenCutFactor(HYPRE_Solver solver,
HYPRE_Int coarsen_cut_factor);
/**
* (Optional) Sets AMG strength threshold. The default is 0.25.
* For 2D Laplace operators, 0.25 is a good value, for 3D Laplace
* operators, 0.5 or 0.6 is a better value. For elasticity problems,
* a large strength threshold, such as 0.9, is often better.
**/
HYPRE_Int HYPRE_BoomerAMGSetStrongThreshold(HYPRE_Solver solver,
HYPRE_Real strong_threshold);
/**
* (Optional) The strong threshold for R is strong connections used
* in building an approximate ideal restriction.
* Default value is 0.25.
**/
HYPRE_Int HYPRE_BoomerAMGSetStrongThresholdR(HYPRE_Solver solver,
HYPRE_Real strong_threshold);
/**
* (Optional) The filter threshold for R is used to eliminate small entries
* of the approximate ideal restriction after building it.
* Default value is 0.0, which disables filtering.
**/
HYPRE_Int HYPRE_BoomerAMGSetFilterThresholdR(HYPRE_Solver solver,
HYPRE_Real filter_threshold);
/**
* (Optional) Deprecated. This routine now has no effect.
**/
HYPRE_Int HYPRE_BoomerAMGSetSCommPkgSwitch(HYPRE_Solver solver,
HYPRE_Real S_commpkg_switch);
/**
* (Optional) Sets a parameter to modify the definition of strength for
* diagonal dominant portions of the matrix. The default is 0.9.
* If \e max_row_sum is 1, no checking for diagonally dominant rows is
* performed.
**/
HYPRE_Int HYPRE_BoomerAMGSetMaxRowSum(HYPRE_Solver solver,
HYPRE_Real max_row_sum);
/**
* (Optional) Defines which parallel coarsening algorithm is used.
* There are the following options for \e coarsen_type:
*
* - 0 : CLJP-coarsening (a parallel coarsening algorithm using independent sets.
* - 1 : classical Ruge-Stueben coarsening on each processor, no boundary treatment
(not recommended!)
* - 3 : classical Ruge-Stueben coarsening on each processor, followed by a third pass,
which adds coarse points on the boundaries
* - 6 : Falgout coarsening (uses 1 first, followed by CLJP using the interior coarse points
* generated by 1 as its first independent set)
* - 7 : CLJP-coarsening (using a fixed random vector, for debugging purposes only)
* - 8 : PMIS-coarsening (a parallel coarsening algorithm using independent sets, generating
* lower complexities than CLJP, might also lead to slower convergence)
* - 9 : PMIS-coarsening (using a fixed random vector, for debugging purposes only)
* - 10 : HMIS-coarsening (uses one pass Ruge-Stueben on each processor independently, followed
* by PMIS using the interior C-points generated as its first independent set)
* - 11 : one-pass Ruge-Stueben coarsening on each processor, no boundary treatment
(not recommended!)
* - 21 : CGC coarsening by M. Griebel, B. Metsch and A. Schweitzer
* - 22 : CGC-E coarsening by M. Griebel, B. Metsch and A.Schweitzer
*
* The default is 10.
**/
HYPRE_Int HYPRE_BoomerAMGSetCoarsenType(HYPRE_Solver solver,
HYPRE_Int coarsen_type);
/**
* (Optional) Defines the non-Galerkin drop-tolerance
* for sparsifying coarse grid operators and thus reducing communication.
* Value specified here is set on all levels.
* This routine should be used before HYPRE_BoomerAMGSetLevelNonGalerkinTol, which
* then can be used to change individual levels if desired
**/
HYPRE_Int HYPRE_BoomerAMGSetNonGalerkinTol (HYPRE_Solver solver,
HYPRE_Real nongalerkin_tol);
/**
* (Optional) Defines the level specific non-Galerkin drop-tolerances
* for sparsifying coarse grid operators and thus reducing communication.
* A drop-tolerance of 0.0 means to skip doing non-Galerkin on that
* level. The maximum drop tolerance for a level is 1.0, although
* much smaller values such as 0.03 or 0.01 are recommended.
*
* Note that if the user wants to set a specific tolerance on all levels,
* HYPRE_BooemrAMGSetNonGalerkinTol should be used. Individual levels
* can then be changed using this routine.
*
* In general, it is safer to drop more aggressively on coarser levels.
* For instance, one could use 0.0 on the finest level, 0.01 on the second level and
* then using 0.05 on all remaining levels. The best way to achieve this is
* to set 0.05 on all levels with HYPRE_BoomerAMGSetNonGalerkinTol and then
* change the tolerance on level 0 to 0.0 and the tolerance on level 1 to 0.01
* with HYPRE_BoomerAMGSetLevelNonGalerkinTol.
* Like many AMG parameters, these drop tolerances can be tuned. It is also common
* to delay the start of the non-Galerkin process further to a later level than
* level 1.
*
* @param solver [IN] solver or preconditioner object to be applied.
* @param nongalerkin_tol [IN] level specific drop tolerance
* @param level [IN] level on which drop tolerance is used
**/
HYPRE_Int HYPRE_BoomerAMGSetLevelNonGalerkinTol (HYPRE_Solver solver,
HYPRE_Real nongalerkin_tol,
HYPRE_Int level);
/**
* (Optional) Defines the non-Galerkin drop-tolerance (old version)
**/
HYPRE_Int HYPRE_BoomerAMGSetNonGalerkTol (HYPRE_Solver solver,
HYPRE_Int nongalerk_num_tol,
HYPRE_Real *nongalerk_tol);
/**
* (Optional) Defines whether local or global measures are used.
**/
HYPRE_Int HYPRE_BoomerAMGSetMeasureType(HYPRE_Solver solver,
HYPRE_Int measure_type);
/**
* (Optional) Defines the number of levels of aggressive coarsening.
* The default is 0, i.e. no aggressive coarsening.
**/
HYPRE_Int HYPRE_BoomerAMGSetAggNumLevels(HYPRE_Solver solver,
HYPRE_Int agg_num_levels);
/**
* (Optional) Defines the degree of aggressive coarsening.
* The default is 1. Larger numbers lead to less aggressive
* coarsening.
**/
HYPRE_Int HYPRE_BoomerAMGSetNumPaths(HYPRE_Solver solver,
HYPRE_Int num_paths);
/**
* (optional) Defines the number of pathes for CGC-coarsening.
**/
HYPRE_Int HYPRE_BoomerAMGSetCGCIts (HYPRE_Solver solver,
HYPRE_Int its);
/**
* (Optional) Sets whether to use the nodal systems coarsening.
* Should be used for linear systems generated from systems of PDEs.
* The default is 0 (unknown-based coarsening,
* only coarsens within same function).
* For the remaining options a nodal matrix is generated by
* applying a norm to the nodal blocks and applying the coarsening
* algorithm to this matrix.
* - 1 : Frobenius norm
* - 2 : sum of absolute values of elements in each block
* - 3 : largest element in each block (not absolute value)
* - 4 : row-sum norm
* - 6 : sum of all values in each block
**/
HYPRE_Int HYPRE_BoomerAMGSetNodal(HYPRE_Solver solver,
HYPRE_Int nodal);
/**
* (Optional) Sets whether to give special treatment to diagonal elements in
* the nodal systems version.
* The default is 0.
* If set to 1, the diagonal entry is set to the negative sum of all off
* diagonal entries.
* If set to 2, the signs of all diagonal entries are inverted.
*/
HYPRE_Int HYPRE_BoomerAMGSetNodalDiag(HYPRE_Solver solver,
HYPRE_Int nodal_diag);
/*
* (Optional) Sets whether to keep same sign in S for nodal > 0
* The default is 0, i.e., discard those elements.
*/
HYPRE_Int HYPRE_BoomerAMGSetKeepSameSign(HYPRE_Solver solver,
HYPRE_Int keep_same_sign);
/**
* (Optional) Defines which parallel interpolation operator is used.
* There are the following options for \e interp_type:
*
* - 0 : classical modified interpolation
* - 1 : LS interpolation (for use with GSMG)
* - 2 : classical modified interpolation for hyperbolic PDEs
* - 3 : direct interpolation (with separation of weights) (also for GPU use)
* - 4 : multipass interpolation
* - 5 : multipass interpolation (with separation of weights)
* - 6 : extended+i interpolation (also for GPU use)
* - 7 : extended+i (if no common C neighbor) interpolation
* - 8 : standard interpolation
* - 9 : standard interpolation (with separation of weights)
* - 10 : classical block interpolation (for use with nodal systems version only)
* - 11 : classical block interpolation (for use with nodal systems version only)
* with diagonalized diagonal blocks
* - 12 : FF interpolation
* - 13 : FF1 interpolation
* - 14 : extended interpolation (also for GPU use)
* - 15 : interpolation with adaptive weights (GPU use only)
* - 16 : extended interpolation in matrix-matrix form
* - 17 : extended+i interpolation in matrix-matrix form
* - 18 : extended+e interpolation in matrix-matrix form
*
* The default is ext+i interpolation (interp_type 6) trunctated to at most 4
* elements per row. (see HYPRE_BoomerAMGSetPMaxElmts).
**/
HYPRE_Int HYPRE_BoomerAMGSetInterpType(HYPRE_Solver solver,
HYPRE_Int interp_type);
/**
* (Optional) Defines a truncation factor for the interpolation. The default is 0.
**/
HYPRE_Int HYPRE_BoomerAMGSetTruncFactor(HYPRE_Solver solver,
HYPRE_Real trunc_factor);
/**
* (Optional) Defines the maximal number of elements per row for the interpolation.
* The default is 4. To turn off truncation, it needs to be set to 0.
**/
HYPRE_Int HYPRE_BoomerAMGSetPMaxElmts(HYPRE_Solver solver,
HYPRE_Int P_max_elmts);
/**
* (Optional) Defines whether separation of weights is used
* when defining strength for standard interpolation or
* multipass interpolation.
* Default: 0, i.e. no separation of weights used.
**/
HYPRE_Int HYPRE_BoomerAMGSetSepWeight(HYPRE_Solver solver,
HYPRE_Int sep_weight);
/**
* (Optional) Defines the interpolation used on levels of aggressive coarsening
* The default is 4, i.e. multipass interpolation.
* The following options exist:
*
* - 1 : 2-stage extended+i interpolation
* - 2 : 2-stage standard interpolation
* - 3 : 2-stage extended interpolation
* - 4 : multipass interpolation
* - 5 : 2-stage extended interpolation in matrix-matrix form
* - 6 : 2-stage extended+i interpolation in matrix-matrix form
* - 7 : 2-stage extended+e interpolation in matrix-matrix form
**/
HYPRE_Int HYPRE_BoomerAMGSetAggInterpType(HYPRE_Solver solver,
HYPRE_Int agg_interp_type);
/**
* (Optional) Defines the truncation factor for the
* interpolation used for aggressive coarsening.
* The default is 0.
**/
HYPRE_Int HYPRE_BoomerAMGSetAggTruncFactor(HYPRE_Solver solver,
HYPRE_Real agg_trunc_factor);
/**
* (Optional) Defines the truncation factor for the
* matrices P1 and P2 which are used to build 2-stage interpolation.
* The default is 0.
**/
HYPRE_Int HYPRE_BoomerAMGSetAggP12TruncFactor(HYPRE_Solver solver,
HYPRE_Real agg_P12_trunc_factor);
/**
* (Optional) Defines the maximal number of elements per row for the
* interpolation used for aggressive coarsening.
* The default is 0.
**/
HYPRE_Int HYPRE_BoomerAMGSetAggPMaxElmts(HYPRE_Solver solver,
HYPRE_Int agg_P_max_elmts);
/**
* (Optional) Defines the maximal number of elements per row for the
* matrices P1 and P2 which are used to build 2-stage interpolation.
* The default is 0.
**/
HYPRE_Int HYPRE_BoomerAMGSetAggP12MaxElmts(HYPRE_Solver solver,
HYPRE_Int agg_P12_max_elmts);
/**
* (Optional) Allows the user to incorporate additional vectors
* into the interpolation for systems AMG, e.g. rigid body modes for
* linear elasticity problems.
* This can only be used in context with nodal coarsening and still
* requires the user to choose an interpolation.
**/
HYPRE_Int HYPRE_BoomerAMGSetInterpVectors (HYPRE_Solver solver,
HYPRE_Int num_vectors,
HYPRE_ParVector *interp_vectors );
/**
* (Optional) Defines the interpolation variant used for
* HYPRE_BoomerAMGSetInterpVectors:
* - 1 : GM approach 1
* - 2 : GM approach 2 (to be preferred over 1)
* - 3 : LN approach
**/
HYPRE_Int HYPRE_BoomerAMGSetInterpVecVariant (HYPRE_Solver solver,
HYPRE_Int var );
/**
* (Optional) Defines the maximal elements per row for Q, the additional
* columns added to the original interpolation matrix P, to reduce complexity.
* The default is no truncation.
**/
HYPRE_Int HYPRE_BoomerAMGSetInterpVecQMax (HYPRE_Solver solver,
HYPRE_Int q_max );
/**
* (Optional) Defines a truncation factor for Q, the additional
* columns added to the original interpolation matrix P, to reduce complexity.
* The default is no truncation.
**/
HYPRE_Int HYPRE_BoomerAMGSetInterpVecAbsQTrunc (HYPRE_Solver solver,
HYPRE_Real q_trunc );
/**
* (Optional) Specifies the use of GSMG - geometrically smooth
* coarsening and interpolation. Currently any nonzero value for
* gsmg will lead to the use of GSMG.
* The default is 0, i.e. (GSMG is not used)
**/
HYPRE_Int HYPRE_BoomerAMGSetGSMG(HYPRE_Solver solver,
HYPRE_Int gsmg);
/**
* (Optional) Defines the number of sample vectors used in GSMG
* or LS interpolation.
**/
HYPRE_Int HYPRE_BoomerAMGSetNumSamples(HYPRE_Solver solver,
HYPRE_Int num_samples);
/**
* (Optional) Defines the type of cycle.
* For a V-cycle, set \e cycle_type to 1, for a W-cycle
* set \e cycle_type to 2. The default is 1.
**/
HYPRE_Int HYPRE_BoomerAMGSetCycleType(HYPRE_Solver solver,
HYPRE_Int cycle_type);
/**
* (Optional) Specifies the use of Full multigrid cycle.
* The default is 0.
**/
HYPRE_Int
HYPRE_BoomerAMGSetFCycle( HYPRE_Solver solver,
HYPRE_Int fcycle );
/**
* (Optional) Defines use of an additive V(1,1)-cycle using the
* classical additive method starting at level 'addlvl'.
* The multiplicative approach is used on levels 0, ...'addlvl+1'.
* 'addlvl' needs to be > -1 for this to have an effect.
* Can only be used with weighted Jacobi and l1-Jacobi(default).
*
* Can only be used when AMG is used as a preconditioner !!!
**/
HYPRE_Int HYPRE_BoomerAMGSetAdditive(HYPRE_Solver solver,
HYPRE_Int addlvl);
/**
* (Optional) Defines use of an additive V(1,1)-cycle using the
* mult-additive method starting at level 'addlvl'.
* The multiplicative approach is used on levels 0, ...'addlvl+1'.
* 'addlvl' needs to be > -1 for this to have an effect.
* Can only be used with weighted Jacobi and l1-Jacobi(default).
*
* Can only be used when AMG is used as a preconditioner !!!
**/
HYPRE_Int HYPRE_BoomerAMGSetMultAdditive(HYPRE_Solver solver,
HYPRE_Int addlvl);
/**
* (Optional) Defines use of an additive V(1,1)-cycle using the
* simplified mult-additive method starting at level 'addlvl'.
* The multiplicative approach is used on levels 0, ...'addlvl+1'.
* 'addlvl' needs to be > -1 for this to have an effect.
* Can only be used with weighted Jacobi and l1-Jacobi(default).
*
* Can only be used when AMG is used as a preconditioner !!!
**/
HYPRE_Int HYPRE_BoomerAMGSetSimple(HYPRE_Solver solver,
HYPRE_Int addlvl);
/**
* (Optional) Defines last level where additive, mult-additive
* or simple cycle is used.
* The multiplicative approach is used on levels > add_last_lvl.
*
* Can only be used when AMG is used as a preconditioner !!!
**/
HYPRE_Int HYPRE_BoomerAMGSetAddLastLvl(HYPRE_Solver solver,
HYPRE_Int add_last_lvl);
/**
* (Optional) Defines the truncation factor for the
* smoothed interpolation used for mult-additive or simple method.
* The default is 0.
**/
HYPRE_Int HYPRE_BoomerAMGSetMultAddTruncFactor(HYPRE_Solver solver,
HYPRE_Real add_trunc_factor);
/**
* (Optional) Defines the maximal number of elements per row for the
* smoothed interpolation used for mult-additive or simple method.
* The default is 0.
**/
HYPRE_Int HYPRE_BoomerAMGSetMultAddPMaxElmts(HYPRE_Solver solver,
HYPRE_Int add_P_max_elmts);
/**
* (Optional) Defines the relaxation type used in the (mult)additive cycle
* portion (also affects simple method.)
* The default is 18 (L1-Jacobi).
* Currently the only other option allowed is 0 (Jacobi) which should be
* used in combination with HYPRE_BoomerAMGSetAddRelaxWt.
**/
HYPRE_Int HYPRE_BoomerAMGSetAddRelaxType(HYPRE_Solver solver,
HYPRE_Int add_rlx_type);
/**
* (Optional) Defines the relaxation weight used for Jacobi within the
* (mult)additive or simple cycle portion.
* The default is 1.
* The weight only affects the Jacobi method, and has no effect on L1-Jacobi
**/
HYPRE_Int HYPRE_BoomerAMGSetAddRelaxWt(HYPRE_Solver solver,
HYPRE_Real add_rlx_wt);
/**
* (Optional) Sets maximal size for agglomeration or redundant coarse grid solve.
* When the system is smaller than this threshold, sequential AMG is used
* on process 0 or on all remaining active processes (if redundant = 1 ).
**/
HYPRE_Int HYPRE_BoomerAMGSetSeqThreshold(HYPRE_Solver solver,
HYPRE_Int seq_threshold);
/**
* (Optional) operates switch for redundancy. Needs to be used with
* HYPRE_BoomerAMGSetSeqThreshold. Default is 0, i.e. no redundancy.
**/
HYPRE_Int HYPRE_BoomerAMGSetRedundant(HYPRE_Solver solver,
HYPRE_Int redundant);
/**
* (Optional) Defines the number of sweeps for the fine and coarse grid,
* the up and down cycle.
*
* Note: This routine will be phased out!!!!
* Use HYPRE_BoomerAMGSetNumSweeps or HYPRE_BoomerAMGSetCycleNumSweeps instead.
**/
HYPRE_Int HYPRE_BoomerAMGSetNumGridSweeps(HYPRE_Solver solver,
HYPRE_Int *num_grid_sweeps);
/**
* (Optional) Sets the number of sweeps. On the finest level, the up and
* the down cycle the number of sweeps are set to \e num_sweeps and on the
* coarsest level to 1. The default is 1.
**/
HYPRE_Int HYPRE_BoomerAMGSetNumSweeps(HYPRE_Solver solver,
HYPRE_Int num_sweeps);
/**
* (Optional) Sets the number of sweeps at a specified cycle.
* There are the following options for \e k:
*
* - 1 : the down cycle
* - 2 : the up cycle
* - 3 : the coarsest level
**/
HYPRE_Int HYPRE_BoomerAMGSetCycleNumSweeps(HYPRE_Solver solver,
HYPRE_Int num_sweeps,
HYPRE_Int k);
/**
* (Optional) Defines which smoother is used on the fine and coarse grid,
* the up and down cycle.
*
* Note: This routine will be phased out!!!!
* Use HYPRE_BoomerAMGSetRelaxType or HYPRE_BoomerAMGSetCycleRelaxType instead.
**/
HYPRE_Int HYPRE_BoomerAMGSetGridRelaxType(HYPRE_Solver solver,
HYPRE_Int *grid_relax_type);
/**
* (Optional) Defines the smoother to be used. It uses the given
* smoother on the fine grid, the up and
* the down cycle and sets the solver on the coarsest level to Gaussian
* elimination (9). The default is \f$\ell_1\f$-Gauss-Seidel, forward solve (13)
* on the down cycle and backward solve (14) on the up cycle.
*
* There are the following options for \e relax_type:
*
* - 0 : Jacobi
* - 1 : Gauss-Seidel, sequential (very slow!)
* - 2 : Gauss-Seidel, interior points in parallel, boundary sequential (slow!)
* - 3 : hybrid Gauss-Seidel or SOR, forward solve
* - 4 : hybrid Gauss-Seidel or SOR, backward solve
* - 5 : hybrid chaotic Gauss-Seidel (works only with OpenMP)
* - 6 : hybrid symmetric Gauss-Seidel or SSOR
* - 7 : Jacobi (uses Matvec)
* - 8 : \f$\ell_1\f$-scaled hybrid symmetric Gauss-Seidel
* - 9 : Gaussian elimination (only on coarsest level)
* - 10 : On-processor direct forward solve for matrices with
* triangular structure
* - 11 : Two Stage approximation to GS. Uses the strict lower
* part of the diagonal matrix
* - 12 : Two Stage approximation to GS. Uses the strict lower
* part of the diagonal matrix and a second iteration
* for additional error approximation
* - 13 : \f$\ell_1\f$ Gauss-Seidel, forward solve
* - 14 : \f$\ell_1\f$ Gauss-Seidel, backward solve
* - 15 : CG (warning - not a fixed smoother - may require FGMRES)
* - 16 : Chebyshev
* - 17 : FCF-Jacobi
* - 18 : \f$\ell_1\f$-scaled jacobi
* - 19 : Gaussian elimination (old version)
* - 21 : The same as 8 except forcing serialization on CPU (#OMP-thread = 1)
* - 29 : Direct solve: use Gaussian elimination & BLAS
* (with pivoting) (old version)
* - 30 : Kaczmarz
* - 88: The same methods as 8 with a convergent l1-term
* - 89: Symmetric l1-hybrid Gauss-Seidel (i.e., 13 followed by 14)
* - 98 : LU with pivoting
* - 99 : LU with pivoting
* -199 : Matvec with the inverse
**/
HYPRE_Int HYPRE_BoomerAMGSetRelaxType(HYPRE_Solver solver,
HYPRE_Int relax_type);
/**
* (Optional) Defines the smoother at a given cycle.
*
* For options of \e relax_type see description of HYPRE_BoomerAMGSetRelaxType.
* In addition, the following options for \e relax_type are available when choosing
* the coarsest level solver (k = 3):
*
* For coarsest level systems formed via a sub-communicator defined with active ranks:
* - 9 : hypre's internal Gaussian elimination (host only).
* - 99 : LU factorization with pivoting.
* - 199 : explicit (dense) inverse.
*
* For coarsest level systems formed via hypre_DataExchangeList:
* - 19 : hypre's internal Gaussian elimination (host only).
* - 98 : LU factorization with pivoting.
* - 198 : explicit (dense) inverse.
*
* Options for \e k are
*
* - 1 : the down cycle
* - 2 : the up cycle
* - 3 : the coarsest level
**/
HYPRE_Int HYPRE_BoomerAMGSetCycleRelaxType(HYPRE_Solver solver,
HYPRE_Int relax_type,
HYPRE_Int k);
/**
* (Optional) Defines in which order the points are relaxed. There are
* the following options for \e relax_order:
*
* - 0 : the points are relaxed in natural or lexicographic order on each processor
* - 1 : CF-relaxation is used, i.e on the fine grid and the down cycle the
* coarse points are relaxed first, followed by the fine points; on the
* up cycle the F-points are relaxed first, followed by the C-points.
* On the coarsest level, if an iterative scheme is used, the points
* are relaxed in lexicographic order.
*
* The default is 0.
**/
HYPRE_Int HYPRE_BoomerAMGSetRelaxOrder(HYPRE_Solver solver,
HYPRE_Int relax_order);
/**
* (Optional) Defines in which order the points are relaxed.
*
* See also HYPRE_BoomerAMGSetRelaxOrder.
**/
HYPRE_Int HYPRE_BoomerAMGSetGridRelaxPoints(HYPRE_Solver solver,
HYPRE_Int **grid_relax_points);
/**
* (Optional) Defines the relaxation weight for smoothed Jacobi and hybrid SOR.
*
* Note: This routine will be phased out!!!!
* Use HYPRE_BoomerAMGSetRelaxWt or HYPRE_BoomerAMGSetLevelRelaxWt instead.
**/
HYPRE_Int HYPRE_BoomerAMGSetRelaxWeight(HYPRE_Solver solver,
HYPRE_Real *relax_weight);
/**
* (Optional) Defines the relaxation weight for smoothed Jacobi and hybrid SOR
* on all levels.
*
* Values for \e relax_weight are
* - > 0 : this assigns the given relaxation weight on all levels
* - = 0 : the weight is determined on each level with the estimate
* \f$3 \over {4\|D^{-1/2}AD^{-1/2}\|}\f$, where \f$D\f$ is the diagonal of \f$A\f$
* (this should only be used with Jacobi)
* - = -k : the relaxation weight is determined with at most k CG steps on each level
* (this should only be used for symmetric positive definite problems)
*
* The default is 1.
**/
HYPRE_Int HYPRE_BoomerAMGSetRelaxWt(HYPRE_Solver solver,
HYPRE_Real relax_weight);
/**
* (Optional) Defines the relaxation weight for smoothed Jacobi and hybrid SOR
* on the user defined level. Note that the finest level is denoted 0, the
* next coarser level 1, etc. For nonpositive \e relax_weight, the parameter is
* determined on the given level as described for HYPRE_BoomerAMGSetRelaxWt.
* The default is 1.
**/
HYPRE_Int HYPRE_BoomerAMGSetLevelRelaxWt(HYPRE_Solver solver,
HYPRE_Real relax_weight,
HYPRE_Int level);
/**
* (Optional) Defines the outer relaxation weight for hybrid SOR.
* Note: This routine will be phased out!!!!
* Use HYPRE_BoomerAMGSetOuterWt or HYPRE_BoomerAMGSetLevelOuterWt instead.
**/
HYPRE_Int HYPRE_BoomerAMGSetOmega(HYPRE_Solver solver,
HYPRE_Real *omega);
/**
* (Optional) Defines the outer relaxation weight for hybrid SOR and SSOR
* on all levels.
*
* Values for \e omega are
* - > 0 : this assigns the same outer relaxation weight omega on each level
* - = -k : an outer relaxation weight is determined with at most k CG steps on each level
* (this only makes sense for symmetric positive definite problems and smoothers
* such as SSOR)
*
* The default is 1.
**/
HYPRE_Int HYPRE_BoomerAMGSetOuterWt(HYPRE_Solver solver,
HYPRE_Real omega);
/**
* (Optional) Defines the outer relaxation weight for hybrid SOR or SSOR
* on the user defined level. Note that the finest level is denoted 0, the
* next coarser level 1, etc. For nonpositive omega, the parameter is
* determined on the given level as described for HYPRE_BoomerAMGSetOuterWt.
* The default is 1.
**/
HYPRE_Int HYPRE_BoomerAMGSetLevelOuterWt(HYPRE_Solver solver,
HYPRE_Real omega,
HYPRE_Int level);
/**
* (Optional) Defines the Order for Chebyshev smoother.
* The default is 2 (valid options are 1-4).
**/
HYPRE_Int HYPRE_BoomerAMGSetChebyOrder(HYPRE_Solver solver,
HYPRE_Int order);
/**
* (Optional) Fraction of the spectrum to use for the Chebyshev smoother.
* The default is .3 (i.e., damp on upper 30% of the spectrum).
**/
HYPRE_Int HYPRE_BoomerAMGSetChebyFraction (HYPRE_Solver solver,
HYPRE_Real ratio);
/**
* (Optional) Defines whether matrix should be scaled.
* The default is 1 (i.e., scaled).
**/
HYPRE_Int HYPRE_BoomerAMGSetChebyScale (HYPRE_Solver solver,
HYPRE_Int scale);
/**
* (Optional) Defines which polynomial variant should be used.
* The default is 0 (i.e., scaled).
**/
HYPRE_Int HYPRE_BoomerAMGSetChebyVariant (HYPRE_Solver solver,
HYPRE_Int variant);
/**
* (Optional) Defines how to estimate eigenvalues.
* The default is 10 (i.e., 10 CG iterations are used to find extreme
* eigenvalues.) If eig_est=0, the largest eigenvalue is estimated
* using Gershgorin, the smallest is set to 0.
* If eig_est is a positive number n, n iterations of CG are used to
* determine the smallest and largest eigenvalue.
**/
HYPRE_Int HYPRE_BoomerAMGSetChebyEigEst (HYPRE_Solver solver,
HYPRE_Int eig_est);
/**
* (Optional) Enables the use of more complex smoothers.
* The following options exist for \e smooth_type:
*
* - 4 : FSAI (routines needed to set: HYPRE_BoomerAMGSetFSAIMaxSteps,
* HYPRE_BoomerAMGSetFSAIMaxStepSize, HYPRE_BoomerAMGSetFSAIEigMaxIters,
* HYPRE_BoomerAMGSetFSAIKapTolerance)
* - 5 : ParILUK (routines needed to set: HYPRE_ILUSetLevelOfFill, HYPRE_ILUSetType)
* - 6 : Schwarz (routines needed to set: HYPRE_BoomerAMGSetDomainType,
* HYPRE_BoomerAMGSetOverlap, HYPRE_BoomerAMGSetVariant,
* HYPRE_BoomerAMGSetSchwarzRlxWeight)
* - 7 : Pilut (routines needed to set: HYPRE_BoomerAMGSetDropTol,
* HYPRE_BoomerAMGSetMaxNzPerRow)
* - 8 : ParaSails (routines needed to set: HYPRE_BoomerAMGSetSym,
* HYPRE_BoomerAMGSetLevel, HYPRE_BoomerAMGSetFilter,
* HYPRE_BoomerAMGSetThreshold)
* - 9 : Euclid (routines needed to set: HYPRE_BoomerAMGSetEuclidFile)
*
* The default is 6. Also, if no smoother parameters are set via the routines
* mentioned in the table above, default values are used.
**/
HYPRE_Int HYPRE_BoomerAMGSetSmoothType(HYPRE_Solver solver,
HYPRE_Int smooth_type);
/**
* (Optional) Sets the number of levels for more complex smoothers.
* The smoothers,
* as defined by HYPRE_BoomerAMGSetSmoothType, will be used
* on level 0 (the finest level) through level \e smooth_num_levels-1.
* The default is 0, i.e. no complex smoothers are used.
**/
HYPRE_Int HYPRE_BoomerAMGSetSmoothNumLevels(HYPRE_Solver solver,
HYPRE_Int smooth_num_levels);
/**
* (Optional) Sets the number of sweeps for more complex smoothers.
* The default is 1.
**/
HYPRE_Int HYPRE_BoomerAMGSetSmoothNumSweeps(HYPRE_Solver solver,
HYPRE_Int smooth_num_sweeps);
/**
* (Optional) Defines which variant of the Schwarz method is used.
* The following options exist for \e variant:
*
* - 0 : hybrid multiplicative Schwarz method (no overlap across processor boundaries)
* - 1 : hybrid additive Schwarz method (no overlap across processor boundaries)
* - 2 : additive Schwarz method
* - 3 : hybrid multiplicative Schwarz method (with overlap across processor boundaries)
*
* The default is 0.
**/
HYPRE_Int HYPRE_BoomerAMGSetVariant(HYPRE_Solver solver,
HYPRE_Int variant);
/**
* (Optional) Defines the overlap for the Schwarz method.
* The following options exist for overlap:
*
* - 0 : no overlap
* - 1 : minimal overlap (default)
* - 2 : overlap generated by including all neighbors of domain boundaries
**/
HYPRE_Int HYPRE_BoomerAMGSetOverlap(HYPRE_Solver solver,
HYPRE_Int overlap);
/**
* (Optional) Defines the type of domain used for the Schwarz method.
* The following options exist for \e domain_type:
*
* - 0 : each point is a domain
* - 1 : each node is a domain (only of interest in "systems" AMG)
* - 2 : each domain is generated by agglomeration (default)