Back to home page

darwin3

 
 

    


Warning, /doc/examples/global_oce_optim/global_oce_optim.rst is written in an unsupported language. File is not indexed.

view on githubraw file Latest commit 4d722833 on 2023-10-06 20:00:24 UTC
1c8cebb321 Jeff*0001 .. _sec_global_oce_optim:
                0002 
d67096e55c Jeff*0003 Global Ocean State Estimation
                0004 =============================
                0005 
                0006 (in directory: :filelink:`verification/tutorial_global_oce_optim/`)
                0007 
                0008 Overview
                0009 --------
                0010 
                0011 This experiment illustrates the optimization capacity of the MITgcm:
                0012 here, a high level description.
                0013 
                0014 In this tutorial, a very simple case is used to illustrate the
                0015 optimization capacity of the MITgcm. Using an ocean configuration with
                0016 realistic geography and bathymetry on a :math:`4\times4^\circ` spherical
                0017 polar grid, we estimate a time-independent surface heat flux adjustment
                0018 :math:`Q_\mathrm{netm}` that attempts to bring the model climatology
                0019 into consistency with observations (Levitus and Boyer (1994a,b)
                0020 :cite:`levitus:94a,levitus:94b`).
                0021 
                0022 This adjustment :math:`Q_\mathrm{netm}` (a 2-D field only function of
                0023 longitude and latitude) is the control variable of an optimization
                

** Warning **

Wide character in print at /usr/local/share/lxr/source line 1030, <$git> line 25.

0024 problem. It is inferred by an iterative procedure using an ‘adjoint

** Warning **

Wide character in print at /usr/local/share/lxr/source line 1030, <$git> line 26.

