Solve interface 3

.github/workflows/actions.yml

@@ -8,7 +8,7 @@ jobs:
     runs-on: ${{ matrix.os }}
     strategy:
       matrix:
-        julia-version: ['1', 'nightly']
+        julia-version: ['1', 'nightly', '1.9']
         julia-arch: [x64]
         os: [ubuntu-latest]
       fail-fast: false

Project.toml

@@ -1,12 +1,13 @@
 name = "Rimu"
 uuid = "c53c40cc-bd84-11e9-2cf4-a9fde2b9386e"
 authors = ["Joachim Brand <j.brand@massey.ac.nz>"]
-version = "0.11.2"
+version = "0.12.0"
 
 [deps]
 Arrow = "69666777-d1a9-59fb-9406-91d4454c9d45"
 BSON = "fbb218c0-5317-5bc6-957e-2ee96dd4b1f0"
 Combinatorics = "861a8166-3701-5b0c-9a16-15d98fcdc6aa"
+CommonSolve = "38540f10-b2f7-11e9-35d8-d573e4eb0ff2"
 ConsoleProgressMonitor = "88cd18e8-d9cc-4ea6-8889-5259c0d15c8b"
 DataFrames = "a93c6f00-e57d-5684-b7b6-d8193f3e46c0"
 DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
@@ -21,6 +22,7 @@ MPI = "da04e1cc-30fd-572f-bb4f-1f8673147195"
 MacroTools = "1914dd2f-81c6-5fcd-8719-6d5c9610ff09"
 Measurements = "eff96d63-e80a-5855-80a2-b1b0885c5ab7"
 MonteCarloMeasurements = "0987c9cc-fe09-11e8-30f0-b96dd679fdca"
+NamedTupleTools = "d9ec5142-1e00-5aa0-9d6a-321866360f50"
 OrderedCollections = "bac558e1-5e72-5ebc-8fee-abe8a469f55d"
 Parameters = "d96e819e-fc66-5662-9728-84c9c7592b0a"
 ProgressLogging = "33c8b6b6-d38a-422a-b730-caa89a2f386c"
@@ -35,32 +37,41 @@ StatsBase = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91"
 StrFormat = "b5087856-efa9-5a6d-8e6f-97303a7af894"
 StrLiterals = "68059f60-971f-57ff-a2d0-18e7de9ccc84"
 TOML = "fa267f1f-6049-4f14-aa54-33bafae1ed76"
+Tables = "bd369af6-aec1-5ad0-b16a-f7cc5008161c"
 TerminalLoggers = "5d786b92-1e48-4d6f-9151-6b4477ca9bed"
 TupleTools = "9d95972d-f1c8-5527-a6e0-b4b365fa01f6"
 VectorInterface = "409d34a3-91d5-4945-b6ec-7529ddf182d8"
 
 [weakdeps]
+Arpack = "7d9fca2a-8960-54d3-9f78-7d1dccf2cb97"
+IterativeSolvers = "42fd0dbc-a981-5370-80f2-aaf504508153"
 KrylovKit = "0b1a1467-8014-51b9-945f-bf0ae24f4b77"
 
 [extensions]
+ArpackExt = "Arpack"
+IterativeSolversExt = "IterativeSolvers"
 KrylovKitExt = "KrylovKit"
 
 [compat]
+Arpack = "0.5"
 Arrow = "1.5, 2"
 BSON = "0.3"
 Combinatorics = "1"
+CommonSolve = "0.2.4"
 ConsoleProgressMonitor = "0.1"
 DataFrames = "1"
 DataStructures = "0.18"
 Distributions = "0.25"
 FFTW = "1"
 Folds = "0.2"
 HypergeometricFunctions = "0.3"
+IterativeSolvers = "0.9"
 KrylovKit = "0.6, 0.7"
 MPI = "0.20"
 MacroTools = "0.5"
 Measurements = "2"
 MonteCarloMeasurements = "1.0.2"
