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 // ************************************************************************* 9 // This is the HYPRE implementation of LinearSystemCore. 10 // ************************************************************************* 11 12 #ifndef _HYPRE_LinSysCore_h_ 13 #define _HYPRE_LinSysCore_h_ 14 15 #define HYPRE_FEI_Version() "FEI/HYPRE 2.7.0R1" 16 17 // ************************************************************************* 18 // system libraries used 19 // ------------------------------------------------------------------------- 20 21 #include <stdio.h> 22 #include <stdlib.h> 23 #include <string.h> 24 #include <math.h> 25 26 #ifdef NOFEI 27 #undef NOFEI 28 #endif 29 30 // ************************************************************************* 31 // FEI-specific include files 32 // ------------------------------------------------------------------------- 33 34 #include "HYPRE_FEI_includes.h" 35 36 // ************************************************************************* 37 // local enumerations and defines 38 // ------------------------------------------------------------------------- 39 40 enum HYsolverID {HYPCG,HYLSICG,HYGMRES,HYFGMRES,HYCGSTAB,HYCGSTABL,HYTFQMR, 41 HYBICGS,HYSYMQMR,HYAMG,HYSUPERLU,HYSUPERLUX,HYDSUPERLU, 42 HYY12M,HYAMGE,HYHYBRID}; 43 enum HYpreconID {HYIDENTITY,HYDIAGONAL,HYPILUT,HYPARASAILS,HYBOOMERAMG,HYML, 44 HYDDILUT,HYPOLY,HYDDICT,HYSCHWARZ,HYEUCLID,HYBLOCK,HYMLI, 45 HYUZAWA,HYMLMAXWELL,HYAMS,HYSYSPDE,HYDSLU}; 46 47 #define HYFEI_HIGHMASK (2147483647-255) 48 #define HYFEI_SPECIALMASK 255 49 #define HYFEI_SLIDEREDUCE1 256 50 #define HYFEI_SLIDEREDUCE2 512 51 #define HYFEI_SLIDEREDUCE3 1024 52 #define HYFEI_PRINTMAT 2048 53 #define HYFEI_PRINTREDMAT 4096 54 #define HYFEI_PRINTSOL 8192 55 #define HYFEI_DDILUT 16384 56 #define HYFEI_SCHURREDUCE1 32768 57 #define HYFEI_SCHURREDUCE2 65536 58 #define HYFEI_SCHURREDUCE3 131072 59 #define HYFEI_PRINTFEINFO 262144 60 #define HYFEI_AMGDEBUG 524288 61 #define HYFEI_STOPAFTERPRINT 1048576 62 #define HYFEI_PRINTPARCSRMAT 2097152 63 #define HYFEI_IMPOSENOBC 4194304 64 65 // ************************************************************************* 66 // substructure definition 67 // ------------------------------------------------------------------------- 68 69 typedef struct 70 { 71 HYPRE_BigInt *EdgeNodeList_; 72 HYPRE_BigInt *NodeNumbers_; 73 HYPRE_Int numEdges_; 74 HYPRE_Int numLocalNodes_; 75 HYPRE_Int numNodes_; 76 HYPRE_Real *NodalCoord_; 77 } HYPRE_FEI_AMSData; 78 79 // ************************************************************************* 80 // class definition 81 // ------------------------------------------------------------------------- 82 83 class HYPRE_LinSysCore 84 #ifndef NOFEI 85 : public LinearSystemCore 86 #endif 87 { 88 public: 89 HYPRE_LinSysCore(MPI_Comm comm); 90 virtual ~HYPRE_LinSysCore(); 91 92 // ---------------------------------------------------------------------- 93 // for creating another one, w/o knowing the run-time type of 'this' one. 94 // ---------------------------------------------------------------------- 95 96 #ifndef NOFEI 97 LinearSystemCore* clone(); 98 #endif 99 100 // ---------------------------------------------------------------------- 101 // parameters : for setting generic argc/argv style parameters. 102 // ---------------------------------------------------------------------- 103 104 int parameters(int numParams, char** params); 105 106 // ---------------------------------------------------------------------- 107 // new functions in FEI 1.5 and above (not implemented here) 108 // ---------------------------------------------------------------------- 109 110 int setLookup(Lookup& lookup); 111 112 int setConnectivities(GlobalID elemBlock, int numElements, 113 int numNodesPerElem, const GlobalID* elemIDs, 114 const int* const* connNodes) ; 115 116 int setStiffnessMatrices(GlobalID elemBlock, int numElems, 117 const GlobalID* elemIDs, const double *const *const *stiff, 118 int numEqnsPerElem, const int *const * eqnIndices); 119 120 int setLoadVectors(GlobalID elemBlock, int numElems, 121 const GlobalID* elemIDs, const double *const * load, 122 int numEqnsPerElem, const int *const * eqnIndices); 123 124 int setMultCREqns(int multCRSetID, int numCRs, int numNodesPerCR, 125 int** nodeNumbers, int** eqnNumbers, int* fieldIDs, 126 int* multiplierEqnNumbers); 127 128 int setPenCREqns(int penCRSetID, int numCRs, int numNodesPerCR, 129 int** nodeNumbers, int** eqnNumbers, int* fieldIDs); 130 131 // ---------------------------------------------------------------------- 132 // setGlobalOffsets : provide info for initial creation of 133 // matrix/vector data, Equation numbers are 1-based, and local sets 134 // of equation numbers are contiguous. 135 // ---------------------------------------------------------------------- 136 137 int setGlobalOffsets(int len, int* nodeOffsets, int* eqnOffsets, 138 int* blkEqnOffsets); 139 140 // ---------------------------------------------------------------------- 141 // setMatrixStructure : provide enough info to allocate the matrix -- 142 // i.e., define the structure. 143 // ---------------------------------------------------------------------- 144 145 int setMatrixStructure(int** ptColIndices, int* ptRowLengths, 146 int** blkColIndices, int* blkRowLengths, int* ptRowsPerBlkRow); 147 148 // ---------------------------------------------------------------------- 149 // resetMatrixAndVector : don't destroy the structure of the matrix, 150 // but set the value 's' throughout the matrix and vectors. 151 // ---------------------------------------------------------------------- 152 153 int resetMatrixAndVector(double s); 154 155 // ---------------------------------------------------------------------- 156 // reset matrix and vector individually 157 // ---------------------------------------------------------------------- 158 159 int resetMatrix(double s); 160 int resetRHSVector(double s); 161 162 // ---------------------------------------------------------------------- 163 // sumIntoSystemMatrix: 164 // this is the primary assembly function. The coefficients 'values' 165 // are to be accumumlated into (added to any values already in place) 166 // global (0-based) equation 'row' of the matrix. 167 // ---------------------------------------------------------------------- 168 169 int sumIntoSystemMatrix(int numPtRows, const int* ptRows, 170 int numPtCols, const int* ptCols, int numBlkRows, 171 const int* blkRows, int numBlkCols, const int* blkCols, 172 const double* const* values); 173 174 // ---------------------------------------------------------------------- 175 // sumIntoSystemMatrix: 176 // this is the primary assembly function. The coefficients 'values' 177 // are to be accumumlated into (added to any values already in place) 178 // global (1-based) [old - 0-based] equation 'row' of the matrix. 179 // ---------------------------------------------------------------------- 180 181 int sumIntoSystemMatrix(int numPtRows, const int* ptRows, 182 int numPtCols, const int* ptCols, 183 const double* const* values); 184 185 // ---------------------------------------------------------------------- 186 // Point-entry matrix data as for 'sumIntoSystemMatrix', but in this case 187 // the data should be "put" into the matrix (i.e., overwrite any coefficients 188 // already present) rather than being "summed" into the matrix. 189 // ---------------------------------------------------------------------- 190 191 int putIntoSystemMatrix(int numPtRows, const int* ptRows, int numPtCols, 192 const int* ptCols, const double* const* values); 193 194 // ---------------------------------------------------------------------- 195 // Get the length of a row of the matrix. 196 // ---------------------------------------------------------------------- 197 198 int getMatrixRowLength(int row, int& length); 199 200 // ---------------------------------------------------------------------- 201 // Obtain the coefficients and indices for a row of the matrix. 202 // ---------------------------------------------------------------------- 203 204 int getMatrixRow(int row, double* coefs, int* indices, int len, 205 int& rowLength); 206 207 // ---------------------------------------------------------------------- 208 // sumIntoRHSVector: 209 // this is the rhs vector equivalent to sumIntoSystemMatrix above. 210 // ---------------------------------------------------------------------- 211 212 int sumIntoRHSVector(int num, const double* values, const int* indices); 213 214 // ---------------------------------------------------------------------- 215 // For putting coefficients into the rhs vector 216 // ---------------------------------------------------------------------- 217 218 int putIntoRHSVector(int num, const double* values, const int* indices); 219 220 // ---------------------------------------------------------------------- 221 // For getting coefficients out of the rhs vector 222 // ---------------------------------------------------------------------- 223 224 int getFromRHSVector(int num, double* values, const int* indices); 225 226 // ---------------------------------------------------------------------- 227 // matrixLoadComplete: 228 // do any internal synchronization/communication. 229 // ---------------------------------------------------------------------- 230 231 int matrixLoadComplete(); 232 233 // ---------------------------------------------------------------------- 234 // Pass nodal data that probably doesn't mean anything to the FEI 235 // implementation, but may mean something to the linear solver. Examples: 236 // geometric coordinates, nullspace data, etc. 237 // ---------------------------------------------------------------------- 238 239 int putNodalFieldData(int fieldID, int fieldSize, int* nodeNumbers, 240 int numNodes, const double* data); 241 242 // ---------------------------------------------------------------------- 243 // functions for enforcing boundary conditions. 244 // ---------------------------------------------------------------------- 245 246 int enforceEssentialBC(int* globalEqn,double* alpha,double* gamma,int len); 247 248 int enforceRemoteEssBCs(int numEqns, int* globalEqns, int** colIndices, 249 int* colIndLen, double** coefs); 250 251 int enforceOtherBC(int* globalEqn, double* alpha, double* beta, 252 double* gamma, int len); 253 254 // ---------------------------------------------------------------------- 255 //functions for getting/setting matrix or vector pointers. 256 // ---------------------------------------------------------------------- 257 258 // ---------------------------------------------------------------------- 259 // getMatrixPtr: 260 // obtain a pointer to the 'A' matrix. This should be considered a 261 // constant pointer -- i.e., this class remains responsible for the 262 // matrix (e.g., de-allocation upon destruction). 263 // ---------------------------------------------------------------------- 264 265 #ifndef NOFEI 266 int getMatrixPtr(Data& data); 267 #endif 268 269 // ---------------------------------------------------------------------- 270 // copyInMatrix: 271 // replaces the internal matrix with a copy of the input argument, scaled 272 // by the coefficient 'scalar'. 273 // ---------------------------------------------------------------------- 274 275 #ifndef NOFEI 276 int copyInMatrix(double scalar, const Data& data); 277 #endif 278 279 // ---------------------------------------------------------------------- 280 // copyOutMatrix: 281 // passes out a copy of the internal matrix, scaled by the coefficient 282 // 'scalar'. 283 // ---------------------------------------------------------------------- 284 285 #ifndef NOFEI 286 int copyOutMatrix(double scalar, Data& data); 287 #endif 288 289 // ---------------------------------------------------------------------- 290 // sumInMatrix: 291 // accumulate (sum) a copy of the input argument into the internal 292 // matrix, scaling the input by the coefficient 'scalar'. 293 // ---------------------------------------------------------------------- 294 295 #ifndef NOFEI 296 int sumInMatrix(double scalar, const Data& data); 297 #endif 298 299 // ---------------------------------------------------------------------- 300 // get/setRHSVectorPtr: 301 // the same semantics apply here as for the matrixPtr functions above. 302 // ---------------------------------------------------------------------- 303 304 #ifndef NOFEI 305 int getRHSVectorPtr(Data& data); 306 #endif 307 308 // ---------------------------------------------------------------------- 309 // copyInRHSVector/copyOutRHSVector/sumInRHSVector: 310 // the same semantics apply here as for the matrix functions above. 311 // ---------------------------------------------------------------------- 312 313 #ifndef NOFEI 314 int copyInRHSVector(double scalar, const Data& data); 315 int copyOutRHSVector(double scalar, Data& data); 316 int sumInRHSVector(double scalar, const Data& data); 317 #endif 318 319 // ---------------------------------------------------------------------- 320 // destroyMatrixData/destroyVectorData: 321 // Utility function for destroying the matrix (or vector) in Data 322 // ---------------------------------------------------------------------- 323 324 #ifndef NOFEI 325 int destroyMatrixData(Data& data); 326 int destroyVectorData(Data& data); 327 #endif 328 329 // ---------------------------------------------------------------------- 330 // functions for managing multiple rhs vectors 331 // ---------------------------------------------------------------------- 332 333 int setNumRHSVectors(int numRHSs, const int* rhsIDs); 334 335 // ---------------------------------------------------------------------- 336 // setRHSID: 337 // set the 'current' rhs context, assuming multiple rhs vectors. 338 // ---------------------------------------------------------------------- 339 340 int setRHSID(int rhsID); 341 342 // ---------------------------------------------------------------------- 343 // putInitialGuess: 344 // function for setting (a subset of) the initial-guess 345 // solution values (i.e., in the 'x' vector). 346 // ---------------------------------------------------------------------- 347 348 int putInitialGuess(const int* eqnNumbers, const double* values,int len); 349 350 // ---------------------------------------------------------------------- 351 // function for getting all of the answers ('x' vector). 352 // ---------------------------------------------------------------------- 353 354 int getSolution(double* answers, int len); 355 356 // ---------------------------------------------------------------------- 357 // function for getting the (single) entry at equation number 'eqnNumber'. 358 // ---------------------------------------------------------------------- 359 360 int getSolnEntry(int eqnNumber, double& answer); 361 362 // ---------------------------------------------------------------------- 363 // This will be called to request that LinearSystemCore form the residual 364 // vector r = b - A*x, and pass the coefficients for r back out in the 365 // 'values' list. 366 // ---------------------------------------------------------------------- 367 368 int formResidual(double* values, int len); 369 double HYPRE_LSC_GetRNorm(); 370 371 // ---------------------------------------------------------------------- 372 // function for launching the linear solver 373 // ---------------------------------------------------------------------- 374 375 int launchSolver(int& solveStatus, int& iterations); 376 377 // ---------------------------------------------------------------------- 378 // other functions 379 // ---------------------------------------------------------------------- 380 381 int writeSystem(const char *); 382 383 // ---------------------------------------------------------------------- 384 // old functions before FEI 1.5 (but still needed here) 385 // ---------------------------------------------------------------------- 386 387 int createMatricesAndVectors(int numGlobalEqns, int firstLocalEqn, 388 int numLocalEqns); 389 390 int allocateMatrix(int** colIndices, int* rowLengths); 391 392 int sumIntoSystemMatrix(int row, int numValues, const double* values, 393 const int* scatterIndices); 394 395 // ---------------------------------------------------------------------- 396 // HYPRE-specific public functions 397 // ---------------------------------------------------------------------- 398 399 void loadConstraintNumbers(int length, int *list); 400 char *getVersion(); 401 void beginCreateMapFromSoln(); 402 void endCreateMapFromSoln(); 403 void putIntoMappedMatrix(int row, int numValues, const double* values, 404 const int* scatterIndices); getFEDataObject(void ** object)405 void getFEDataObject(void **object) { (*object) = feData_; } 406 int HYPRE_LSC_Matvec(void *x, void *y); 407 int HYPRE_LSC_Axpby(double, void *, double, void *); 408 void *HYPRE_LSC_GetRHSVector(); 409 void *HYPRE_LSC_GetSolVector(); 410 void *HYPRE_LSC_GetMatrix(); 411 void *HYPRE_LSC_SetColMap(int, int); 412 void *HYPRE_LSC_MatMatMult(void *); 413 414 // ---------------------------------------------------------------------- 415 // MLI-specific public functions 416 // ---------------------------------------------------------------------- 417 418 void FE_initFields(int nFields, int *fieldSizes, int *fieldIDs); 419 void FE_initElemBlock(int nElems, int nNodesPerElem, int numNodeFields, 420 int *nodeFieldIDs); 421 void FE_initElemNodeList(int elemID, int nNodesPerElem, int *nodeIDs); 422 void FE_initSharedNodes(int nShared, int *sharedIDs, int *sharedPLengs, 423 int **sharedProcs); 424 void FE_initComplete(); 425 void FE_loadElemMatrix(int elemID, int nNodes, int *elemNodeList, 426 int matDim, double **elemMat); 427 428 private: //functions 429 430 // ---------------------------------------------------------------------- 431 // HYPRE specific private functions 432 // ---------------------------------------------------------------------- 433 434 void setupPCGPrecon(); 435 void setupLSICGPrecon(); 436 void setupGMRESPrecon(); 437 void setupFGMRESPrecon(); 438 void setupBiCGSTABPrecon(); 439 void setupBiCGSTABLPrecon(); 440 void setupTFQmrPrecon(); 441 void setupBiCGSPrecon(); 442 void setupSymQMRPrecon(); 443 void setupPreconBoomerAMG(); 444 void setupPreconParaSails(); 445 void setupPreconDDICT(); 446 void setupPreconDDILUT(); 447 void setupPreconPILUT(); 448 void setupPreconPoly(); 449 void setupPreconSchwarz(); 450 void setupPreconML(); 451 void setupPreconMLMaxwell(); 452 void setupPreconAMS(); 453 void setupPreconBlock(); 454 void setupPreconEuclid(); 455 void setupPreconSysPDE(); 456 void solveUsingBoomeramg(int&); 457 double solveUsingSuperLU(int&); 458 double solveUsingSuperLUX(int&); 459 double solveUsingDSuperLU(int&); 460 void solveUsingY12M(int&); 461 void solveUsingAMGe(int&); 462 void buildSlideReducedSystem(); 463 void buildSlideReducedSystem2(); 464 double buildSlideReducedSoln(); 465 double buildSlideReducedSoln2(); 466 void buildSchurReducedSystem(); 467 void buildSchurReducedSystem2(); 468 void buildSlideReducedSystemPartA(int*,int*,int,int,int*,int*); 469 void buildSlideReducedSystemPartB(int*,int*,int,int,int*,int*, 470 HYPRE_ParCSRMatrix *); 471 void buildSlideReducedSystemPartC(int*,int*,int,int,int*,int*, 472 HYPRE_ParCSRMatrix); 473 void buildSchurReducedRHS(); 474 void buildSchurInitialGuess(); 475 double buildSchurReducedSoln(); 476 void computeAConjProjection(HYPRE_ParCSRMatrix A_csr, HYPRE_ParVector x, 477 HYPRE_ParVector b); 478 void computeMinResProjection(HYPRE_ParCSRMatrix A_csr, HYPRE_ParVector x, 479 HYPRE_ParVector b); 480 void addToAConjProjectionSpace(HYPRE_IJVector x, HYPRE_IJVector b); 481 void addToMinResProjectionSpace(HYPRE_IJVector x, HYPRE_IJVector b); 482 int HYPRE_Schur_Search(int,int,int*,int*,int,int); 483 void HYPRE_LSI_BuildNodalCoordinates(); 484 485 // ---------------------------------------------------------------------- 486 // private functions for selecting solver/preconditioner 487 // ---------------------------------------------------------------------- 488 489 void selectSolver(char* name); 490 void selectPreconditioner(char* name); 491 492 private: //variables 493 494 // ---------------------------------------------------------------------- 495 // parallel communication information and output levels 496 // ---------------------------------------------------------------------- 497 498 MPI_Comm comm_; // MPI communicator 499 int numProcs_; // number of processors 500 int mypid_; // my processor ID 501 int HYOutputLevel_; // control print information 502 int memOptimizerFlag_; // turn on memory optimizer 503 504 // ---------------------------------------------------------------------- 505 // for storing information about how to load matrix directly (bypass FEI) 506 // ---------------------------------------------------------------------- 507 508 int mapFromSolnFlag_; 509 int mapFromSolnLeng_; 510 int mapFromSolnLengMax_; 511 int *mapFromSolnList_; 512 int *mapFromSolnList2_; 513 514 // ---------------------------------------------------------------------- 515 // matrix and vectors 516 // ---------------------------------------------------------------------- 517 518 HYPRE_IJMatrix HYA_; // the system matrix 519 HYPRE_IJMatrix HYnormalA_; // normalized system matrix 520 HYPRE_IJVector HYb_; // the current RHS 521 HYPRE_IJVector HYnormalB_; // normalized system rhs 522 HYPRE_IJVector *HYbs_; // an array of RHSs 523 HYPRE_IJVector HYx_; // the solution vector 524 HYPRE_IJVector HYr_; // temporary vector for residual 525 HYPRE_IJVector *HYpxs_; // an array of previous solutions 526 HYPRE_IJVector *HYpbs_; // an array of previous rhs 527 int numGlobalRows_; 528 int localStartRow_; 529 int localEndRow_; 530 int localStartCol_; 531 int localEndCol_; 532 int *rowLengths_; 533 int **colIndices_; 534 double **colValues_; 535 double truncThresh_; 536 double rnorm_; 537 538 // ---------------------------------------------------------------------- 539 // matrix and vectors for reduction 540 // ---------------------------------------------------------------------- 541 542 HYPRE_IJMatrix reducedA_; // matrix for reduction 543 HYPRE_IJVector reducedB_; // RHS vector for reduction 544 HYPRE_IJVector reducedX_; // solution vector for reduction 545 HYPRE_IJVector reducedR_; // temporary vector for reduction 546 HYPRE_IJMatrix HYA21_; // (2,1) block in reduction 547 HYPRE_IJMatrix HYA12_; // (1,2) block in reduction 548 int A21NRows_; // number of rows in (2,1) block 549 int A21NCols_; // number of cols in (2,1) block 550 int reducedAStartRow_; // Nrows in reduced system 551 HYPRE_IJMatrix HYinvA22_; // inv(A22) in slide reduction 552 553 // ---------------------------------------------------------------------- 554 // pointers to current matrix and vectors for solver 555 // ---------------------------------------------------------------------- 556 557 HYPRE_IJMatrix currA_; 558 HYPRE_IJVector currB_; 559 HYPRE_IJVector currX_; 560 HYPRE_IJVector currR_; 561 int currentRHS_; 562 int *rhsIDs_; 563 int numRHSs_; 564 int nStored_; 565 int *storedIndices_; 566 int *auxStoredIndices_; 567 int mRHSFlag_; 568 int mRHSNumGEqns_; 569 int *mRHSGEqnIDs_; 570 int *mRHSNEntries_; 571 int *mRHSBCType_; 572 int **mRHSRowInds_; 573 double **mRHSRowVals_; 574 575 // ---------------------------------------------------------------------- 576 // flags for matrix assembly, various reductions, and projections 577 // ---------------------------------------------------------------------- 578 579 int matrixVectorsCreated_; 580 int systemAssembled_; 581 int slideReduction_; 582 double slideReductionMinNorm_; 583 int slideReductionScaleMatrix_; 584 int schurReduction_; 585 int schurReductionCreated_; 586 int projectionScheme_; 587 int projectSize_; 588 int projectCurrSize_; 589 double **projectionMatrix_; 590 int normalEqnFlag_; 591 void *slideObj_; 592 593 // ---------------------------------------------------------------------- 594 // variables for slide and Schur reduction 595 // ---------------------------------------------------------------------- 596 597 int *selectedList_; 598 int *selectedListAux_; 599 int nConstraints_; 600 int *constrList_; 601 int matrixPartition_; 602 603 // ---------------------------------------------------------------------- 604 // variables for the selected solver and preconditioner 605 // ---------------------------------------------------------------------- 606 607 char *HYSolverName_; 608 HYPRE_Solver HYSolver_; 609 HYsolverID HYSolverID_; 610 int gmresDim_; 611 int fgmresUpdateTol_; 612 int maxIterations_; 613 double tolerance_; 614 int normAbsRel_; 615 int pcgRecomputeRes_; 616 617 char *HYPreconName_; 618 HYPRE_Solver HYPrecon_; 619 HYpreconID HYPreconID_; 620 int HYPreconReuse_; 621 int HYPreconSetup_; 622 623 // ---------------------------------------------------------------------- 624 // preconditioner specific variables 625 // ---------------------------------------------------------------------- 626 627 int amgMaxLevels_; 628 int amgCoarsenType_; 629 int amgMaxIter_; 630 int amgMeasureType_; 631 int amgNumSweeps_[4]; 632 int amgRelaxType_[4]; 633 int amgGridRlxType_; 634 double amgRelaxWeight_[25]; 635 double amgRelaxOmega_[25]; 636 double amgStrongThreshold_; 637 int amgSystemSize_; 638 int amgSmoothType_; 639 int amgSmoothNumLevels_; 640 int amgSmoothNumSweeps_; 641 int amgCGSmoothNumSweeps_; 642 double amgSchwarzRelaxWt_; 643 int amgSchwarzVariant_; 644 int amgSchwarzOverlap_; 645 int amgSchwarzDomainType_; 646 int amgUseGSMG_; 647 int amgGSMGNSamples_; 648 int amgAggLevels_; 649 int amgInterpType_; 650 int amgPmax_; 651 int pilutFillin_; 652 double pilutDropTol_; 653 int pilutMaxNnzPerRow_; 654 int parasailsSym_; 655 double parasailsThreshold_; 656 int parasailsNlevels_; 657 double parasailsFilter_; 658 double parasailsLoadbal_; 659 int parasailsReuse_; 660 int mlMethod_; 661 int mlNumPreSweeps_; 662 int mlNumPostSweeps_; 663 int mlPresmootherType_; 664 int mlPostsmootherType_; 665 double mlRelaxWeight_; 666 double mlStrongThreshold_; 667 int mlCoarseSolver_; 668 int mlCoarsenScheme_; 669 int mlNumPDEs_; 670 int superluOrdering_; 671 char superluScale_[1]; 672 double ddilutFillin_; 673 double ddilutDropTol_; 674 int ddilutOverlap_; 675 int ddilutReorder_; 676 double ddictFillin_; 677 double ddictDropTol_; 678 double schwarzFillin_; 679 int schwarzNblocks_; 680 int schwarzBlksize_; 681 int polyOrder_; 682 int euclidargc_; 683 char **euclidargv_; 684 HYPRE_IJVector amsX_; 685 HYPRE_IJVector amsY_; 686 HYPRE_IJVector amsZ_; 687 int localStartRowAMSV_; 688 int localEndRowAMSV_; 689 HYPRE_IJMatrix amsG_; 690 HYPRE_IJMatrix amsD0_; 691 HYPRE_IJMatrix amsD1_; 692 int localStartRowAMSG_; 693 int localEndRowAMSG_; 694 int localStartColAMSG_; 695 int localEndColAMSG_; 696 HYPRE_ParCSRMatrix amsBetaPoisson_; 697 int amsNumPDEs_; 698 int amsMaxIter_; 699 double amsTol_; 700 int amsCycleType_; 701 int amsRelaxType_; 702 int amsRelaxTimes_; 703 double amsRelaxWt_; 704 double amsRelaxOmega_; 705 int amsPrintLevel_; 706 int amsAlphaCoarsenType_; 707 int amsAlphaAggLevels_; 708 int amsAlphaRelaxType_; 709 double amsAlphaStrengthThresh_; 710 int amsAlphaInterpType_; 711 int amsAlphaPmax_; 712 int amsBetaCoarsenType_; 713 int amsBetaAggLevels_; 714 int amsBetaRelaxType_; 715 double amsBetaStrengthThresh_; 716 int amsBetaInterpType_; 717 int amsBetaPmax_; 718 int sysPDEMethod_; 719 int sysPDEFormat_; 720 double sysPDETol_; 721 int sysPDEMaxIter_; 722 int sysPDENumPre_; 723 int sysPDENumPost_; 724 int sysPDENVars_; 725 726 // ---------------------------------------------------------------------- 727 // FEI and MLI variables 728 // ---------------------------------------------------------------------- 729 730 void *feData_; 731 int haveFEData_; 732 Lookup *lookup_; 733 int haveLookup_; 734 int MLI_NumNodes_; 735 int MLI_FieldSize_; 736 int *MLI_EqnNumbers_; 737 double *MLI_NodalCoord_; 738 int MLI_Hybrid_NSIncr_; 739 int MLI_Hybrid_GSA_; 740 int MLI_Hybrid_MaxIter_; 741 double MLI_Hybrid_ConvRate_; 742 int MLI_Hybrid_NTrials_; 743 HYPRE_FEI_AMSData AMSData_; 744 int FEI_mixedDiagFlag_; 745 double *FEI_mixedDiag_; 746 747 // ---------------------------------------------------------------------- 748 // ML Maxwell variables 749 // ---------------------------------------------------------------------- 750 751 HYPRE_ParCSRMatrix maxwellANN_; // Maxwell nodal matrix 752 HYPRE_ParCSRMatrix maxwellGEN_; // Maxwell gradient matrix 753 754 // ---------------------------------------------------------------------- 755 // temporary functions for testing purposes 756 // ---------------------------------------------------------------------- 757 758 friend void fei_hypre_test(int argc, char *argv[]); 759 friend void fei_hypre_domaindecomposition(int argc, char *argv[]); 760 761 }; 762 763 #endif 764 765