boomerAMG on GPU for SolidMechanicsLagrangianSSLE

examples/GPU/beamBending.xml

@@ -0,0 +1,161 @@
+<?xml version="1.0" ?>
+
+<Problem>
+  <!-- SPHINX_Gravity -->
+  <Solvers
+    gravityVector="0.0, 0.0, 0.0">
+    <!-- SPHINX_GravityEnd -->
+    <!-- SPHINX_SolidMechanicsSolver -->
+    <SolidMechanicsLagrangianSSLE
+      name="lagsolve"
+      timeIntegrationOption="QuasiStatic"
+      discretization="FE1"
+      logLevel="0"
+      targetRegions="{ Region2 }"
+      solidMaterialNames="{ shale }">
+      <!-- SPHINX_SolidMechanicsSolverEnd -->
+      <NonlinearSolverParameters
+        newtonTol="1.0e-6"
+        newtonMaxIter="8"/>
+      <LinearSolverParameters
+        logLevel="1"
+        solverType="gmres"
+        krylovTol="1.0e-12"
+        preconditionerType="amg"/>
+    </SolidMechanicsLagrangianSSLE>
+  </Solvers>
+
+  <!-- SPHINX_BeamBendingMesh -->
+  <Mesh>
+    <InternalMesh
+      name="mesh1"
+      elementTypes="{ C3D8 }"
+      xCoords="{ 0, 80 }"
+      yCoords="{ 0, 8 }"
+      zCoords="{ 0, 4 }"
+      nx="{ 320 }"
+      ny="{ 32 }"
+      nz="{ 16 }"
+      cellBlockNames="{ cb1 }"/>
+  </Mesh>
+
+  <!-- SPHINX_BeamBendingMeshEnd -->
+  <Events
+    maxTime="2.0">
+    <!-- This event is applied every cycle, and overrides the
+    solver time-step request -->
+    <!-- SPHINX_Execution -->
+    <PeriodicEvent
+      name="solverApplications"
+      forceDt="1.0"
+      target="/Solvers/lagsolve"/>
+
+    <!-- SPHINX_ExecutionEnd -->
+    <!-- This event is applied every 5.0e-5s.  The targetExactTimestep
+    flag allows this event to request a dt modification to match an
+    integer multiple of the timeFrequency. -->
+    <!--PeriodicEvent
+      name="outputs"
+      timeFrequency="1.0"
+      targetExactTimestep="1"
+      target="/Outputs/siloOutput"/-->
+
+    <!--PeriodicEvent
+      name="restarts"
+      timeFrequency="1e99"
+      targetExactTimestep="0"
+      target="/Outputs/restartOutput"/-->
+  </Events>
+
+  <!-- SPHINX_FEDiscretization -->
+  <NumericalMethods>
+    <FiniteElements>
+      <FiniteElementSpace
+        name="FE1"
+        order="1"/>
+    </FiniteElements>
+  </NumericalMethods>
+
+  <!-- SPHINX_FEDiscretizationEnd -->
+  <ElementRegions>
+    <CellElementRegion
+      name="Region2"
+      cellBlocks="{ cb1 }"
+      materialList="{ shale }"/>
+  </ElementRegions>
+
+  <Constitutive>
+    <LinearElasticIsotropic
+      name="granite"
+      defaultDensity="2700"
+      defaultBulkModulus="5.5556e9"
+      defaultShearModulus="4.16667e9"/>
+
+    <!-- SPHINX_Material -->
+    <LinearElasticIsotropic
+      name="shale"
+      defaultDensity="2700"
+      defaultBulkModulus="5.5556e9"
+      defaultShearModulus="4.16667e9"/>
+
+    <!-- SPHINX_MaterialEnd -->
+  </Constitutive>
+
+  <!-- SPHINX_BoundaryConditions -->
+  <FieldSpecifications>
+    <FieldSpecification
+      name="xnegconstraint"
+      objectPath="nodeManager"
+      fieldName="TotalDisplacement"
+      component="0"
+      scale="0.0"
+      setNames="{ xneg }"/>
+
+    <FieldSpecification
+      name="yconstraint"
+      objectPath="nodeManager"
+      fieldName="TotalDisplacement"
+      component="1"
+      scale="0.0"
+      setNames="{ xneg }"/>
+
+    <FieldSpecification
+      name="zconstraint"
+      objectPath="nodeManager"
+      fieldName="TotalDisplacement"
+      component="2"
+      scale="0.0"
+      setNames="{ zneg, zpos }"/>
+
+    <FieldSpecification
+      name="xposconstraint"
+      objectPath="faceManager"
+      fieldName="Traction"
+      component="1"
+      scale="1.0e6"
+      functionName="timeFunction"
+      setNames="{ xpos }"/>
+  </FieldSpecifications>
+
+  <!-- SPHINX_BoundaryConditionsEnd -->
+  <!-- SPHINX_TableFunction -->
+  <Functions>
+    <TableFunction
+      name="timeFunction"
+      inputVarNames="{ time }"
+      coordinates="{ 0.0, 10.0 }"
+      values="{ 0.0, 10.0 }"/>
+  </Functions>
+
+  <!-- SPHINX_TableFunctionEnd -->
+  <Outputs>
+    <Silo
+      name="siloOutput"
+      parallelThreads="32"
+      plotFileRoot="plot"
+      childDirectory="sub"/>
+
+    <Restart
+      name="restartOutput"/>
+  </Outputs>
+</Problem>

