Demo CUDSS solver

.clang-format

@@ -0,0 +1,95 @@
+---
+Language:        Cpp
+# BasedOnStyle:  LLVM
+AccessModifierOffset: -2
+AlignAfterOpenBracket: Align
+AlignConsecutiveAssignments: false
+AlignConsecutiveDeclarations: false
+AlignEscapedNewlinesLeft: false
+AlignOperands:   true
+AlignTrailingComments: true
+AllowAllParametersOfDeclarationOnNextLine: true
+AllowShortBlocksOnASingleLine: false
+AllowShortCaseLabelsOnASingleLine: false
+AllowShortFunctionsOnASingleLine: All
+AllowShortIfStatementsOnASingleLine: false
+AllowShortLoopsOnASingleLine: false
+AlwaysBreakAfterDefinitionReturnType: None
+AlwaysBreakAfterReturnType: None
+AlwaysBreakBeforeMultilineStrings: false
+AlwaysBreakTemplateDeclarations: true
+BinPackArguments: true
+BinPackParameters: true
+BraceWrapping:
+  AfterClass:      false
+  AfterControlStatement: false
+  AfterEnum:       false
+  AfterFunction:   false
+  AfterNamespace:  false
+  AfterObjCDeclaration: false
+  AfterStruct:     false
+  AfterUnion:      false
+  BeforeCatch:     false
+  BeforeElse:      false
+  IndentBraces:    false
+BreakBeforeBinaryOperators: All
+BreakBeforeBraces: Allman
+BreakBeforeTernaryOperators: true
+BreakConstructorInitializersBeforeComma: false
+BreakAfterJavaFieldAnnotations: false
+BreakStringLiterals: true
+ColumnLimit:     80
+CommentPragmas:  '^ IWYU pragma:'
+ConstructorInitializerAllOnOneLineOrOnePerLine: false
+ConstructorInitializerIndentWidth: 4
+ContinuationIndentWidth: 4
+Cpp11BracedListStyle: true
+DerivePointerAlignment: false
+DisableFormat:   false
+ExperimentalAutoDetectBinPacking: false
+ForEachMacros:   [ foreach, Q_FOREACH, BOOST_FOREACH ]
+IncludeCategories:
+  - Regex:           '^"(llvm|llvm-c|clang|clang-c)/'
+    Priority:        2
+  - Regex:           '^(<|"(gtest|isl|json)/)'
+    Priority:        3
+  - Regex:           '.*'
+    Priority:        1
+IncludeIsMainRegex: '$'
+IndentCaseLabels: false
+IndentWidth:     2
+IndentWrappedFunctionNames: false
+JavaScriptQuotes: Leave
+JavaScriptWrapImports: true
+KeepEmptyLinesAtTheStartOfBlocks: true
+MacroBlockBegin: ''
+MacroBlockEnd:   ''
+MaxEmptyLinesToKeep: 1
+NamespaceIndentation: None
+ObjCBlockIndentWidth: 2
+ObjCSpaceAfterProperty: false
+ObjCSpaceBeforeProtocolList: true
+PenaltyBreakBeforeFirstCallParameter: 19
+PenaltyBreakComment: 300
+PenaltyBreakFirstLessLess: 120
+PenaltyBreakString: 1000
+PenaltyExcessCharacter: 1000000
+PenaltyReturnTypeOnItsOwnLine: 60
+PointerAlignment: Left
+ReflowComments:  true
+SortIncludes:    true
+SpaceAfterCStyleCast: false
+SpaceAfterTemplateKeyword: true
+SpaceBeforeAssignmentOperators: true
+SpaceBeforeParens: ControlStatements
+SpaceInEmptyParentheses: false
+SpacesBeforeTrailingComments: 1
+SpacesInAngles:  false
+SpacesInContainerLiterals: true
+SpacesInCStyleCastParentheses: false
+SpacesInParentheses: false
+SpacesInSquareBrackets: false
+Standard:        Cpp11
+TabWidth:        8
+UseTab:          Never
+...

dolfinx-gpu/LICENSE

@@ -0,0 +1,21 @@
+MIT License
+
+Copyright (c) 2025 FEniCS Project
+
+Permission is hereby granted, free of charge, to any person obtaining a copy
+of this software and associated documentation files (the "Software"), to deal
+in the Software without restriction, including without limitation the rights
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
+copies of the Software, and to permit persons to whom the Software is
+furnished to do so, subject to the following conditions:
+
+The above copyright notice and this permission notice shall be included in all
+copies or substantial portions of the Software.
+
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE
+SOFTWARE.