+NamedTupleTools = "v0.14"
 OrderedCollections = "1"
 Parameters = "0.12"
 ProgressLogging = "0.1.3"
@@ -73,12 +84,15 @@ StatsBase = "0.33, 0.34"
 StrFormat = "1"
 StrLiterals = "1"
 TOML = "1"
+Tables = "1.9"
 TerminalLoggers = "0.1.4"
 TupleTools = "1"
 VectorInterface = "0.2, 0.3, 0.4"
-julia = "1.7"
+julia = "1.9"
 
 [extras]
+Arpack = "7d9fca2a-8960-54d3-9f78-7d1dccf2cb97"
+IterativeSolvers = "42fd0dbc-a981-5370-80f2-aaf504508153"
 KrylovKit = "0b1a1467-8014-51b9-945f-bf0ae24f4b77"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 

benchmark/benchmarks.jl

@@ -68,9 +68,9 @@ const SUITE = @benchmarkset "Rimu" begin
             ham = BoseHubbardMom1D2C(addr, v=0.1)
             dv = PDVec(addr => 1.0f0; style=IsDynamicSemistochastic{Float32}())
             s_strat = DoubleLogUpdate(targetwalkers=10_000)
-            replica = AllOverlaps(2, ntuple(i -> G2Correlator(i - 1), 7))
+            replica_strategy = AllOverlaps(2; operator = ntuple(i -> G2Correlator(i - 1), 7))
 
-            lomc!(ham, dv; s_strat, replica, laststep=2000)
+            lomc!(ham, dv; s_strat, replica_strategy, laststep=2000)
         end seconds=150
 
         @case "(50, 50) Real space" begin

docs/make.jl