host-configs/LLNL/lassen-clang@upstream.cmake

@@ -34,4 +34,4 @@ set(ESSL_LIBRARIES /usr/tcetmp/packages/essl/essl-6.2/lib64/libesslsmpcuda.so
 # PETSc doesn't seem to work correctly with clang.
 set(ENABLE_PETSC OFF CACHE BOOL "" FORCE )
 #set(PETSC_OMP_DIR ${GEOSX_TPL_ROOT_DIR}/omp-links-for-petsc CACHE STRING "")
-set(ENABLE_HYPRE OFF CACHE BOOL "" FORCE )
+set(ENABLE_HYPRE ON CACHE BOOL "" FORCE )

src/coreComponents/linearAlgebra/interfaces/hypre/HypreInterface.cpp

@@ -18,16 +18,27 @@
 
 #include "HypreInterface.hpp"
 #include "linearAlgebra/interfaces/hypre/HyprePreconditioner.hpp"
+#include "HYPRE_utilities.h"
+#include "_hypre_utilities.h"
 
 namespace geosx
 {
 
 void HypreInterface::initialize( int & GEOSX_UNUSED_PARAM( argc ),
                                  char * * & GEOSX_UNUSED_PARAM( argv ) )
-{}
+{
+  HYPRE_Init();
+#if defined(USE_CUDA)
+  HYPRE_ExecutionPolicy default_exec_policy = HYPRE_EXEC_DEVICE;
+  hypre_HandleDefaultExecPolicy(hypre_handle()) = HYPRE_EXEC_DEVICE;
+  hypre_HandleSpgemmUseCusparse(hypre_handle()) = 1;
+#endif
+}
 
 void HypreInterface::finalize()
-{}
+{
+  HYPRE_Finalize();
+}
 
 std::unique_ptr< PreconditionerBase< HypreInterface > >
 geosx::HypreInterface::createPreconditioner( LinearSolverParameters params )

src/coreComponents/linearAlgebra/interfaces/hypre/HypreMatrix.cpp

@@ -267,6 +267,85 @@ void HypreMatrix::close()
   m_assembled = true;
 }
 
+//void HypreMatrix::create( CRSMatrixView< real64 const, globalIndex const > const & localMatrix,
+//                          MPI_Comm const & comm )
+//{
+//
+//
+//
+//  HYPRE_BigInt *row_starts;
+//  HYPRE_BigInt *col_starts;
+//  HYPRE_Int num_cols_offd;
+//  HYPRE_Int num_nonzeros_diag;
+//  HYPRE_Int num_nonzeros_offd;
+//
+//  HYPRE_Int *diag_i;
+//  HYPRE_Int *diag_j;
+//  double *diag_data;
+//  HYPRE_Int *offd_i;
+//  HYPRE_Int *offd_j;
+//  double *offd_data;
+//  HYPRE_Int offd_num_cols;
+//  HYPRE_Int *offd_col_map;
+//
+//
+//
+//  m_parcsr_mat = hypre_ParCSRMatrixCreate( comm,
+//                                           localMatrix.numRows(),
+//                                           localMatrix.numRows(),
+//                                           row_starts,
+//                                           col_starts,
+//                                           num_cols_offd,
+//                                           num_nonzeros_diag,
+//                                           num_nonzeros_offd);
+//
+//
+//
+//
+//
+//     hypre_ParCSRMatrixSetDataOwner(m_parcsr_mat,1);
+//     hypre_ParCSRMatrixSetRowStartsOwner(m_parcsr_mat,0);
+//     hypre_ParCSRMatrixSetColStartsOwner(m_parcsr_mat,0);
+//
+//     HYPRE_Int local_num_rows = hypre_CSRMatrixNumRows(m_parcsr_mat->diag);
+//
+//     hypre_CSRMatrixSetDataOwner(m_parcsr_mat->diag,0);
+//     hypre_CSRMatrixI(m_parcsr_mat->diag) = diag_i;
+//     hypre_CSRMatrixJ(m_parcsr_mat->diag) = diag_j;
+//     hypre_CSRMatrixData(m_parcsr_mat->diag) = diag_data;
+//     hypre_CSRMatrixNumNonzeros(m_parcsr_mat->diag) = diag_i[local_num_rows];
+//     hypre_CSRMatrixSetRownnz(m_parcsr_mat->diag);
+//     // Prevent hypre from destroying m_parcsr_mat->diag->{i,j,data}, own m_parcsr_mat->diag->{i,j,data}
+//     //diagOwner = 3;
+//
+//     hypre_CSRMatrixSetDataOwner(m_parcsr_mat->offd,0);
+//     hypre_CSRMatrixI(m_parcsr_mat->offd) = offd_i;
+//     hypre_CSRMatrixJ(m_parcsr_mat->offd) = offd_j;
+//     hypre_CSRMatrixData(m_parcsr_mat->offd) = offd_data;
+//     hypre_CSRMatrixNumNonzeros(m_parcsr_mat->offd) = offd_i[local_num_rows];
+//     hypre_CSRMatrixSetRownnz(m_parcsr_mat->offd);
+//     // Prevent hypre from destroying m_parcsr_mat->offd->{i,j,data}, own m_parcsr_mat->offd->{i,j,data}
+//     //offdOwner = 3;
+//
+//     hypre_ParCSRMatrixColMapOffd(m_parcsr_mat) = offd_col_map;
+//     // Prevent hypre from destroying m_parcsr_mat->col_map_offd, own m_parcsr_mat->col_map_offd
+//     //colMapOwner = 1;
+//
+//     hypre_ParCSRMatrixSetNumNonzeros(m_parcsr_mat);
+//
+//     /* Make sure that the first entry in each row is the diagonal one. */
+//     if (row_starts == col_starts)
+//     {
+//        hypre_CSRMatrixReorder(hypre_ParCSRMatrixDiag(m_parcsr_mat));
+//     }
+//
+//     hypre_MatvecCommPkgCreate(m_parcsr_mat);
+//
+//     height = GetNumRows();
+//     width = GetNumCols();
+//}
+
+
 bool HypreMatrix::created() const
 {
   return m_ij_mat != nullptr;

src/coreComponents/linearAlgebra/interfaces/hypre/HypreMatrix.hpp

@@ -118,6 +118,9 @@ class HypreMatrix final : public virtual LinearOperator< HypreVector >,
 
   virtual void close() override;
 
+//  virtual void create( CRSMatrixView< real64 const, globalIndex const > const & localMatrix,
+//                       MPI_Comm const & comm ) override;
+
   virtual bool created() const override;
 
   virtual void reset() override;
@@ -336,6 +339,9 @@ class HypreMatrix final : public virtual LinearOperator< HypreVector >,
    */
   HYPRE_ParCSRMatrix m_parcsr_mat = nullptr;
 
+//  CRSMatrix< real64, globalIndex > & m_diagCRS;
+//  CRSMatrix< real64, globalIndex > & m_offDiagCRS;
+
 };
 
 } // namespace geosx

