Merge branch 'develop' into feature/save-state

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.

src/StatsTools/growth_witness.jl

@@ -51,7 +51,7 @@ function growth_witness(
     shift=:shift, norm=:norm, time_step=nothing, kwargs...
 )
     df = DataFrame(sim)
-    time_step = determine_constant_time_step(df)
+    time_step = isnothing(time_step) ? determine_constant_time_step(df) : time_step
 
     shift_vec = getproperty(df, Symbol(shift))
     norm_vec = getproperty(df, Symbol(norm))

src/StatsTools/reweighting.jl

@@ -20,7 +20,7 @@ function determine_constant_time_step(df)
         elseif hasproperty(df, "dτ_1")
             return df.dτ_1[end]
         else
-            throw(ArgumentError("Time step not found in `df`"))
+            throw(ArgumentError("key `\"time_step\"` not found in `df` metadata"))
         end
     else
         throw(ArgumentError("Time step not constant"))
@@ -219,6 +219,7 @@ Returns a `NamedTuple` with the fields
 * `threading = Threads.nthreads() > 1`: if `false` a progress meter is displayed
 * `shift_name = :shift` name of column in `df` with shift data
 * `norm_name = :norm` name of column in `df` with walkernumber data
+* `time_step = determine_constant_time_step(df)` the time step
 * `warn = true` whether to log warning messages when blocking fails or denominators are
   small
 
@@ -242,14 +243,15 @@ function growth_estimator_analysis(
     threading=Threads.nthreads() > 1,
     shift_name=:shift,
     norm_name=:norm,
+    time_step=nothing,
     warn=true,
     kwargs...
 )
     df = DataFrame(sim)
     shift_v = Vector(getproperty(df, Symbol(shift_name))) # casting to `Vector` to make SIMD loops efficient
     norm_v = Vector(getproperty(df, Symbol(norm_name)))
     num_reps = length(filter(startswith("norm"), names(df)))
-    time_step = determine_constant_time_step(df)
+    time_step = isnothing(time_step) ? determine_constant_time_step(df) : time_step
     se = blocking_analysis(shift_v; skip)
     E_r = se.mean
     correlation_estimate = 2^(se.k - 1)
@@ -388,6 +390,7 @@ Returns a `NamedTuple` with the fields
 * `shift_name = :shift` name of column in `df` with shift data
 * `hproj_name = :hproj` name of column in `df` with operator overlap data
 * `vproj_name = :vproj` name of column in `df` with projector overlap data
+* `time_step = determine_constant_time_step(df)` the time step
 * `warn = true` whether to log warning messages when blocking fails or denominators are small
 
 ## Example
@@ -412,14 +415,15 @@ function mixed_estimator_analysis(
     hproj_name=:hproj,
     vproj_name=:vproj,
     warn=true,
+    time_step=nothing,
     kwargs...
 )
     shift_v = Vector(getproperty(df, Symbol(shift_name))) # casting to `Vector` to make SIMD loops efficient
     hproj_v = Vector(getproperty(df, Symbol(hproj_name)))
     vproj_v = Vector(getproperty(df, Symbol(vproj_name)))
     num_reps = length(filter(startswith("norm"), names(df)))
 
-    time_step = determine_constant_time_step(df)
+    time_step = isnothing(time_step) ? determine_constant_time_step(df) : time_step
     se = blocking_analysis(shift_v; skip)
     E_r = se.mean
     correlation_estimate = 2^(se.k - 1)
@@ -555,11 +559,12 @@ function rayleigh_replica_estimator(
     h=0,
     skip=0,
     Anorm=1,
+    time_step=nothing,
     kwargs...
 )
     df = DataFrame(sim)
     num_reps = length(filter(startswith("norm"), names(df)))
-    time_step = determine_constant_time_step(df)
+    time_step = isnothing(time_step) ? determine_constant_time_step(df) : time_step
     T = eltype(df[!, Symbol(shift_name, "_1")])
     shift_v = Vector{T}[]
     for a in 1:num_reps
@@ -625,11 +630,12 @@ function rayleigh_replica_estimator_analysis(
     vec_name="dot",
     Anorm=1,
     warn=true,
+    time_step=nothing,
     kwargs...
 )
     df = DataFrame(sim)
     num_reps = length(filter(startswith("norm"), names(df)))
-    time_step = determine_constant_time_step(df)
+    time_step = isnothing(time_step) ? determine_constant_time_step(df) : time_step
     # estimate the correlation time by blocking the shift data
     T = eltype(df[!, Symbol(shift_name, "_1")])
     shift_v = Vector{T}[]

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.maxlength[]
-        if length(state.spectral_states) > 1
-            @error "`maxlength` reached in single state $(s_state.id). Aborting."
-        else
-            @error "`maxlength` reached. Aborting."
-        end
+    if len > state.max_length[]
+        @error "`max_length` reached in state $(s_state.id). Aborting."
         return false
     end
     return proceed # Bool

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)
@@ -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
@@ -156,7 +158,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
@@ -193,15 +195,18 @@ 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
-    @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 +218,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,

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).
 """
@@ -80,7 +90,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)
 
@@ -114,30 +124,33 @@ 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}
 
     # set up the initial state
     state = ReplicaState(
         spectral_states,
-        Ref(maxlength),
+        Ref(max_length),
         Ref(simulation_plan.starting_step),
         simulation_plan,
         reporting_strategy,
@@ -242,7 +255,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 +313,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 +323,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 +349,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 +365,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 +433,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