Merge pull request NCAR#322 from NCAR/tracer_advection_siparcs

Tracer advection siparcs
Benjamin-Gunn · Mar 10, 2022 · 90b6dcc · 90b6dcc
2 parents 3076d43 + d4d0192
commit 90b6dcc
Show file tree

Hide file tree

Showing 62 changed files with 4,878 additions and 19 deletions.
diff --git a/CHANGELOG.rst b/CHANGELOG.rst
@@ -22,6 +22,14 @@ individual files.
 
 The changes are now listed with the most recent at the top.
 
+**March 10th 2022 :: Tracer advection model. Tag v9.14.0**
+
+- New model for tracer advection based on Lorenz_96 using a
+  Semi-Lagrangian scheme. 
+- Forward operator documentation updated to Manhattan
+
+*lorenz_96_tracer_advection contributed by Fairuz Ishraque (SIParCS)*
+
 **February 11 2022 :: Bug-fix for nag compiler. Tag v9.13.2**
 
 - Bug-fix for fixsystem preprocessing for NAG compiler

diff --git a/README.rst b/README.rst
@@ -390,6 +390,7 @@ References
    models/lorenz_84/readme
    models/lorenz_96/readme
    models/lorenz_96_2scale/readme
+   models/lorenz_96_tracer_advection/readme
    models/forced_lorenz_96/readme
    models/MITgcm_ocean/readme
    models/mpas_atm/readme

diff --git a/assimilation_code/modules/observations/oned_quantities_mod.f90 b/assimilation_code/modules/observations/oned_quantities_mod.f90
@@ -27,5 +27,7 @@
 !   QTY_LARGE_SCALE_STATE            desc="state varies with large time/space scale"
 !   QTY_SMALL_SCALE_STATE            desc="state varies with small time/space scale"
 !   QTY_1D_PARAMETER
+!   QTY_TRACER_CONCENTRATION
+!   QTY_TRACER_SOURCE
 ! 
 ! END DART PREPROCESS QUANTITY DEFINITIONS
diff --git a/conf.py b/conf.py
@@ -21,7 +21,7 @@
 author = 'Data Assimilation Research Section'
 
 # The full version, including alpha/beta/rc tags
-release = '9.13.2'
+release = '9.14.0'
 master_doc = 'README'
 
 # -- General configuration ---------------------------------------------------

diff --git a/diagnostics/matlab/plot_bins.m b/diagnostics/matlab/plot_bins.m
@@ -70,7 +70,8 @@ function plot_bins(input_file)
 switch lower(pinfo.model)
 
    case {'9var','lorenz_63','lorenz_84','lorenz_96','lorenz_96_2scale', ...
-	 'lorenz_04', 'forced_lorenz_96','ikeda','simple_advection', 'null'}
+	 'lorenz_04', 'forced_lorenz_96','ikeda','simple_advection', ...
+         'lorenz_96_tracer_advection', 'null'}
 
       varid = SetVariableID(pinfo);      % queries for variable IDs
       pinfo.var        = varid.var;

diff --git a/diagnostics/matlab/plot_correl.m b/diagnostics/matlab/plot_correl.m
@@ -91,7 +91,7 @@ function plot_correl(diagn_file)
       fprintf('Correlating variable %s index %d at time %d.\n', ...
            pinfo.base_var,pinfo.base_var_index, pinfo.base_time)
 
-   case {'simple_advection'}
+   case {'simple_advection', 'lorenz_96_tracer_advection'}
 
       disp('Your choice of variables are:')
       disp(pinfo.vars)

diff --git a/diagnostics/matlab/plot_ens_err_spread.m b/diagnostics/matlab/plot_ens_err_spread.m
@@ -71,7 +71,7 @@ function plot_ens_err_spread(input_file)
 switch lower(pinfo.model)
 
    case {'9var','lorenz_63','lorenz_84','lorenz_96','lorenz_96_2scale', ...
-	 'forced_lorenz_96','lorenz_04','ikeda','simple_advection', 'null'}
+	 'forced_lorenz_96','lorenz_04','ikeda','simple_advection','lorenz_96_tracer_advection', 'null'}
 
       varid          = SetVariableID(pinfo);
       pinfo.var      = varid.var;