dolfinx-gpu/README.md

@@ -0,0 +1,9 @@
+# dolfinx-gpu
+
+In this folder:
+
+* a working C++ example, which demonstrates the use of `cuDSS` and `cuSPARSE` - only available for NVIDIA
+
+* Python wrappers, using the same C++ headers, which can be used alongside `dolfinx` and has an interface to `cupy`
+
+The Python wrappers will compile for both CUDA and HIP, but cuDSS (sparse solver) is only available in CUDA.

dolfinx-gpu/cpp/demo/heat_equation/CMakeLists.txt

@@ -0,0 +1,36 @@
+# This file was generated by running
+#
+# python cmake/scripts/generate-cmakefiles from dolfinx/cpp
+#
+cmake_minimum_required(VERSION 3.21)
+
+set(PROJECT_NAME demo_heat)
+project(${PROJECT_NAME} LANGUAGES C CXX CUDA)
+
+if(NOT TARGET dolfinx)
+  find_package(DOLFINX REQUIRED)
+endif()
+
+add_custom_command(
+  OUTPUT heat.c
+  COMMAND ffcx ${CMAKE_CURRENT_SOURCE_DIR}/heat.py
+  VERBATIM
+  DEPENDS heat.py
+  COMMENT "Compile heat.py using FFCx"
+)
+
+set(CMAKE_INCLUDE_CURRENT_DIR ON)
+set(CMAKE_CUDA_ARCHITECTURES 89)
+set(CMAKE_CUDA_FLAGS "${CMAKE_CUDA_FLAGS} -g -march=native -std=c++20 -O3 -DNDEBUG --extended-lambda --expt-relaxed-constexpr -lineinfo")
+
+set_source_files_properties(main.cpp PROPERTIES LANGUAGE CUDA)
+
+find_package(cudss)
+find_package(CUDAToolkit)
+
+add_executable(${PROJECT_NAME} main.cpp ${CMAKE_CURRENT_BINARY_DIR}/heat.c)
+target_link_libraries(${PROJECT_NAME} dolfinx cudss ${CUDA_cusparse_LIBRARY})
+
+# Set C++20 standard
+set(CMAKE_CXX_EXTENSIONS OFF)
+target_compile_features(${PROJECT_NAME} PUBLIC cxx_std_20)

dolfinx-gpu/cpp/demo/heat_equation/README.md

@@ -0,0 +1,23 @@
+# Solver for the heat equation using cuDSS from NVIDIA
+
+This code solves the time-dependent heat equation on a 100x100x100
+cube, using NVIDIA's cuDSS sparse direct solver, and cuSPARSE to apply
+a mass matrix.
+
+## Prerequisites
+
+* NVIDIA CUDA 12.9+
+* NVIDIA cuSPARSE
+* NVIDIA cuDSS library
+* dolfinx
+* ffcx
+
+## Building
+
+Create a `build` directory, then:
+
+```
+cd build
+cmake ..
+make
+```

dolfinx-gpu/cpp/demo/heat_equation/heat.py

@@ -0,0 +1,46 @@
+# The first step is to define the variational problem at hand. We define
+# the variational problem in UFL terms in a separate form file
+# {download}`demo_poisson/poisson.py`.  We begin by defining the finite
+# element:
+
+from basix.ufl import element
+from ufl import (
+    Coefficient,
+    Constant,
+    FunctionSpace,
+    Mesh,
+    TestFunction,
+    TrialFunction,
+    ds,
+    dx,
+    grad,
+    inner,
+)
+
+e = element("Lagrange", "tetrahedron", 1)
+
+# The first argument to :py:class:`FiniteElement` is the finite element
+# family, the second argument specifies the domain, while the third
+# argument specifies the polynomial degree. Thus, in this case, our
+# element `element` consists of first-order, continuous Lagrange basis
+# functions on triangles (or in order words, continuous piecewise linear
+# polynomials on triangles).
+#
+# Next, we use this element to initialize the trial and test functions
+# ($u$ and $v$) and the coefficient functions ($f$ and $g$):
+
+coord_element = element("Lagrange", "tetrahedron", 1, shape=(3,))
+mesh = Mesh(coord_element)
+
+V = FunctionSpace(mesh, e)
+
+u = TrialFunction(V)
+v = TestFunction(V)
+f = Coefficient(V)
+kappa = Constant(mesh)
+
+# Finally, we define the bilinear and linear forms according to the
+# variational formulation of the equations:
+
+a = kappa * inner(grad(u), grad(v)) * dx + inner(u, v) * dx
+L = inner(f, v) * dx