Changes between Initial Version and Version 1 of Lorenz_63_model

Timestamp:

Jul 4, 2019, 3:19:56 PM (5 years ago)

Author:

lnerger

Comment:

--

Lorenz_63_model

@@ +1 @@
+= Lorenz-63 model with PDAF = 
+
+[[PageOutline(2-3,Contents of this page)]]
+
+The implementation of the Lorenz-63 model coupled to PDAF is in the directory `testsuite/src/lorenz63/` of the PDAF package.
+Provided is a full implementation of PDAF with the nonlinear Lorenz63 model (E. N. Lorenz (1963) Deterministic non-periodic flows. J. Atmos. Sci. 20, 130-141) providing various filter and smoother methods. The model has state dimenion 3. So it's not usable for
+localized filters, but it's a good test case for the particle filter.
+
+Next to the implementation of the Lorenz-63 model with PDAF, the test case provides tool programs and scripts that
+allow to run a test case and to display the outputs. Note, that this implementation runs without parallelization.
+
+The data assimilation with the Lorenz-63 model has been added in version 1.14 of PDAF.
+
+
+== Running the test case  ==
+
+Runnning a data assimilation experiment with the Lorenz-63 model is a two step process: First one runs the model without PDAF to generate a file holding the trajectory of a forward run. Then one generates files with observations and a covariance matrix for the initialization of the initial ensemble. In the second step, one compiles the Lorenz-63 model with activated coupling to PDAF and runs the data assimilation experiments. 
+
+=== 1. Compile and run the forward model without assimilation ===
+
+First change in the file `make.arch/linux_gfortran.h` the line
+{{{
+CPP_DEFS = -DUSE_PDAF
+}}}
+to 
+{{{
+CPP_DEFS = #-DUSE_PDAF
+}}}
+to deactivate the coupling to PDAF in the model code. Now build the forward model program with
+{{{
+  cd testsuite/src
+  make lorenz_63 PDAF_ARCH=linux_gfortran
+}}}
+in the directory `testsuite/src/` of the PDAF package. You have to ensure
+that in the machine-specific make include file `linux_gfortran.h` `-DUSE_PDAF` is not defined
+for CPP_DEFS (such that calls to PDAF are not active). You can replace `linux_gfortran` by any other make include file from `make.arch/`, e.g. specify `macos_gfortran` for compiling on MacOS. The executable
+is generated in `testsuite/bin/`.
+
+'''Note''': The implementation uses the NetCDF library for file outputs. If the compilation above fails, please ensure the netcdf-library ist installed. On computers running Linux, it is usually available as a package of the operating system. On MacOS one can install the netcdf library e.g. using Fink or !MacPorts. NetCDF is a self-describing binary output format, but here it is not required that you know details about it. Anyway, if you like to look 'into' a NetCDF file, please try to use `ncdump FILENAME | less`.  
+
+To run the forward model use
+{{{
+  cd ../bin
+  ./lorenz_63 -total_steps 10000
+}}}
+This runs the Lorenz-63 model for 10000 time steps and the trajectory is written into a file `state_l63.nc`. 
+
+=== 2. Generate observations ===
+
+To build the executable to generate observations from the state trajectore use
+{{{
+  cd ../src/lorenz63/tools
+  make generate_obs PDAF_ARCH=linux_gfortran
+}}}
+
+Now run
+{{{
+  ./generate_obs
+}}}
+to generate a file holding observations (`obs_l63.nc` in testsuite/bin/).
+
+
+=== 3. Build and run the assimilation program ===
+
+Now do
+{{{
+  cd ../../../../make.arch
+}}}
+and change in the file `linux_gfortran.h` the line
+{{{
+CPP_DEFS = #-DUSE_PDAF
+}}}
+back to 
+{{{
+CPP_DEFS = -DUSE_PDAF
+}}}
+to activate the coupling to PDAF in the model code. 
+
+Now compile the Lorenz-63 model with activated
+coupling to PDAF. First clean the directories for the main driver and the Lorenz-63 model using
+{{{
+  cd ../testsuite/src
+  make cleandriver PDAF_ARCH=linux_gfortran
+  make cleanlorenz_63 PDAF_ARCH=linux_gfortran
+}}}
+(This removes object files that were compiled without support for PDAF)
+Then build the executable using
+{{{
+  make pdaf_lorenz_63 PDAF_ARCH=linux_gfortran
+}}} 
+The program `pdaf_lorenz_63` is generated in testsuite/bin.
+
+
+To run the assimilation program, do
+{{{
+  cd ../bin
+  ../src/lorenz63/tools/run_ESTKF.sh 
+}}}
+The script run_ESTKF.sh runs an experiment with the ESTKF filter method for ensmeble size 20 in which only the variable X os the model is observed at each 10th time step. The execution should only take seconds. 
+Further you can run
+{{{
+  cd ../bin
+  ../src/lorenz63/tools/run_PF.sh 
+}}}
+to run a similar experiment with the particle filter.
+
+
+
+=== 4. Plot output from the assimilation experiments ===
+
+To display the output of the assimilation experiments we provide plotting scripts for Matlab, Octave and Python. To use them do
+{{{
+cd ../src/lorenz63/plotting/
+}}}
+
+The following plot functions are available:
+||= Script =||= Description =||
+|| plot_eofs || Plot the mean state and the eigenvectors (EOFs) from the file holdig the covariance matrix that is generated with tools/generate_covar.F90. ||
+|| plot_obs || Plot the observations for a specified time step. The observation file is generated using tools/generate_obs.F90. ||
+|| plot_obs_series || Plot a time series of the observations for a specified variable of the model. The observation file is generated using tools/generate_obs.F90. ||
+|| plot_rms || Plot the true and estimated RMS errors over time from an assimilation experiment. ||
+|| plot_sigma || Plot the singular values from the file holdig the covariance matrix that is generated with tools/generate_covar.F90. ||
+|| plot_state || Plot the state from the file holding a true state trajectory (state.nc) or the estimated forecast or analysis states from an assimilation experiment. ||
+|| plot_state_series || Plot a time series of the state from the file holding a true state trajectory (state._l63nc) or the estimated forecast or analysis states from an assimilation experiment for a selected state variable. ||
+
+The scripts need the specification of the input file and the time step or variable to plot. The syntax of the functions is identical for Matlab and Python (except for `plot_rms`):
+
+||= Matlab =||= Python =||= Comment =||
+|| plot_eofs(filename, index) || plot_eofs(filename, index) || index=0 for mean state, index>1 for EOFs ||
+|| plot_obs(filename, timestep) || plot_obs(filename, timestep) || ||
+|| plot_obs_series(filename, variable) || plot_obs_series(filename, variable) || ||
+|| plot_rms(filename [, plot_forecast=1, plot_analysis=1]) || plot_rms(filename [, plot_forecast=True, plot_analysis=True]) || plot_forecast/plot_analysis are optional and switch on/off a plot ||
+|| plot_sigma(filename) || plot_sigma(filename) || ||
+|| plot_state(filename, timestep [, choice='t']) || plot_state(filename, timestep [, choice='t']) || `choice` is optional; options for `choice`: 't' true (use with truth file `state.nc`); for use with assimilation experiment output files: 'f' forecast, 'a' analysis, 'i' initial ||
+|| plot_state_series(filename, variable [, choice='t']) || plot_state(filename, variable [, choice='t']) || `choice` is optional; options for `choice`: 't' true (use with truth file `state_l63.nc`); for use with assimilation experiment output files: 'f' forecast, 'a' analysis, 'i' initial ||
+
+Here `filename` is the name of the file including its path. 
+In Matlab use 'help' to display the information about required input. 
+
+'''!Matlab/Octave plotting examples:'''[[BR]]
+`plot_obs('../../../bin/obs.nc',100)` plots the observation at time step 100[[BR]]
+`plot_obs_series('../../../bin/obs.nc',1)` plots the time series of observations for variable X (likely set 2=Y, 3=Z) [[BR]]
+`plot_state('../../../bin/ESTKF_N20.nc',100,'f')` plots the forecast state estimate from the ESTKF run example at the 100th analysis step[[BR]]
+`plot_rms('../../../bin/t1_N30_f0.97.nc')` plots the true and estimated RMS errors over time for the chosen experiment[[BR]]
+`plot_state('../../../bin/state_l63.nc',1101)` plots the true state at model time step 1101 (= analysis step 100)[[BR]] 
+`plot_state_series('../../../bin/ESTKF_N20.nc',1,'f')` plots the time series for the ESTKF run example of the forecast state for variable X[[BR]] 
+
+'''Python plotting'''[[BR]]
+For Python the scripts are provided by `plot_l63.py`.
+
+The module can either be imported, e.g. to use its functions interactively:[[BR]]
+`>>> import plot_l63`[[BR]]
+`>>> plot_l63.plot_obs('../../../bin/obs_l63.nc', 4)`[[BR]]
+
+Alternatively the script can be run from the command line, providing the function
+name and its argument as command line parameters:[[BR]]
+`./plot_l63.py plot_obs ../../../bin/obs_l63.nc 4`[[BR]]
+(If this fails you can also try to run the script as
+`python plot_l63.py plot_obs ../../../bin/obs_l63.nc 4`). To plot only the Analysis RMS error one can use[[BR]]
+`plot_rms ../../../bin/ESTKF_N20.nc False True`
+
+
+== Run options ==
+
+The implementation for the Lorenz-63 model has a wide range of options. For the full set, we recommend to check the list in the file init_pdaf.F90 in the directory lorenz63. Essentially all of the options can be specified on the command line in the format '-VARIABLE VALUE', where VARIABLE is the name of the variable in the program and VALUE is its value. The experiment in `run_ESTKF.sh` is executed like
+{{{
+./pdaf_lorenz_63 -total_steps 5000 -step_null 1000 -dim_ens 20 -forget 0.8 -filtertype 6 -file_asml ESTKF_N20.nc
+}}}
+and the experiment in `run_PF.sh` is executed like
+{{{
+./pdaf_lorenz_63 -total_steps 5000 -step_null 1000 -dim_ens 20 -pf_res_type 2 -pf_noise_type 2 -pf_noise_amp 0.2 -filtertype 12 -file_asml PF_N20.nc
+}}}
+The meaning of the different options is the following:
+
+||= Variable =||= Description =||= Default value =||
+|| dim_ens || ensemble size || 20
+|| file_asml || the name of the output file || ESTKF_N20.nc ||
+|| filtertype || Specifies the filter to use; 6 sets the ESTKF filter, 12 the particle filter || 6 || 
+|| forget  || value of the forgetting factor for multiplicative variance inflation || 1.0 ||
+|| step_null || initial time step of experiment. This specifies the offset of the observations which are generated from a model forward run. A value 1000 avoids the spin-up phase of the model || 0 ||
+|| total_steps || number of time step in the experiment || 5000 ||
+|| pf_res_type || type of resmpling scheme in the particle filter (2 set stochastic unversal resampling) || 1 ||
+|| pf_noise_type || type of noise to perturb resampled particles (2 sets noise with amplitude relative to ensemble standard deviation) || 0 ||
+|| pf_noise_amp || amplitude of noise to perturb resampled particles || 0.0 ||
+
+Some further options:
+
+||= Variable =||= Description =||= Default value =||
+|| delt_obs  || time interval between observations, i.e. the forecast length || 10 ||
+|| dim_lag   || set the time lag for the smoother (the smoother is active for dim_lag>0) || 0 ||
+|| dx_obs    || distance between observations on the grid; let's you specify an incompletely observated state (only used when `use_obs_mask=.true.`) || 3 ||
+|| numobs || number of observed grid points; the points 1 to numobs are observed (can be combined with dx_obs; only used when `use_obs_mask=.true.`) || 1 ||
+|| use_obs_mask || whether to use incomplete observations (need to be .true. for any other settings on the observation density to be used) || .true. ||
+
+For the full set of options, please see init_pdaf.F90 in the lorenz63 directory. There, also the possible settings e.g. for `filtertype` are described.
+
+