From 0b7d7ae4da131906cec70954cad964ea0a48b6c6 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Sat, 21 Dec 2024 23:07:50 +1300
Subject: [PATCH 01/17] change maxlength to max_length in
 ProjectorMonteCarloProblem

---
 src/Hamiltonians/HubbardRealSpace.jl  |  2 +-
 src/fciqmc.jl                         |  6 ++--
 src/lomc.jl                           |  6 ++--
 src/pmc_simulation.jl                 |  4 +--
 src/projector_monte_carlo_problem.jl  | 49 ++++++++++++++++-----------
 src/qmc_states.jl                     |  4 +--
 test/lomc.jl                          |  2 +-
 test/projector_monte_carlo_problem.jl |  2 +-
 8 files changed, 42 insertions(+), 33 deletions(-)

diff --git a/src/Hamiltonians/HubbardRealSpace.jl b/src/Hamiltonians/HubbardRealSpace.jl
index 533204da2..8b7488d5b 100644
--- a/src/Hamiltonians/HubbardRealSpace.jl
+++ b/src/Hamiltonians/HubbardRealSpace.jl
@@ -266,7 +266,7 @@ function Base.show(io::IO, h::HubbardRealSpace{C}) where C
         println(io, "  u = ", Float64.(h.u), ",")
     end
     !isnothing(h.v) && println(io, "  v = ", Float64.(h.v), ",")
-    println(io, ")")
+    print(io, ")")
 end
 
 # Overload equality due to stored potential energy arrays.