diff --git a/diagnostics/matlab/plot_ens_mean_time_series.m b/diagnostics/matlab/plot_ens_mean_time_series.m
@@ -80,7 +80,7 @@ function plot_ens_mean_time_series(input_file)
 switch lower(pinfo.model)
 
    case {'9var','lorenz_63','lorenz_84','lorenz_96','lorenz_96_2scale', ...
-	 'forced_lorenz_96','lorenz_04','ikeda','simple_advection', 'null'}
+	 'forced_lorenz_96','lorenz_04','ikeda','simple_advection','lorenz_96_tracer_advection', 'null'}
 
       varid          = SetVariableID(pinfo);
       pinfo.var      = varid.var;

diff --git a/diagnostics/matlab/plot_ens_time_series.m b/diagnostics/matlab/plot_ens_time_series.m
@@ -81,7 +81,7 @@ function plot_ens_time_series(input_file)
 
 switch lower(pinfo.model)
 
-   case {'9var','lorenz_63','lorenz_84','lorenz_96','lorenz_96_2scale', ...
+   case {'9var','lorenz_63','lorenz_84','lorenz_96','lorenz_96_2scale','lorenz_96_tracer_advection' ...
 	 'forced_lorenz_96','lorenz_04','ikeda','simple_advection', 'null'}
 
       varid          = SetVariableID(pinfo);

diff --git a/diagnostics/matlab/plot_jeff_correl.m b/diagnostics/matlab/plot_jeff_correl.m
@@ -75,7 +75,7 @@ function plot_jeff_correl(fname)
       inputstring = input('Input variable and index for correlation \n','s');
       [pinfo.state_var, pinfo.state_var_index] = ParseAlphaNumerics(inputstring);
 
-   case {'simple_advection'}
+   case {'simple_advection', 'lorenz_96_tracer_advection'}
 
       disp('Your choice of variables are:')
       disp(pinfo.vars)

diff --git a/diagnostics/matlab/plot_phase_space.m b/diagnostics/matlab/plot_phase_space.m
@@ -159,7 +159,7 @@ function plot_phase_space(fname)
                      'ens_mem' , ens_mem  , ...
                      'ltype'   , ltype   );
 
-   case {'simple_advection'}
+   case {'simple_advection','lorenz_96_tracer_advection'}
 
       if (ishold), clear var1 var2 var3 ens_mem ltype; end
 

diff --git a/diagnostics/matlab/plot_sawtooth.m b/diagnostics/matlab/plot_sawtooth.m
@@ -60,7 +60,7 @@ function plot_sawtooth(posterior_file,prior_file)
 switch lower(pinfo.model)
 
    case {'9var','lorenz_63','lorenz_84','lorenz_96','lorenz_96_2scale', ...
-	 'lorenz_04','forced_lorenz_96','ikeda','simple_advection', 'null'}
+	 'lorenz_04','forced_lorenz_96','ikeda','simple_advection','lorenz_96_tracer_advection','null'}
 
       % determine which variable ['state','X','Y' ...] and
       % determine which variable IDs (locations), as well as

diff --git a/diagnostics/matlab/plot_total_err.m b/diagnostics/matlab/plot_total_err.m
@@ -62,7 +62,7 @@ function plot_total_err(input_file)
 pinfo   = CheckModelCompatibility(truth_file,input_file);
 
 switch lower(pinfo.model)
