Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

feat: add PETScSNES #482

Merged
merged 6 commits into from
Oct 27, 2024
Merged

feat: add PETScSNES #482

merged 6 commits into from
Oct 27, 2024

Conversation

avik-pal
Copy link
Member

@avik-pal avik-pal commented Oct 25, 2024

PETSc 0.2 resolves old compatibility issues we were having with wrapping this library.

Warning

I am still seeing occasional segfaults with PETSc. Very common in 1.11. And for some reason PETSc hijacks Ctrl + C and trying to use it leads to a segfault even when no code is running

fixes #160

TODOs

  • Sparsity Handling?
  • Testing
  • Rewrite one of the examples from PETSc.jl

@avik-pal
Copy link
Member Author

This 150 lines of code https://github.com/JuliaParallel/PETSc.jl/blob/main/examples/SNES_ex2.jl translates to 30 lines in nonlinearsolve with automatic sparsity detection

# This implements src/snes/examples/tutorials/ex2.c from PETSc using the PETSc.jl package, using SNES
#
# This solves the equations sequentially
# 
# Newton method to solve u'' + u^{2} = f, sequentially.
using NonlinearSolve, PETSc, LinearAlgebra, SparseConnectivityTracer

n = 21
u0 = fill(0.5, n)

function form_residual!(resid, x, _)
    n = length(x)
    xp = LinRange(0.0, 1.0, n)
    F = 6xp .+ (xp .+ 1e-12) .^ 6

    dx = 1 / (n - 1)
    resid[1] = x[1]
    for i in 2:(n - 1)
        resid[i] = (x[i - 1] - 2x[i] + x[i + 1]) / dx^2 + x[i] * x[i] - F[i]
    end
    resid[n] = x[n] - 1

    return
end

nlfunc = NonlinearFunction{true}(form_residual!; sparsity = TracerSparsityDetector())
nlprob = NonlinearProblem(nlfunc, u0)

solve(nlprob, NewtonRaphson())
solve(nlprob, PETScSNES())

@avik-pal
Copy link
Member Author

PETSc is faster than us in the sparse case at least on my local computer. Once all the refactoring is done we will have to revisit the benchmarks

@avik-pal
Copy link
Member Author

Looking at the profile all of the time goes in the UMFPACK

Comment on lines 61 to 84
## Runtimes

### Dense Jacobian

```@example snes_ex2
@benchmark solve($(nlprob_dense), $(NewtonRaphson()); abstol = 1e-8)
nothing # hide
```

```@example snes_ex2
@benchmark solve($(nlprob_dense), $(PETScSNES()); abstol = 1e-8)
nothing # hide
```

### Sparse Jacobian

```@example snes_ex2
@benchmark solve($(nlprob_sparse), $(NewtonRaphson()); abstol = 1e-8)
nothing # hide
```

```@example snes_ex2
@benchmark solve($(nlprob_sparse), $(PETScSNES()); abstol = 1e-8)
nothing # hide
Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

What's the result? Just curious

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

How optimized is our sparse solver selection for LinearSolve?

Copy link
Member

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Relatively optimal? What sparse linear solver is PETSc using here? Try all of UMFPACK, KLU, and MKLPardiso. Other thing to try is using MKL with Pardiso. Then note that there's a new sparse solver ParU that @rayegun is getting wrapped SciML/LinearSolve.jl#394.

One thing to double check is what the extra time is spent in. Is it only doing one symbolic factorization?

Copy link
Member Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Let's get this in, since it is non-blocking. I anyways need to redo benchmarks for the paper, so I can revisit it after the splitting is done.

@ChrisRackauckas ChrisRackauckas merged commit 28c0189 into master Oct 27, 2024
16 of 39 checks passed
@ChrisRackauckas ChrisRackauckas deleted the ap/petsc branch October 27, 2024 10:02
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

Successfully merging this pull request may close these issues.

An interface to PETSc SNES nonlinear solve would be nice
2 participants