diff --git a/src/fciqmc.jl b/src/fciqmc.jl
index 8b630edfd..f2d22361c 100644
--- a/src/fciqmc.jl
+++ b/src/fciqmc.jl
@@ -177,11 +177,11 @@ function advance!(algorithm::FCIQMC, report, state::ReplicaState, s_state::Singl
         end
         return false
     end
-    if len > state.maxlength[]
+    if len > state.max_length[]
         if length(state.spectral_states) > 1
-            @error "`maxlength` reached in single state $(s_state.id). Aborting."
+            @error "`max_length` reached in single state $(s_state.id). Aborting."
         else
-            @error "`maxlength` reached. Aborting."
+            @error "`max_length` reached. Aborting."
         end
         return false
     end
diff --git a/src/lomc.jl b/src/lomc.jl
index e88eb66b5..cc3a693c0 100644
--- a/src/lomc.jl
+++ b/src/lomc.jl
@@ -156,7 +156,7 @@ function lomc!(
         replica_strategy = replica,
         reporting_strategy = r_strat,
         post_step_strategy = post_step,
-        maxlength,
+        max_length = maxlength,
         metadata,
         display_name = name,
         random_seed = false
@@ -201,7 +201,7 @@ function lomc!(state::ReplicaState, df=DataFrame(); laststep=0, name="lomc!", me
     if !iszero(laststep)
         state = @set state.simulation_plan.last_step = laststep
     end
-    @unpack spectral_states, maxlength, step, simulation_plan,
+    @unpack spectral_states, max_length, step, simulation_plan,
     reporting_strategy, post_step_strategy, replica_strategy = state
     first_replica = first(state) # SingleState
     @unpack hamiltonian = first_replica
@@ -213,7 +213,7 @@ function lomc!(state::ReplicaState, df=DataFrame(); laststep=0, name="lomc!", me
         replica_strategy,
         reporting_strategy,
         post_step_strategy,
-        maxlength=maxlength[],
+        max_length=max_length[],
         simulation_plan,
         metadata,
         display_name=name,
diff --git a/src/pmc_simulation.jl b/src/pmc_simulation.jl
index 60178e1ab..a5ee8da1f 100644
--- a/src/pmc_simulation.jl
+++ b/src/pmc_simulation.jl
@@ -80,7 +80,7 @@ function PMCSimulation(problem::ProjectorMonteCarloProblem; copy_vectors=true)
     @unpack algorithm, hamiltonian, start_at, style, threading, simulation_plan,
         replica_strategy, initial_shift_parameters,
         reporting_strategy, post_step_strategy,
-        maxlength, metadata, initiator, random_seed, spectral_strategy, minimum_size = problem
+        max_length, metadata, initiator, random_seed, spectral_strategy, minimum_size = problem
 
     reporting_strategy = refine_reporting_strategy(reporting_strategy)
 
@@ -137,7 +137,7 @@ function PMCSimulation(problem::ProjectorMonteCarloProblem; copy_vectors=true)
     # set up the initial state
     state = ReplicaState(
         spectral_states,
-        Ref(maxlength),
+        Ref(max_length),
         Ref(simulation_plan.starting_step),
         simulation_plan,
         reporting_strategy,
diff --git a/src/projector_monte_carlo_problem.jl b/src/projector_monte_carlo_problem.jl
index c88de80ce..28d4531c1 100644
--- a/src/projector_monte_carlo_problem.jl
+++ b/src/projector_monte_carlo_problem.jl
@@ -54,7 +54,7 @@ the [`FCIQMC`](@ref) algorithm.
   simulations, see [`ReplicaStrategy`](@ref).
 - `n_spectral = 1`: Number of targeted spectral states. Set `n_spectral > 1` to find excited
   states.
-- `spectral_strategy = GramSchmidt(n_spectral)`: The [`SpectralStrategy`](@ref) used for 
+- `spectral_strategy = GramSchmidt(n_spectral)`: The [`SpectralStrategy`](@ref) used for
   orthogonalizing spectral states.
 
 # Example
@@ -90,7 +90,7 @@ julia> size(DataFrame(simulation))
     Hamiltonian and the starting vectors.
 - `initial_shift_parameters`: Initial shift parameters or collection of initial shift
     parameters. Overrides `shift` if provided.
-- `maxlength = 2 * target_walkers + 100`: Maximum length of the vectors.
+- `max_length = 2 * target_walkers + 100`: Maximum length of the vectors.
 - `display_name = "PMCSimulation"`: Name displayed in progress bar (via `ProgressLogging`).
 - `metadata`: User-supplied metadata to be added to the report. Must be an iterable of
   pairs or a `NamedTuple`, e.g. `metadata = ("key1" => "value1", "key2" => "value2")`.
@@ -119,7 +119,7 @@ struct ProjectorMonteCarloProblem{N,S} # is not type stable but does not matter
     reporting_strategy::ReportingStrategy
     post_step_strategy::Tuple
     spectral_strategy::SpectralStrategy{S}
-    maxlength::Int
+    max_length::Int
     metadata::LittleDict{String,String} # user-supplied metadata + display_name
     random_seed::Union{Nothing,UInt64}
     minimum_size::Int
@@ -129,19 +129,19 @@ function Base.show(io::IO, p::ProjectorMonteCarloProblem)
     nr = num_replicas(p)
     ns = num_spectral_states(p)
     println(io, "ProjectorMonteCarloProblem with $nr replica(s) and $ns spectral state(s):")
-    isnothing(p.algorithm) || println(io, "  Algorithm: ", p.algorithm)
-    println(io, "  Hamiltonian: ", p.hamiltonian)
-    println(io, "  Style: ", p.style)
-    println(io, "  Initiator: ", p.initiator)
-    println(io, "  Threading: ", p.threading)
-    println(io, "  Simulation Plan: ", p.simulation_plan)
-    println(io, "  Replica Strategy: ", p.replica_strategy)
-    print(io, "  Reporting Strategy: ", p.reporting_strategy)
-    println(io, "  Post Step Strategy: ", p.post_step_strategy)
-    println(io, "  Spectral Strategy: ", p.spectral_strategy)
-    println(io, "  Maxlength: ", p.maxlength)
-    println(io, "  Metadata: ", p.metadata)
-    print(io, "  Random Seed: ", p.random_seed)
+    isnothing(p.algorithm) || println(io, "  algorithm = ", p.algorithm)
+    println(io, "  hamiltonian = ", p.hamiltonian)
+    println(io, "  style = ", p.style)
+    println(io, "  initiator = ", p.initiator)
+    println(io, "  threading = ", p.threading)
+    println(io, "  simulation_plan = ", p.simulation_plan)
+    println(io, "  replica_strategy = ", p.replica_strategy)
+    print(io, "  reporting_strategy = ", p.reporting_strategy)
+    println(io, "  post_step_strategy = ", p.post_step_strategy)
+    println(io, "  spectral_strategy = ", p.spectral_strategy)
+    println(io, "  max_length = ", p.max_length)
+    println(io, "  metadata = ", p.metadata)
+    print(io, "  random_seed = ", p.random_seed)
 end
 
 
@@ -172,7 +172,8 @@ function ProjectorMonteCarloProblem(
     n_spectral = 1,
     spectral_strategy = GramSchmidt(n_spectral),
     minimum_size = 2*num_spectral_states(spectral_strategy),
-    maxlength = nothing,
+    max_length = nothing,
+    maxlength = nothing, # deprecated
     metadata = nothing,
     display_name = "PMCSimulation",
     random_seed = true
@@ -219,10 +220,18 @@ function ProjectorMonteCarloProblem(
     end
 
     shift_strategy = algorithm.shift_strategy
-    if isnothing(maxlength)
-        maxlength = round(Int, 2 * abs(shift_strategy.target_walkers) + 100)
+
+    if isnothing(maxlength) # deprecated
+        if isdefined(shift_strategy, :target_walkers)
+            maxlength = round(Int, 2 * abs(shift_strategy.target_walkers) + 100)
+        else
+            maxlength = round(Int, 2 * abs(target_walkers) + 100)
+        end
         # padding for small walkernumbers
+    else
+        @warn "The keyword argument `maxlength` is deprecated. Use `max_length` instead."
     end
+    max_length = isnothing(max_length) ? maxlength : max_length
 
     # convert metadata to LittleDict
     report = Report()
@@ -254,7 +263,7 @@ function ProjectorMonteCarloProblem(
         reporting_strategy,
         post_step_strategy,
         spectral_strategy,
-        maxlength,
+        max_length,
         metadata,
         random_seed,
         minimum_size
diff --git a/src/qmc_states.jl b/src/qmc_states.jl
index c115ec1f5..2c280dc84 100644
--- a/src/qmc_states.jl
+++ b/src/qmc_states.jl
@@ -100,7 +100,7 @@ struct ReplicaState{
     PS<:NTuple{<:Any,PostStepStrategy},
 } <: AbstractMatrix{SingleState}
     spectral_states::R
-    maxlength::Ref{Int}
+    max_length::Ref{Int}
     step::Ref{Int}
     simulation_plan::SimulationPlan
     reporting_strategy::RS
@@ -172,7 +172,7 @@ function report_default_metadata!(report::Report, state::ReplicaState)
     report_metadata!(report, "time_step", shift_parameters.time_step)
     report_metadata!(report, "step", state.step[])
     report_metadata!(report, "shift", shift_parameters.shift)
-    report_metadata!(report, "maxlength", state.maxlength[])
+    report_metadata!(report, "max_length", state.max_length[])
     report_metadata!(report, "post_step_strategy", state.post_step_strategy)
     report_metadata!(report, "v_summary", summary(s_state.v))
     report_metadata!(report, "v_type", typeof(s_state.v))
diff --git a/test/lomc.jl b/test/lomc.jl
index 41dcd5d36..987df06a2 100644
--- a/test/lomc.jl
+++ b/test/lomc.jl
@@ -287,7 +287,7 @@ Random.seed!(1234)
         @test all(df.len_5[1:end-1] .≤ 10)
         @test all(df.len_6[1:end-1] .≤ 10)
 
-        state.maxlength[] += 1000
+        state.max_length[] += 1000
         df_cont = lomc!(state).df
         @test size(df_cont, 1) == 100 - size(df, 1)
     end
diff --git a/test/projector_monte_carlo_problem.jl b/test/projector_monte_carlo_problem.jl
index a2aebebb5..f8a492c60 100644
--- a/test/projector_monte_carlo_problem.jl
+++ b/test/projector_monte_carlo_problem.jl
@@ -22,7 +22,7 @@ using OrderedCollections: freeze
     @test sp.shift == diagonal_element(h, starting_address(h))
     @test sp.pnorm == walkernumber(only(state_vectors(simulation)))
     @test sp.pnorm isa Float64
-    @test p.maxlength == 2 * p.algorithm.shift_strategy.target_walkers + 100
+    @test p.max_length == 2 * p.algorithm.shift_strategy.target_walkers + 100
 
     ps = ProjectorMonteCarloProblem(h; initial_shift_parameters=sp, threading=false)
     @test ps.initial_shift_parameters === sp

From ad5fd23461529f1ff88e431d97e80ee43405b9f7 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Sat, 21 Dec 2024 23:36:24 +1300
Subject: [PATCH 02/17] eliminate lomc! from allocation tests

---
 test/allocation_tests.jl | 16 +++++++---------
 1 file changed, 7 insertions(+), 9 deletions(-)

diff --git a/test/allocation_tests.jl b/test/allocation_tests.jl
index 026da2824..a75131017 100644
--- a/test/allocation_tests.jl
+++ b/test/allocation_tests.jl
@@ -90,13 +90,13 @@ end
                     dv = dv_type(addr => 1.0, style=IsDynamicSemistochastic())
                     sizehint!(dv, 500_000)
 
-                    # Warmup for lomc!
-                    _, st = lomc!(
-                        H, dv;
-                        params=RunTillLastStep(;shift=float(diagonal_element(H, addr)), dτ),
-                        maxlength=10_000,
-                        laststep=1,
+                    p = ProjectorMonteCarloProblem(
+                        H;
+                        start_at=dv, last_step=200, time_step=dτ, max_length=10_000
                     )
+                    # Warmup for solve
+                    res = solve!(init(p); last_step=1)
+                    st = res.state
 
                     r = only(only(st.spectral_states).single_states)
 
@@ -111,9 +111,7 @@ end
                     @test allocs_step ≤ 512
 
                     dv = dv_type(addr => 1.0, style=IsDynamicSemistochastic())
-                    allocs_full = @allocated lomc!(
-                        H, dv; dτ, laststep=200, maxlength=10_000
-                    )
+                    allocs_full = @allocated solve(p)
                     @test allocs_full ≤ 1e8 # 100MiB
 
                     # Print out the results to make it easier to find problems.

From 53d4e8d044b111ac29fff8ae0c4bf1f7c6740681 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Sun, 22 Dec 2024 00:00:43 +1300
Subject: [PATCH 03/17] rename walltime to wall_time

---
 src/pmc_simulation.jl                 | 25 +++++++++++++------------
 src/projector_monte_carlo_problem.jl  | 18 +++++++++---------
 test/projector_monte_carlo_problem.jl |  8 ++++----
 3 files changed, 26 insertions(+), 25 deletions(-)

diff --git a/src/pmc_simulation.jl b/src/pmc_simulation.jl
index a5ee8da1f..63a476868 100644
--- a/src/pmc_simulation.jl
+++ b/src/pmc_simulation.jl
@@ -242,7 +242,7 @@ end
 
 Advance the simulation by one step.
 
-Calling [`solve!`](@ref) will advance the simulation until the last step or the walltime is
+Calling [`solve!`](@ref) will advance the simulation until the last step or the wall time is
 exceeded. When completing the simulation without calling [`solve!`](@ref), the simulation
 report needs to be finalised by calling [`Rimu.finalize_report!`](@ref).
 
@@ -300,7 +300,7 @@ end
     solve(::ProjectorMonteCarloProblem)::PMCSimulation
 
 Initialize and solve a [`ProjectorMonteCarloProblem`](@ref) until the last step is completed
-or the walltime limit is reached.
+or the wall time limit is reached.
 
 See also [`init`](@ref), [`solve!`](@ref), [`step!`](@ref), [`Rimu.PMCSimulation`](@ref),
 and [`solve(::ExactDiagonalizationProblem)`](@ref).
@@ -310,16 +310,17 @@ CommonSolve.solve
 """
     solve!(sm::PMCSimulation; kwargs...)::PMCSimulation
 
-Solve a [`Rimu.PMCSimulation`](@ref) until the last step is completed or the walltime limit
+Solve a [`Rimu.PMCSimulation`](@ref) until the last step is completed or the wall time limit
 is reached.
 
-To continue a previously completed simulation, set a new `last_step` or `walltime` using the
-keyword arguments. Optionally, changes can be made to the `replica_strategy`, the
+To continue a previously completed simulation, set a new `last_step` or `wall_time` using
+the keyword arguments. Optionally, changes can be made to the `replica_strategy`, the
 `post_step_strategy`, or the `reporting_strategy`.
 
 # Optional keyword arguments:
 * `last_step = nothing`: Set the last step to a new value and continue the simulation.
-* `walltime = nothing`: Set the allowed walltime to a new value and continue the simulation.
+* `wall_time = nothing`: Set the allowed wall time to a new value and continue the
+    simulation.
 * `reset_time = false`: Reset the `elapsed_time` counter and continue the simulation.
 * `empty_report = false`: Empty the report before continuing the simulation.
 * `replica_strategy = nothing`: Change the replica strategy. Requires the number of replicas
@@ -335,7 +336,7 @@ See also [`ProjectorMonteCarloProblem`](@ref), [`init`](@ref), [`solve`](@ref),
 """
 function CommonSolve.solve!(sm::PMCSimulation;
     last_step = nothing,
-    walltime = nothing,
+    wall_time = nothing,
     reset_time = false,
     replica_strategy=nothing,
     post_step_strategy=nothing,
@@ -351,9 +352,9 @@ function CommonSolve.solve!(sm::PMCSimulation;
         report_metadata!(sm.report, "laststep", last_step)
         reset_flags = true
     end
-    if !isnothing(walltime)
+    if !isnothing(wall_time)
         state = sm.state
-        sm.state = @set state.simulation_plan.walltime = walltime
+        sm.state = @set state.simulation_plan.wall_time = wall_time
         reset_flags = true
     end
     if !isnothing(replica_strategy)
@@ -419,10 +420,10 @@ function CommonSolve.solve!(sm::PMCSimulation;
     name = get_metadata(sm.report, "display_name")
 
     @withprogress name = while !sm.aborted && !sm.success
-        if time() - starting_time > simulation_plan.walltime
+        if time() - starting_time > simulation_plan.wall_time
             sm.aborted = true
-            sm.message = "Walltime limit reached."
-            @warn "Walltime limit reached. Aborting simulation."
+            sm.message = "Wall time limit reached."
+            @warn "Wall time limit reached. Aborting simulation."
         else
             step!(sm)
         end
diff --git a/src/projector_monte_carlo_problem.jl b/src/projector_monte_carlo_problem.jl
index 28d4531c1..bbd7f8e46 100644
--- a/src/projector_monte_carlo_problem.jl
+++ b/src/projector_monte_carlo_problem.jl
@@ -7,21 +7,21 @@ See [`ProjectorMonteCarloProblem`](@ref), [`FCIQMC`](@ref).
 abstract type PMCAlgorithm end
 
 """
-    SimulationPlan(; starting_step = 1, last_step = 100, walltime = Inf)
+    SimulationPlan(; starting_step = 1, last_step = 100, wall_time = Inf)
 Defines the duration of the simulation. The simulation ends when the `last_step` is reached
-or the `walltime` is exceeded.
+or the `wall_time` is exceeded.
 
 See [`ProjectorMonteCarloProblem`](@ref), [`PMCSimulation`](@ref).
 """
 Base.@kwdef struct SimulationPlan
     starting_step::Int = 0
     last_step::Int = 100
-    walltime::Float64 = Inf
+    wall_time::Float64 = Inf
 end
 function Base.show(io::IO, plan::SimulationPlan)
     print(
         io, "SimulationPlan(starting_step=", plan.starting_step,
-        ", last_step=", plan.last_step, ", walltime=", plan.walltime, ")"
+        ", last_step=", plan.last_step, ", wall_time=", plan.wall_time, ")"
     )
 end
 
@@ -75,9 +75,9 @@ julia> size(DataFrame(simulation))
 
 # Further keyword arguments:
 - `starting_step = 1`: Starting step of the simulation.
-- `walltime = Inf`: Maximum time allowed for the simulation.
-- `simulation_plan = SimulationPlan(; starting_step, last_step, walltime)`: Defines the
-    duration of the simulation. Takes precedence over `last_step` and `walltime`.
+- `wall_time = Inf`: Maximum time allowed for the simulation.
+- `simulation_plan = SimulationPlan(; starting_step, last_step, wall_time)`: Defines the
+    duration of the simulation. Takes precedence over `last_step` and `wall_time`.
 - `ζ = 0.08`: Damping parameter for the shift update.
 - `ξ = ζ^2/4`: Forcing parameter for the shift update.
 - `shift_strategy = DoubleLogUpdate(; target_walkers, ζ, ξ)`: How to update the `shift`,
@@ -156,8 +156,8 @@ function ProjectorMonteCarloProblem(
     time_step = 0.01,
     starting_step = 0,
     last_step = 100,
-    walltime = Inf,
-    simulation_plan = SimulationPlan(starting_step, last_step, walltime),
+    wall_time = Inf,
+    simulation_plan = SimulationPlan(starting_step, last_step, wall_time),
     replica_strategy = NoStats(n_replicas),
     targetwalkers = nothing, # deprecated
     target_walkers = 1_000,
diff --git a/test/projector_monte_carlo_problem.jl b/test/projector_monte_carlo_problem.jl
index f8a492c60..5cdde872c 100644
--- a/test/projector_monte_carlo_problem.jl
+++ b/test/projector_monte_carlo_problem.jl
@@ -181,14 +181,14 @@ using Rimu: num_replicas, num_spectral_states
     @test sm.success == true == parse(Bool, (Rimu.get_metadata(sm.report, "success")))
 
     # time out
-    p = ProjectorMonteCarloProblem(h; last_step=500, walltime=1e-3)
+    p = ProjectorMonteCarloProblem(h; last_step=500, wall_time=1e-3)
     sm = init(p)
-    @test_logs (:warn, Regex("(Walltime)")) solve!(sm)
+    @test_logs (:warn, Regex("(Wall time)")) solve!(sm)
     @test sm.success == false
     @test sm.aborted == true
-    @test sm.message == "Walltime limit reached."
+    @test sm.message == "Wall time limit reached."
 
-    sm2 = solve!(sm; walltime=1.0)
+    sm2 = solve!(sm; wall_time=1.0)
     @test sm2 === sm
     @test sm.success == true
     @test sm.state.step[] == 500 == size(sm.df)[1]

From 717baedf443686634427fc436990b60f39b0738b Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Sun, 22 Dec 2024 12:34:26 +1300
Subject: [PATCH 04/17] add deprecation warning for lomc!

---
 src/lomc.jl  | 5 +++++
 test/lomc.jl | 4 ++--
 2 files changed, 7 insertions(+), 2 deletions(-)

diff --git a/src/lomc.jl b/src/lomc.jl
index cc3a693c0..0ab87e568 100644
--- a/src/lomc.jl
+++ b/src/lomc.jl
@@ -118,6 +118,8 @@ function lomc!(
     maxlength = nothing,
     wm = nothing
 )
+    @warn "The use of `lomc!` is deprecated. Use `ProjectorMonteCarloProblem` and `solve` instead."
+
     if !isnothing(wm)
         @warn "The `wm` argument has been removed and will be ignored."
     end
@@ -193,11 +195,14 @@ function lomc!(
 end
 
 function lomc!(::AbstractMatrix, v=nothing; kwargs...)
+    @warn "The use of `lomc!` is deprecated. Use `ProjectorMonteCarloProblem` and `solve` instead."
     throw(ArgumentError("Using lomc! with a matrix is no longer supported. Use `MatrixHamiltonian` instead."))
 end
 
 # methods for backward compatibility
 function lomc!(state::ReplicaState, df=DataFrame(); laststep=0, name="lomc!", metadata=nothing)
+    @warn "The use of `lomc!` is deprecated. Use `ProjectorMonteCarloProblem` and `solve` instead."
+
     if !iszero(laststep)
         state = @set state.simulation_plan.last_step = laststep
     end
diff --git a/test/lomc.jl b/test/lomc.jl
index 987df06a2..9cfd74dc1 100644
--- a/test/lomc.jl
+++ b/test/lomc.jl
@@ -53,7 +53,7 @@ Random.seed!(1234)
         address = BoseFS{5,2}((2,3))
         H = HubbardReal1D(address; u=0.1)
         dv = DVec(address => 1; style=IsStochasticInteger())
-        df, state = @test_logs (:warn, Regex("(Simulation)")) lomc!(H, dv; laststep=0, shift=23.1, dτ=0.002)
+        df, state = @test_logs (:warn, Regex("(Simulation)")) match_mode = :any lomc!(H, dv; laststep=0, shift=23.1, dτ=0.002)
         @test state.spectral_states[1].single_states[1].shift_parameters.time_step  == 0.002
         @test state.spectral_states[1].single_states[1].shift_parameters.shift == 23.1
         @test state.replica_strategy == NoStats{1}() # uses getfield method
@@ -91,7 +91,7 @@ Random.seed!(1234)
         walkers = lomc!(H, copy(dv); s_strat, laststep=1000).df.norm
         @test median(walkers) ≈ 1000 rtol=0.1
 
-        _, state = @test_logs (:warn, Regex("(Simulation)")) lomc!(H, copy(dv); targetwalkers=500, laststep=0)
+        _, state = @test_logs (:warn, Regex("(Simulation)")) match_mode = :any lomc!(H, copy(dv); targetwalkers=500, laststep=0)
         @test only(state).algorithm.shift_strategy.target_walkers == 500
     end
 

From 703116a9ecb1b6ec611e46e1721aeb250044a4d1 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Sun, 22 Dec 2024 12:35:24 +1300
Subject: [PATCH 05/17] bypass lomc! tests

---
 test/runtests.jl | 6 +++---
 1 file changed, 3 insertions(+), 3 deletions(-)

diff --git a/test/runtests.jl b/test/runtests.jl
index a566bc46b..d9a7ce896 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -42,9 +42,9 @@ end
     include("projector_monte_carlo_problem.jl")
 end
 
-@safetestset "lomc!" begin
-    include("lomc.jl")
-end
+# @safetestset "lomc!" begin
+#     include("lomc.jl")
+# end
 
 @safetestset "RimuIO" begin
     include("RimuIO.jl")

From e3e3d19ce31f8c925400bc98a5dcc25387fc30e7 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Sun, 22 Dec 2024 13:03:16 +1300
Subject: [PATCH 06/17] remove lomc! from tests

---
 test/Hamiltonians.jl     | 20 ++++++++++++--------
 test/Interfaces.jl       |  6 +-----
 test/StochasticStyles.jl |  3 ++-
 3 files changed, 15 insertions(+), 14 deletions(-)

diff --git a/test/Hamiltonians.jl b/test/Hamiltonians.jl
index 03451a2b2..61174082a 100644
--- a/test/Hamiltonians.jl
+++ b/test/Hamiltonians.jl
@@ -675,8 +675,6 @@ end
 end
 
 @testset "MatrixHamiltonian" begin
-    # using lomc! with an `AbstractMatrix` was removed in Rimu.jl v0.12.0
-
     # generate matrix
     ham = HubbardReal1D(BoseFS((1, 1, 1, 1)))
     dim = dimension(ham)
@@ -685,9 +683,9 @@ end
     sparse_matrix, basis = sparse(bsr), bsr.basis
     @test dim == length(basis)
 
-    # run lomc! in deterministic mode with Hamiltonian and DVec
+    # run ProjectorMonteCarloProblem in deterministic mode with Hamiltonian and DVec
     v = DVec(k=>1.0 for k in basis; style=IsDeterministic()) # corresponds to `ones(dim)`
-    a = lomc!(ham, v).df
+    a = solve(ProjectorMonteCarloProblem(ham; start_at=v)).df
 
     # MatrixHamiltonian
     @test_throws ArgumentError MatrixHamiltonian([1 2 3; 4 5 6])
@@ -706,17 +704,23 @@ end
     @test starting_address(mh) == 1
     @test dimension(mh) == dim
 
-    # lomc!
+    # ProjectorMonteCarloProblem with MatrixHamiltonian
     # float walkernumber triggers IsDeterministic algorithm
-    d = lomc!(mh, DVec(pairs(ones(dim)))).df
+    sim = solve(ProjectorMonteCarloProblem(mh; start_at=DVec(pairs(ones(dim)))))
+    @test StochasticStyle(only(state_vectors(sim))) == IsDeterministic()
+    d = DataFrame(sim)
     @test d.shift ≈ a.shift
     # integer walkernumber triggers IsStochasticInteger algorithm
     Random.seed!(15)
-    e = lomc!(mh, DVec(pairs(ones(Int,dim)))).df
+    sim = solve(ProjectorMonteCarloProblem(mh; start_at=DVec(pairs(ones(Int, dim)))))
+    @test StochasticStyle(only(state_vectors(sim))) == IsStochasticInteger()
+    e = DataFrame(sim)
     @test ≈(e.shift[end], a.shift[end], atol=0.3)
     # wrap full matrix as MatrixHamiltonian
     fmh =  MatrixHamiltonian(Matrix(sparse_matrix))
-    f = lomc!(fmh, DVec(pairs(ones(dim)))).df
+    sim = solve(ProjectorMonteCarloProblem(fmh; start_at=DVec(pairs(ones(dim)))))
+    @test StochasticStyle(only(state_vectors(sim))) == IsDeterministic()
+    f = DataFrame(sim)
     @test f.shift ≈ a.shift
 end
 
diff --git a/test/Interfaces.jl b/test/Interfaces.jl
index 6b1905ad3..18fa74f50 100644
--- a/test/Interfaces.jl
+++ b/test/Interfaces.jl
@@ -44,7 +44,6 @@ end
            0 0 1 1 2;
            0 1 0 1 0]
     vector = ones(5)
-    @test_throws ArgumentError lomc!(ham, vector; laststep=10_000)
 
     # rephrase with MatrixHamiltonian
     mh = MatrixHamiltonian(ham)
@@ -54,13 +53,10 @@ end
     # solve with new API
     p = ProjectorMonteCarloProblem(mh; start_at=sv, last_step=10_000, post_step_strategy)
     sm = solve(p)
-    last_shift = DataFrame(sm).shift[end]
+    df = DataFrame(sm)
 
-    # solve with old API
-    df, _ = lomc!(mh, sv; laststep=10_000, post_step_strategy)
     eigs = eigen(ham)
 
-    @test eigs.values[1] ≈ last_shift rtol = 0.01
     @test df.shift[end] ≈ eigs.values[1] rtol=0.01
     @test df.hproj[end] / df.vproj[end] ≈ eigs.values[1] rtol=0.01
     @test normalize(state_vectors(sm)[1]) ≈ DVec(pairs(eigs.vectors[:, 1])) rtol = 0.01
diff --git a/test/StochasticStyles.jl b/test/StochasticStyles.jl
index 52b5a8f38..8b3eb5747 100644
--- a/test/StochasticStyles.jl
+++ b/test/StochasticStyles.jl
@@ -241,7 +241,8 @@ end
     add = BoseFS((0,0,0,10,0,0,0))
     H = HubbardMom1D(add)
     dv = DVec(add => 1.0, style=IsDeterministic())
-    lomc!(H, dv)
+    sim = solve(ProjectorMonteCarloProblem(H; start_at=dv))
+    dv = only(state_vectors(sim))
     normalize!(dv)
 
     for compression in (

From a5d7a399ed78d63aad1ee67213ed68978dcb1fab Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Sun, 22 Dec 2024 13:52:55 +1300
Subject: [PATCH 07/17] suppress warning in doctest

---
 src/lomc.jl      | 4 ++--
 test/doctests.jl | 2 +-
 2 files changed, 3 insertions(+), 3 deletions(-)

diff --git a/src/lomc.jl b/src/lomc.jl
index 0ab87e568..e47d80ea1 100644
--- a/src/lomc.jl
+++ b/src/lomc.jl
@@ -55,9 +55,9 @@ julia> address = BoseFS(1,2,3);
 
 julia> hamiltonian = HubbardReal1D(address);
 
-julia> df1, state = lomc!(hamiltonian; targetwalkers=500, laststep=100);
+julia> df1, state = @suppress lomc!(hamiltonian; targetwalkers=500, laststep=100);
 
-julia> df2, _ = lomc!(state, df1; laststep=200, metadata=(;info="cont")); # Continuation run
+julia> df2, _ = @suppress lomc!(state, df1; laststep=200, metadata=(;info="cont")); # Continuation run
 
 julia> size(df1)
 (100, 9)
diff --git a/test/doctests.jl b/test/doctests.jl
index 11e61cda3..2b3b6a07e 100644
--- a/test/doctests.jl
+++ b/test/doctests.jl
@@ -5,7 +5,7 @@ DocMeta.setdocmeta!(
     Rimu,
     :DocTestSetup,
     :(using Rimu; using Rimu.StatsTools; using DataFrames; using Random;
-      using LinearAlgebra);
+      using LinearAlgebra; using Suppressor);
     recursive=true
 )
 # Run with fix=true to fix docstrings

From 5ed6497ed31d5be4be8af1f603f5c8a4b341b8d7 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Wed, 25 Dec 2024 18:51:06 +1300
Subject: [PATCH 08/17] eliminate lomc! and consistenly apply rng seeding with
 ProjectorMonteCarloProblem

---
 src/projector_monte_carlo_problem.jl |  8 +++--
 test/Hamiltonians.jl                 |  9 ++++--
 test/StatsTools.jl                   | 46 ++++++++++++++++++----------
 test/mpi_runtests.jl                 |  8 ++---
 test/runtests.jl                     | 36 +++++++++++-----------
 5 files changed, 62 insertions(+), 45 deletions(-)

diff --git a/src/projector_monte_carlo_problem.jl b/src/projector_monte_carlo_problem.jl
index bbd7f8e46..93a5aed68 100644
--- a/src/projector_monte_carlo_problem.jl
+++ b/src/projector_monte_carlo_problem.jl
@@ -96,9 +96,11 @@ julia> size(DataFrame(simulation))
   pairs or a `NamedTuple`, e.g. `metadata = ("key1" => "value1", "key2" => "value2")`.
   All metadata is converted to strings.
 - `random_seed = true`: Provide and store a seed for the random number generator. If set to
-    `true`, a random seed is generated. If set to number, this number is used as the seed.
-    The seed is used by `solve` such that `solve`ing the problem twice will yield identical
-    results. If set to `false`, no seed is used and results are not reproducible.
+    `true`, a new random seed is generated from `RandomDevice()`. If set to number, this
+    number is used as the seed. This seed is used by `solve` (and `init`) to re-seed the
+    default random number generator (consistently on each MPI rank) such that
+    `solve`ing the same `ProjectorMonteCarloProblem` twice will yield identical results. If
+    set to `false`, no seed is used and consecutive random numbers are used.
 - `minimum_size = 2*num_spectral_states(spectral_strategy)`: The minimum size of the basis
     used to construct starting vectors for simulations of spectral states, if `start_at`
     is not provided.
diff --git a/test/Hamiltonians.jl b/test/Hamiltonians.jl
index 61174082a..61f480414 100644
--- a/test/Hamiltonians.jl
+++ b/test/Hamiltonians.jl
@@ -711,14 +711,17 @@ end
     d = DataFrame(sim)
     @test d.shift ≈ a.shift
     # integer walkernumber triggers IsStochasticInteger algorithm
-    Random.seed!(15)
-    sim = solve(ProjectorMonteCarloProblem(mh; start_at=DVec(pairs(ones(Int, dim)))))
+    sim = solve(
+        ProjectorMonteCarloProblem(mh; start_at=DVec(pairs(ones(Int, dim))), random_seed=18)
+    )
     @test StochasticStyle(only(state_vectors(sim))) == IsStochasticInteger()
     e = DataFrame(sim)
     @test ≈(e.shift[end], a.shift[end], atol=0.3)
     # wrap full matrix as MatrixHamiltonian
     fmh =  MatrixHamiltonian(Matrix(sparse_matrix))
-    sim = solve(ProjectorMonteCarloProblem(fmh; start_at=DVec(pairs(ones(dim)))))
+    sim = solve(
+        ProjectorMonteCarloProblem(fmh; start_at=DVec(pairs(ones(dim))), random_seed=15)
+    )
     @test StochasticStyle(only(state_vectors(sim))) == IsDeterministic()
     f = DataFrame(sim)
     @test f.shift ≈ a.shift
diff --git a/test/StatsTools.jl b/test/StatsTools.jl
index 0f15e2e54..9478b6821 100644
--- a/test/StatsTools.jl
+++ b/test/StatsTools.jl
@@ -193,25 +193,29 @@ using Rimu.StatsTools: x_by_y_linear, ratio_estimators, particles
 end
 
 @testset "Reweighting" begin
-    Random.seed!(133)
     ham = HubbardReal1D(BoseFS((1, 1, 1, 1)), u=6.0, t=1.0)
     # using KrylovKit
     # fv = DVec(starting_address(ham)=>1.0; capacity=dimension(ham))
     # kkresults = eigsolve(ham, fv, 1, :SR; issymmetric = true)
-    # exact_energy = kkresults[1][1]
-    exact_energy = -2.869739978337469
+    # exact_energy_value = kkresults[1][1]
+    exact_energy_value = -2.869739978337469
     # run integer walker FCIQMC to get significant bias
     v = DVec(starting_address(ham) => 2; capacity=dimension(ham))
     steps_equi = 200
     steps_meas = 2^10
-    p = RunTillLastStep(laststep=steps_equi + steps_meas)
-    post_step = ProjectedEnergy(ham, v)
-    s_strat = DoubleLogUpdate(target_walkers=10)
-    df = lomc!(ham, v; params=p, s_strat, post_step).df
+    # p = RunTillLastStep(laststep=steps_equi + steps_meas)
+    post_step_strategy = ProjectedEnergy(ham, v)
+    shift_strategy = DoubleLogUpdate(target_walkers=10)
+    p = ProjectorMonteCarloProblem(
+        ham;
+        start_at=v, last_step=steps_equi + steps_meas, shift_strategy, post_step_strategy,
+        random_seed=1234
+    )
+    df = DataFrame(solve(p))
     @test_throws ArgumentError variational_energy_estimator(df) # see next testset
     bs = shift_estimator(df; skip=steps_equi)
     @test bs == blocking_analysis(df.shift[steps_equi+1:end])
-    pcb = bs.mean - exact_energy
+    pcb = bs.mean - exact_energy_value
 
     @test pcb > 0.0 # the shift has a large population control bias
     # test growth_estimator
@@ -264,12 +268,16 @@ end
     num_reps = 2
     tw = 10
 
-    params = RunTillLastStep(laststep = skipsteps + runsteps)
-    s_strat = DoubleLogUpdate(target_walkers = tw)
+    shift_strategy = DoubleLogUpdate(target_walkers = tw)
     G2list = ([G2RealCorrelator(i) for i in dvals]...,)
 
-    Random.seed!(174)
-    df = lomc!(ham, dv; params, s_strat, replica_strategy = AllOverlaps(num_reps; operator = G2list)).df
+    p = ProjectorMonteCarloProblem(
+        ham;
+        start_at=dv, last_step=skipsteps+runsteps, shift_strategy,
+        replica_strategy = AllOverlaps(num_reps; operator = G2list),
+        random_seed=179
+    )
+    df = DataFrame(solve(p))
 
     for d in dvals
         # without reweighting
@@ -306,13 +314,17 @@ using Rimu.StatsTools: replica_fidelity
     v = DVec(starting_address(ham) => 2; capacity=dimension(ham))
     steps_equi = 200
     steps_meas = 2^10
-    p = RunTillLastStep(laststep=steps_equi + steps_meas)
-    post_step = (Projector(vproj=gs), Projector(hproj=os))
-    s_strat = DoubleLogUpdate(target_walkers=10)
+    post_step_strategy = (Projector(vproj=gs), Projector(hproj=os))
+    shift_strategy = DoubleLogUpdate(target_walkers=10)
+    replica_strategy = AllOverlaps(2)
 
     # run replica fciqmc
-    Random.seed!(170)
-    rr = lomc!(ham, v; params=p, s_strat, post_step, replica_strategy=AllOverlaps()).df
+    p = ProjectorMonteCarloProblem(
+        ham;
+        start_at=v, last_step=steps_equi+steps_meas, shift_strategy, post_step_strategy,
+        random_seed=170, replica_strategy
+    )
+    rr = DataFrame(solve(p))
 
     # check fidelity with ground state
     fid_gs = replica_fidelity(rr; p_field=:vproj, skip=steps_equi)
diff --git a/test/mpi_runtests.jl b/test/mpi_runtests.jl
index fbbc47f05..82261916f 100644
--- a/test/mpi_runtests.jl
+++ b/test/mpi_runtests.jl
@@ -171,7 +171,7 @@ end
                 Projector(proj_1=Norm2Projector()),
             )
             prob = ProjectorMonteCarloProblem(
-                H; start_at=dv, post_step_strategy, last_step=5000
+                H; start_at=dv, post_step_strategy, last_step=5000, random_seed=13
             )
             df = DataFrame(solve(prob))
 
@@ -194,7 +194,7 @@ end
             dv = PDVec(addr => 3; initiator_threshold=1)
             shift_strategy = DoubleLogUpdate(target_walkers=100)
             prob = ProjectorMonteCarloProblem(
-                H; start_at=dv, last_step=5000, shift_strategy
+                H; start_at=dv, last_step=5000, shift_strategy, random_seed=13
             )
             df = DataFrame(solve(prob))
 
@@ -214,7 +214,7 @@ end
                 # Diagonal
                 replica_strategy = AllOverlaps(2; operator=ntuple(DensityMatrixDiagonal, M))
                 prob = ProjectorMonteCarloProblem(
-                    H; start_at=dv, replica_strategy, last_step=10_000
+                    H; start_at=dv, replica_strategy, last_step=10_000, random_seed=13
                 )
                 df = DataFrame(solve(prob))
 
@@ -229,7 +229,7 @@ end
                 ops = ntuple(x -> G2MomCorrelator(x - cld(M, 2)), M)
                 replica_strategy = AllOverlaps(2; operator=ops)
                 prob = ProjectorMonteCarloProblem(
-                    H; start_at=dv, replica_strategy, last_step=10_000
+                    H; start_at=dv, replica_strategy, last_step=10_000, random_seed=13
                 )
                 df = DataFrame(solve(prob))
 
diff --git a/test/runtests.jl b/test/runtests.jl
index d9a7ce896..d5e48658c 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -14,6 +14,24 @@ using ExplicitImports: check_no_implicit_imports
 
 @test Rimu.PACKAGE_VERSION == VersionNumber(TOML.parsefile(pkgdir(Rimu, "Project.toml"))["version"])
 
+@safetestset "ExplicitImports" begin
+    using Rimu
+    using ExplicitImports
+    # Check that no implicit imports are used in the Rimu module.
+    # See https://ericphanson.github.io/ExplicitImports.jl/stable/
+    @test check_no_implicit_imports(Rimu; skip=(Rimu, Base, Core, VectorInterface)) === nothing
+    # If this test fails, make your import statements explicit.
+    # For example, replace `using Foo` with `using Foo: bar, baz`.
+end
+
+@safetestset "doctests" begin
+    include("doctests.jl")
+end
+
+@safetestset "excited states" begin
+    include("excited_states_tests.jl")
+end
+
 @safetestset "Interfaces" begin
     include("Interfaces.jl")
 end
@@ -120,22 +138,4 @@ end
     @test default_logger() isa Logging.ConsoleLogger
 end
 
-@safetestset "doctests" begin
-    include("doctests.jl")
-end
-
-@safetestset "ExplicitImports" begin
-    using Rimu
-    using ExplicitImports
-    # Check that no implicit imports are used in the Rimu module.
-    # See https://ericphanson.github.io/ExplicitImports.jl/stable/
-    @test check_no_implicit_imports(Rimu; skip=(Rimu, Base, Core, VectorInterface)) === nothing
-    # If this test fails, make your import statements explicit.
-    # For example, replace `using Foo` with `using Foo: bar, baz`.
-end
-
-@safetestset "excited states" begin
-    include("excited_states_tests.jl")
-end
-
 # Note: Running Rimu with several MPI ranks is tested seperately on GitHub CI and not here.

From be591d2440ae5ddbe0ee5d144c60d3958d9dd41d Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Wed, 25 Dec 2024 20:11:59 +1300
Subject: [PATCH 09/17] more tests for ProjectorMonteCarloProblem

---
 src/fciqmc.jl                         | 12 +----
 test/projector_monte_carlo_problem.jl | 64 +++++++++++++++++++++++++++
 2 files changed, 66 insertions(+), 10 deletions(-)

diff --git a/src/fciqmc.jl b/src/fciqmc.jl
index f2d22361c..661932d18 100644
--- a/src/fciqmc.jl
+++ b/src/fciqmc.jl
@@ -170,19 +170,11 @@ function advance!(algorithm::FCIQMC, report, state::ReplicaState, s_state::Singl
     end
 
     if len == 0
-        if length(state.spectral_states) > 1
-            @error "population in single state $(s_state.id) is dead. Aborting."
-        else
-            @error "population is dead. Aborting."
-        end
+        @error "Population in state $(s_state.id) is dead. Aborting."
         return false
     end
     if len > state.max_length[]
-        if length(state.spectral_states) > 1
-            @error "`max_length` reached in single state $(s_state.id). Aborting."
-        else
-            @error "`max_length` reached. Aborting."
-        end
+        @error "`max_length` reached in state $(s_state.id). Aborting."
         return false
     end
     return proceed # Bool
diff --git a/test/projector_monte_carlo_problem.jl b/test/projector_monte_carlo_problem.jl
index 5cdde872c..9ffbf402e 100644
--- a/test/projector_monte_carlo_problem.jl
+++ b/test/projector_monte_carlo_problem.jl
@@ -172,6 +172,8 @@ using Rimu: num_replicas, num_spectral_states
 
     # Tables.jl interface
     @test Tables.istable(sm)
+    @test Tables.columnaccess(sm)
+    @test Tables.schema(sm) == Tables.schema(DataFrame(sm))
     @test map(NamedTuple, Tables.rows(sm)) == map(NamedTuple, Tables.rows(df))
 
     # continue simulation
@@ -220,3 +222,65 @@ using Rimu: num_replicas, num_spectral_states
     @test solve!(sm; last_step=300, reporting_strategy=ReportDFAndInfo()) === sm
     @test size(sm.df)[1] == 200 # the report was not emptied
 end
+
+@testset "Dead population" begin
+    address = BoseFS{5,2}((2, 3))
+    H = HubbardReal1D(address; u=20)
+    dv = DVec(address => 10; style=IsStochasticInteger())
+
+    # Only population is dead.
+    p = ProjectorMonteCarloProblem(H; start_at=dv, last_step=100, shift=0.0, random_seed=7)
+    sim = @suppress_err solve(p)
+    @test sim.aborted == true
+    @test sim.success == false
+    @test sim.modified == true
+    @test sim.message == "Aborted in step 5."
+    @test size(sim.df, 1) < 100
+
+    # population does not die with sensible default shift
+    p = ProjectorMonteCarloProblem(H; start_at=dv, last_step=100, random_seed=7)
+    sim = solve(p)
+    @test sim.aborted == false
+    @test sim.success == true
+    @test sim.modified == true
+    @test size(sim.df, 1) == 100
+
+    # Populations in replicas are dead.
+    p = ProjectorMonteCarloProblem(
+        H;
+        start_at=dv, n_replicas=3, last_step=100, shift=0.0, random_seed=7
+    )
+    sim = @suppress_err solve(p)
+    @test sim.aborted == true
+    @test sim.success == false
+    @test sim.modified == true
+    @test sim.message == "Aborted in step 3."
+    @test size(sim.df, 1) < 100
+end
+
+@testset "max_length" begin
+    # walker number blows up when time_step is too large
+    h = HubbardReal1D(BoseFS(1, 3, 5, 2, 1))
+    p = ProjectorMonteCarloProblem(h; time_step=0.1, target_walkers=100, random_seed=7)
+    sm = init(p)
+    @test size(DataFrame(sm)) == (0, 0)
+    @test size(sm.df) == (0, 0)
+    @test sm.state.max_length[] > 100 # default max_length
+    @test sm.state.step[] == 0
+    @test @suppress_err solve!(sm) === sm
+    @test sm.modified == true
+    @test sm.success == false
+    @test sm.aborted == true
+    @test sm.message == "Aborted in step 6."
+    @test is_finalized(sm.report) == true
+    @test @suppress_err step!(sm) === sm # no effect, aborted
+
+
+    # runs fine with a smaller time_step
+    p = ProjectorMonteCarloProblem(h; time_step=0.01, target_walkers=100, random_seed=7)
+    sm = solve!(init(p))
+    @test sm.success == true
+    @test sm.aborted == false
+    @test size(sm.df, 1) == 100
+    @test @suppress_err step!(sm) === sm # no effect, already finalized
+end

From e24b5d55beaa31ab8db54d9f6cd7d005d148972b Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Wed, 25 Dec 2024 22:43:01 +1300
Subject: [PATCH 10/17] fix tests

---
 test/projector_monte_carlo_problem.jl | 2 +-
 1 file changed, 1 insertion(+), 1 deletion(-)

diff --git a/test/projector_monte_carlo_problem.jl b/test/projector_monte_carlo_problem.jl
index 9ffbf402e..9ddda143e 100644
--- a/test/projector_monte_carlo_problem.jl
+++ b/test/projector_monte_carlo_problem.jl
@@ -1,5 +1,5 @@
 using Rimu
-using Test
+using Test, Suppressor
 import Random
 
 using Rimu: is_finalized

From 504032d35255b94c51a0b4b25886eee5465ce519 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Wed, 25 Dec 2024 23:10:15 +1300
Subject: [PATCH 11/17] remove _n_walkers and improve tests

---
 src/qmc_states.jl                     | 13 -------------
 test/projector_monte_carlo_problem.jl |  4 ++++
 2 files changed, 4 insertions(+), 13 deletions(-)

diff --git a/src/qmc_states.jl b/src/qmc_states.jl
index 2c280dc84..56937ad45 100644
--- a/src/qmc_states.jl
+++ b/src/qmc_states.jl
@@ -67,19 +67,6 @@ function state_vectors(state::SpectralState)
     return SMatrix{1, num_spectral_states(state)}(s.v for s in state.single_states)
 end
 
-"""
-    _n_walkers(v, shift_strategy)
-Returns an estimate of the expected number of walkers as an integer.
-"""
-function _n_walkers(v, shift_strategy)
-    n = if hasfield(typeof(shift_strategy), :target_walkers)
-        shift_strategy.target_walkers
-    else # e.g. for LogUpdate()
-        walkernumber(v)
-    end
-    return ceil(Int, max(real(n), imag(n)))
-end
-
 """
     ReplicaState <: AbstractMatrix{SingleState}
 
diff --git a/test/projector_monte_carlo_problem.jl b/test/projector_monte_carlo_problem.jl
index 9ddda143e..d7104f18b 100644
--- a/test/projector_monte_carlo_problem.jl
+++ b/test/projector_monte_carlo_problem.jl
@@ -95,6 +95,10 @@ end
         @test sm.state[2].shift_parameters.shift ≡ 2.0
         @test state_vectors(sm.state)[1][BoseFS(1, 3)] == 10
         @test state_vectors(sm.state)[2][BoseFS(3, 1)] == 10
+        @test startswith(sprint(show, sm.state), "ReplicaState")
+        @test startswith(sprint(show, sm.state[1]), "SingleState")
+        @test startswith(sprint(show, sm.state.spectral_states), "(SpectralState")
+
     end
 
     @testset "random seeds" begin

From 9e3ba1f505152d2950865d0d37d96628320ef88c Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Wed, 25 Dec 2024 23:33:48 +1300
Subject: [PATCH 12/17] test deprecated keyword arguments

---
 src/projector_monte_carlo_problem.jl  | 23 ++++++++++++++++++++---
 test/projector_monte_carlo_problem.jl | 10 ++++++++++
 2 files changed, 30 insertions(+), 3 deletions(-)

diff --git a/src/projector_monte_carlo_problem.jl b/src/projector_monte_carlo_problem.jl
index 93a5aed68..4128fb0b5 100644
--- a/src/projector_monte_carlo_problem.jl
+++ b/src/projector_monte_carlo_problem.jl
@@ -159,15 +159,16 @@ function ProjectorMonteCarloProblem(
     starting_step = 0,
     last_step = 100,
     wall_time = Inf,
-    simulation_plan = SimulationPlan(starting_step, last_step, wall_time),
+    walltime = nothing, # deprecated
+    simulation_plan = nothing,
     replica_strategy = NoStats(n_replicas),
     targetwalkers = nothing, # deprecated
     target_walkers = 1_000,
     ζ = 0.08,
     ξ = ζ^2/4,
-    shift_strategy = DoubleLogUpdate(; target_walkers, ζ, ξ),
+    shift_strategy = nothing,
     time_step_strategy=ConstantTimeStep(),
-    algorithm=FCIQMC(; shift_strategy, time_step_strategy),
+    algorithm=nothing,
     initial_shift_parameters=nothing,
     reporting_strategy = ReportDFAndInfo(),
     post_step_strategy = (),
@@ -180,11 +181,27 @@ function ProjectorMonteCarloProblem(
     display_name = "PMCSimulation",
     random_seed = true
 )
+    if !isnothing(walltime)
+        @warn "The keyword argument `walltime` is deprecated. Use `wall_time` instead."
+        wall_time = walltime
+    end
+    if isnothing(simulation_plan)
+        simulation_plan = SimulationPlan(starting_step, last_step, wall_time)
+    end
+
     if !isnothing(targetwalkers)
         @warn "The keyword argument `targetwalkers` is deprecated. Use `target_walkers` instead."
         target_walkers = targetwalkers
     end
 
+    if isnothing(shift_strategy)
+        shift_strategy = DoubleLogUpdate(; target_walkers, ζ, ξ)
+    end
+
+    if isnothing(algorithm)
+        algorithm = FCIQMC(; shift_strategy, time_step_strategy)
+    end
+
     n_replicas = num_replicas(replica_strategy) # replica_strategy may override n_replicas
 
     n_spectral = num_spectral_states(spectral_strategy) # spectral_strategy may override n_spectral
diff --git a/test/projector_monte_carlo_problem.jl b/test/projector_monte_carlo_problem.jl
index d7104f18b..dfc793c91 100644
--- a/test/projector_monte_carlo_problem.jl
+++ b/test/projector_monte_carlo_problem.jl
@@ -288,3 +288,13 @@ end
     @test size(sm.df, 1) == 100
     @test @suppress_err step!(sm) === sm # no effect, already finalized
 end
+
+@testset "deprecated keyword arguments" begin
+    h = HubbardReal1D(BoseFS(1, 3))
+    p = @suppress_err ProjectorMonteCarloProblem(
+        h; shift=1.0, targetwalkers=100, maxlength=200, walltime=23
+    )
+    @test p.algorithm.shift_strategy.target_walkers == 100
+    @test p.max_length == 200
+    @test p.simulation_plan.wall_time == 23
+end

From 7e7c236de4aab26ef1da848e4ca0d85c89fdd025 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Wed, 25 Dec 2024 23:59:41 +1300
Subject: [PATCH 13/17] more tests

---
 test/projector_monte_carlo_problem.jl | 16 ++++++++++++++++
 test/runtests.jl                      |  6 +++---
 2 files changed, 19 insertions(+), 3 deletions(-)

diff --git a/test/projector_monte_carlo_problem.jl b/test/projector_monte_carlo_problem.jl
index dfc793c91..5010c49cb 100644
--- a/test/projector_monte_carlo_problem.jl
+++ b/test/projector_monte_carlo_problem.jl
@@ -98,9 +98,25 @@ end
         @test startswith(sprint(show, sm.state), "ReplicaState")
         @test startswith(sprint(show, sm.state[1]), "SingleState")
         @test startswith(sprint(show, sm.state.spectral_states), "(SpectralState")
+    end
+
+    @testset "Default DVec" begin
+        address = BoseFS(2, 3)
+        H = HubbardReal1D(address; u=20)
+        sm = init(ProjectorMonteCarloProblem(H; threading=false))
+        @test only(state_vectors(sm)) isa DVec
+        @test StochasticStyle(only(state_vectors(sm))) isa IsDynamicSemistochastic
 
+        sm = init(ProjectorMonteCarloProblem(H; threading=true))
+        @test only(state_vectors(sm)) isa PDVec
+        @test StochasticStyle(only(state_vectors(sm))) isa IsDynamicSemistochastic
+
+        sm = init(ProjectorMonteCarloProblem(H; threading=false, initiator=true))
+        @test only(state_vectors(sm)) isa InitiatorDVec
+        @test StochasticStyle(only(state_vectors(sm))) isa IsDynamicSemistochastic
     end
 
+
     @testset "random seeds" begin
         p = ProjectorMonteCarloProblem(h) # generates random_seed
         @test p.random_seed isa UInt64
diff --git a/test/runtests.jl b/test/runtests.jl
index d5e48658c..37f46602e 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -60,9 +60,9 @@ end
     include("projector_monte_carlo_problem.jl")
 end
 
-# @safetestset "lomc!" begin
-#     include("lomc.jl")
-# end
+@suppress_err @safetestset "lomc!" begin
+    include("lomc.jl")
+end
 
 @safetestset "RimuIO" begin
     include("RimuIO.jl")

From d0ee8ea4cb4c32e7dd7563823ab535456b3d273f Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Thu, 26 Dec 2024 13:31:52 +1300
Subject: [PATCH 14/17] docstring improvements

---
 src/pmc_simulation.jl                         | 10 +++
 src/qmc_states.jl                             | 73 ++++++++++++++-----
 .../reportingstrategy.jl                      | 46 +++++++-----
 test/projector_monte_carlo_problem.jl         | 11 ++-
 4 files changed, 99 insertions(+), 41 deletions(-)

diff --git a/src/pmc_simulation.jl b/src/pmc_simulation.jl
index 63a476868..f45b0de0e 100644
--- a/src/pmc_simulation.jl
+++ b/src/pmc_simulation.jl
@@ -6,6 +6,16 @@ Is returned by [`init(::ProjectorMonteCarloProblem)`](@ref) and solved with
 
 Obtain the results of a simulation `sm` as a DataFrame with `DataFrame(sm)`.
 
+## Fields
+- `problem::ProjectorMonteCarloProblem`: The problem that was solved
+- `state::Rimu.ReplicaState`: The current state of the simulation
+- `report::Rimu.Report`: The report of the simulation
+- `modified::Bool`: Whether the simulation has been modified
+- `aborted::Bool`: Whether the simulation has been aborted
+- `success::Bool`: Whether the simulation has been completed successfully
+- `message::String`: A message about the simulation status
+- `elapsed_time::Float64`: The time elapsed during the simulation
+
 See also [`state_vectors`](@ref),
 [`ProjectorMonteCarloProblem`](@ref), [`init`](@ref), [`solve!`](@ref).
 """
diff --git a/src/qmc_states.jl b/src/qmc_states.jl
index 56937ad45..216431679 100644
--- a/src/qmc_states.jl
+++ b/src/qmc_states.jl
@@ -2,11 +2,20 @@
     SingleState(hamiltonian, algorithm, v, wm, pnorm, params, id)
 
 Struct that holds a single state vector and all information needed for an independent run
-of the algorithm. Can be advanced a step forward with [`advance!`](@ref).
+of the algorithm. Can be advanced a step forward with [`Rimu.advance!`](@ref).
 
-See also [`SpectralState`](@ref), [`SpectralStrategy`](@ref),
-[`ReplicaState`](@ref), [`ReplicaStrategy`](@ref), [`replica_stats`](@ref),
-[`PMCSimulation`](@ref).
+## Fields
+- `hamiltonian`: Hamiltonian
+- `algorithm`: Algorithm
+- `v`: Vector
+- `pv`: Previous vector
+- `wm`: Working memory
+- `shift_parameters`: Shift parameters
+- `id::String`: id is appended to column names
+
+See also [`SpectralStrategy`](@ref), [`ReplicaStrategy`](@ref),
+[`Rimu.SpectralState`](@ref), [`Rimu.ReplicaState`](@ref), [`Rimu.replica_stats`](@ref),
+[`Rimu.PMCSimulation`](@ref).
 """
 mutable struct SingleState{H,A,V,W,SP}
     # Future TODO: rename these fields, add interface for accessing them.
@@ -19,11 +28,12 @@ mutable struct SingleState{H,A,V,W,SP}
     id::String # id is appended to column names
 end
 
-function Base.show(io::IO, r::SingleState)
+Base.show(io::IO, r::SingleState) = show(io, MIME("text/plain"), r)
+function Base.show(io::IO, ::MIME"text/plain", r::SingleState)
     print(
         io,
-        "SingleState(v: ", length(r.v), "-element ", nameof(typeof(r.v)),
-        ", wm: ", length(r.wm), "-element ", nameof(typeof(r.wm)), ")"
+        "Rimu.SingleState with ", length(r.v), "-element ", nameof(typeof(r.v)),
+        " and id \"", r.id, "\""
     )
 end
 
@@ -31,9 +41,13 @@ end
     SpectralState <: AbstractVector{SingleState}
 Holds one or several [`SingleState`](@ref)s representing the ground state and excited
 states of a single replica.
-
 Indexing the `SpectralState` `state[i]` returns the `i`th `SingleState`.
 
+## Fields
+- `single_states`: Tuple of `SingleState`s
+- `spectral_strategy`: Strategy for computing the spectral states
+- `id::String`: id is appended to column names
+
 See also [`SpectralStrategy`](@ref), [`ReplicaState`](@ref), [`SingleState`](@ref),
 [`PMCSimulation`](@ref).
 """
@@ -54,9 +68,11 @@ num_spectral_states(::SpectralState{N}) where {N} = N
 Base.size(s::SpectralState) = (num_spectral_states(s),)
 Base.getindex(s::SpectralState, i::Int) = s.single_states[i]
 
-function Base.show(io::IO, s::SpectralState)
-    print(io, "SpectralState")
-    print(io, " with ", num_spectral_states(s), " spectral states")
+Base.show(io::IO, s::SpectralState) = show(io, MIME("text/plain"), s)
+function Base.show(io::IO, ::MIME"text/plain", s::SpectralState)
+    ns = num_spectral_states(s)
+    print(io, "$ns-element Rimu.SpectralState")
+    print(io, " with ", ns, " spectral state(s) of type ", nameof(typeof(s[1])))
     print(io, "\n    spectral_strategy: ", s.spectral_strategy)
     for (i, r) in enumerate(s.single_states)
         print(io, "\n      $i: ", r)
@@ -71,10 +87,18 @@ end
     ReplicaState <: AbstractMatrix{SingleState}
 
 Holds information about multiple replicas of [`SpectralState`](@ref)s.
-
 Indexing the `ReplicaState` `state[i, j]` returns a `SingleState` from the `i`th replica
 and `j`th spectral state.
 
+## Fields
+- `spectral_states`: Tuple of `SpectralState`s
+- `max_length::Ref{Int}`: Maximum length of the simulation
+- `step::Ref{Int}`: Current step of the simulation
+- `simulation_plan`: Simulation plan
+- `reporting_strategy`: Reporting strategy
+- `post_step_strategy`: Post-step strategy
+- `replica_strategy`: Replica strategy
+
 See also [`ReplicaStrategy`](@ref), [`SpectralState`](@ref), [`SingleState`](@ref),
 [`PMCSimulation`](@ref).
 """
@@ -99,10 +123,13 @@ end
 num_replicas(::ReplicaState{N}) where {N} = N
 num_spectral_states(::ReplicaState{<:Any, S}) where {S} = S
 
-function Base.show(io::IO, st::ReplicaState)
-    print(io, "ReplicaState")
-    print(io, " with ", num_replicas(st), " replicas  and ")
-    print(io, num_spectral_states(st), " spectral states")
+Base.show(io::IO, r::ReplicaState) = show(io, MIME("text/plain"), r)
+function Base.show(io::IO, ::MIME"text/plain", st::ReplicaState)
+    r = num_replicas(st)
+    s = num_spectral_states(st)
+    print(io, "$r×$s Rimu.ReplicaState")
+    print(io, " with ", r, " replica(s) and ")
+    print(io, s, " spectral state(s) of type ", nameof(typeof(st[1])))
     print(io, "\n  H:    ", first(st).hamiltonian)
     print(io, "\n  step: ", st.step[], " / ", st.simulation_plan.last_step)
     print(io, "\n  replicas: ")
@@ -125,17 +152,25 @@ struct StateVectors{V,R} <: AbstractMatrix{V}
 end
 Base.size(sv::StateVectors) = size(sv.state)
 Base.getindex(sv::StateVectors, i::Int, j::Int) = sv.state[i, j].v
+function Base.show(io::IO, ::MIME"text/plain", sv::StateVectors)
+    r = num_replicas(sv.state)
+    s = num_spectral_states(sv.state)
+    print(io, "$r×$s Rimu.StateVectors")
+    print(io, " with $r replica(s) and $s spectral state(s)")
+    print(io, " of type ", nameof(typeof(first(sv.state).v)))
+end
 
 """
     state_vectors(state::ReplicaState)
     state_vectors(sim::PMCSimulation)
-Return an `AbstractMatrix` of configuration vectors from the `state`.
+Return an `r×s` `AbstractMatrix` of configuration vectors from the `state`, or the result of
+[`solve(::ProjectorMonteCarloProblem)`](@ref).
 The vectors can be accessed by indexing the resulting collection, where the row index
 corresponds to the replica index and the column index corresponds to the spectral state
 index.
 
-See also [`SingleState`](@ref), [`ReplicaState`](@ref), [`SpectralState`](@ref),
-[`PMCSimulation`](@ref).
+See also [`ProjectorMonteCarloProblem`](@ref), [`Rimu.PMCSimulation`](@ref),
+[`Rimu.SingleState`](@ref), [`Rimu.ReplicaState`](@ref), [`Rimu.SpectralState`](@ref).
 """
 function state_vectors(state::R) where {R<:ReplicaState}
     V = typeof(first(state).v)
diff --git a/src/strategies_and_params/reportingstrategy.jl b/src/strategies_and_params/reportingstrategy.jl
index cdbea4ed3..c30eaae13 100644
--- a/src/strategies_and_params/reportingstrategy.jl
+++ b/src/strategies_and_params/reportingstrategy.jl
@@ -175,22 +175,25 @@ end
 """
     ReportingStrategy
 
-Abstract type for strategies for reporting data in a DataFrame with [`report!()`](@ref).
+Abstract type for strategies for reporting data during a simulation of a
+[`ProjectorMonteCarloProblem`](@ref).
+
 
 # Implemented strategies:
 
 * [`ReportDFAndInfo`](@ref)
 * [`ReportToFile`](@ref)
 
-# Interface:
+# Extended help
+## Interface:
 
 A `ReportingStrategy` can define any of the following:
 
-* [`refine_reporting_strategy`](@ref)
-* [`report!`](@ref)
-* [`report_after_step!`](@ref)
-* [`finalize_report!`](@ref)
-* [`reporting_interval`](@ref)
+* [`Rimu.refine_reporting_strategy`](@ref)
+* [`Rimu.report!`](@ref)
+* [`Rimu.report_after_step!`](@ref)
+* [`Rimu.finalize_report!`](@ref)
+* [`Rimu.reporting_interval`](@ref)
 
 """
 abstract type ReportingStrategy end
@@ -208,9 +211,9 @@ refine_reporting_strategy(reporting_strategy::ReportingStrategy) = reporting_str
      report!(::ReportingStrategy, step, report::Report, nt, id="")
 
 Report `keys` and `values` to `report`, which will be converted to a `DataFrame` before
-[`ProjectorMonteCarloProblem`](@ref) exits. Alternatively, a `nt::NamedTuple` can be passed in place of `keys`
-and `values`. If `id` is specified, it is appended to all `keys`. This is used to
-differentiate between values reported by different replicas.
+[`ProjectorMonteCarloProblem`](@ref) exits. Alternatively, a `nt::NamedTuple` can be passed
+in place of `keys` and `values`. If `id` is specified, it is appended to all `keys`. This is
+used to differentiate between values reported by different replicas.
 
 To overload this function for a new `ReportingStrategy`, overload
 `report!(::ReportingStrategy, step, args...)` and apply the report by calling
@@ -224,8 +227,8 @@ end
 """
     report_after_step!(::ReportingStrategy, step, report, state) -> report
 
-This function is called at the very end of a step, after [`reporting_interval`](@ref) steps.
-It may modify the `report`.
+This function is called at the very end of a step, after [`Rimu.reporting_interval`](@ref)
+steps. It may modify the `report`.
 """
 function report_after_step!(::ReportingStrategy, _, report, args...)
     return report
@@ -261,10 +264,13 @@ end
 """
     ReportDFAndInfo(; reporting_interval=1, info_interval=100, io=stdout, writeinfo=false) <: ReportingStrategy
 
-The default [`ReportingStrategy`](@ref). Report every `reporting_interval`th step to a `DataFrame`
-and write info message to `io` every `info_interval`th reported step (unless `writeinfo == false`). The flag
-`writeinfo` is useful for controlling info messages in MPI codes, e.g. by setting
+The default [`ReportingStrategy`](@ref) for [`ProjectorMonteCarloProblem`](@ref). Report
+every `reporting_interval`th step to a `DataFrame` and write info message to `io` every
+`info_interval`th reported step (unless `writeinfo == false`). The flag `writeinfo` is
+useful for controlling info messages in MPI codes, e.g. by setting
 `writeinfo = `[`is_mpi_root()`](@ref Rimu.is_mpi_root).
+
+See also [`ProjectorMonteCarloProblem`](@ref), [`ReportToFile`](@ref).
 """
 @with_kw struct ReportDFAndInfo <: ReportingStrategy
     reporting_interval::Int = 1
@@ -288,9 +294,10 @@ end
 """
     ReportToFile(; kwargs...) <: ReportingStrategy
 
-[`ReportingStrategy`](@ref) that writes the report directly to a file in the
-[`Arrow`](https://arrow.apache.org/julia/dev/) format. Useful when dealing with long
-jobs or large numbers of replicas, when the report can incur a significant memory cost.
+[`ReportingStrategy`](@ref) for [`ProjectorMonteCarloProblem`](@ref) that writes the report
+directly to a file in the [`Arrow`](https://arrow.apache.org/julia/dev/) format. Useful when
+dealing with long jobs or large numbers of replicas, when the report can incur a significant
+memory cost.
 
 The arrow file can be read back in with [`load_df(filename)`](@ref) or
 `using Arrow; Arrow.Table(filename)`.
@@ -311,7 +318,8 @@ The arrow file can be read back in with [`load_df(filename)`](@ref) or
   messages printed out.
 * `compress = :zstd`: compression algorithm to use. Can be `:zstd`, `:lz4` or `nothing`.
 
-See also [`load_df`](@ref) and [`save_df`](@ref).
+See also [`load_df`](@ref), [`save_df`](@ref), [`ReportDFAndInfo`](@ref), and
+[`ProjectorMonteCarloProblem`](@ref).
 """
 mutable struct ReportToFile{C} <: ReportingStrategy
     filename::String
diff --git a/test/projector_monte_carlo_problem.jl b/test/projector_monte_carlo_problem.jl
index 5010c49cb..7c2e0fe90 100644
--- a/test/projector_monte_carlo_problem.jl
+++ b/test/projector_monte_carlo_problem.jl
@@ -95,9 +95,14 @@ end
         @test sm.state[2].shift_parameters.shift ≡ 2.0
         @test state_vectors(sm.state)[1][BoseFS(1, 3)] == 10
         @test state_vectors(sm.state)[2][BoseFS(3, 1)] == 10
-        @test startswith(sprint(show, sm.state), "ReplicaState")
-        @test startswith(sprint(show, sm.state[1]), "SingleState")
-        @test startswith(sprint(show, sm.state.spectral_states), "(SpectralState")
+        @test startswith(sprint(show, sm.state), "2×1 Rimu.ReplicaState")
+        @test startswith(sprint(show, sm.state[1]), "Rimu.SingleState")
+        @test startswith(sprint(show, sm.state.spectral_states),
+            "(1-element Rimu.SpectralState"
+        )
+        @test startswith(sprint(show, sm.state.spectral_states[1]),
+            "1-element Rimu.SpectralState"
+        )
     end
 
     @testset "Default DVec" begin

From d3dbef42eed2eb652050be5a7d650feea0567ef9 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Thu, 26 Dec 2024 17:59:17 +1300
Subject: [PATCH 15/17] append replica id for row index, then spectral id for
 column index

---
 src/pmc_simulation.jl | 17 ++++++++++-------
 src/qmc_states.jl     | 18 ++++++++----------
 2 files changed, 18 insertions(+), 17 deletions(-)

diff --git a/src/pmc_simulation.jl b/src/pmc_simulation.jl
index f45b0de0e..4a545b9a3 100644
--- a/src/pmc_simulation.jl
+++ b/src/pmc_simulation.jl
@@ -124,23 +124,26 @@ function PMCSimulation(problem::ProjectorMonteCarloProblem; copy_vectors=true)
     # set up the spectral_states
     wm = working_memory(vectors[1, 1])
     spectral_states = ntuple(n_replicas) do i
+        replica_id = if n_replicas == 1
+            ""
+        else
+            "_$(i)"
+        end
         SpectralState(
             ntuple(n_spectral) do j
                 v = vectors[i, j]
                 sp = shift_parameters[i, j]
-                id = if n_replicas * n_spectral == 1
+                spectral_id = if n_spectral == 1
                     ""
-                elseif n_spectral == 1
-                    "_$(i)"
                 else
-                    "_s$(j)_$(i)" # j is the spectral state index, i is the replica index
-                end # we have to think about how to label the spectral states
+                    "_s$(j)" # j is the spectral state index, i is the replica index
+                end
                 SingleState(
                     hamiltonian, algorithm, v, zerovector(v),
                     wm isa PDWorkingMemory ? wm : working_memory(v), # reuse for PDVec
-                    sp, id
+                    sp, replica_id * spectral_id
                 )
-            end, spectral_strategy)
+            end, spectral_strategy, replica_id)
     end
     @assert spectral_states isa NTuple{n_replicas, <:SpectralState}
 
diff --git a/src/qmc_states.jl b/src/qmc_states.jl
index 216431679..acfc04860 100644
--- a/src/qmc_states.jl
+++ b/src/qmc_states.jl
@@ -39,17 +39,17 @@ end
 
 """
     SpectralState <: AbstractVector{SingleState}
-Holds one or several [`SingleState`](@ref)s representing the ground state and excited
+Holds one or several [`Rimu.SingleState`](@ref)s representing the ground state and excited
 states of a single replica.
 Indexing the `SpectralState` `state[i]` returns the `i`th `SingleState`.
 
 ## Fields
 - `single_states`: Tuple of `SingleState`s
 - `spectral_strategy`: Strategy for computing the spectral states
-- `id::String`: id is appended to column names
+- `id::String`: Identifies the replica
 
-See also [`SpectralStrategy`](@ref), [`ReplicaState`](@ref), [`SingleState`](@ref),
-[`PMCSimulation`](@ref).
+See also [`SpectralStrategy`](@ref), [`Rimu.ReplicaState`](@ref), [`Rimu.SingleState`](@ref),
+[`Rimu.PMCSimulation`](@ref).
 """
 struct SpectralState{
     N,
@@ -58,10 +58,7 @@ struct SpectralState{
 } <: AbstractVector{SingleState}
     single_states::NS # Tuple of SingleState
     spectral_strategy::SS # SpectralStrategy
-    id::String # id is appended to column names
-end
-function SpectralState(t::Tuple, ss::SpectralStrategy, id="")
-    return SpectralState(t, ss, id)
+    id::String # identifies the replica and is appended to row names
 end
 num_spectral_states(::SpectralState{N}) where {N} = N
 
@@ -73,6 +70,7 @@ function Base.show(io::IO, ::MIME"text/plain", s::SpectralState)
     ns = num_spectral_states(s)
     print(io, "$ns-element Rimu.SpectralState")
     print(io, " with ", ns, " spectral state(s) of type ", nameof(typeof(s[1])))
+    print(io, " and id \"", s.id, "\"")
     print(io, "\n    spectral_strategy: ", s.spectral_strategy)
     for (i, r) in enumerate(s.single_states)
         print(io, "\n      $i: ", r)
@@ -99,8 +97,8 @@ and `j`th spectral state.
 - `post_step_strategy`: Post-step strategy
 - `replica_strategy`: Replica strategy
 
-See also [`ReplicaStrategy`](@ref), [`SpectralState`](@ref), [`SingleState`](@ref),
-[`PMCSimulation`](@ref).
+See also [`ReplicaStrategy`](@ref), [`Rimu.SpectralState`](@ref), [`Rimu.SingleState`](@ref),
+[`Rimu.PMCSimulation`](@ref).
 """
 struct ReplicaState{
     N, # number of replicas

From 4e3f531ebff2aa54116f3d866cb380d07240b844 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Thu, 26 Dec 2024 18:02:51 +1300
Subject: [PATCH 16/17] fix tests

---
 test/excited_states_tests.jl | 20 ++++++++++----------
 1 file changed, 10 insertions(+), 10 deletions(-)

diff --git a/test/excited_states_tests.jl b/test/excited_states_tests.jl
index dc097bac0..5d2922cc1 100644
--- a/test/excited_states_tests.jl
+++ b/test/excited_states_tests.jl
@@ -13,9 +13,9 @@ using Test
     p = ProjectorMonteCarloProblem(ham; spectral_strategy, last_step, style)
     s = solve(p)
     df = DataFrame(s)
-    energy1 = shift_estimator(df, shift="shift_s1_1", skip=1000)
-    energy2 = shift_estimator(df, shift="shift_s2_1", skip=1000)
-    energy3 = shift_estimator(df, shift="shift_s3_1", skip=1000)
+    energy1 = shift_estimator(df, shift="shift_s1", skip=1000)
+    energy2 = shift_estimator(df, shift="shift_s2", skip=1000)
+    energy3 = shift_estimator(df, shift="shift_s3", skip=1000)
 
     @test energy1.mean ≈ vals[1]
     @test energy2.mean ≈ vals[2]
@@ -25,13 +25,13 @@ using Test
     p = ProjectorMonteCarloProblem(ham; spectral_strategy, last_step, style, n_replicas)
     s = solve(p)
     df = DataFrame(s)
-    energy1 = shift_estimator(df, shift="shift_s1_1", skip=1000)
-    energy2 = shift_estimator(df, shift="shift_s2_1", skip=1000)
-    energy3 = shift_estimator(df, shift="shift_s3_1", skip=1000)
-    energy4 = shift_estimator(df, shift="shift_s1_2", skip=1000)
-    energy5 = shift_estimator(df, shift="shift_s2_2", skip=1000)
-    energy6 = shift_estimator(df, shift="shift_s3_2", skip=1000)
-    
+    energy1 = shift_estimator(df, shift="shift_1_s1", skip=1000)
+    energy2 = shift_estimator(df, shift="shift_1_s2", skip=1000)
+    energy3 = shift_estimator(df, shift="shift_1_s3", skip=1000)
+    energy4 = shift_estimator(df, shift="shift_2_s1", skip=1000)
+    energy5 = shift_estimator(df, shift="shift_2_s2", skip=1000)
+    energy6 = shift_estimator(df, shift="shift_2_s3", skip=1000)
+
     @test energy1.mean ≈ vals[1]
     @test energy2.mean ≈ vals[2]
     @test energy3.mean ≈ vals[3]

From 1a2270ef9ed21d008c6bf99b1c9b71e2a4335a39 Mon Sep 17 00:00:00 2001
From: Joachim Brand <joachim.brand@gmail.com>
Date: Thu, 26 Dec 2024 22:28:38 +1300
Subject: [PATCH 17/17] test show for StateVectors

---
 src/qmc_states.jl                     | 2 ++
 test/projector_monte_carlo_problem.jl | 1 +
 2 files changed, 3 insertions(+)

diff --git a/src/qmc_states.jl b/src/qmc_states.jl
index acfc04860..9e6fe17bf 100644
--- a/src/qmc_states.jl
+++ b/src/qmc_states.jl
@@ -150,6 +150,8 @@ struct StateVectors{V,R} <: AbstractMatrix{V}
 end
 Base.size(sv::StateVectors) = size(sv.state)
 Base.getindex(sv::StateVectors, i::Int, j::Int) = sv.state[i, j].v
+
+Base.show(io::IO, sv::StateVectors) = show(io, MIME("text/plain"), sv)
 function Base.show(io::IO, ::MIME"text/plain", sv::StateVectors)
     r = num_replicas(sv.state)
     s = num_spectral_states(sv.state)
diff --git a/test/projector_monte_carlo_problem.jl b/test/projector_monte_carlo_problem.jl
index 7c2e0fe90..7ef4d9465 100644
--- a/test/projector_monte_carlo_problem.jl
+++ b/test/projector_monte_carlo_problem.jl
@@ -103,6 +103,7 @@ end
         @test startswith(sprint(show, sm.state.spectral_states[1]),
             "1-element Rimu.SpectralState"
         )
+        @test startswith(sprint(show, state_vectors(sm)), "2×1 Rimu.StateVectors")
     end
 
     @testset "Default DVec" begin