@@ -53,6 +53,7 @@ makedocs(;
         "Guide" => "index.md",
         "Examples" => EXAMPLES_PAIRS[sortperm(EXAMPLES_NUMS)],
         "User documentation" => [
+            "Exact Diagonalization" => "exactdiagonalization.md",
             "StatsTools" => "statstools.md",
             "Using MPI" => "mpi.md",
         ],
@@ -74,6 +75,7 @@ makedocs(;
     authors="Joachim Brand <j.brand@massey.ac.nz>",
     checkdocs=:exports,
     doctest=false, # Doctests are done while testing.
+    # warnonly = true, # should be diabled for a release
 )
 
 deploydocs(

docs/src/API.md

@@ -6,6 +6,10 @@ Modules = [Rimu]
 ```
 # Reexported Submodules
 
+## ExactDiagonalization
+
+See [Exact Diagonalization](@ref)
+
 ## Interfaces
 
 See [Module `Interfaces`](@ref)
@@ -28,7 +32,7 @@ See [Module `DictVectors`](@ref)
 
 ## StatsTools
 
-See [Module `Rimu/StatsTools`](@ref)
+See [Module `StatsTools`](@ref)
 
 ## RMPI
 

docs/src/exactdiagonalization.md

@@ -0,0 +1,40 @@
+# Exact Diagonalization
+
+The main functionality of Rimu for exact diagonalization is contained in the module `ExactDiagonalization`.
+
+```@docs
+ExactDiagonalization
+```
+
+## `ExactDiagonalizationProblem`
+
+```@docs
+ExactDiagonalizationProblem
+solve(::ExactDiagonalizationProblem)
+init(::ExactDiagonalizationProblem)
+```
+
+## Solver algorithms
+
+```@docs
+KrylovKitSolver
+LinearAlgebraSolver
+ArpackSolver
+LOBPCGSolver
+```
+
+## Converting a Hamiltonian in to a matrix
+
+```@docs
+BasisSetRepresentation
+build_basis
+Matrix
+sparse
+```
+
+## Deprecated
+```@docs
+BasisSetRep
+```
+
+

docs/src/hamiltonians.md

@@ -5,41 +5,12 @@ physical models of interest. These are organised by means of an interface
 around the abstract type [`AbstractHamiltonian`](@ref), in the spirit of the
 `AbstractArray` interface as discussed in the [Julia Documentation](https://docs.julialang.org/en/v1/manual/interfaces/).
 
+The Hamiltonians can be used for projector quantum Monte Carlo with [`ProjectorMonteCarloProblem`](@ref) or for exact diagonalization with [`ExactDiagonalizationProblem`](@ref), see [Exact Diagonalization](@ref).
+
 ```@docs
 Hamiltonians
 ```
 
-## Usage with FCIQMC and exact diagonalisation
-
-In order to define a specific model Hamiltonian with relevant parameters
-for the model, instantiate the model like this in the input file:
-
-```julia-repl
-hubb = HubbardReal1D(BoseFS((1,2,0,3)); u=1.0, t=1.0)
-```
-
-The Hamiltonian `hubb` is now ready to be used for FCIQMC in [`lomc!`](@ref)
-and for exact diagonalisation with [`KrylovKit.jl`](https://github.com/Jutho/KrylovKit.jl) directly, or after
-transforming into a sparse matrix first with
-```julia-repl
-using SparseArrays
-sh = sparse(hubb)
-```
-or into a full matrix with
-```julia-repl
-using LinearAlgebra
-fh = Matrix(hubb)
-```
-This functionality relies on
-```@docs
-Hamiltonians.BasisSetRep
-sparse
-Matrix
-```
-If only the basis is required and not the matrix representation it is more efficient to use
-```@docs
-Hamiltonians.build_basis
-```
 
 ## Model Hamiltonians
 

docs/src/index.md

@@ -37,7 +37,7 @@ Didier Adrien, Chris Scott (NeSI), Alexander Pletzer (NeSI).
 
 `Rimu` is a registered package and can be installed with the package manager.
 Hit the `]` key at the Julia REPL to get into `Pkg` mode and type
-```julia-repl
+```
 pkg> add Rimu
 ```
 Alternatively, use
@@ -53,7 +53,7 @@ etc.,
 clone the git repository with `git clone` to a convenient location, e.g.
 `~/mygitpackagefolder/`. Then
 hit the `]` key at the Julia REPL to get into `Pkg` mode and type
-```julia-repl
+```
 pkg> develop ~/mygitpackagefolder/rimu.jl
 ```
 where the file path has to be adjusted to the location of the cloned git
@@ -116,6 +116,10 @@ be stored in a distributed fashion among the MPI ranks and only communicated bet
 when necessary. Additional MPI-related functionality is provided by the module [`RMPI`](@ref
 Rimu.RMPI).
 
+## Compatibility
+
+We recommend using `Rimu` with the latest Julia release version. Rimu requires at least julia `v1.9`.
+
 ## References
 The code implements the FCIQMC algorithm originally described in
 - "Fermion Monte Carlo without fixed nodes: A game of life, death, and annihilation in Slater determinant space", G. H. Booth, A. J. W. Thom, A. Alavi, [*J. Chem. Phys.* **131**, 054106 (2009)](https://doi.org/10.1063/1.3193710).

docs/src/statstools.md

@@ -1,6 +1,6 @@
-# Module `Rimu/StatsTools`
+# Module `StatsTools`
 
-The  module `Rimu/StatsTools` contains helper function for analysis and post
+The  module `StatsTools` contains helper function for analysis and post
 processing of Monte Carlo data.
 
 ## Blocking analysis

docs/src/testing.md

@@ -12,17 +12,3 @@ More tests should be added over time to test core functionality of the code. To
 
 GitHub Actions are set up to run the test script automatically on the GitHub cloud server every time a new commit to the master branch is pushed to the server. The setup for this to happen is configured in the file
 `actions.yml` in the `Rimu/.github/workflows` folder.
-
-## Modifying the `Project.toml` file
-
-In order for the testing code to be able to run on the cloud server, external packages that are accessed in the
-code with `using` or `import` need to be installed first.
-This is done in the script `actions.yml` via the package manager, based on the information contained in the file
-`test/Project.toml`. More packages can be added to this file using the package manager in the following way: Say we want to install the package `DelimitedFiles`. At the Julia REPL, type the following:
-```julia-repl
-julia> cd("test")
-julia> ]
-(v1.0) pkg> activate .
-(test) pkg> add DelimitedFiles
-```
-This will a new line to the file `Project.toml` with the name of the package and the corresponding uuid. When Pipelines now runs the commands in yml script, it will install the package `DelimitedFiles` before running the `runtest.jl` script.

ext/ArpackExt.jl

@@ -0,0 +1,77 @@
+module ArpackExt
+
+using Arpack: Arpack, eigs
+using CommonSolve: CommonSolve, solve
+using NamedTupleTools: delete
+using LinearAlgebra: norm
+
+using Rimu: Rimu, DVec, replace_keys
+using Rimu.ExactDiagonalization: ArpackSolver, MatrixEDSolver,
+    LazyDVecs, EDResult
+
+struct ArpackConvergenceInfo
+    converged::Int
+    numiter::Int
+    numops::Int
+    residual::Vector{Float64}
+end
+function Base.show(io::IO, info::ArpackConvergenceInfo)
+    print(io, "converged = $(info.converged), ")
+    print(io, "numiter = $(info.numiter), ")
+    print(io, "numops = $(info.numops), ")
+    print(io, "residual norm = ")
+    show(io, norm(info.residual))
+end
+
+function CommonSolve.solve(s::S; kwargs...
+) where {S<:MatrixEDSolver{<:ArpackSolver}}
+    # combine keyword arguments and set defaults for `howmany` and `which`
+    kw_nt = (; howmany=1, which=:SR, s.kw_nt..., kwargs...)
+    # check if universal keyword arguments are present
+    kw_nt = replace_keys(kw_nt, (:abstol=>:tol, :maxiters=>:maxiter))
+    verbose = get(kw_nt, :verbose, false)
+    kw_nt = delete(kw_nt, (:verbose,))
+
+    # Remove the `howmany` key from the kwargs.
+    kw_nt = (; nev=kw_nt.howmany, kw_nt..., ritzvec=true)
+    kw_nt = delete(kw_nt, (:howmany,))
+    howmany = kw_nt.nev
+
+    # set up the starting vector
+    v0 = if isnothing(s.v0)
+        zeros((0,))
+    else
+        # convert v0 to a DVec to use it like a dictionary
+        dvec = DVec(s.v0)
+        [dvec[a] for a in s.basissetrep.basis]
+    end
+    # solve the problem
+    vals, vec_matrix, nconv, niter, nmult, resid = eigs(s.basissetrep.sparse_matrix; v0, kw_nt...)
+
+    verbose && @info "Arpack.eigs: $nconv converged out of $howmany requested eigenvalues,"*
+        " $niter iterations," *
+        " $nmult matrix vector multiplications, norm of residual = $(norm(resid))"
+    success = nconv ≥ howmany
+    # vecs = [view(vec_matrix, :, i) for i in 1:length(vals)] # convert to array of vectors
+    coefficient_vectors = eachcol(vec_matrix)
+    vectors = LazyDVecs(coefficient_vectors, s.basissetrep.basis)
+    info = ArpackConvergenceInfo(nconv, niter, nmult, resid)
+    if !success
+        @warn "Arpack.eigs did not converge for all requested eigenvalues:" *
+              " $nconv converged out of $howmany requested value(s)."
+    end
+    return EDResult(
+        s.algorithm,
+        s.problem,
+        vals,
+        vectors,
+        coefficient_vectors,
+        s.basissetrep.basis,
+        info,
+        howmany,
+        vec_matrix,
+        success
+    )
+end
+
+end # module