-    case{'9var','lorenz_63','lorenz_84','lorenz_96','lorenz_96_2scale', ...
+    case{'9var','lorenz_63','lorenz_84','lorenz_96','lorenz_96_2scale','lorenz_96_tracer_advection', ...
          'lorenz_04', 'forced_lorenz_96','ikeda','simple_advection', 'null'}
 
     case{'fms_bgrid'}

diff --git a/diagnostics/matlab/plot_var_var_correl.m b/diagnostics/matlab/plot_var_var_correl.m
@@ -85,7 +85,7 @@ function plot_var_var_correl(fname)
       inputstring = input('Input variable and index for correlation \n','s');
       [pinfo.state_var, pinfo.state_var_index] = ParseAlphaNumerics(inputstring);
 
-   case {'simple_advection'}
+   case {'simple_advection', 'lorenz_96_tracer_advection'}
 
       disp('Your choice of variables are:')
       disp(pinfo.vars)

diff --git a/diagnostics/matlab/private/CheckModel.m b/diagnostics/matlab/private/CheckModel.m
@@ -139,7 +139,7 @@
 
       vars.fname = fname;
 
-   case {'simple_advection'}
+   case {'simple_advection', 'lorenz_96_tracer_advection'}
 
       num_locs = dim_length(fname,'location'); % # of X variables
       loc1d    =     ncread(fname,'location');

diff --git a/diagnostics/matlab/private/PlotBins.m b/diagnostics/matlab/private/PlotBins.m
@@ -102,7 +102,7 @@ function PlotBins(pinfo)
          axis tight
       end
 
-   case {'lorenz_96_2scale','simple_advection'}
+   case {'lorenz_96_2scale','simple_advection', 'lorenz_96_tracer_advection'}
 
       clf; iplot = 0;
       for ivar = pinfo.var_inds,

diff --git a/diagnostics/matlab/private/PlotCorrel.m b/diagnostics/matlab/private/PlotCorrel.m
@@ -38,7 +38,7 @@ function PlotCorrel( pinfo )
 switch(lower(pinfo.model))
 
     case {'9var','lorenz_63','lorenz_84','lorenz_96', 'lorenz_96_2scale', ...
-            'lorenz_04','forced_lorenz_96','ikeda','simple_advection', 'null'}
+            'lorenz_04','forced_lorenz_96','ikeda','simple_advection', 'lorenz_96_tracer_advection', 'null'}
 
         % The Base Variable Index must be a valid state variable
         if ( pinfo.base_var_index > pinfo.num_state_vars )

diff --git a/diagnostics/matlab/private/PlotEnsErrSpread.m b/diagnostics/matlab/private/PlotEnsErrSpread.m
@@ -104,7 +104,7 @@ function PlotEnsErrSpread( pinfo )
         end
 
     case {'lorenz_63','lorenz_84','lorenz_96', 'lorenz_96_2scale', ...
-          'lorenz_04','forced_lorenz_96','ikeda','simple_advection', 'null'}
+          'lorenz_04','forced_lorenz_96','ikeda','simple_advection','lorenz_96_tracer_advection','null'}
 
         truth      = get_hyperslab('fname',pinfo.truth_file, ...
             'varname',pinfo.var, ...

diff --git a/diagnostics/matlab/private/PlotEnsMeanTimeSeries.m b/diagnostics/matlab/private/PlotEnsMeanTimeSeries.m
@@ -167,7 +167,7 @@ function PlotEnsMeanTimeSeries( pinfo )
       legend boxoff
 
    case {'lorenz_96', 'lorenz_96_2scale', 'forced_lorenz_96', 'lorenz_04', ...
-         'ikeda', 'simple_advection', 'null'}
+         'ikeda', 'simple_advection', 'lorenz_96_tracer_advection','null'}
 
       %% Plot all variables in own subplot ... might get cluttered.
       figure(1); clf; iplot = 0;

diff --git a/diagnostics/matlab/private/PlotEnsTimeSeries.m b/diagnostics/matlab/private/PlotEnsTimeSeries.m
@@ -181,7 +181,7 @@ function PlotEnsTimeSeries( pinfo )
       legend boxoff
 
    case {'lorenz_96', 'lorenz_96_2scale', 'forced_lorenz_96', 'lorenz_04', ...
-         'ikeda', 'simple_advection', 'null'}
+         'ikeda', 'simple_advection', 'lorenz_96_tracer_advection', 'null'}
 
       % Use one figure with subplots
       figure(1); clf; iplot = 0;

diff --git a/diagnostics/matlab/private/PlotPhaseSpace.m b/diagnostics/matlab/private/PlotPhaseSpace.m
@@ -113,7 +113,7 @@ function PlotPhaseSpace( pinfo )
       legend boxoff
 
    case {'9var', 'lorenz_63', 'lorenz_84', 'lorenz_96', 'lorenz_96_2scale', ...
-	 'lorenz_04', 'forced_lorenz_96', 'simple_advection', 'null'}
+	 'lorenz_04', 'forced_lorenz_96', 'simple_advection', 'lorenz_96_tracer_advection', 'null'}
 
       BulletProof(pinfo, X, Y, Z)          % rudimentary bulletproofing
 

