diff --git a/gen/solver_options.jl b/gen/solver_options.jl
index b30fe32..7cdf98f 100644
--- a/gen/solver_options.jl
+++ b/gen/solver_options.jl
@@ -57,6 +57,7 @@ open(joinpath(@__DIR__, "..", "src", "solver_options.jl"), "w") do io
     generate_options(io, "BoomerAMG", "HYPRE_BoomerAMGSet")
     # generate_options(io, "FSAI", "HYPRE_FSAISet")
     generate_options(io, "GMRES", "HYPRE_ParCSRGMRESSet", "HYPRE_GMRESSet")
+    generate_options(io, "ILU", "HYPRE_ILUSet")
     generate_options(io, "ParaSails", "HYPRE_ParCSRParaSailsSet", "HYPRE_ParaSailsSet")
     generate_options(io, "PCG", "HYPRE_ParCSRPCGSet", "HYPRE_PCGSet")
 end
diff --git a/src/solver_options.jl b/src/solver_options.jl
index 789026d..50d6604 100644
--- a/src/solver_options.jl
+++ b/src/solver_options.jl
@@ -323,6 +323,41 @@ function Internals.set_options(s::GMRES, kwargs)
     end
 end
 
+function Internals.set_options(s::ILU, kwargs)
+    solver = s.solver
+    for (k, v) in kwargs
+        if k === :DropThreshold
+            @check HYPRE_ILUSetDropThreshold(solver, v)
+        elseif k === :DropThresholdArray
+            @check HYPRE_ILUSetDropThresholdArray(solver, v)
+        elseif k === :LevelOfFill
+            @check HYPRE_ILUSetLevelOfFill(solver, v)
+        elseif k === :LocalReordering
+            @check HYPRE_ILUSetLocalReordering(solver, v)
+        elseif k === :Logging
+            @check HYPRE_ILUSetLogging(solver, v)
+        elseif k === :MaxIter
+            @check HYPRE_ILUSetMaxIter(solver, v)
+        elseif k === :MaxNnzPerRow
+            @check HYPRE_ILUSetMaxNnzPerRow(solver, v)
+        elseif k === :NSHDropThreshold
+            @check HYPRE_ILUSetNSHDropThreshold(solver, v)
+        elseif k === :NSHDropThresholdArray
+            @check HYPRE_ILUSetNSHDropThresholdArray(solver, v)
+        elseif k === :PrintLevel
+            @check HYPRE_ILUSetPrintLevel(solver, v)
+        elseif k === :SchurMaxIter
+            @check HYPRE_ILUSetSchurMaxIter(solver, v)
+        elseif k === :Tol
+            @check HYPRE_ILUSetTol(solver, v)
+        elseif k === :Type
+            @check HYPRE_ILUSetType(solver, v)
+        else
+            throw(ArgumentError("unknown option $k for HYPRE.ILU"))
+        end
+    end
+end
+
 function Internals.set_options(s::ParaSails, kwargs)
     solver = s.solver
     for (k, v) in kwargs
diff --git a/src/solvers.jl b/src/solvers.jl
index b5ceeee..a341958 100644
--- a/src/solvers.jl
+++ b/src/solvers.jl
@@ -227,6 +227,41 @@ function Internals.set_precond(gmres::GMRES, p::HYPRESolver)
 end
 
 
+#######
+# ILU #
+#######
+
+mutable struct ILU <: HYPRESolver
+    solver::HYPRE_Solver
+    function ILU(; kwargs...)
+        solver = new(C_NULL)
+        solver_ref = Ref{HYPRE_Solver}(C_NULL)
+        @check HYPRE_ILUCreate(solver_ref)
+        solver.solver = solver_ref[]
+        # Attach a finalizer
+        finalizer(x -> HYPRE_ILUDestroy(x.solver), solver)
+        # Set the options
+        Internals.set_options(solver, kwargs)
+        return solver
+    end
+end
+
+function solve!(ilu::ILU, x::HYPREVector, A::HYPREMatrix, b::HYPREVector)
+    @check HYPRE_ILUSetup(ilu.solver, A.parmatrix, b.parvector, x.parvector)
+    @check HYPRE_ILUSolve(ilu.solver, A.parmatrix, b.parvector, x.parvector)
+    return x
+end
+
+Internals.setup_func(::ILU) = HYPRE_ILUSetup
+Internals.solve_func(::ILU) = HYPRE_ILUSolve
+
+function Internals.set_precond_defaults(ilu::ILU)
+    defaults = (; Tol = 0.0, MaxIter = 1)
+    Internals.set_options(ilu, pairs(defaults))
+    return nothing
+end
+
+
 #####################
 # (ParCSR)ParaSails #
 #####################
diff --git a/test/runtests.jl b/test/runtests.jl
index 3d7bde8..b65e0c0 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -409,6 +409,46 @@ end
     @test x ≈ A \ b atol=tol
 end
 
+@testset "ILU" begin
+    # Solver constructor and options
+    @test_throws(
+        ArgumentError("unknown option UnknownOption for HYPRE.ILU"),
+        HYPRE.ILU(; UnknownOption = 1)
+    )
+    # Setup
+    A = sprand(100, 100, 0.05); A = A'A + 5I
+    b = rand(100)
+    x = zeros(100)
+    A_h = HYPREMatrix(A)
+    b_h = HYPREVector(b)
+    x_h = HYPREVector(x)
+    # Solve
+    tol = 1e-9
+    ilu = HYPRE.ILU(; Tol = tol)
+    HYPRE.solve!(ilu, x_h, A_h, b_h)
+    copy!(x, x_h)
+    # Test result with direct solver
+    @test x ≈ A \ b atol=tol
+    # Test without passing initial guess
+    x_h = HYPRE.solve(ilu, A_h, b_h)
+    copy!(x, x_h)
+    @test x ≈ A \ b atol=tol
+
+    # Use as preconditioner to PCG
+    precond = HYPRE.ILU()
+    pcg = HYPRE.PCG(; Tol = tol, Precond = precond)
+    x_h = HYPREVector(zeros(100))
+    HYPRE.solve!(pcg, x_h, A_h, b_h)
+    copy!(x, x_h)
+    # Test result with direct solver
+    @test x ≈ A \ b atol=tol
+    # Test without passing initial guess
+    x_h = HYPRE.solve(pcg, A_h, b_h)
+    copy!(x, x_h)
+    @test x ≈ A \ b atol=tol
+end
+
+
 @testset "(ParCSR)ParaSails" begin
     # Solver constructor and options
     @test_throws(