1 /****************************************************************************** 2 * Copyright 1998-2019 Lawrence Livermore National Security, LLC and other 3 * HYPRE Project Developers. See the top-level COPYRIGHT file for details. 4 * 5 * SPDX-License-Identifier: (Apache-2.0 OR MIT) 6 ******************************************************************************/ 7 8 #ifndef HYPRE_STRUCT_LS_HEADER 9 #define HYPRE_STRUCT_LS_HEADER 10 11 #include "HYPRE_utilities.h" 12 #include "HYPRE_struct_mv.h" 13 #include "HYPRE_lobpcg.h" 14 15 #ifdef __cplusplus 16 extern "C" { 17 #endif 18 19 /*-------------------------------------------------------------------------- 20 *--------------------------------------------------------------------------*/ 21 22 /** 23 * @defgroup StructSolvers Struct Solvers 24 * 25 * These solvers use matrix/vector storage schemes that are tailored 26 * to structured grid problems. 27 * 28 * @memo Linear solvers for structured grids 29 * 30 * @{ 31 **/ 32 33 /*-------------------------------------------------------------------------- 34 *--------------------------------------------------------------------------*/ 35 36 /** 37 * @name Struct Solvers 38 * 39 * @{ 40 **/ 41 42 struct hypre_StructSolver_struct; 43 /** 44 * The solver object. 45 **/ 46 typedef struct hypre_StructSolver_struct *HYPRE_StructSolver; 47 48 typedef HYPRE_Int (*HYPRE_PtrToStructSolverFcn)(HYPRE_StructSolver, 49 HYPRE_StructMatrix, 50 HYPRE_StructVector, 51 HYPRE_StructVector); 52 53 #ifndef HYPRE_MODIFYPC 54 #define HYPRE_MODIFYPC 55 /* if pc not defined, then may need HYPRE_SOLVER also */ 56 57 #ifndef HYPRE_SOLVER_STRUCT 58 #define HYPRE_SOLVER_STRUCT 59 struct hypre_Solver_struct; 60 typedef struct hypre_Solver_struct *HYPRE_Solver; 61 #endif 62 63 typedef HYPRE_Int (*HYPRE_PtrToModifyPCFcn)(HYPRE_Solver, 64 HYPRE_Int, 65 HYPRE_Real); 66 #endif 67 68 /**@}*/ 69 70 /*-------------------------------------------------------------------------- 71 *--------------------------------------------------------------------------*/ 72 73 /** 74 * @name Struct Jacobi Solver 75 * 76 * @{ 77 **/ 78 79 /** 80 * Create a solver object. 81 **/ 82 HYPRE_Int HYPRE_StructJacobiCreate(MPI_Comm comm, 83 HYPRE_StructSolver *solver); 84 85 /** 86 * Destroy a solver object. An object should be explicitly destroyed 87 * using this destructor when the user's code no longer needs direct 88 * access to it. Once destroyed, the object must not be referenced 89 * again. Note that the object may not be deallocated at the 90 * completion of this call, since there may be internal package 91 * references to the object. The object will then be destroyed when 92 * all internal reference counts go to zero. 93 **/ 94 HYPRE_Int HYPRE_StructJacobiDestroy(HYPRE_StructSolver solver); 95 96 /** 97 * Prepare to solve the system. The coefficient data in \e b and \e x is 98 * ignored here, but information about the layout of the data may be used. 99 **/ 100 HYPRE_Int HYPRE_StructJacobiSetup(HYPRE_StructSolver solver, 101 HYPRE_StructMatrix A, 102 HYPRE_StructVector b, 103 HYPRE_StructVector x); 104 105 /** 106 * Solve the system. 107 **/ 108 HYPRE_Int HYPRE_StructJacobiSolve(HYPRE_StructSolver solver, 109 HYPRE_StructMatrix A, 110 HYPRE_StructVector b, 111 HYPRE_StructVector x); 112 113 /** 114 * (Optional) Set the convergence tolerance. 115 **/ 116 HYPRE_Int HYPRE_StructJacobiSetTol(HYPRE_StructSolver solver, 117 HYPRE_Real tol); 118 119 120 HYPRE_Int HYPRE_StructJacobiGetTol(HYPRE_StructSolver solver, 121 HYPRE_Real *tol ); 122 123 /** 124 * (Optional) Set maximum number of iterations. 125 **/ 126 HYPRE_Int HYPRE_StructJacobiSetMaxIter(HYPRE_StructSolver solver, 127 HYPRE_Int max_iter); 128 129 HYPRE_Int HYPRE_StructJacobiGetMaxIter(HYPRE_StructSolver solver, 130 HYPRE_Int *max_iter ); 131 132 /** 133 * (Optional) Use a zero initial guess. This allows the solver to cut corners 134 * in the case where a zero initial guess is needed (e.g., for preconditioning) 135 * to reduce compuational cost. 136 **/ 137 HYPRE_Int HYPRE_StructJacobiSetZeroGuess(HYPRE_StructSolver solver); 138 139 HYPRE_Int HYPRE_StructJacobiGetZeroGuess(HYPRE_StructSolver solver, 140 HYPRE_Int *zeroguess ); 141 142 /** 143 * (Optional) Use a nonzero initial guess. This is the default behavior, but 144 * this routine allows the user to switch back after using \e SetZeroGuess. 145 **/ 146 HYPRE_Int HYPRE_StructJacobiSetNonZeroGuess(HYPRE_StructSolver solver); 147 148 /** 149 * Return the number of iterations taken. 150 **/ 151 HYPRE_Int HYPRE_StructJacobiGetNumIterations(HYPRE_StructSolver solver, 152 HYPRE_Int *num_iterations); 153 154 /** 155 * Return the norm of the final relative residual. 156 **/ 157 HYPRE_Int HYPRE_StructJacobiGetFinalRelativeResidualNorm(HYPRE_StructSolver solver, 158 HYPRE_Real *norm); 159 160 /**@}*/ 161 162 /*-------------------------------------------------------------------------- 163 *--------------------------------------------------------------------------*/ 164 165 /** 166 * @name Struct PFMG Solver 167 * 168 * PFMG is a semicoarsening multigrid solver that uses pointwise relaxation. 169 * For periodic problems, users should try to set the grid size in periodic 170 * dimensions to be as close to a power-of-two as possible. That is, if the 171 * grid size in a periodic dimension is given by \f$N = 2^m * M\f$ where \f$M\f$ 172 * is not a power-of-two, then \f$M\f$ should be as small as possible. Large 173 * values of \f$M\f$ will generally result in slower convergence rates. 174 * 175 * @{ 176 **/ 177 178 /** 179 * Create a solver object. 180 **/ 181 HYPRE_Int HYPRE_StructPFMGCreate(MPI_Comm comm, 182 HYPRE_StructSolver *solver); 183 184 /** 185 * Destroy a solver object. 186 **/ 187 HYPRE_Int HYPRE_StructPFMGDestroy(HYPRE_StructSolver solver); 188 189 /** 190 * Prepare to solve the system. The coefficient data in \e b and \e x is 191 * ignored here, but information about the layout of the data may be used. 192 **/ 193 HYPRE_Int HYPRE_StructPFMGSetup(HYPRE_StructSolver solver, 194 HYPRE_StructMatrix A, 195 HYPRE_StructVector b, 196 HYPRE_StructVector x); 197 198 /** 199 * Solve the system. 200 **/ 201 HYPRE_Int HYPRE_StructPFMGSolve(HYPRE_StructSolver solver, 202 HYPRE_StructMatrix A, 203 HYPRE_StructVector b, 204 HYPRE_StructVector x); 205 206 /** 207 * (Optional) Set the convergence tolerance. 208 **/ 209 HYPRE_Int HYPRE_StructPFMGSetTol(HYPRE_StructSolver solver, 210 HYPRE_Real tol); 211 212 HYPRE_Int HYPRE_StructPFMGGetTol (HYPRE_StructSolver solver, 213 HYPRE_Real *tol); 214 215 /** 216 * (Optional) Set maximum number of iterations. 217 **/ 218 HYPRE_Int HYPRE_StructPFMGSetMaxIter(HYPRE_StructSolver solver, 219 HYPRE_Int max_iter); 220 221 HYPRE_Int HYPRE_StructPFMGGetMaxIter(HYPRE_StructSolver solver, 222 HYPRE_Int *max_iter); 223 224 /** 225 * (Optional) Set maximum number of multigrid grid levels. 226 **/ 227 HYPRE_Int HYPRE_StructPFMGSetMaxLevels(HYPRE_StructSolver solver, 228 HYPRE_Int max_levels); 229 230 HYPRE_Int HYPRE_StructPFMGGetMaxLevels (HYPRE_StructSolver solver, 231 HYPRE_Int *max_levels ); 232 233 /** 234 * (Optional) Additionally require that the relative difference in 235 * successive iterates be small. 236 **/ 237 HYPRE_Int HYPRE_StructPFMGSetRelChange(HYPRE_StructSolver solver, 238 HYPRE_Int rel_change); 239 240 HYPRE_Int HYPRE_StructPFMGGetRelChange (HYPRE_StructSolver solver, 241 HYPRE_Int *rel_change); 242 243 /** 244 * (Optional) Use a zero initial guess. This allows the solver to cut corners 245 * in the case where a zero initial guess is needed (e.g., for preconditioning) 246 * to reduce compuational cost. 247 **/ 248 HYPRE_Int HYPRE_StructPFMGSetZeroGuess(HYPRE_StructSolver solver); 249 250 HYPRE_Int HYPRE_StructPFMGGetZeroGuess(HYPRE_StructSolver solver, 251 HYPRE_Int *zeroguess); 252 253 /** 254 * (Optional) Use a nonzero initial guess. This is the default behavior, but 255 * this routine allows the user to switch back after using \e SetZeroGuess. 256 **/ 257 HYPRE_Int HYPRE_StructPFMGSetNonZeroGuess(HYPRE_StructSolver solver); 258 259 /** 260 * (Optional) Set relaxation type. 261 * 262 * Current relaxation methods set by \e relax_type are: 263 * 264 * - 0 : Jacobi 265 * - 1 : Weighted Jacobi (default) 266 * - 2 : Red/Black Gauss-Seidel (symmetric: RB pre-relaxation, BR post-relaxation) 267 * - 3 : Red/Black Gauss-Seidel (nonsymmetric: RB pre- and post-relaxation) 268 **/ 269 HYPRE_Int HYPRE_StructPFMGSetRelaxType(HYPRE_StructSolver solver, 270 HYPRE_Int relax_type); 271 272 HYPRE_Int HYPRE_StructPFMGGetRelaxType(HYPRE_StructSolver solver, 273 HYPRE_Int *relax_type); 274 275 /* 276 * (Optional) Set Jacobi weight (this is purposely not documented) 277 */ 278 HYPRE_Int HYPRE_StructPFMGSetJacobiWeight(HYPRE_StructSolver solver, 279 HYPRE_Real weight); 280 HYPRE_Int HYPRE_StructPFMGGetJacobiWeight(HYPRE_StructSolver solver, 281 HYPRE_Real *weight); 282 283 284 /** 285 * (Optional) Set type of coarse-grid operator to use. 286 * 287 * Current operators set by \e rap_type are: 288 * 289 * - 0 : Galerkin (default) 290 * - 1 : non-Galerkin 5-pt or 7-pt stencils 291 * 292 * Both operators are constructed algebraically. The non-Galerkin option 293 * maintains a 5-pt stencil in 2D and a 7-pt stencil in 3D on all grid levels. 294 * The stencil coefficients are computed by averaging techniques. 295 **/ 296 HYPRE_Int HYPRE_StructPFMGSetRAPType(HYPRE_StructSolver solver, 297 HYPRE_Int rap_type); 298 299 HYPRE_Int HYPRE_StructPFMGGetRAPType(HYPRE_StructSolver solver, 300 HYPRE_Int *rap_type ); 301 302 /** 303 * (Optional) Set number of relaxation sweeps before coarse-grid correction. 304 **/ 305 HYPRE_Int HYPRE_StructPFMGSetNumPreRelax(HYPRE_StructSolver solver, 306 HYPRE_Int num_pre_relax); 307 308 HYPRE_Int HYPRE_StructPFMGGetNumPreRelax(HYPRE_StructSolver solver, 309 HYPRE_Int *num_pre_relax); 310 311 /** 312 * (Optional) Set number of relaxation sweeps after coarse-grid correction. 313 **/ 314 HYPRE_Int HYPRE_StructPFMGSetNumPostRelax(HYPRE_StructSolver solver, 315 HYPRE_Int num_post_relax); 316 317 HYPRE_Int HYPRE_StructPFMGGetNumPostRelax(HYPRE_StructSolver solver, 318 HYPRE_Int *num_post_relax); 319 320 /** 321 * (Optional) Skip relaxation on certain grids for isotropic problems. This can 322 * greatly improve efficiency by eliminating unnecessary relaxations when the 323 * underlying problem is isotropic. 324 **/ 325 HYPRE_Int HYPRE_StructPFMGSetSkipRelax(HYPRE_StructSolver solver, 326 HYPRE_Int skip_relax); 327 328 HYPRE_Int HYPRE_StructPFMGGetSkipRelax(HYPRE_StructSolver solver, 329 HYPRE_Int *skip_relax); 330 331 /* 332 * RE-VISIT 333 **/ 334 HYPRE_Int HYPRE_StructPFMGSetDxyz(HYPRE_StructSolver solver, 335 HYPRE_Real *dxyz); 336 337 /** 338 * (Optional) Set the amount of logging to do. 339 **/ 340 HYPRE_Int HYPRE_StructPFMGSetLogging(HYPRE_StructSolver solver, 341 HYPRE_Int logging); 342 343 HYPRE_Int HYPRE_StructPFMGGetLogging(HYPRE_StructSolver solver, 344 HYPRE_Int *logging); 345 346 /** 347 * (Optional) Set the amount of printing to do to the screen. 348 **/ 349 HYPRE_Int HYPRE_StructPFMGSetPrintLevel(HYPRE_StructSolver solver, 350 HYPRE_Int print_level); 351 352 HYPRE_Int HYPRE_StructPFMGGetPrintLevel(HYPRE_StructSolver solver, 353 HYPRE_Int *print_level); 354 355 /** 356 * Return the number of iterations taken. 357 **/ 358 HYPRE_Int HYPRE_StructPFMGGetNumIterations(HYPRE_StructSolver solver, 359 HYPRE_Int *num_iterations); 360 361 /** 362 * Return the norm of the final relative residual. 363 **/ 364 HYPRE_Int HYPRE_StructPFMGGetFinalRelativeResidualNorm(HYPRE_StructSolver solver, 365 HYPRE_Real *norm); 366 367 #if 0 //defined(HYPRE_USING_CUDA) || defined(HYPRE_USING_HIP) 368 HYPRE_Int 369 HYPRE_StructPFMGSetDeviceLevel( HYPRE_StructSolver solver, 370 HYPRE_Int device_level ); 371 #endif 372 /**@}*/ 373 374 /*-------------------------------------------------------------------------- 375 *--------------------------------------------------------------------------*/ 376 377 /** 378 * @name Struct SMG Solver 379 * 380 * SMG is a semicoarsening multigrid solver that uses plane smoothing (in 3D). 381 * The plane smoother calls a 2D SMG algorithm with line smoothing, and the line 382 * smoother is cyclic reduction (1D SMG). For periodic problems, the grid size 383 * in periodic dimensions currently must be a power-of-two. 384 * 385 * @{ 386 **/ 387 388 /** 389 * Create a solver object. 390 **/ 391 HYPRE_Int HYPRE_StructSMGCreate(MPI_Comm comm, 392 HYPRE_StructSolver *solver); 393 394 /** 395 * Destroy a solver object. 396 **/ 397 HYPRE_Int HYPRE_StructSMGDestroy(HYPRE_StructSolver solver); 398 399 /** 400 * Prepare to solve the system. The coefficient data in \e b and \e x is 401 * ignored here, but information about the layout of the data may be used. 402 **/ 403 HYPRE_Int HYPRE_StructSMGSetup(HYPRE_StructSolver solver, 404 HYPRE_StructMatrix A, 405 HYPRE_StructVector b, 406 HYPRE_StructVector x); 407 408 /** 409 * Solve the system. 410 **/ 411 HYPRE_Int HYPRE_StructSMGSolve(HYPRE_StructSolver solver, 412 HYPRE_StructMatrix A, 413 HYPRE_StructVector b, 414 HYPRE_StructVector x); 415 416 /* 417 * RE-VISIT 418 **/ 419 HYPRE_Int HYPRE_StructSMGSetMemoryUse(HYPRE_StructSolver solver, 420 HYPRE_Int memory_use); 421 422 HYPRE_Int HYPRE_StructSMGGetMemoryUse(HYPRE_StructSolver solver, 423 HYPRE_Int *memory_use); 424 425 /** 426 * (Optional) Set the convergence tolerance. 427 **/ 428 HYPRE_Int HYPRE_StructSMGSetTol(HYPRE_StructSolver solver, 429 HYPRE_Real tol); 430 431 HYPRE_Int HYPRE_StructSMGGetTol(HYPRE_StructSolver solver, 432 HYPRE_Real *tol); 433 434 /** 435 * (Optional) Set maximum number of iterations. 436 **/ 437 HYPRE_Int HYPRE_StructSMGSetMaxIter(HYPRE_StructSolver solver, 438 HYPRE_Int max_iter); 439 440 HYPRE_Int HYPRE_StructSMGGetMaxIter(HYPRE_StructSolver solver, 441 HYPRE_Int *max_iter); 442 443 /** 444 * (Optional) Additionally require that the relative difference in 445 * successive iterates be small. 446 **/ 447 HYPRE_Int HYPRE_StructSMGSetRelChange(HYPRE_StructSolver solver, 448 HYPRE_Int rel_change); 449 450 HYPRE_Int HYPRE_StructSMGGetRelChange(HYPRE_StructSolver solver, 451 HYPRE_Int *rel_change); 452 453 /** 454 * (Optional) Use a zero initial guess. This allows the solver to cut corners 455 * in the case where a zero initial guess is needed (e.g., for preconditioning) 456 * to reduce compuational cost. 457 **/ 458 HYPRE_Int HYPRE_StructSMGSetZeroGuess(HYPRE_StructSolver solver); 459 460 HYPRE_Int HYPRE_StructSMGGetZeroGuess(HYPRE_StructSolver solver, 461 HYPRE_Int *zeroguess); 462 463 /** 464 * (Optional) Use a nonzero initial guess. This is the default behavior, but 465 * this routine allows the user to switch back after using \e SetZeroGuess. 466 **/ 467 HYPRE_Int HYPRE_StructSMGSetNonZeroGuess(HYPRE_StructSolver solver); 468 469 /** 470 * (Optional) Set number of relaxation sweeps before coarse-grid correction. 471 **/ 472 HYPRE_Int HYPRE_StructSMGSetNumPreRelax(HYPRE_StructSolver solver, 473 HYPRE_Int num_pre_relax); 474 475 HYPRE_Int HYPRE_StructSMGGetNumPreRelax(HYPRE_StructSolver solver, 476 HYPRE_Int *num_pre_relax); 477 478 /** 479 * (Optional) Set number of relaxation sweeps after coarse-grid correction. 480 **/ 481 HYPRE_Int HYPRE_StructSMGSetNumPostRelax(HYPRE_StructSolver solver, 482 HYPRE_Int num_post_relax); 483 484 HYPRE_Int HYPRE_StructSMGGetNumPostRelax(HYPRE_StructSolver solver, 485 HYPRE_Int *num_post_relax); 486 487 /** 488 * (Optional) Set the amount of logging to do. 489 **/ 490 HYPRE_Int HYPRE_StructSMGSetLogging(HYPRE_StructSolver solver, 491 HYPRE_Int logging); 492 493 HYPRE_Int HYPRE_StructSMGGetLogging(HYPRE_StructSolver solver, 494 HYPRE_Int *logging); 495 496 /** 497 * (Optional) Set the amount of printing to do to the screen. 498 **/ 499 HYPRE_Int HYPRE_StructSMGSetPrintLevel(HYPRE_StructSolver solver, 500 HYPRE_Int print_level); 501 502 HYPRE_Int HYPRE_StructSMGGetPrintLevel(HYPRE_StructSolver solver, 503 HYPRE_Int *print_level); 504 505 /** 506 * Return the number of iterations taken. 507 **/ 508 HYPRE_Int HYPRE_StructSMGGetNumIterations(HYPRE_StructSolver solver, 509 HYPRE_Int *num_iterations); 510 511 /** 512 * Return the norm of the final relative residual. 513 **/ 514 HYPRE_Int HYPRE_StructSMGGetFinalRelativeResidualNorm(HYPRE_StructSolver solver, 515 HYPRE_Real *norm); 516 517 #if 0 //defined(HYPRE_USING_CUDA) || defined(HYPRE_USING_HIP) 518 HYPRE_Int 519 HYPRE_StructSMGSetDeviceLevel( HYPRE_StructSolver solver, 520 HYPRE_Int device_level ); 521 #endif 522 523 /**@}*/ 524 525 /*-------------------------------------------------------------------------- 526 *--------------------------------------------------------------------------*/ 527 528 /** 529 * @name Struct CycRed Solver 530 * 531 * CycRed is a cyclic reduction solver that simultaneously solves a collection 532 * of 1D tridiagonal systems embedded in a d-dimensional grid. 533 * 534 * @{ 535 **/ 536 537 /** 538 * Create a solver object. 539 **/ 540 HYPRE_Int HYPRE_StructCycRedCreate(MPI_Comm comm, 541 HYPRE_StructSolver *solver); 542 543 /** 544 * Destroy a solver object. 545 **/ 546 HYPRE_Int HYPRE_StructCycRedDestroy(HYPRE_StructSolver solver); 547 548 /** 549 * Prepare to solve the system. The coefficient data in \e b and \e x is 550 * ignored here, but information about the layout of the data may be used. 551 **/ 552 HYPRE_Int HYPRE_StructCycRedSetup(HYPRE_StructSolver solver, 553 HYPRE_StructMatrix A, 554 HYPRE_StructVector b, 555 HYPRE_StructVector x); 556 557 /** 558 * Solve the system. 559 **/ 560 HYPRE_Int HYPRE_StructCycRedSolve(HYPRE_StructSolver solver, 561 HYPRE_StructMatrix A, 562 HYPRE_StructVector b, 563 HYPRE_StructVector x); 564 565 /** 566 * 567 * (Optional) Set the dimension number for the embedded 1D tridiagonal systems. 568 * The default is \e tdim = 0. 569 **/ 570 HYPRE_Int HYPRE_StructCycRedSetTDim(HYPRE_StructSolver solver, 571 HYPRE_Int tdim); 572 573 /** 574 * (Optional) Set the base index and stride for the embedded 1D systems. The 575 * stride must be equal one in the dimension corresponding to the 1D systems 576 * (see \ref HYPRE_StructCycRedSetTDim). 577 **/ 578 HYPRE_Int HYPRE_StructCycRedSetBase(HYPRE_StructSolver solver, 579 HYPRE_Int ndim, 580 HYPRE_Int *base_index, 581 HYPRE_Int *base_stride); 582 583 /**@}*/ 584 585 /*-------------------------------------------------------------------------- 586 *--------------------------------------------------------------------------*/ 587 588 /** 589 * @name Struct PCG Solver 590 * 591 * These routines should be used in conjunction with the generic interface in 592 * \ref KrylovSolvers. 593 * 594 * @{ 595 **/ 596 597 /** 598 * Create a solver object. 599 **/ 600 HYPRE_Int HYPRE_StructPCGCreate(MPI_Comm comm, 601 HYPRE_StructSolver *solver); 602 603 /** 604 * Destroy a solver object. 605 **/ 606 HYPRE_Int HYPRE_StructPCGDestroy(HYPRE_StructSolver solver); 607 608 HYPRE_Int HYPRE_StructPCGSetup(HYPRE_StructSolver solver, 609 HYPRE_StructMatrix A, 610 HYPRE_StructVector b, 611 HYPRE_StructVector x); 612 613 HYPRE_Int HYPRE_StructPCGSolve(HYPRE_StructSolver solver, 614 HYPRE_StructMatrix A, 615 HYPRE_StructVector b, 616 HYPRE_StructVector x); 617 618 HYPRE_Int HYPRE_StructPCGSetTol(HYPRE_StructSolver solver, 619 HYPRE_Real tol); 620 621 HYPRE_Int HYPRE_StructPCGSetAbsoluteTol(HYPRE_StructSolver solver, 622 HYPRE_Real tol); 623 624 HYPRE_Int HYPRE_StructPCGSetMaxIter(HYPRE_StructSolver solver, 625 HYPRE_Int max_iter); 626 627 HYPRE_Int HYPRE_StructPCGSetTwoNorm(HYPRE_StructSolver solver, 628 HYPRE_Int two_norm); 629 630 HYPRE_Int HYPRE_StructPCGSetRelChange(HYPRE_StructSolver solver, 631 HYPRE_Int rel_change); 632 633 HYPRE_Int HYPRE_StructPCGSetPrecond(HYPRE_StructSolver solver, 634 HYPRE_PtrToStructSolverFcn precond, 635 HYPRE_PtrToStructSolverFcn precond_setup, 636 HYPRE_StructSolver precond_solver); 637 638 HYPRE_Int HYPRE_StructPCGSetLogging(HYPRE_StructSolver solver, 639 HYPRE_Int logging); 640 641 HYPRE_Int HYPRE_StructPCGSetPrintLevel(HYPRE_StructSolver solver, 642 HYPRE_Int level); 643 644 HYPRE_Int HYPRE_StructPCGGetNumIterations(HYPRE_StructSolver solver, 645 HYPRE_Int *num_iterations); 646 647 HYPRE_Int HYPRE_StructPCGGetFinalRelativeResidualNorm(HYPRE_StructSolver solver, 648 HYPRE_Real *norm); 649 650 HYPRE_Int HYPRE_StructPCGGetResidual(HYPRE_StructSolver solver, 651 void **residual); 652 653 /** 654 * Setup routine for diagonal preconditioning. 655 **/ 656 HYPRE_Int HYPRE_StructDiagScaleSetup(HYPRE_StructSolver solver, 657 HYPRE_StructMatrix A, 658 HYPRE_StructVector y, 659 HYPRE_StructVector x); 660 661 /** 662 * Solve routine for diagonal preconditioning. 663 **/ 664 HYPRE_Int HYPRE_StructDiagScale(HYPRE_StructSolver solver, 665 HYPRE_StructMatrix HA, 666 HYPRE_StructVector Hy, 667 HYPRE_StructVector Hx); 668 669 /**@}*/ 670 671 /*-------------------------------------------------------------------------- 672 *--------------------------------------------------------------------------*/ 673 674 /** 675 * @name Struct GMRES Solver 676 * 677 * These routines should be used in conjunction with the generic interface in 678 * \ref KrylovSolvers. 679 * 680 * @{ 681 **/ 682 683 /** 684 * Create a solver object. 685 **/ 686 HYPRE_Int HYPRE_StructGMRESCreate(MPI_Comm comm, 687 HYPRE_StructSolver *solver); 688 689 690 /** 691 * Destroy a solver object. 692 **/ 693 HYPRE_Int HYPRE_StructGMRESDestroy(HYPRE_StructSolver solver); 694 695 HYPRE_Int HYPRE_StructGMRESSetup(HYPRE_StructSolver solver, 696 HYPRE_StructMatrix A, 697 HYPRE_StructVector b, 698 HYPRE_StructVector x); 699 700 HYPRE_Int HYPRE_StructGMRESSolve(HYPRE_StructSolver solver, 701 HYPRE_StructMatrix A, 702 HYPRE_StructVector b, 703 HYPRE_StructVector x); 704 705 HYPRE_Int HYPRE_StructGMRESSetTol(HYPRE_StructSolver solver, 706 HYPRE_Real tol); 707 708 HYPRE_Int HYPRE_StructGMRESSetAbsoluteTol(HYPRE_StructSolver solver, 709 HYPRE_Real tol); 710 711 HYPRE_Int HYPRE_StructGMRESSetMaxIter(HYPRE_StructSolver solver, 712 HYPRE_Int max_iter); 713 714 HYPRE_Int HYPRE_StructGMRESSetKDim(HYPRE_StructSolver solver, 715 HYPRE_Int k_dim); 716 717 HYPRE_Int HYPRE_StructGMRESSetPrecond(HYPRE_StructSolver solver, 718 HYPRE_PtrToStructSolverFcn precond, 719 HYPRE_PtrToStructSolverFcn precond_setup, 720 HYPRE_StructSolver precond_solver); 721 722 HYPRE_Int HYPRE_StructGMRESSetLogging(HYPRE_StructSolver solver, 723 HYPRE_Int logging); 724 725 HYPRE_Int HYPRE_StructGMRESSetPrintLevel(HYPRE_StructSolver solver, 726 HYPRE_Int level); 727 728 HYPRE_Int HYPRE_StructGMRESGetNumIterations(HYPRE_StructSolver solver, 729 HYPRE_Int *num_iterations); 730 731 HYPRE_Int HYPRE_StructGMRESGetFinalRelativeResidualNorm(HYPRE_StructSolver solver, 732 HYPRE_Real *norm); 733 734 HYPRE_Int HYPRE_StructGMRESGetResidual(HYPRE_StructSolver solver, 735 void **residual); 736 /**@}*/ 737 738 /*-------------------------------------------------------------------------- 739 *--------------------------------------------------------------------------*/ 740 741 /** 742 * @name Struct FlexGMRES Solver 743 * 744 * These routines should be used in conjunction with the generic interface in 745 * \ref KrylovSolvers. 746 * 747 * @{ 748 **/ 749 750 /** 751 * Create a solver object. 752 **/ 753 HYPRE_Int HYPRE_StructFlexGMRESCreate(MPI_Comm comm, 754 HYPRE_StructSolver *solver); 755 756 /** 757 * Destroy a solver object. 758 **/ 759 HYPRE_Int HYPRE_StructFlexGMRESDestroy(HYPRE_StructSolver solver); 760 761 HYPRE_Int HYPRE_StructFlexGMRESSetup(HYPRE_StructSolver solver, 762 HYPRE_StructMatrix A, 763 HYPRE_StructVector b, 764 HYPRE_StructVector x); 765 766 HYPRE_Int HYPRE_StructFlexGMRESSolve(HYPRE_StructSolver solver, 767 HYPRE_StructMatrix A, 768 HYPRE_StructVector b, 769 HYPRE_StructVector x); 770 771 HYPRE_Int HYPRE_StructFlexGMRESSetTol(HYPRE_StructSolver solver, 772 HYPRE_Real tol); 773 774 HYPRE_Int HYPRE_StructFlexGMRESSetAbsoluteTol(HYPRE_StructSolver solver, 775 HYPRE_Real tol); 776 777 HYPRE_Int HYPRE_StructFlexGMRESSetMaxIter(HYPRE_StructSolver solver, 778 HYPRE_Int max_iter); 779 780 HYPRE_Int HYPRE_StructFlexGMRESSetKDim(HYPRE_StructSolver solver, 781 HYPRE_Int k_dim); 782 783 HYPRE_Int HYPRE_StructFlexGMRESSetPrecond(HYPRE_StructSolver solver, 784 HYPRE_PtrToStructSolverFcn precond, 785 HYPRE_PtrToStructSolverFcn precond_setup, 786 HYPRE_StructSolver precond_solver); 787 788 HYPRE_Int HYPRE_StructFlexGMRESSetLogging(HYPRE_StructSolver solver, 789 HYPRE_Int logging); 790 791 HYPRE_Int HYPRE_StructFlexGMRESSetPrintLevel(HYPRE_StructSolver solver, 792 HYPRE_Int level); 793 794 HYPRE_Int HYPRE_StructFlexGMRESGetNumIterations(HYPRE_StructSolver solver, 795 HYPRE_Int *num_iterations); 796 797 HYPRE_Int HYPRE_StructFlexGMRESGetFinalRelativeResidualNorm(HYPRE_StructSolver solver, 798 HYPRE_Real *norm); 799 800 HYPRE_Int HYPRE_StructFlexGMRESGetResidual(HYPRE_StructSolver solver, 801 void **residual); 802 803 HYPRE_Int HYPRE_StructFlexGMRESSetModifyPC(HYPRE_StructSolver solver, 804 HYPRE_PtrToModifyPCFcn modify_pc); 805 806 /**@}*/ 807 808 /*-------------------------------------------------------------------------- 809 *--------------------------------------------------------------------------*/ 810 811 /** 812 * @name Struct LGMRES Solver 813 * 814 * These routines should be used in conjunction with the generic interface in 815 * \ref KrylovSolvers. 816 * 817 * @{ 818 **/ 819 820 /** 821 * Create a solver object. 822 **/ 823 HYPRE_Int HYPRE_StructLGMRESCreate(MPI_Comm comm, 824 HYPRE_StructSolver *solver); 825 826 /** 827 * Destroy a solver object. 828 **/ 829 HYPRE_Int HYPRE_StructLGMRESDestroy(HYPRE_StructSolver solver); 830 831 HYPRE_Int HYPRE_StructLGMRESSetup(HYPRE_StructSolver solver, 832 HYPRE_StructMatrix A, 833 HYPRE_StructVector b, 834 HYPRE_StructVector x); 835 836 HYPRE_Int HYPRE_StructLGMRESSolve(HYPRE_StructSolver solver, 837 HYPRE_StructMatrix A, 838 HYPRE_StructVector b, 839 HYPRE_StructVector x); 840 841 HYPRE_Int HYPRE_StructLGMRESSetTol(HYPRE_StructSolver solver, 842 HYPRE_Real tol); 843 844 HYPRE_Int HYPRE_StructLGMRESSetAbsoluteTol(HYPRE_StructSolver solver, 845 HYPRE_Real tol); 846 847 HYPRE_Int HYPRE_StructLGMRESSetMaxIter(HYPRE_StructSolver solver, 848 HYPRE_Int max_iter); 849 850 HYPRE_Int HYPRE_StructLGMRESSetKDim(HYPRE_StructSolver solver, 851 HYPRE_Int k_dim); 852 853 HYPRE_Int HYPRE_StructLGMRESSetAugDim(HYPRE_StructSolver solver, 854 HYPRE_Int aug_dim); 855 856 HYPRE_Int HYPRE_StructLGMRESSetPrecond(HYPRE_StructSolver solver, 857 HYPRE_PtrToStructSolverFcn precond, 858 HYPRE_PtrToStructSolverFcn precond_setup, 859 HYPRE_StructSolver precond_solver); 860 861 HYPRE_Int HYPRE_StructLGMRESSetLogging(HYPRE_StructSolver solver, 862 HYPRE_Int logging); 863 864 HYPRE_Int HYPRE_StructLGMRESSetPrintLevel(HYPRE_StructSolver solver, 865 HYPRE_Int level); 866 867 HYPRE_Int HYPRE_StructLGMRESGetNumIterations(HYPRE_StructSolver solver, 868 HYPRE_Int *num_iterations); 869 870 HYPRE_Int HYPRE_StructLGMRESGetFinalRelativeResidualNorm(HYPRE_StructSolver solver, 871 HYPRE_Real *norm); 872 873 HYPRE_Int HYPRE_StructLGMRESGetResidual(HYPRE_StructSolver solver, 874 void **residual); 875 /**@}*/ 876 877 /*-------------------------------------------------------------------------- 878 *--------------------------------------------------------------------------*/ 879 880 /** 881 * @name Struct BiCGSTAB Solver 882 * 883 * These routines should be used in conjunction with the generic interface in 884 * \ref KrylovSolvers. 885 * 886 * @{ 887 **/ 888 889 /** 890 * Create a solver object. 891 **/ 892 HYPRE_Int HYPRE_StructBiCGSTABCreate(MPI_Comm comm, 893 HYPRE_StructSolver *solver); 894 895 /** 896 * Destroy a solver object. 897 **/ 898 HYPRE_Int HYPRE_StructBiCGSTABDestroy(HYPRE_StructSolver solver); 899 900 HYPRE_Int HYPRE_StructBiCGSTABSetup(HYPRE_StructSolver solver, 901 HYPRE_StructMatrix A, 902 HYPRE_StructVector b, 903 HYPRE_StructVector x); 904 905 HYPRE_Int HYPRE_StructBiCGSTABSolve(HYPRE_StructSolver solver, 906 HYPRE_StructMatrix A, 907 HYPRE_StructVector b, 908 HYPRE_StructVector x); 909 910 HYPRE_Int HYPRE_StructBiCGSTABSetTol(HYPRE_StructSolver solver, 911 HYPRE_Real tol); 912 913 HYPRE_Int HYPRE_StructBiCGSTABSetAbsoluteTol(HYPRE_StructSolver solver, 914 HYPRE_Real tol); 915 916 HYPRE_Int HYPRE_StructBiCGSTABSetMaxIter(HYPRE_StructSolver solver, 917 HYPRE_Int max_iter); 918 919 HYPRE_Int HYPRE_StructBiCGSTABSetPrecond(HYPRE_StructSolver solver, 920 HYPRE_PtrToStructSolverFcn precond, 921 HYPRE_PtrToStructSolverFcn precond_setup, 922 HYPRE_StructSolver precond_solver); 923 924 HYPRE_Int HYPRE_StructBiCGSTABSetLogging(HYPRE_StructSolver solver, 925 HYPRE_Int logging); 926 927 HYPRE_Int HYPRE_StructBiCGSTABSetPrintLevel(HYPRE_StructSolver solver, 928 HYPRE_Int level); 929 930 HYPRE_Int HYPRE_StructBiCGSTABGetNumIterations(HYPRE_StructSolver solver, 931 HYPRE_Int *num_iterations); 932 933 HYPRE_Int HYPRE_StructBiCGSTABGetFinalRelativeResidualNorm(HYPRE_StructSolver solver, 934 HYPRE_Real *norm); 935 936 HYPRE_Int HYPRE_StructBiCGSTABGetResidual( HYPRE_StructSolver solver, 937 void **residual); 938 /**@}*/ 939 940 /*-------------------------------------------------------------------------- 941 *--------------------------------------------------------------------------*/ 942 943 /** 944 * @name Struct Hybrid Solver 945 * 946 * @{ 947 **/ 948 949 /** 950 * Create a solver object. 951 **/ 952 HYPRE_Int HYPRE_StructHybridCreate(MPI_Comm comm, 953 HYPRE_StructSolver *solver); 954 955 /** 956 * Destroy a solver object. 957 **/ 958 HYPRE_Int HYPRE_StructHybridDestroy(HYPRE_StructSolver solver); 959 960 /** 961 * Prepare to solve the system. The coefficient data in \e b and \e x is 962 * ignored here, but information about the layout of the data may be used. 963 **/ 964 HYPRE_Int HYPRE_StructHybridSetup(HYPRE_StructSolver solver, 965 HYPRE_StructMatrix A, 966 HYPRE_StructVector b, 967 HYPRE_StructVector x); 968 969 /** 970 * Solve the system. 971 **/ 972 HYPRE_Int HYPRE_StructHybridSolve(HYPRE_StructSolver solver, 973 HYPRE_StructMatrix A, 974 HYPRE_StructVector b, 975 HYPRE_StructVector x); 976 977 /** 978 * (Optional) Set the convergence tolerance. 979 **/ 980 HYPRE_Int HYPRE_StructHybridSetTol(HYPRE_StructSolver solver, 981 HYPRE_Real tol); 982 983 /** 984 * (Optional) Set an accepted convergence tolerance for diagonal scaling (DS). 985 * The solver will switch preconditioners if the convergence of DS is slower 986 * than \e cf_tol. 987 **/ 988 HYPRE_Int HYPRE_StructHybridSetConvergenceTol(HYPRE_StructSolver solver, 989 HYPRE_Real cf_tol); 990 991 /** 992 * (Optional) Set maximum number of iterations for diagonal scaling (DS). The 993 * solver will switch preconditioners if DS reaches \e ds_max_its. 994 **/ 995 HYPRE_Int HYPRE_StructHybridSetDSCGMaxIter(HYPRE_StructSolver solver, 996 HYPRE_Int ds_max_its); 997 998 /** 999 * (Optional) Set maximum number of iterations for general preconditioner (PRE). 1000 * The solver will stop if PRE reaches \e pre_max_its. 1001 **/ 1002 HYPRE_Int HYPRE_StructHybridSetPCGMaxIter(HYPRE_StructSolver solver, 1003 HYPRE_Int pre_max_its); 1004 1005 /** 1006 * (Optional) Use the two-norm in stopping criteria. 1007 **/ 1008 HYPRE_Int HYPRE_StructHybridSetTwoNorm(HYPRE_StructSolver solver, 1009 HYPRE_Int two_norm); 1010 1011 HYPRE_Int HYPRE_StructHybridSetStopCrit(HYPRE_StructSolver solver, 1012 HYPRE_Int stop_crit); 1013 1014 /** 1015 * (Optional) Additionally require that the relative difference in 1016 * successive iterates be small. 1017 **/ 1018 HYPRE_Int HYPRE_StructHybridSetRelChange(HYPRE_StructSolver solver, 1019 HYPRE_Int rel_change); 1020 1021 /** 1022 * (Optional) Set the type of Krylov solver to use. 1023 * 1024 * Current krylov methods set by \e solver_type are: 1025 * 1026 * - 0 : PCG (default) 1027 * - 1 : GMRES 1028 * - 2 : BiCGSTAB 1029 **/ 1030 HYPRE_Int HYPRE_StructHybridSetSolverType(HYPRE_StructSolver solver, 1031 HYPRE_Int solver_type); 1032 1033 /** 1034 * (Optional) Set recompute residual (don't rely on 3-term recurrence). 1035 **/ 1036 HYPRE_Int 1037 HYPRE_StructHybridSetRecomputeResidual( HYPRE_StructSolver solver, 1038 HYPRE_Int recompute_residual ); 1039 1040 /** 1041 * (Optional) Get recompute residual option. 1042 **/ 1043 HYPRE_Int 1044 HYPRE_StructHybridGetRecomputeResidual( HYPRE_StructSolver solver, 1045 HYPRE_Int *recompute_residual ); 1046 1047 /** 1048 * (Optional) Set recompute residual period (don't rely on 3-term recurrence). 1049 * 1050 * Recomputes residual after every specified number of iterations. 1051 **/ 1052 HYPRE_Int 1053 HYPRE_StructHybridSetRecomputeResidualP( HYPRE_StructSolver solver, 1054 HYPRE_Int recompute_residual_p ); 1055 1056 /** 1057 * (Optional) Get recompute residual period option. 1058 **/ 1059 HYPRE_Int 1060 HYPRE_StructHybridGetRecomputeResidualP( HYPRE_StructSolver solver, 1061 HYPRE_Int *recompute_residual_p ); 1062 1063 /** 1064 * (Optional) Set the maximum size of the Krylov space when using GMRES. 1065 **/ 1066 HYPRE_Int HYPRE_StructHybridSetKDim(HYPRE_StructSolver solver, 1067 HYPRE_Int k_dim); 1068 1069 /** 1070 * (Optional) Set the preconditioner to use. 1071 **/ 1072 HYPRE_Int HYPRE_StructHybridSetPrecond(HYPRE_StructSolver solver, 1073 HYPRE_PtrToStructSolverFcn precond, 1074 HYPRE_PtrToStructSolverFcn precond_setup, 1075 HYPRE_StructSolver precond_solver); 1076 1077 /** 1078 * (Optional) Set the amount of logging to do. 1079 **/ 1080 HYPRE_Int HYPRE_StructHybridSetLogging(HYPRE_StructSolver solver, 1081 HYPRE_Int logging); 1082 1083 /** 1084 * (Optional) Set the amount of printing to do to the screen. 1085 **/ 1086 HYPRE_Int HYPRE_StructHybridSetPrintLevel(HYPRE_StructSolver solver, 1087 HYPRE_Int print_level); 1088 1089 /** 1090 * Return the number of iterations taken. 1091 **/ 1092 HYPRE_Int HYPRE_StructHybridGetNumIterations(HYPRE_StructSolver solver, 1093 HYPRE_Int *num_its); 1094 1095 /** 1096 * Return the number of diagonal scaling iterations taken. 1097 **/ 1098 HYPRE_Int HYPRE_StructHybridGetDSCGNumIterations(HYPRE_StructSolver solver, 1099 HYPRE_Int *ds_num_its); 1100 1101 /** 1102 * Return the number of general preconditioning iterations taken. 1103 **/ 1104 HYPRE_Int HYPRE_StructHybridGetPCGNumIterations(HYPRE_StructSolver solver, 1105 HYPRE_Int *pre_num_its); 1106 1107 /** 1108 * Return the norm of the final relative residual. 1109 **/ 1110 HYPRE_Int HYPRE_StructHybridGetFinalRelativeResidualNorm(HYPRE_StructSolver solver, 1111 HYPRE_Real *norm); 1112 1113 HYPRE_Int HYPRE_StructHybridSetPCGAbsoluteTolFactor(HYPRE_StructSolver solver, 1114 HYPRE_Real pcg_atolf ); 1115 1116 /**@}*/ 1117 1118 /*-------------------------------------------------------------------------- 1119 *--------------------------------------------------------------------------*/ 1120 1121 /* 1122 * @name Struct SparseMSG Solver 1123 **/ 1124 1125 HYPRE_Int HYPRE_StructSparseMSGCreate(MPI_Comm comm, 1126 HYPRE_StructSolver *solver); 1127 1128 HYPRE_Int HYPRE_StructSparseMSGDestroy(HYPRE_StructSolver solver); 1129 1130 HYPRE_Int HYPRE_StructSparseMSGSetup(HYPRE_StructSolver solver, 1131 HYPRE_StructMatrix A, 1132 HYPRE_StructVector b, 1133 HYPRE_StructVector x); 1134 1135 HYPRE_Int HYPRE_StructSparseMSGSolve(HYPRE_StructSolver solver, 1136 HYPRE_StructMatrix A, 1137 HYPRE_StructVector b, 1138 HYPRE_StructVector x); 1139 1140 HYPRE_Int HYPRE_StructSparseMSGSetTol(HYPRE_StructSolver solver, 1141 HYPRE_Real tol); 1142 1143 HYPRE_Int HYPRE_StructSparseMSGSetMaxIter(HYPRE_StructSolver solver, 1144 HYPRE_Int max_iter); 1145 1146 HYPRE_Int HYPRE_StructSparseMSGSetJump(HYPRE_StructSolver solver, 1147 HYPRE_Int jump); 1148 1149 HYPRE_Int HYPRE_StructSparseMSGSetRelChange(HYPRE_StructSolver solver, 1150 HYPRE_Int rel_change); 1151 1152 HYPRE_Int HYPRE_StructSparseMSGSetZeroGuess(HYPRE_StructSolver solver); 1153 1154 HYPRE_Int HYPRE_StructSparseMSGSetNonZeroGuess(HYPRE_StructSolver solver); 1155 1156 HYPRE_Int HYPRE_StructSparseMSGSetRelaxType(HYPRE_StructSolver solver, 1157 HYPRE_Int relax_type); 1158 1159 HYPRE_Int HYPRE_StructSparseMSGSetJacobiWeight(HYPRE_StructSolver solver, 1160 HYPRE_Real weight); 1161 1162 HYPRE_Int HYPRE_StructSparseMSGSetNumPreRelax(HYPRE_StructSolver solver, 1163 HYPRE_Int num_pre_relax); 1164 1165 HYPRE_Int HYPRE_StructSparseMSGSetNumPostRelax(HYPRE_StructSolver solver, 1166 HYPRE_Int num_post_relax); 1167 1168 HYPRE_Int HYPRE_StructSparseMSGSetNumFineRelax(HYPRE_StructSolver solver, 1169 HYPRE_Int num_fine_relax); 1170 1171 HYPRE_Int HYPRE_StructSparseMSGSetLogging(HYPRE_StructSolver solver, 1172 HYPRE_Int logging); 1173 1174 HYPRE_Int HYPRE_StructSparseMSGSetPrintLevel(HYPRE_StructSolver solver, 1175 HYPRE_Int print_level); 1176 1177 1178 HYPRE_Int HYPRE_StructSparseMSGGetNumIterations(HYPRE_StructSolver solver, 1179 HYPRE_Int *num_iterations); 1180 1181 HYPRE_Int HYPRE_StructSparseMSGGetFinalRelativeResidualNorm(HYPRE_StructSolver solver, 1182 HYPRE_Real *norm); 1183 1184 /*-------------------------------------------------------------------------- 1185 *--------------------------------------------------------------------------*/ 1186 1187 /** 1188 * @name Struct LOBPCG Eigensolver 1189 * 1190 * These routines should be used in conjunction with the generic interface in 1191 * \ref Eigensolvers. 1192 * 1193 * @{ 1194 **/ 1195 1196 /** 1197 * Load interface interpreter. Vector part loaded with hypre_StructKrylov 1198 * functions and multivector part loaded with mv_TempMultiVector functions. 1199 **/ 1200 HYPRE_Int 1201 HYPRE_StructSetupInterpreter(mv_InterfaceInterpreter *i); 1202 1203 /** 1204 * Load Matvec interpreter with hypre_StructKrylov functions. 1205 **/ 1206 HYPRE_Int 1207 HYPRE_StructSetupMatvec(HYPRE_MatvecFunctions *mv); 1208 1209 /**@}*/ 1210 1211 /*-------------------------------------------------------------------------- 1212 *--------------------------------------------------------------------------*/ 1213 1214 /**@}*/ 1215 1216 #ifdef __cplusplus 1217 } 1218 #endif 1219 1220 #endif 1221 1222