0025 technique’ and a least-squares method (see, for example, 0026 Stammer et al. (2002) :cite:`stammer:02` and Ferriera et a. (2005) :cite:`ferriera:05`. 0027 0028 The ocean model is run forward in time and the quality of the solution 0029 is determined by a cost function, :math:`J_1`, a measure of the 0030 departure of the model climatology from observations: 0031 0032 .. math:: 0033 J_1=\frac{1}{N}\sum_{i=1}^N \left[ \frac{\overline{T}_i-\overline{T}_i^{lev}}{\sigma_i^T}\right]^2 0034 :label: cost-temp 0035 0036 where :math:`\overline{T}_i` and :math:`\overline{T}_i^{lev}` are, 0037 respectively, the model and observed potential temperature at each grid 0038 point :math:`i`. The differences are weighted by an *a priori* 0039 uncertainty :math:`\sigma_i^T` on observations (as provided by 0040 Levitus and Boyer (1994a) 0041 :cite:`levitus:94a`). The error :math:`\sigma_i^T` is only a 0042 function of depth and varies from 0.5 K at the surface to 0.05K at the 0043 bottom of the ocean, mainly reflecting the decreasing temperature 0044 variance with depth (see :numref:`tut_global_optim_errors`\ a). A value of :math:`J_1` of order 1 0045 means that the model is, on average, within observational uncertainties. 0046 0047 .. figure:: figs/Error.png 0048 :width: 80% 0049 :align: center 0050 :alt: A priori errors on potential temperature and surface hf 0051 :name: tut_global_optim_errors 0052 0053 *A priori* errors on potential temperature (left, in :sup:`o`\ C) and surface heat flux 0054 (right, in Wm\ :sup:`-2`) used to compute the cost terms :math:`J_1` and :math:`J_2`, respectively. 0055 0056 The cost function also places constraints on the adjustment to insure it

** Warning **

Wide character in print at /usr/local/share/lxr/source line 1030, <$git> line 58.

0057 is “reasonable”, i.e., of order of the uncertainties on the observed 0058 surface heat flux: 0059 0060 .. math:: J_2 = \frac{1}{N} \sum_{i=1}^N \left[\frac{Q_\mathrm{netm}}{\sigma^Q_i} \right]^2 0061 0062 where :math:`\sigma^Q_i` are the *a priori* errors on the observed heat 0063 flux as estimated by Stammer et al. (2002) :cite:`stammer:02` from 30% of local 0064 root-mean-square variability of the NCEP forcing field (see :numref:`tut_global_optim_errors`\ b). 0065 0066 The total cost function is defined as 0067 :math:`J=\lambda_1 J_1+ \lambda_2 J_2` where :math:`\lambda_1` and 0068 :math:`\lambda_2` are weights controlling the relative contribution of 0069 the two components. The adjoint model then yields the sensitivities 0070 :math:`\partial J/\partial Q_\mathrm{netm}` of :math:`J` relative to the 0071 2-D fields :math:`Q_\mathrm{netm}`. Using a line-searching algorithm 0072 (Gilbert and Lemarchal 1989 :cite:`gil-lem:89`), :math:`Q_\mathrm{netm}` is adjusted

** Warning **

Wide character in print at /usr/local/share/lxr/source line 1030, <$git> line 74.

0073 then in the sense to reduce :math:`J` — the procedure is repeated until 0074 convergence. 0075 0076 :numref:`tut_global_optim_tutfig` shows the results of such an optimization. The model is 0077 started from rest and from January-mean temperature and salinity initial 0078 conditions taken from the Levitus dataset. The experiment is run a year 0079 and the averaged temperature over the whole run (i.e., annual mean) is 0080 used in the cost function :eq:`cost-temp` to evaluate the model [1]_. 0081 Only the top 2 levels are used. The first guess :math:`Q_\mathrm{netm}` 0082 is chosen to be zero. The weights :math:`\lambda_1` and 0083 :math:`\lambda_2` are set to 1 and 2, respectively. The total cost 0084 function converges after 15 iterations, decreasing from 6.1 to 2.7 (the 0085 temperature contribution decreases from 6.1 to 1.8 while the heat flux 0086 one increases from 0 to 0.42). The right panels of :numref:`tut_global_optim_tutfig` 0087 illustrate the evolution of the temperature error at the surface from 0088 iteration 0 to iteration 15. Unsurprisingly, the largest errors at 0089 iteration 0 (up to 6 :sup:`o`\ C, top left panels) are found in 0090 the Western boundary currents. After optimization, the departure of the 0091 model temperature from observations is reduced to 1 :sup:`o`\ C 0092 or less almost everywhere except in the Pacific equatorial cold tongue. 0093 Comparison of the initial temperature error (top, right) and heat flux 0094 adjustment (bottom, left) shows that the system basically increased the 0095 heat flux out of the ocean where temperatures were too warm and 0096 vice-versa. Obviously, heat flux uncertainties are not solely 0097 responsible for temperature errors, and the heat flux adjustment partly 0098 compensates the poor representation of narrow currents (Western boundary 0099 currents, equatorial currents) at :math:`4\times4^\circ` resolution. 0100 This is allowed by the large *a priori* error on the heat flux :numref:`tut_global_optim_errors`. 0101 The Pacific cold tongue is a counter example: there, heat 0102 fluxes uncertainties are fairly small (about 20W m\ :sup:`-2`), and a 0103 large temperature errors remains after optimization. 0104 0105 .. figure:: figs/Tutorial_fig.png 0106 :width: 100% 0107 :align: center 0108 :alt: Surface HF and Temp Iter 0 v 15 0109 :name: tut_global_optim_tutfig 0110 0111 Initial annual mean surface heat flux (top right in W m\ :sup:`-2`) and adjustment obtained at iteration 15 (bottom right). 0112 Averaged difference between model and observed potential temperatures at the surface (in :math:`^\circ`\ C) 0113 before optimization (iteration 0, top right) and after optimization (iteration 15, bottom right). 0114 Contour intervals for heat flux and temperature are 25W m\ :sup:`-2` and 1 :sup:`o`\ C, respectively. A positive flux is out of the ocean. 0115 0116 Implementation of the control variable and the cost function 0117 ------------------------------------------------------------ 0118 0119 One of the goals of this tutorial is to illustrate how to implement a new 0120 control variable. Most of this is fairly generic and is done in :filelink:`pkg/ctrl` 0121 and :filelink:`pkg/cost`. The modifications can be 0122 tracked by the CPP option :varlink:`ALLOW_HFLUXM_CONTROL` or the comment 0123 ``cHFLUXM_CONTROL``. The more specific modifications required for the 0124 experiment are found in 0125 :filelink:`verification/tutorial_global_oce_optim/code_ad`. Here follows a brief 0126 description of the implementation. 0127 0128 The control variable 0129 ~~~~~~~~~~~~~~~~~~~~ 0130 0131 The adjustment :math:`Q_\mathrm{netm}` is activated by setting ``#define`` 0132 :varlink:`ALLOW_HFLUXM_CONTROL` in :filelink:`code_ad/CTRL_OPTIONS.h <verification/tutorial_global_oce_optim/code_ad//CTRL_OPTIONS.h>`. 0133

** Warning **

Wide character in print at /usr/local/share/lxr/source line 1030, <$git> line 135.

0134 It is first implemented as a “normal” forcing variable. It is defined in 0135 :filelink:`model/inc/FFIELDS.h`, initialized to zero in :filelink:`model/src/ini_forcing.F`, and then used in 0136 :filelink:`model/src/external_forcing_surf.F`. :math:`Q_\mathrm{netm}` is made a control 0137 variable in :filelink:`pkg/ctrl` by modifying the following subroutines: 0138 0139 - :filelink:`pkg/ctrl/ctrl_init.F` where :math:`Q_\mathrm{netm}` is defined as the control 0140 variable number 24, 0141 0142 - :filelink:`pkg/ctrl/ctrl_pack.F` which writes, at the end of each iteration, the 0143 sensitivity of the cost function 0144 :math:`\partial J/\partial Q_\mathrm{netm}` in to a file to be used 0145 by the line-search algorithm, 0146 0147 - :filelink:`pkg/ctrl/ctrl_unpack.F` which reads, at the start of each iteration, the 0148 updated adjustment as provided by the line-search algorithm, 0149 0150 - :filelink:`pkg/ctrl/ctrl_map_forcing.F` in which the updated adjustment is added to the 0151 first guess :math:`Q_\mathrm{netm}`. 0152 0153 Cost functions 0154 ~~~~~~~~~~~~~~ 0155 0156 The cost functions are implemented using :filelink:`pkg/cost`. 0157 0158 - The temperature cost function :math:`J_1` which measures the drift of 0159 the mean model temperature from the Levitus climatology is 0160 implemented in :filelink:`/verification/tutorial_global_oce_optim/code_ad/cost_temp.F`. 0161 It is activated by ``#define`` :varlink:`ALLOW_COST_TEMP` in 0162 :filelink:`code_ad/COST_OPTIONS.h <verification/tutorial_global_oce_optim/code_ad//COST_OPTIONS.h>`. 0163 It requires the mean temperature of the model which 0164 is obtained by accumulating the temperature in :filelink:`pkg/cost/cost_tile.F` (called 0165 at each time step). The value of the cost function is stored in 0166 :varlink:`objf_temp` and its weight :math:`\lambda_1` in :varlink:`mult_temp`. 0167 0168 - The heat flux cost function, penalizing the departure of the surface 0169 heat flux from observations is implemented in :filelink:`/verification/tutorial_global_oce_optim/code_ad/cost_hflux.F`, and 0170 activated by ``#define`` :varlink:`ALLOW_COST_HFLUXM` in 0171 :filelink:`code_ad/COST_OPTIONS.h <verification/tutorial_global_oce_optim/code_ad//COST_OPTIONS.h>`. The 0172 value of the cost function is stored in :varlink:`objf_hfluxm` and its weight 0173 :math:`\lambda_2` in :varlink:`mult_hflux`. 0174 0175 - The subroutine :filelink:`pkg/cost/cost_final.F` calls the cost function subroutines 0176 and makes the (weighted) sum of the various contributions. 0177 0178 - The various weights used in the cost functions are read in 0179 :filelink:`/verification/tutorial_global_oce_optim/code_ad/cost_weights.F`. The weight of the cost functions are read in 0180 :filelink:`pkg/cost/cost_readparms.F` from the input file :filelink:`verification/tutorial_global_oce_optim/input_ad/data.cost`. 0181 0182 Code Configuration 0183 ------------------ 0184 0185 The experiment files in :filelink:`verification/tutorial_global_oce_optim/code_ad/` 0186 and :filelink:`verification/tutorial_global_oce_optim/input_ad/` contain the code customizations and parameter 0187 settings. Most of them are identical to those used in the Global Ocean ( 0188 experiment :filelink:`verification/tutorial_global_oce_latlon/`). Below, we describe some of 0189 the customizations required for this experiment. 0190 0191 Compilation-time customizations in :filelink:`code_ad <verification/tutorial_global_oce_optim/code_ad/>` 0192 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 0193 0194 In :filelink:`code_ad/CTRL_OPTIONS.h <verification/tutorial_global_oce_optim/code_ad//CTRL_OPTIONS.h>`: 0195 0196 - ``#define`` :varlink:`ALLOW_ECCO_OPTIMIZATION` 0197 0198 .. _tut_global_oce_runsect: 0199 0200 Running-time customizations in :filelink:`input_ad <verification/tutorial_global_oce_optim/input_ad/>` 0201 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 0202 0203 - :filelink:`input_ad/data <verification/tutorial_global_oce_optim/input_ad/data>`: note the smaller :varlink:`cg2dTargetResidual` than in the 0204 forward-only experiment, 0205 0206 - :filelink:`input_ad/data.optim <verification/tutorial_global_oce_optim/input_ad/data.optim>` specifies the iteration number, 0207 0208 - :filelink:`input_ad/data.ctrl <verification/tutorial_global_oce_optim/input_ad/data.ctrl>` is used, in particular, to specify the name of the 0209 sensitivity and adjustment files associated to a control variable, 0210 0211 - :filelink:`input_ad/data.cost <verification/tutorial_global_oce_optim/input_ad/data.cost>`: parameters of the cost functions, in particular 0212 :varlink:`lastinterval` specifies the length of time-averaging for the model 0213 temperature to be used in the cost function :eq:`cost-temp`, 0214 0215 - :filelink:`input_ad/data.pkg <verification/tutorial_global_oce_optim/input_ad/data.pkg>`: note that the Gradient Check package is turned on by 0216 default (:varlink:`useGrdchk` ``=.TRUE.``), 0217 0218 - ``Err_hflux.bin`` and ``Err_levitus_15layer.bin`` are the files 0219 containing the heat flux and potential temperature uncertainties, 0220 respectively. 0221 0222 Compiling 0223 --------- 0224 0225 The optimization experiment requires two executables: 1) the MITgcm and 0226 its adjoint (``mitgcmuv_ad``) and 2) the line-search algorithm 0227 (``optim.x``). 0228 0229 Compilation of MITgcm and its adjoint: ``mitcgmuv_ad`` 0230 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 0231 0232 Before compiling, first note that in the directory :filelink:`code_ad <verification/tutorial_global_oce_optim/code_ad/>`, two files 0233 must be updated: 0234 0235 - :filelink:`code_ad/code_ad_diff.list <verification/tutorial_global_oce_optim/code_ad/code_ad_diff.list>` which lists new subroutines to be compiled by the 0236 TAF software (:filelink:`code_ad/cost_temp.F <verification/tutorial_global_oce_optim/code_ad/cost_temp.F>` 0237 and :filelink:`code_ad/cost_hflux.F <verification/tutorial_global_oce_optim/code_ad/cost_hflux.F>`), 0238 0239 - the file :filelink:`code_ad/ad_optfile.local <verification/tutorial_global_oce_optim/code_ad/ad_optfile.local>` provides a list of the control 0240 variables and the name of cost function to the TAF software. 0241 0242 Then, in the directory :filelink:`build <verification/tutorial_global_oce_optim/build/>`, type: 0243 0244 :: 0245 0246 % ../../../tools/genmake2 -mods=../code_ad -adof=../code_ad/ad_optfile.local 0247 % make depend 0248 % make adall 0249 0250 to generate the MITgcm executable ``mitgcmuv_ad``. 0251 0252 Compilation of the line-search algorithm: ``optim.x`` 0253 ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ 0254 0255 This is done from the directories :filelink:`lsopt/` and :filelink:`optim/` (found in the top MITgcm directory). In 0256 :filelink:`lsopt/`, unzip the ``blash1`` library adapted to your platform (see :filelink:`lsopt/README`), and change 0257 the ``Makefile`` accordingly. Compile with: 0258 0259 :: 0260 0261 % make all 0262 0263 (more details in :filelink:`lsopt/lsopt_doc.txt`) 0264 0265 In :filelink:`optim/`, the path of the directory where ``mitgcm_ad`` was compiled 0266 must be specified in the ``Makefile`` in the variable :varlink:`INCLUDEDIRS`. The file 0267 name of the control variable (here, :varlink:`xx_hfluxm_file`) must be added to 0268 the namelist read by :filelink:`optim/optim_numbmod.F`. Then use 0269 0270 :: 0271 0272 % make depend 0273 0274 and 0275 0276 :: 0277 0278 % make 0279 0280 to generate the line-search executable ``optim.x``. 0281 0282 Running the estimation 0283 ---------------------- 0284 0285 Make a new subdirectory ``input_ad/OPTIM``. 0286 Copy the ``mitgcmuv_ad`` executable to :filelink:`input_ad <verification/tutorial_global_oce_optim/input_ad/>` 0287 and ``optim.x`` to this subdirectory. 0288 ``cd`` into :filelink:`input_ad/<verification/tutorial_global_oce_optim/input_ad/>`. The first iteration

** Warning **

Wide character in print at /usr/local/share/lxr/source line 1030, <$git> line 290.

0289 is somewhat particular and is best done “by hand” while the following 0290 iterations can be run automatically (see below). Check that the 0291 iteration number is set to 0 in :filelink:`input_ad/data.optim <verification/tutorial_global_oce_optim/input_ad/data.optim>` and run MITgcm: 0292 0293 :: 0294 0295 % ./mitgcmuv_ad 0296 0297 The output files ``adxx_hfluxm.0000000000.*`` and ``xx_hfluxm.0000000000.*`` 0298 contain the sensitivity of the cost function to :math:`Q_\mathrm{netm}` 0299 and the adjustment to :math:`Q_\mathrm{netm}` (zero at the first 0300 iteration), respectively. Two other files called 0301 ``costhflux_tut_MITgcm.opt0000`` and ``ctrlhflux_tut_MITgcm.opt0000`` are 0302 also generated. They essentially contain the same information as the 0303 ``adxx_.hfluxm*`` and ``xx_hfluxm*`` files, but in a compressed format. 0304 These two files are the only ones involved in the communication between 0305 the adjoint model ``mitgcmuv_ad`` and the line-search algorithm 0306 ``optim.x``. Only at the first iteration, are they both generated by 0307 ``mitgcmuv_ad``. Subsequently, ``costhflux_tut_MITgcm.opt`` :math:`n` is 0308 an output of the adjoint model at iteration :math:`n` and an input of 0309 the line-search. The latter returns an updated adjustment in 0310 ``ctrlhflux_tut_MITgcm.opt`` :math:`n+1` to be used as an input of the 0311 adjoint model at iteration :math:`n+1`. 0312 0313 At the first iteration, move ``costhflux_tut_MITgcm.opt0000`` and 0314 ``ctrlhflux_tut_MITgcm.opt0000`` to ``input_ad/OPTIM``, 0315 move into this directory and link :filelink:`input_ad/data.optim <verification/tutorial_global_oce_optim/input_ad/data.optim>` 0316 and :filelink:`input_ad/data.ctrl <verification/tutorial_global_oce_optim/input_ad/data.ctrl>` locally: 0317 0318 :: 0319 0320 % cd OPTIM/ 0321 % ln -s ../data.optim . 0322 % ln -s ../data.ctrl . 0323 0324 The target cost function :varlink:`fmin` needs to be specified 0325 in :filelink:`input_ad/data.optim <verification/tutorial_global_oce_optim/input_ad/data.optim>`: as a rule of 0326 thumb, it should be about 0.95-0.90 times the value of the cost function 0327 at the first iteration. This value is only used at the first iteration 0328 and does not need to be updated afterward. However, it implicitly

** Warning **

Wide character in print at /usr/local/share/lxr/source line 1030, <$git> line 330.

0329 specifies the “pace” at which the cost function is going down (if you 0330 are lucky and it does indeed diminish!). 0331 0332 Once this is done, run the line-search algorithm: 0333 0334 :: 0335 0336 % ./optim.x 0337 0338 which computes the updated adjustment for iteration 1, 0339 ``ctrlhflux_tut_MITgcm.opt0001``. 0340 0341 The following iterations can be executed automatically using the shell 0342 script :filelink:`input_ad/cycsh <verification/tutorial_global_oce_optim/input_ad/cycsh>`. This script will take care of 0343 changing the iteration numbers in :filelink:`input_ad/data.optim <verification/tutorial_global_oce_optim/input_ad/data.optim>`, launch the adjoint 0344 model, clean and store the outputs, move the ``costhflux*`` and ``ctrlhflux*`` 0345 files, and run the line-search algorithm. Edit :filelink:`input_ad/cycsh <verification/tutorial_global_oce_optim/input_ad/cycsh>` to specify the 0346 prefix of the directories used to store the outputs and the maximum 0347 number of iteration. 0348 0349 .. [1] 0350 Because of the daily automatic testing, the experiment as found in 0351 the repository is set-up with a very small number of time-steps. To 0352 reproduce the results shown here, one needs to set :varlink:`nTimeSteps` = 360 0353 and :varlink:`lastinterval` =31104000 (both corresponding to a year, see :numref:`tut_global_oce_runsect` for further details).