src/coreComponents/linearAlgebra/interfaces/hypre/HypreSolver.cpp

@@ -311,6 +311,12 @@ void HypreSolver::solve_krylov( HypreMatrix & mat,
 {
   Stopwatch watch;
 
+#if defined(USE_CUDA)
+  HYPRE_ExecutionPolicy default_exec_policy = HYPRE_EXEC_DEVICE;
+  hypre_HandleDefaultExecPolicy(hypre_handle()) = HYPRE_EXEC_DEVICE;
+  hypre_HandleSpgemmUseCusparse(hypre_handle()) = 1;
+#endif
+
   // Create the preconditioner, but don't compute (this is done by solver setup)
   HyprePreconditioner precond( m_parameters, dofManager );
 
@@ -340,6 +346,13 @@ void HypreSolver::solve_krylov( HypreMatrix & mat,
     HYPRE_BoomerAMGSetAggNumLevels( uu_amg_solver, 1 );
     HYPRE_BoomerAMGSetNumFunctions( uu_amg_solver, 3 );
 
+    // Relaxation options: Use options 18 or 7.
+    // Coarsening options: Only PMIS is supported
+    // Interpolation options: Use options 3, 6, 14 or 15.
+    HYPRE_BoomerAMGSetRelaxType( uu_amg_solver, 7 );
+    HYPRE_BoomerAMGSetCoarsenType( uu_amg_solver, 8 );
+    HYPRE_BoomerAMGSetInterpType( uu_amg_solver, 3 );
+
     HYPRE_BoomerAMGSetup( uu_amg_solver, separateComponentMatrix.unwrapped(), nullptr, nullptr );
 
     HYPRE_MGRSetFSolver( precond.unwrapped(), HYPRE_BoomerAMGSolve, HYPRE_BoomerAMGSetup, uu_amg_solver );

src/coreComponents/managers/initialization.cpp

@@ -445,6 +445,8 @@ void basicSetup( int argc, char * argv[], bool const parseCommandLine )
   setupMKL();
   setupLAI( argc, argv );
 
+
+
   if( parseCommandLine )
   {
     internal::parseCommandLineOptions( argc, argv );

src/coreComponents/physicsSolvers/solidMechanics/integratedTests/SSLE-QS-beamBending.xml

@@ -18,6 +18,7 @@
         newtonTol="1.0e-6"
         newtonMaxIter="8"/>
       <LinearSolverParameters
+        logLevel="1"
         solverType="gmres"
         krylovTol="1.0e-12"
         preconditionerType="amg"/>
@@ -32,9 +33,9 @@
       xCoords="{ 0, 80 }"
       yCoords="{ 0, 8 }"
       zCoords="{ 0, 4 }"
-      nx="{ 80 }"
-      ny="{ 8 }"
-      nz="{ 4 }"
+      nx="{ 160 }"
+      ny="{ 16 }"
+      nz="{ 8 }"
       cellBlockNames="{ cb1 }"/>
   </Mesh>