-
Notifications
You must be signed in to change notification settings - Fork 145
Inflation Testing
You need a case where we expect the inflation to do something. Since all the models use the same inflation code, it doesn't matter which model. The test should be fast. Since some of these algorithms adapt, the test should have multiple assimilation cycles. All locations are treated similarly, so any location will do. If the localization is large enough to impact more than one location - this is a sufficient test to ensure that spatially-varying inflations vary spatially and spatially-constant inflations are spatially constant.
The following scripts/files can be used to exhaustively test all the inflation
algorithms that we support - as both prior and posterior inflation. At present,
all tests are repeated on the Manhattan branch and the feature branch - but this
requires two separate clones, and the mkmf.template file should be identical on
both clones (each clone has one of the branches). The test_inflation.csh and
CheckInflation.m will require customization to reflect the location of your
DART branches.
There is a duality between the inflation flavor specified as an integer and that specified by a string.
| character string | integer |
|---|---|
| 'NO_INFLATION' | 0 |
| 'OBS_INFLATION' | 1 |
| 'VARYING_SS_INFLATION' | 2 |
| 'SINGLE_SS_INFLATION' | 3 |
| 'RELAXATION_TO_PRIOR_SPREAD' | 4 |
| 'ENHANCED_SS_INFLATION' | 5 |
I chose to use the Lorenz 96 model for the test. I created a single identity observation (-10) and repeated it for 3 cycles - since time_step_seconds is 3600 seconds (i.e. times 0,0 ... 3600,0 ... 7200,0) ran it through perfect_model_obs, and then manually biased the observation to force the separation of the observation value and the ensemble. I also changed the model forcing (pmo was run with 8.0, filter with 9.0) for the same reason. I left the cutoff at 0.02 which results in impacting 3 spatial locations (#9,#10,#11). I used the model states that result from running workshop_setup.csh or quickbuild.csh.
To recap the input.nml I was using for the inflation testing:
inf_initial = 1.1, 1.2
inf_sd_initial = 0.6, 0.6
inf_initial = 1.1, 1.2
inf_lower_bound = 0.1, 0.1
inf_upper_bound = 1000000.0, 1000000.0
inf_damping = 0.9, 0.9
inf_sd_initial = 0.6, 0.6
inf_sd_lower_bound = 0.5, 0.5
inf_sd_max_change = 1.05, 1.05
...
cutoff = 0.02
...
forcing = 9.00
Put another way:
diff inf_3_4/input.nml ../work/input.nml
55c55
< obs_sequence_in_name = 'obs_seq.out.inf_test'
---
> obs_sequence_in_name = "obs_seq.out",
65c65
< inf_flavor = 3, 4
---
> inf_flavor = 2, 0,
69,75c69,75
< inf_initial = 1.1, 1.2
< inf_lower_bound = 0.0, 0.0
< inf_upper_bound = 1000000.0, 1000000.0
< inf_damping = 0.9, 0.9
< inf_sd_initial = 0.6, 0.6
< inf_sd_lower_bound = 0.5, 0.5
< inf_sd_max_change = 1.05, 1.05
---
> inf_initial = 1.0, 1.0,
> inf_lower_bound = 1.0, 1.0,
> inf_upper_bound = 1000000.0, 1000000.0,
> inf_damping = 0.9, 1.0,
> inf_sd_initial = 0.6, 0.0,
> inf_sd_lower_bound = 0.6, 0.0,
> inf_sd_max_change = 1.05, 1.05,
97c97
< cutoff = 0.02
---
> cutoff = 0.02,
133c133
< forcing = 9.0
---
> forcing = 8.00,
I thought there was a way to attach a file to the wiki page, but I guess not ... so here are the three files in question.
0[1539]flex_filter_input ~/<3>lorenz_96/work % cat obs_seq.out.inf_test
obs_sequence
obs_type_definitions
0
num_copies: 2 num_qc: 1
num_obs: 3 max_num_obs: 3
observations
truth
Quality Control
first: 1 last: 3
OBS 1
5.25
5.0897032269825484
0.0000000000000000
-1 2 -1
obdef
loc1d
0.2250000000000000
kind
-10
0 0
1.0000000000000000E-002
OBS 2
5.26
5.1082190492496267
0.0000000000000000
1 3 -1
obdef
loc1d
0.2250000000000000
kind
-10
3600 0
1.0000000000000000E-002
OBS 3
4.20
4.4661680535183397
0.0000000000000000
2 -1 -1
obdef
loc1d
0.2250000000000000
kind
-10
7200 0
1.0000000000000000E-002
and what I've been calling test_inflation.csh:
#!/bin/csh
set Manhattan = /Users/thoar/git/DART_Manhattan/models/lorenz_96
set Test = /Users/thoar/git/DART_development/models/lorenz_96
sed -e "s/obs_type_definitions/obs_kind_definitions/" \
$Test/work/obs_seq.out.inf_test >! $Manhattan/work/obs_seq.out.inf_test
foreach REPO ( $Manhattan $Test )
echo "=============================================================================="
echo "Testing inflation for $REPO"
cd $REPO/work
# Make sure we are starting from the same states
ncgen -o filter_input.nc filter_input.cdl
# Make sure we are starting from current code
# Should make sure that the mkmf.template has same compiler settings, etc.
if (! -e filter ) then
./quickbuild.csh
endif
# Test all possible combinarions of inflation - even impossible ones
# Note ... obs_inflation (1) should gracefully die.
# Note ... prior RTPS (4) should gracefully die.
# Note ... inflation (6) does not exist, should gracefully die.
foreach PRIOR ( 0 1 2 3 4 5 )
foreach POSTERIOR ( 0 1 2 3 4 5 6)
set CASE = "inf_${PRIOR}_${POSTERIOR}"
echo "------------------------------------------------------------------------------"
echo "Testing inflation case $CASE"
mkdir $CASE
cd $CASE
\rm -f input.nml filter_log.txt dart_log.*
cp ../obs_seq.out.inf_test .
cp ../filter_input.nc .
echo "filter_input.nc" >! filter_input_list.txt
echo "filter_output.nc" >! filter_output_list.txt
sed -e "s/obs_sequence_in_name .*/obs_sequence_in_name = 'obs_seq.out.inf_test'/g" \
-e "s/inf_flavor .*/inf_flavor = ${PRIOR}, ${POSTERIOR}/" \
-e "s/inf_initial .*/inf_initial = 1.1, 1.2/" \
-e "s/inf_lower_bound .*/inf_lower_bound = 0.0, 0.0/" \
-e "s/inf_upper_bound .*/inf_upper_bound = 1000000.0, 1000000.0/" \
-e "s/inf_damping .*/inf_damping = 0.9, 0.9/" \
-e "s/inf_sd_initial .*/inf_sd_initial = 0.6, 0.6/" \
-e "s/inf_sd_lower_bound .*/inf_sd_lower_bound = 0.5, 0.5/" \
-e "s/inf_sd_max_change .*/inf_sd_max_change = 1.05, 1.05/" \
-e "s/cutoff .*/cutoff = 0.02/" \
-e "s/forcing .*/forcing = 9.0/" \
../input.nml >! input.nml || exit 1
../filter |& tee filter_log.txt
cd ..
end
end
end
echo "=============================================================================="
echo "Comparing observation-space results"
foreach OUTPUT ( $Test/work/inf_?_?/obs_seq.final )
set GOLD = `echo $OUTPUT | sed -e "s#$Test#$Manhattan#"`
echo "Comparing "
echo $OUTPUT
echo $GOLD
echo
end
and the Matlab script to compare the netCDF files that have the states as well
as the inflation values. It relies on the DART matlab function Compare_netCDF_files
so the DART diagnostics/matlab directory must be part of your matlabpath.
0[1218]flex_filter_input ~/<3>lorenz_96/matlab % cat CheckInflation.m
% Easy way to compare the results of two inflation experiments.
%
% The assumption is that you have configured and run the 'test_inflation.csh'
% script in two directories (repeated here).
% This only tests the state-space results since test_inflation.csh compares
% the obs_seq.final files.
%
% I suppose we could run the test_inflation.csh script from inside this
% script by using the 'system' function. left for a later date ...
% or we could modify the test_inflation.csh script to call the
gold = '/Users/thoar/git/DART_Manhattan/models/lorenz_96/work';
test = '/Users/thoar/git/DART_development/models/lorenz_96/work';
for prior = 0:5
for posterior = 0:6
directory = sprintf('inf_%d_%d',prior,posterior);
fprintf('Comparing inflation case %s \n',directory)
fnamegold = sprintf('%s/%s/analysis.nc',gold,directory);
fnametest = sprintf('%s/%s/analysis.nc',test,directory);
if ( exist(fnamegold,'file') == 2 & exist(fnametest,'file') == 2 )
Compare_netCDF_files(fnamegold,fnametest)
else
fprintf('WARNING: %s did not create output.\n',directory)
end
fprintf('\n')
end
end