diff --git a/diagnostics/matlab/private/PlotTotalErr.m b/diagnostics/matlab/private/PlotTotalErr.m
@@ -282,7 +282,59 @@ function PlotTotalErr( pinfo )
          title({string1,string2},'interpreter','none','fontweight','bold')
          xlabel(sprintf('model "days" (%d timesteps)',pinfo.time_series_length))
          ylabel('Total Error')
+      end
+
+
+   case {'lorenz_96_tracer_advection'}
+
+      varlist = {'state_variable','tracer_concentration','source'};
+
+      for ivar = 1:length(varlist)
+
+         % Get the appropriate netcdf variables
+
+         truth  = get_hyperslab('fname', pinfo.truth_file, ...
+                   'varname', varlist{ivar}, ...
+                   'permute', 'T', ...
+                   'squeeze', 'T', ...
+                   'tindex1', pinfo.truth_time(1), ...
+                   'tcount', pinfo.truth_time(2));
+
+         varname = sprintf('%s_mean',varlist{ivar});
+         ens    = get_hyperslab('fname', pinfo.diagn_file, ...
+                   'varname', varname, ...
+                   'permute', 'T', ...
+                   'tindex1', pinfo.diagn_time(1), ...
+                   'tcount', pinfo.diagn_time(2));
+
+         varname = sprintf('%s_sd',varlist{ivar});
+         spread = get_hyperslab('fname', pinfo.diagn_file, ...
+                   'varname', varname, ...
+                   'permute', 'T', ...
+                   'tindex1', pinfo.diagn_time(1), ...
+                   'tcount', pinfo.diagn_time(2));
+
+         [~,num_vars] = size(spread);
+
+         % Also need to compute the spread; zero truth for this and
+         % compute distance from 0
+         err        = total_err(truth, ens);
+         err_spread = total_err(zeros(size(spread)), spread);
+         errTotal   = sum(err)/pinfo.time_series_length;
+         spreadTotal= sum(err_spread)/pinfo.time_series_length;
+         string1 = ['time-mean Ensemble Mean Total Error = ' num2str(errTotal)];
+         string2 = ['time-mean Ensemble Spread = ' num2str(spreadTotal)];
 
+         figure(ivar); clf(ivar);
+         plot(pinfo.time,err, 'b', pinfo.time,err_spread, 'r','LineWidth',figdata.linewidth);
+         set(gca,'FontSize',figdata.fontsize)
+         legend(string1,string2,'Location','NorthEast','FontSize',figdata.fontsize)
+         legend boxoff
+         string1 = sprintf('%s Total Error over all %d variables', pinfo.model, num_vars);
+         string2 = sprintf('''%s'' %s', varlist{ivar}, pinfo.diagn_file);
+         title({string1,string2},'interpreter','none','fontweight','bold')
+         xlabel(sprintf('model "days" (%d timesteps)',pinfo.time_series_length))
+         ylabel('Total Error')
       end
 
    case 'fms_bgrid'

diff --git a/diagnostics/matlab/private/SetVariableID.m b/diagnostics/matlab/private/SetVariableID.m
@@ -73,7 +73,7 @@
          varid.var_inds = vrbl_inds;
       end
 
-   case 'simple_advection'
+   case {'simple_advection','lorenz_96_tracer_advection'}
 
       % query to see if the defaults are OK ...
 

diff --git a/models/README.rst b/models/README.rst
@@ -23,6 +23,7 @@ DART supported models:
 - :doc:`lorenz_84/readme`
 - :doc:`lorenz_96/readme`
 - :doc:`lorenz_96_2scale/readme`
+- :doc:`lorenz_96_tracer_advection/readme`
 - :doc:`forced_lorenz_96/readme`
 - :doc:`MITgcm_ocean/readme`
 - :doc:`mpas_atm/readme`

diff --git a/models/lorenz_96_tracer_advection/images/Semi_lag.png b/models/lorenz_96_tracer_advection/images/Semi_lag.png