1 | % ================================================================ |
---|
2 | % Chapter stochastic parametrization of EOS (STO) |
---|
3 | % ================================================================ |
---|
4 | \chapter{Stochastic parametrization of EOS (STO)} |
---|
5 | \label{STO} |
---|
6 | |
---|
7 | Authors: P.-A. Bouttier |
---|
8 | |
---|
9 | \minitoc |
---|
10 | |
---|
11 | \newpage |
---|
12 | $\ $\newline % force a new line |
---|
13 | |
---|
14 | The stochastic parametrization module aims to explicitly simulate uncertainties in the model. More particularly, \cite{Brankart_OM2013} has shown that, because of the nonlinearity of the seawater equation of state, unresolved scales represent a major source of uncertainties in the computation of the large scale horizontal density gradient (from T/S large scale fields), and that the impact of these uncertainties can be simulated by random processes representing unresolved T/S fluctuations. |
---|
15 | |
---|
16 | The stochastic formulation of the equation of state can be written as: |
---|
17 | \begin{equation} |
---|
18 | \label{eq:eos_sto} |
---|
19 | \rho = \frac{1}{2} \sum_{i=1}^m\{ \rho[T+\Delta T_i,S+\Delta S_i,p_o(z)] + \rho[T-\Delta T_i,S-\Delta S_i,p_o(z)] \} |
---|
20 | \end{equation} |
---|
21 | where $p_o(z)$ is the reference pressure depending on the depth and $\Delta T_i$ and $\Delta S_i$ are a set of T/S perturbations defined as the scalar product of the respective local T/S gradients with random walks $\mathbf{\xi}$: |
---|
22 | \begin{equation} |
---|
23 | \label{eq:sto_pert} |
---|
24 | \Delta T_i = \mathbf{\xi}_i \cdot \nabla T \qquad \hbox{and} \qquad \Delta S_i = \mathbf{\xi}_i \cdot \nabla S |
---|
25 | \end{equation} |
---|
26 | $\mathbf{\xi}_i$ are produced by a first-order autoregressive processes (AR-1) with a parametrized decorrelation time scale, and horizontal and vertical standard deviations $\sigma_s$. $\mathbf{\xi}$ are uncorrelated over the horizontal and fully correlated along the vertical. |
---|
27 | |
---|
28 | |
---|
29 | \section{Stochastic processes} |
---|
30 | \label{STO_the_details} |
---|
31 | |
---|
32 | The starting point of our implementation of stochastic parameterizations |
---|
33 | in NEMO is to observe that many existing parameterizations are based |
---|
34 | on autoregressive processes, which are used as a basic source of randomness |
---|
35 | to transform a deterministic model into a probabilistic model. |
---|
36 | A generic approach is thus to add one single new module in NEMO, |
---|
37 | generating processes with appropriate statistics |
---|
38 | to simulate each kind of uncertainty in the model |
---|
39 | (see \cite{Brankart_al_GMD2015} for more details). |
---|
40 | |
---|
41 | In practice, at every model grid point, independent Gaussian autoregressive |
---|
42 | processes~$\xi^{(i)},\,i=1,\ldots,m$ are first generated |
---|
43 | using the same basic equation: |
---|
44 | |
---|
45 | \begin{equation} |
---|
46 | \label{eq:autoreg} |
---|
47 | \xi^{(i)}_{k+1} = a^{(i)} \xi^{(i)}_k + b^{(i)} w^{(i)} + c^{(i)} |
---|
48 | \end{equation} |
---|
49 | |
---|
50 | \noindent |
---|
51 | where $k$ is the index of the model timestep; and |
---|
52 | $a^{(i)}$, $b^{(i)}$, $c^{(i)}$ are parameters defining |
---|
53 | the mean ($\mu^{(i)}$) standard deviation ($\sigma^{(i)}$) |
---|
54 | and correlation timescale ($\tau^{(i)}$) of each process: |
---|
55 | |
---|
56 | \begin{itemize} |
---|
57 | \item for order~1 processes, $w^{(i)}$ is a Gaussian white noise, |
---|
58 | with zero mean and standard deviation equal to~1, and the parameters |
---|
59 | $a^{(i)}$, $b^{(i)}$, $c^{(i)}$ are given by: |
---|
60 | |
---|
61 | \begin{equation} |
---|
62 | \label{eq:ord1} |
---|
63 | \left\{ |
---|
64 | \begin{array}{l} |
---|
65 | a^{(i)} = \varphi \\ |
---|
66 | b^{(i)} = \sigma^{(i)} \sqrt{ 1 - \varphi^2 } |
---|
67 | \qquad\qquad\mbox{with}\qquad\qquad |
---|
68 | \varphi = \exp \left( - 1 / \tau^{(i)} \right) \\ |
---|
69 | c^{(i)} = \mu^{(i)} \left( 1 - \varphi \right) \\ |
---|
70 | \end{array} |
---|
71 | \right. |
---|
72 | \end{equation} |
---|
73 | |
---|
74 | \item for order~$n>1$ processes, $w^{(i)}$ is an order~$n-1$ autoregressive process, |
---|
75 | with zero mean, standard deviation equal to~$\sigma^{(i)}$; correlation timescale |
---|
76 | equal to~$\tau^{(i)}$; and the parameters |
---|
77 | $a^{(i)}$, $b^{(i)}$, $c^{(i)}$ are given by: |
---|
78 | |
---|
79 | \begin{equation} |
---|
80 | \label{eq:ord2} |
---|
81 | \left\{ |
---|
82 | \begin{array}{l} |
---|
83 | a^{(i)} = \varphi \\ |
---|
84 | b^{(i)} = \frac{n-1}{2(4n-3)} \sqrt{ 1 - \varphi^2 } |
---|
85 | \qquad\qquad\mbox{with}\qquad\qquad |
---|
86 | \varphi = \exp \left( - 1 / \tau^{(i)} \right) \\ |
---|
87 | c^{(i)} = \mu^{(i)} \left( 1 - \varphi \right) \\ |
---|
88 | \end{array} |
---|
89 | \right. |
---|
90 | \end{equation} |
---|
91 | |
---|
92 | \end{itemize} |
---|
93 | |
---|
94 | \noindent |
---|
95 | In this way, higher order processes can be easily generated recursively using the same piece of code implementing Eq.~(\ref{eq:autoreg}), and using succesively processes from order $0$ to~$n-1$ as~$w^{(i)}$. |
---|
96 | The parameters in Eq.~(\ref{eq:ord2}) are computed so that this recursive application |
---|
97 | of Eq.~(\ref{eq:autoreg}) leads to processes with the required standard deviation |
---|
98 | and correlation timescale, with the additional condition that |
---|
99 | the $n-1$ first derivatives of the autocorrelation function |
---|
100 | are equal to zero at~$t=0$, so that the resulting processes |
---|
101 | become smoother and smoother as $n$ is increased. |
---|
102 | |
---|
103 | Overall, this method provides quite a simple and generic way of generating a wide class of stochastic processes. However, this also means that new model parameters are needed to specify each of these stochastic processes. As in any parameterization of lacking physics, a very important issues then to tune these new parameters using either first principles, model simulations, or real-world observations. |
---|
104 | |
---|
105 | \section{Implementation details} |
---|
106 | \label{STO_thech_details} |
---|
107 | The computer code implementing stochastic parametrisations is made of one single FORTRAN module, |
---|
108 | with 3 public routines to be called by the model (in our case, NEMO): |
---|
109 | |
---|
110 | The first routine ({sto\_par}) is a direct implementation of Eq.~(\ref{eq:autoreg}), |
---|
111 | applied at each model grid point (in 2D or 3D), |
---|
112 | and called at each model time step ($k$) to update |
---|
113 | every autoregressive process ($i=1,\ldots,m$). |
---|
114 | This routine also includes a filtering operator, applied to $w^{(i)}$, |
---|
115 | to introduce a spatial correlation between the stochastic processes. |
---|
116 | |
---|
117 | The second routine ({sto\_par\_init}) |
---|
118 | is an initialization routine mainly dedicated |
---|
119 | to the computation of parameters $a^{(i)}, b^{(i)}, c^{(i)}$ |
---|
120 | for each autoregressive process, as a function of the statistical properties |
---|
121 | required by the model user (mean, standard deviation, time correlation, |
---|
122 | order of the process,\ldots). Parameters for the processes can be specified through the following namelist parameters: |
---|
123 | \begin{alltt} |
---|
124 | \tiny |
---|
125 | \begin{verbatim} |
---|
126 | nn_sto_eos = 1 ! number of independent random walks |
---|
127 | rn_eos_stdxy = 1.4 ! random walk horz. standard deviation (in grid points) |
---|
128 | rn_eos_stdz = 0.7 ! random walk vert. standard deviation (in grid points) |
---|
129 | rn_eos_tcor = 1440.0 ! random walk time correlation (in timesteps) |
---|
130 | nn_eos_ord = 1 ! order of autoregressive processes |
---|
131 | nn_eos_flt = 0 ! passes of Laplacian filter |
---|
132 | rn_eos_lim = 2.0 ! limitation factor (default = 3.0) |
---|
133 | \end{verbatim} |
---|
134 | \end{alltt} |
---|
135 | This routine also includes the initialization (seeding) |
---|
136 | of the random number generator. |
---|
137 | |
---|
138 | The third routine ({sto\_rst\_write}) writes a ``restart file'' |
---|
139 | with the current value of all autoregressive processes |
---|
140 | to allow restarting a simulation from where it has been interrupted. |
---|
141 | This file also contains the current state of the random number generator. |
---|
142 | In case of a restart, this file is then read by the initialization routine |
---|
143 | ({sto\_par\_init}), so that the simulation can continue exactly |
---|
144 | as if it was not interrupted. |
---|
145 | Restart capabilities of the module are driven by the following namelist parameters: |
---|
146 | \begin{alltt} |
---|
147 | \tiny |
---|
148 | \begin{verbatim} |
---|
149 | ln_rststo = .false. ! start from mean parameter (F) or from restart file (T) |
---|
150 | ln_rstseed = .true. ! read seed of RNG from restart file |
---|
151 | cn_storst_in = "restart_sto" ! suffix of stochastic parameter restart file (input) |
---|
152 | cn_storst_out = "restart_sto" ! suffix of stochastic parameter restart file (output) |
---|
153 | \end{verbatim} |
---|
154 | \end{alltt} |
---|
155 | |
---|
156 | In the particular case of the stochastic equation of state, there is also an additional module ({sto\_pts}) implementing Eq~\ref{eq:sto_pert} and specific piece of code in the equation of state implementing Eq~\ref{eq:eos_sto}. |
---|
157 | |
---|
158 | |
---|