11using ParticleDA
22using ParticleDA. Kalman
3- using HDF5, LinearAlgebra, MPI, PDMats, Random, StableRNGs, Statistics, Test, YAML
3+ using HDF5, LinearAlgebra, MPI, PDMats, Random, StableRNGs, Statistics, Test, YAML
44
55include(joinpath(@__DIR__, " models" , " llw2d.jl" ))
66include(joinpath(@__DIR__, " models" , " lorenz63.jl" ))
102102 check_hdf5_group_valid(file, group, group_name)
103103 @test read(group, dataset_name) == test_array
104104 @test all([
105- read_attribute(group[dataset_name], k) == test_attributes[k]
105+ read_attribute(group[dataset_name], k) == test_attributes[k]
106106 for k in keys(test_attributes)
107107 ])
108108 # test writing to existing dataset name results in warning and does not update
117117 # test writing timer data
118118 timer_strings = [" ab" , " cde" , " fg" , " hij" ]
119119 ParticleDA. write_timers(
120- map(length, timer_strings),
121- length(timer_strings),
122- codeunits(join(timer_strings)),
120+ map(length, timer_strings),
121+ length(timer_strings),
122+ codeunits(join(timer_strings)),
123123 output_filename
124124 )
125125 h5open(output_filename, " cw" ) do file
134134
135135
136136@testset (
137- " Generic model interface unit tests - $(parentmodule(typeof(model))) "
137+ " Generic model interface unit tests - $(parentmodule(typeof(model))) "
138138) for model in (
139139 LLW2d. init(Dict()),
140140 Lorenz63. init(Dict()),
161161 " x_length" => 100e3 ,
162162 " y_length" => 100e3 ,
163163 " station_boundary_x" => 30e3 ,
164- " station_boundary_y" => 30e3 ,
164+ " station_boundary_y" => 30e3 ,
165165 )
166166 )
167167 ),
184184 config. model,
185185 seed,
186186 estimate_bound_constant,
187- config. estimate_n_samples;
187+ config. estimate_n_samples;
188188 RNGType= StableRNG
189189 )
190190end
225225 avg= mean(states; dims= 2 ), var= var(states, corrected= true ; dims= 2 )
226226 )
227227 for statistics_type in (
228- ParticleDA. NaiveMeanSummaryStat,
228+ ParticleDA. NaiveMeanSummaryStat,
229229 ParticleDA. NaiveMeanAndVarSummaryStat,
230230 ParticleDA. MeanSummaryStat,
231231 ParticleDA. MeanAndVarSummaryStat
@@ -281,12 +281,12 @@ end
281281 new_states = copy(states)
282282 Random. seed!(rng, seed)
283283 ParticleDA. sample_proposal_and_compute_log_weights!(
284- new_states,
285- log_weights,
286- observation,
287- time_index,
288- model,
289- filter_data,
284+ new_states,
285+ log_weights,
286+ observation,
287+ time_index,
288+ model,
289+ filter_data,
290290 filter_type,
291291 rng
292292 )
@@ -301,12 +301,12 @@ end
301301 # Check filter update gives deterministic updates when rng state is fixed
302302 Random. seed!(rng, seed)
303303 ParticleDA. sample_proposal_and_compute_log_weights!(
304- new_states_2,
305- log_weights_2,
306- observation,
307- time_index,
308- model,
309- filter_data,
304+ new_states_2,
305+ log_weights_2,
306+ observation,
307+ time_index,
308+ model,
309+ filter_data,
310310 filter_type,
311311 rng
312312 )
@@ -320,10 +320,10 @@ end
320320 n_task = 1
321321 model_params_dict = Dict(
322322 " llw2d" => Dict(
323- " nx" => 32 ,
324- " ny" => 32 ,
325- " n_stations_x" => 4 ,
326- " n_stations_y" => 4 ,
323+ " nx" => 32 ,
324+ " ny" => 32 ,
325+ " n_stations_x" => 4 ,
326+ " n_stations_y" => 4 ,
327327 " padding" => 0
328328 )
329329 )
342342 online_matrices = ParticleDA. init_online_matrices(model, 1 )
343343 for f in nfields(online_matrices)
344344 matrix = getfield(online_matrices, f)
345- @test isa(matrix, AbstractMatrix)
346- end
345+ @test isa(matrix, AbstractMatrix)
346+ end
347347 state_dimension = ParticleDA. get_state_dimension(model)
348348 updated_indices = ParticleDA. get_state_indices_correlated_to_observations(
349349 model
@@ -383,10 +383,10 @@ end
383383 # X = F(x) + U and Y = HX + V where x is the previous state value, F the forward
384384 # operator for the deterministic state dynamics, U ~ Normal(0, Q) the additive
385385 # state noise, X the state at the next time step, H the linear observation
386- # operator, V ~ Normal(0, R) the additive observation noise and Y the modelled
387- # observations, is Normal(m, C) where
386+ # operator, V ~ Normal(0, R) the additive observation noise and Y the modelled
387+ # observations, is Normal(m, C) where
388388 # m = F(x) + QHᵀ(HQHᵀ + R)⁻¹(y − HF(x))
389- # = F(x) + cov(X, Y) @ cov(Y, Y)⁻¹ (y − HF(x))
389+ # = F(x) + cov(X, Y) @ cov(Y, Y)⁻¹ (y − HF(x))
390390 # and C = Q − QHᵀ(HQHᵀ + R)⁻¹HQ = cov(X, X) - cov(X, Y) cov(Y, Y)⁻¹ cov(X, Y)ᵀ
391391 analytic_mean = copy(state)
392392 @view(analytic_mean[updated_indices]) .+ = (
406406 # Create set of state ' particles' all equal to propagated state
407407 states = Matrix{ParticleDA.get_state_eltype(model)}(
408408 undef, (state_dimension, nprt)
409- )
409+ )
410410 states .= state
411411 updated_states = copy(states)
412412 for state in eachcol(updated_states)
@@ -416,21 +416,21 @@ end
416416 # Mean of noise added by update_particle_noise! should be zero in all components
417417 # and empirical mean should therefore be zero to within Monte Carlo error. The
418418 # constant in the tolerance below was set by looking at scale of typical
419- # deviation, the point of check is that errors scale at expected O(1/√N) rate.
420- @test maximum(abs.(mean(noise, dims=2))) < (10. / sqrt(nprt))
419+ # deviation, the point of check is that errors scale at expected O(1/√N) rate.
420+ @test maximum(abs.(mean(noise, dims=2))) < (10. / sqrt(nprt))
421421 # Covariance of noise added by update_particle_noise! to observed state
422422 # components should be cov_X_X as computed above and empirical covariance of
423423 # these components should therefore be within Monte Carlo error of cov_X_X. The
424424 # constant in tolerance below was set by looking at scale of typical deviations,
425- # the point of check is that errors scale at expected O(1/√N) rate.
425+ # the point of check is that errors scale at expected O(1/√N) rate.
426426 noise_cov = cov(noise, dims=2)
427- @test maximum(abs.(noise_cov .- cov_X_X)) < (10. / sqrt(nprt))
427+ @test maximum(abs.(noise_cov .- cov_X_X)) < (10. / sqrt(nprt))
428428 ParticleDA.update_states_given_observations!(
429429 updated_states, observation, model, filter_data, rng
430430 )
431431 updated_mean = mean(updated_states, dims=2)
432432 updated_cov = cov(updated_states, dims=2)
433- # Mean and covariance of updated particles should be within O(1/√N) Monte Carlo
433+ # Mean and covariance of updated particles should be within O(1/√N) Monte Carlo
434434 # error of analytic values - constants in tolerances were set by looking at
435435 # scale of typical deviations, main point of checks are that errors scale at
436436 # expected O(1/√N) rate.
@@ -475,24 +475,24 @@ end
475475end
476476
477477
478- @testset "Integration test -- $(input_file) with $(filter_type) and $(stat_type)" for
478+ @testset "Integration test -- $(input_file) with $(filter_type) and $(stat_type)" for
479479 filter_type in (ParticleDA.BootstrapFilter, ParticleDA.OptimalFilter),
480480 stat_type in (ParticleDA.MeanSummaryStat, ParticleDA.MeanAndVarSummaryStat),
481481 input_file in ["integration_test_$i.yaml" for i in 1:6]
482482 observation_file_path = tempname()
483483 ParticleDA.simulate_observations_from_model(
484- LLW2d.init,
485- joinpath(@__DIR__, input_file),
484+ LLW2d.init,
485+ joinpath(@__DIR__, input_file),
486486 observation_file_path
487487 )
488488 observation_sequence = h5open(
489489 ParticleDA.read_observation_sequence, observation_file_path, "r"
490490 )
491491 @test !any(isnan.(observation_sequence))
492492 states, statistics = ParticleDA.run_particle_filter(
493- LLW2d.init,
494- joinpath(@__DIR__, input_file),
495- observation_file_path,
493+ LLW2d.init,
494+ joinpath(@__DIR__, input_file),
495+ observation_file_path,
496496 filter_type,
497497 stat_type,
498498 )
0 commit comments