Excited states #269
Excited states #269
Conversation
Benchmark resultJudge resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsA ratio greater than
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoTargetBaselineTarget resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoBaseline resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoRuntime information
Architecture: x86_64
Benchmark resultJudge resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsA ratio greater than
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoTargetBaselineTarget resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoBaseline resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoRuntime information
Architecture: x86_64
Benchmark resultJudge resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsA ratio greater than
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoTargetBaselineTarget resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoBaseline resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoRuntime information
Architecture: x86_64
Benchmark resultJudge resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsA ratio greater than
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoTargetBaselineTarget resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoBaseline resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoRuntime information
Architecture: x86_64
Benchmark resultJudge resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsA ratio greater than
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoTargetBaselineTarget resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoBaseline resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoRuntime information
Architecture: x86_64
Benchmark resultJudge resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsA ratio greater than
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoTargetBaselineTarget resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoBaseline resultBenchmark Report for /home/runner/work/Rimu.jl/Rimu.jlJob Properties
ResultsBelow is a table of this job's results, obtained by running the benchmarks.
Benchmark Group ListHere's a list of all the benchmark groups executed by this job:
Julia versioninfoRuntime information
Architecture: x86_64 |
There was a problem hiding this comment.
Thanks for the nice PR!
I like the idea of orthogonalising the state before the spawning step. It makes life easier for the norm control machinery. It does mean, however, that we cannot rely on the instantaneous states being exactly orthogonal.
It would be great to have this work smoothly with ProjectorMonteCarloProblem and init! to create an appropriate set of starting vectors. I think currently the number of permissible spectral states is still hardcoded to be one (see line 39 of "pmc_simulation.jl").
It would further be good to make it easy to pass suitable initial state. It would also be super convenient to have suitable initial state generated automatically if no starting vectors (or an insufficient number of them) are provided.
A possible strategy is the following: Say only a single address is available (e.g. from starting_address(ham). Generate further addresses with offdiagonals until you have at least num_spectral_states of them and then sort them by the value of diagonal_element. Set as starting vectors vectors with one of these addresses each. This should work well for diagonally dominated Hamiltonians.
Yes it would be good to be able to automatically generate starting vectors, thanks, I'll look into doing that. |
Rewrite orthogonalization of vectors, make orthogonalization_interval a keyword argument, fix spelling Co-authored-by: Joachim Brand <j.brand@massey.ac.nz> Co-authored-by: mtsch <matijacufar@gmail.com>
Pull Request Test Coverage Report for Build 10866553544Details
💛 - Coveralls |
There was a problem hiding this comment.
Looks great! The next step would be to add some tests. It doesn't have to be anything complicated, just something that makes sure the code works correctly. I suggest something like setting up a small Hamiltonian (maybe HubbardReal1D(BoseFS(1,1,1,1,1))), solving it with FCIQMC with the style set to IsDeterministic(), then comparing the resulting energies and vectors with the result that comes out of solveing an ExactDiagonalizationProblem.
If you run into issues getting this to work feel free to ask here or on Zulip!
| basis = build_basis(ham; cutoff=diagonal_element(ham, starting_address(ham)) + 5) | ||
| basis = sort!(basis; by=x -> diagonal_element(ham, x))[1:2*n_spectral] | ||
| trunc = BasisSetRepresentation(ham, basis; filter=Returns(false), sizelim=Inf) | ||
| vecs = eigvecs(Matrix(trunc)) | ||
| if threading | ||
| eigs = [PDVec(zip(basis, 10*vecs[:,i]); style, initiator) for i in 1:n_spectral] | ||
| else | ||
| if initiator == NonInitiator() | ||
| eigs = [DVec(zip(basis, 10*vecs[:,i]); style) for i in 1:n_spectral] | ||
| else | ||
| eigs = [InitiatorDVec(zip(basis,10*vecs[:,i]);style,initiator) for i in 1:n_spectral] | ||
| end | ||
| end |
There was a problem hiding this comment.
I like this idea a lot, but we should be careful. It's possible to setup a Hamiltonian that's huge, but still has all diagonal elements below diagonal_element(ham, starting_address(ham)) + 5. Running this on it will not work.
There was a problem hiding this comment.
As discussed, I have the same concern. I would suggest to not use an energy cutoff, but instead implement a new option for build_basis where you specify a lower bound, e.g. called mininum_size, that allows build_basis to stop generating new addresses once the mininum_size has been reached.
This will require modifying the code of build_basis, but should be fairly easy to implement.
There was a problem hiding this comment.
Very nice PR! I've found two issues:
- The new functionality is not yet well documented. I suggest documenting the option of passing a
spectral_strategyin the docstring ofProjectorMonteCarloProblem, and adding a reference toGramSchmidtin the docstring ofSpectralStrategy. I'm not completely sure whether the latter should be exported. The used would not need to use it directly in most cases, but I'd like to properly document the available options (which currently is onlyGramSchmidt. - I've checked that this works nicely when generating excited state with a single replica. However, combining replicas with excited states throws an error:
julia> h = HubbardReal1D(near_uniform(BoseFS{10,10}))
HubbardReal1D(fs"|1 1 1 1 1 1 1 1 1 1⟩"; u=1.0, t=1.0)
julia> spectral_strategy = GramSchmidt(3)
GramSchmidt(3; orthogonalization_interval=1)
julia> p = ProjectorMonteCarloProblem(h; spectral_strategy, replica_strategy=AllOverlaps(2), last_step=2_000)
ProjectorMonteCarloProblem with 2 replica(s) and 3 spectral state(s):
Algorithm: FCIQMC(DoubleLogUpdate{Int64}(1000, 0.08, 0.0016), ConstantTimeStep())
Hamiltonian: HubbardReal1D(fs"|1 1 1 1 1 1 1 1 1 1⟩"; u=1.0, t=1.0)
Style: IsDynamicSemistochastic{Float64,ThresholdCompression,DynamicSemistochastic}()
Initiator: NonInitiator()
Threading: false
Simulation Plan: SimulationPlan(starting_step=0, last_step=2000, walltime=Inf)
Replica Strategy: AllOverlaps{2, 0, Tuple{}, true}(())
Reporting Strategy: ReportDFAndInfo
reporting_interval: Int64 1
info_interval: Int64 100
io: Base.TTY
writeinfo: Bool false
Post Step Strategy: ()
Spectral Strategy: GramSchmidt(3; orthogonalization_interval=1)
Maxlength: 2100
Metadata: OrderedCollections.LittleDict("display_name" => "PMCSimulation")
Random Seed: 1943348664599472552
julia> sim = solve(p)
ERROR: ArgumentError: Collection has multiple elements, must contain exactly 1 element
Stacktrace:
[1] only
@ ./iterators.jl:1527 [inlined]
[2] #34
@ ~/git/code/Rimu/src/strategies_and_params/replicastrategy.jl:136 [inlined]
[3] ntuple
@ ./ntuple.jl:49 [inlined]
[4] replica_stats(rs::AllOverlaps{…}, spectral_states::Tuple{…})
@ Rimu ~/git/code/Rimu/src/strategies_and_params/replicastrategy.jl:136
[5] step!(sm::Rimu.PMCSimulation)
@ Rimu ~/git/code/Rimu/src/pmc_simulation.jl:283
[6] macro expansion
@ ~/git/code/Rimu/src/pmc_simulation.jl:428 [inlined]
[7] macro expansion
@ ~/.julia/packages/ProgressLogging/6KXlp/src/ProgressLogging.jl:328 [inlined]
[8]
@ Rimu ~/git/code/Rimu/src/pmc_simulation.jl:422
[9] solve!
@ ~/git/code/Rimu/src/pmc_simulation.jl:337 [inlined]
[10] solve(args::ProjectorMonteCarloProblem{2, 3})
@ CommonSolve ~/.julia/packages/CommonSolve/JfpfI/src/CommonSolve.jl:23
[11] top-level scope
@ REPL[77]:1
Some type information was truncated. Use `show(err)` to see complete types.
julia> Co-authored-by: Joachim Brand <joachim.brand@gmail.com>
Added a new method for the
advancefunction that takes multiple states and orthogonalises them. Theorthogonalisation_intervalparameter of theGramSchmidtSpectralStrategydetermines how often the states are orthogonalised.