Skip to main content
Log in

Improving Prediction Skill of Imperfect Turbulent Models Through Statistical Response and Information Theory

  • Published:
Journal of Nonlinear Science Aims and scope Submit manuscript

Abstract

Turbulent dynamical systems with a large phase space and a high degree of instabilities are ubiquitous in climate science and engineering applications. Statistical uncertainty quantification (UQ) to the response to the change in forcing or uncertain initial data in such complex turbulent systems requires the use of imperfect models due to the lack of both physical understanding and the overwhelming computational demands of Monte Carlo simulation with a large-dimensional phase space. Thus, the systematic development of reduced low-order imperfect statistical models for UQ in turbulent dynamical systems is a grand challenge. This paper applies a recent mathematical strategy for calibrating imperfect models in a training phase and accurately predicting the response by combining information theory and linear statistical response theory in a systematic fashion. A systematic hierarchy of simple statistical imperfect closure schemes for UQ for these problems is designed and tested which are built through new local and global statistical energy conservation principles combined with statistical equilibrium fidelity. The forty mode Lorenz 96 (L-96) model which mimics forced baroclinic turbulence is utilized as a test bed for the calibration and predicting phases for the hierarchy of computationally cheap imperfect closure models both in the full phase space and in a reduced three-dimensional subspace containing the most energetic modes. In all of phase spaces, the nonlinear response of the true model is captured accurately for the mean and variance by the systematic closure model, while alternative methods based on the fluctuation-dissipation theorem alone are much less accurate. For reduced-order model for UQ in the three-dimensional subspace for L-96, the systematic low-order imperfect closure models coupled with the training strategy provide the highest predictive skill over other existing methods for general forced response yet have simple design principles based on a statistical global energy equation. The systematic imperfect closure models and the calibration strategies for UQ for the L-96 model serve as a new template for similar strategies for UQ with model error in vastly more complex realistic turbulent dynamical systems.

This is a preview of subscription content, log in via an institution to check access.

Access this article

Price excludes VAT (USA)
Tax calculation will be finalised during checkout.

Instant access to the full article PDF.

Institutional subscriptions

Fig. 1
Fig. 2
Fig. 3
Fig. 4
Fig. 5
Fig. 6
Fig. 7
Fig. 8
Fig. 9
Fig. 10
Fig. 11
Fig. 12
Fig. 13
Fig. 14
Fig. 15

Similar content being viewed by others

References

  • Abramov, R.V., Majda, A.J.: Blended response algorithms for linear fluctuation-dissipation for complex nonlinear dynamical systems. Nonlinearity 20(12), 2793 (2007)

    Article  MathSciNet  MATH  Google Scholar 

  • Abramov, R.V., Majda, A.J.: New approximations and tests of linear fluctuation-response for chaotic nonlinear forced-dissipative dynamical systems. J. Nonlinear Sci. 18(3), 303–341 (2008)

    Article  MathSciNet  MATH  Google Scholar 

  • Abramov, R.V., Majda, A.J.: A new algorithm for low-frequency climate response. J. Atmos. Sci. 66(2), 286–309 (2009)

    Article  Google Scholar 

  • Bell, T.L.: Climate sensitivity from fluctuation dissipation: some simple model tests. J. Atmos. Sci. 37(8), 1700–1707 (1980)

    Article  Google Scholar 

  • Bengtsson, T., Bickel, P., Li, B., et al.: Curse-of-dimensionality revisited: collapse of the particle filter in very large scale systems. In: Probability and Statistics: Essays in Honor of David A. Freedman, pp. 316–334. Institute of Mathematical Statistics (2008)

  • Branicki, M., Majda, A.J.: Quantifying uncertainty for predictions with model error in non-gaussian systems with intermittency. Nonlinearity 25(9), 2543 (2012)

    Article  MathSciNet  MATH  Google Scholar 

  • Branicki, M., Chen, N., Majda, A.J.: Non-Gaussian test models for prediction and state estimation with model errors. Chin. Ann. Math. Ser. B 34(1), 29–64 (2013)

    Article  MathSciNet  MATH  Google Scholar 

  • Carnevale, G., Falcioni, M., Isola, S., Purini, R., Vulpiani, A.: Fluctuation–response relations in systems with chaotic behavior. Phys. Fluids A Fluid Dyn. (1989–1993) 3(9), 2247–2254 (1991)

    Article  MATH  Google Scholar 

  • Csiszar, I., Körner, J.: Information Theory: Coding Theorems for Discrete Memoryless Systems. Cambridge University Press, Cambridge (2011)

    Book  Google Scholar 

  • DelSole, T.: Predictability and information theory. Part II: imperfect forecasts. J. Atmos. Sci. 62(9), 3368–3381 (2005)

    Article  Google Scholar 

  • DelSole, T., Shukla, J.: Model fidelity versus skill in seasonal forecasting. J. Clim. 23(18), 4794–4806 (2010)

    Article  Google Scholar 

  • Emanuel, K., Wyngaard, J., McWilliams, J., Randall, D., Yung, Y.: Improving the scientific foundation for atmosphere-land-ocean simulations. National Academic Press, Washington, DC (2005). 72

    Google Scholar 

  • Epstein, E.S.: Stochastic dynamic prediction1. Tellus 21(6), 739–759 (1969)

    Article  Google Scholar 

  • Gershgorin, B., Majda, A.J.: Quantifying uncertainty for climate change and long-range forecasting scenarios with model errors. part i: Gaussian models. J. Clim. 25(13), 4523–4548 (2012)

    Article  Google Scholar 

  • Gritsun, A., Branstator, G.: Climate response using a three-dimensional operator based on the fluctuation–dissipation theorem. J. Atmos. Sci. 64(7), 2558–2575 (2007)

    Article  Google Scholar 

  • Gritsun, A., Branstator, G., Majda, A.: Climate response of linear and quadratic functionals using the fluctuation–dissipation theorem. J. Atmos. Sci. 65(9), 2824–2841 (2008)

    Article  Google Scholar 

  • Hairer, M., Majda, A.J.: A simple framework to justify linear response theory. Nonlinearity 23(4), 909 (2010)

    Article  MathSciNet  MATH  Google Scholar 

  • Kullback, S., Leibler, R.A.: On information and sufficiency. Ann. Math. Stat. 22(1), 79–86 (1951)

    Article  MathSciNet  MATH  Google Scholar 

  • Leith, C.: Climate response and fluctuation dissipation. J. Atmos. Sci. 32(10), 2022–2026 (1975)

    Article  Google Scholar 

  • Lorenz, E.N.: Predictability: a problem partly solved. In: Proceedings of Seminar on Predictability, vol. 1 (1996)

  • Majda, A.: Introduction to PDEs and Waves for the Atmosphere and Ocean, vol. 9. American Mathematical Society, Providence (2003)

    MATH  Google Scholar 

  • Majda, A., Wang, X.: Nonlinear Dynamics and Statistical Theories for Basic Geophysical Flows. Cambridge University Press, Cambridge (2006)

    Book  MATH  Google Scholar 

  • Majda, A., Abramov, R.V., Grote, M.J.: Information Theory and Stochastics for Multiscale Nonlinear Systems, vol. 25. American Mathematical Society, Providence (2005)

    MATH  Google Scholar 

  • Majda, A., Kleeman, R., Cai, D.: A mathematical framework for quantifying predictability through relative entropy. Methods Appl. Anal. 9(3), 425–444 (2002)

    MathSciNet  MATH  Google Scholar 

  • Majda, A., Wang, X., et al.: Linear response theory for statistical ensembles in complex systems with time-periodic forcing. Commun. Math. Sci. 8(1), 145–172 (2010)

    Article  MathSciNet  MATH  Google Scholar 

  • Majda, A.J.: Statistical energy conservation principle for inhomogeneous turbulent dynamical systems. Proc. Natl. Acad. Sci. USA 112(29), 8937–8941 (2015)

    Article  Google Scholar 

  • Majda, A.J., Branicki, M.: Lessons in uncertainty quantification for turbulent dynamical systems. Discrete Cont. Dyn. Syst. 32(9), 3133–3221 (2012)

    Article  MathSciNet  MATH  Google Scholar 

  • Majda, A.J., Gershgorin, B.: Quantifying uncertainty in climate change science through empirical information theory. Proc. Natl. Acad. Sci. USA 107(34), 14958–14963 (2010)

    Article  MathSciNet  Google Scholar 

  • Majda, A.J., Gershgorin, B.: Improving model fidelity and sensitivity for complex systems through empirical information theory. Proc. Natl. Acad. Sci. USA 108(25), 10044–10049 (2011)

    Article  MathSciNet  MATH  Google Scholar 

  • Majda, A.J., Gershgorin, B.: Link between statistical equilibrium fidelity and forecasting skill for complex systems with model error. Proc. Natl. Acad. Sci. USA 108(31), 12599–12604 (2011b)

    Article  MathSciNet  MATH  Google Scholar 

  • Majda, A.J., Abramov, R., Gershgorin, B.: High skill in low-frequency climate response through fluctuation dissipation theorems despite structural instability. Proc. Natl. Acad. Sci. USA 107(2), 581–586 (2010a)

    Article  MathSciNet  MATH  Google Scholar 

  • Majda, A.J., Gershgorin, B., Yuan, Y.: Low-frequency climate response and fluctuation-dissipation theorems: theory and practice. J. Atmos. Sci. 67(4), 1186–1201 (2010b)

    Article  Google Scholar 

  • Majda, A.J., Qi, D., Sapsis, T.P.: Blended particle filters for large-dimensional chaotic dynamical systems. Proc. Natl. Acad. Sci. USA 111(21), 7511–7516 (2014)

    Article  MathSciNet  Google Scholar 

  • Neelin, J., Münnich, M., Su, H., Meyerson, J., Holloway, C.: Tropical drying trends in global warming models and observations. Proc. Natl. Acad. Sci. USA 103(16), 6110–6115 (2006)

    Article  Google Scholar 

  • Qi, D., Majda, A.J.: Blended particle methods with adaptive subspaces for filtering turbulent dynamical systems. Phys. D 298, 21–41 (2015)

    Article  MathSciNet  Google Scholar 

  • Sapsis, T.P., Majda, A.J.: Blended reduced subspace algorithms for uncertainty quantification of quadratic systems with a stable mean state. Phys. D 258, 61–76 (2013a)

    Article  MathSciNet  MATH  Google Scholar 

  • Sapsis, T.P., Majda, A.J.: Blending modified Gaussian closure and non-gaussian reduced subspace methods for turbulent dynamical systems. J. Nonlinear Sci. 23(6), 1039–1071 (2013b)

    Article  MathSciNet  MATH  Google Scholar 

  • Sapsis, T.P., Majda, A.J.: Statistically accurate low-order models for uncertainty quantification in turbulent dynamical systems. Proc. Natl. Acad. Sci. USA 110(34), 13705–13710 (2013c)

    Article  MathSciNet  MATH  Google Scholar 

  • Sapsis, T.P., Majda, A.J.: A statistically accurate modified quasilinear Gaussian closure for uncertainty quantification in turbulent dynamical systems. Phys. D 252, 34–45 (2013d)

    Article  MathSciNet  MATH  Google Scholar 

  • Tung, S.: On lower and upper bounds of the difference between the arithmetic and the geometric mean. Math. Comput. 29(131), 834–836 (1975)

    Article  MathSciNet  MATH  Google Scholar 

  • Yaglom, A.M.: An introduction to the theory of stationary random functions. Courier Corporation, New York (2004)

    Google Scholar 

Download references

Acknowledgments

The research of Andrew Majda is partially supported by Office of Naval Research Grant ONR MURI N00014-12-1-0912. Di Qi is supported as a graduate research assistant on this grant.

Author information

Authors and Affiliations

Authors

Corresponding author

Correspondence to Di Qi.

Additional information

Communicated by Paul Newton.

Appendices

Appendix 1: Derivation of the Moment Equations for L-96 System

Here we display the detailed derivations about the moment equations and the properties under homogeneous assumption in Sect. 3.1. Under the Fourier representation of the basis \(\left\{ \mathbf {v}_{j}\right\} \), we can write explicit formulas for each part of the L-96 system (16) as follows

  1. 1.

    The quadratic interaction between Fourier mode \(\mathbf {v}_{i}\) and \(\mathbf {v}_{j}\) can be written explicitly from the definition

    $$\begin{aligned} \mathbf {B}\left( \mathbf {v}_{i},\mathbf {v}_{j}\right)= & {} \left\{ v_{i}^{l-1*}\left( v_{j}^{l+1}-v_{j}^{l-2}\right) \right\} _{l=0}^{J-1}\\= & {} \left\{ \frac{1}{J}e^{-2\pi \imath i\frac{l-1}{J}}\left( e^{2\pi \imath j\frac{l+1}{J}}-e^{-2\pi \imath j\frac{l-2}{J}}\right) \right\} _{l=0}^{J-1}\\= & {} \frac{1}{\sqrt{J}}e^{2\pi \imath \frac{i+j}{J}}\mathbf {v}_{j-i}-\frac{1}{\sqrt{J}}e^{2\pi \imath \frac{i-2j}{J}}\mathbf {v}_{j-i}\\= & {} \frac{1}{\sqrt{J}}e^{2\pi \imath \frac{i}{J}}\left( e^{2\pi \imath \frac{j}{J}}-e^{2\pi \imath \frac{-2j}{J}}\right) \mathbf {v}_{j-i}. \end{aligned}$$

    Particularly, we can find that the nonlinear interaction term vanishes if and only if its second component is only under the zero base mode \(\mathbf {v}_{0}=\frac{1}{\sqrt{J}}\left( 1,\ldots ,1\right) ^{T}\), that is,

    $$\begin{aligned} \mathbf {B}\left( \mathbf {v}_{i},\mathbf {v}_{j}\right) =\mathbf {0}\Leftrightarrow e^{2\pi \imath \frac{3j}{J}}=1\Leftrightarrow j=\frac{nJ}{3},\; n\in \mathbb {N},\, J=40\Leftrightarrow j=0. \end{aligned}$$
  2. 2.

    Using the expression above, we can calculate the explicit formulas for the second-order interaction form in the mean dynamics

    $$\begin{aligned} \sum _{i,j}R_{ij}\mathbf {B}\left( \mathbf {v}_{i},\mathbf {v}_{j}\right)= & {} \frac{1}{\sqrt{J}}\sum _{i,j}r_{ij}e^{2\pi \imath \frac{i}{J}}\left( e^{2\pi \imath \frac{j}{J}}-e^{2\pi \imath \frac{-2j}{J}}\right) \mathbf {v}_{j-i}\\= & {} \left\{ \frac{1}{J}\sum _{i,j}r_{ij}e^{-2\pi \imath i\frac{l-1}{J}}\left( e^{2\pi \imath j\frac{l+1}{J}}-e^{-2\pi \imath j\frac{l-2}{J}}\right) \right\} _{l=0}^{J-1}. \end{aligned}$$

    Then to calculate the explicit formula for the covariance dynamics, we need to get the linear interaction part \(L_{v}\) first

    $$\begin{aligned} \mathbf {B}\left( \bar{\mathbf {u}},\mathbf {v}_{j}\right) _{l}= & {} \frac{1}{\sqrt{J}}\bar{u}_{l-1}\left( e^{2\pi \imath j\frac{l+1}{J}}-e^{2\pi \imath j\frac{l-2}{J}}\right) ,\\ \mathbf {B}\left( \mathbf {v}_{j},\bar{\mathbf {u}}\right) _{l}= & {} \frac{1}{\sqrt{J}}e^{-2\pi \imath j\frac{l-1}{J}}\left( \bar{u}_{l+1}-\bar{u}_{l-2}\right) . \end{aligned}$$

    Therefore, by definition, we get

    $$\begin{aligned} L_{v,ij}= & {} \left( \mathbf {L}\left( t\right) \mathbf {v}_{j}+\mathbf {B}\left( \bar{\mathbf {u}},\mathbf {v}_{j}\right) +\mathbf {B}\left( \mathbf {v}_{j},\bar{\mathbf {u}}\right) \right) \cdot \mathbf {v}_{i}^{*}\\= & {} -d_{j}\left( t\right) \delta _{ij}+\frac{1}{J}\left( e^{2\pi \imath \frac{j}{J}}-e^{-2\pi \imath \frac{2j}{J}}\right) \sum _{l}\bar{u}_{l-1}e^{2\pi \imath l\frac{j-i}{J}}\\&+\frac{1}{J}e^{2\pi \imath \frac{j}{J}}\sum _{l}\left( \bar{u}_{l+1}-\bar{u}_{l-2}\right) e^{-2\pi \imath \frac{j+i}{J}}. \end{aligned}$$

    where the index l represents the l-th component of a vector.

  3. 3.

    Finally we need to calculate the formula for the nonlinear flux term \(Q_{F}\) in the covariance dynamics. Again use the explicit formula for the nonlinear interaction term between modes \(\mathbf {B}\left( \mathbf {v}_{i},\mathbf {v}_{j}\right) \) and the definition of the nonlinear flux, we have

    $$\begin{aligned} Q_{F,ij}= & {} \sum _{m,n}\left\langle Z_{m}Z_{n}^{*}Z_{j} \right\rangle \mathbf {B}\left( \mathbf {v}_{m},\mathbf {v}_{n}\right) ^{*}\cdot \mathbf {v}_{i} +\left\langle Z_{m}^{*}Z_{n}Z_{i}^{*} \right\rangle \mathbf {B}\left( \mathbf {v}_{m},\mathbf {v}_{n}\right) \cdot \mathbf {v}_{j}^{*}\\= & {} \frac{1}{\sqrt{J}}\sum _{m,n}\left\langle Z_{m}Z_{n}^{*}Z_{j}\right\rangle e^{-2\pi \imath \frac{m}{J}}\left( e^{-2\pi \imath \frac{n}{J}}-e^{2\pi \imath \frac{2n}{J}}\right) \delta _{n-m,i}\\&+\frac{1}{\sqrt{J}}\sum _{m,n}\left\langle Z_{m}^{*}Z_{n}Z_{i}^{*}\right\rangle e^{2\pi \imath \frac{m}{J}}\left( e^{2\pi \imath \frac{n}{J}}-e^{-2\pi \imath \frac{2n}{J}}\right) \delta _{n-m,j}. \end{aligned}$$

Therefore, we get the explicit forms for each part of the moment equations for the mean and covariance matrix.

Another important issue is due to the simplification of each order of moments under the homogeneous assumption of the L-96 system. Specifically, it tells that if each order of moments is invariant under shifting in grid points

$$\begin{aligned} \left\langle u_{i_{1}}u_{i_{2}}\ldots u_{i_{n}}\right\rangle =\left\langle u_{i_{1}+l}u_{i_{2}+l}\ldots u_{i_{n}+l}\right\rangle ,\quad \forall \, l, \end{aligned}$$

then the first three moments under the Fourier basis become

$$\begin{aligned} \bar{\mathbf {u}}\left( t\right)= & {} \bar{u}\left( t\right) \left( 1,1,\dots ,1\right) ^{T},\\ R\left( t\right)= & {} \mathrm {diag}\left( r_{-J/2+1}\left( t\right) ,\ldots ,r_{0}\left( t\right) ,\dots ,r_{J/2}\left( t\right) \right) ,\\ \left\langle Z_{i}Z_{j}Z_{k}\right\rangle\ne & {} 0,\quad \mathrm {only}\,\mathrm {if}\; i+j+k=0. \end{aligned}$$

The first equation for the mean state \(\bar{\mathbf {u}}\) is direct from the definition. To get the simplified forms for the second and third moments R and \(\left\langle Z_{i}Z_{j}Z_{k}\right\rangle \), first we write the Fourier transform as a change of basis as

$$\begin{aligned} \mathbf {u}^{\prime }= & {} V\mathbf {Z},\\ Z_{j}= & {} \sum _{l}\mathbf {v}_{j}^{l*}u_{l}^{\prime }. \end{aligned}$$

where \(\mathbf {u}^{\prime }=\mathbf {u}-\bar{\mathbf {u}}=\left( u_{1}^{\prime },u_{2}^{\prime },\ldots ,u_{J}^{\prime }\right) ^{T}\), \(\mathbf {Z}=\left( Z_{-{J}/{2}+1},Z_{-{J}/{2}+2},\ldots ,Z_{{J}/{2}}\right) ^{T}\) are the coefficients under natural basis and Fourier basis separately, and \(V=\left( \mathbf {v}_{-{J}/{2}+1},\mathbf {v}_{-{J}/{2}+2},\ldots ,\mathbf {v}_{{J}/{2}}\right) \) is the transformation matrix formed by Fourier basis. Therefore, the components of the second-order moments under the transform invariant property become

$$\begin{aligned} R_{ij}= & {} \left\langle Z_{i}Z_{j}^{*}\right\rangle =\sum _{m,n}v_{i}^{m*}\left\langle u_{m}^{\prime }u_{n}^{\prime }\right\rangle v_{j}^{n}\\= & {} \sum _{n}\sum _{l}\left\langle u_{0}^{\prime }u_{n}^{\prime }\right\rangle v_{i}^{l*}v_{j}^{n+l}\\= & {} \sum _{n}\left\langle u_{0}^{\prime }u_{n}^{\prime }\right\rangle e^{2\pi \imath \frac{nj}{J}}\sum _{l}e^{2\pi \imath l\frac{j-i}{J}}\\= & {} J\sum _{n}\left\langle u_{0}^{\prime }u_{n}^{\prime }\right\rangle e^{2\pi \imath \frac{nj}{J}}\delta _{ij}. \end{aligned}$$

Here the second equality is using the homogeneity in physical space \(\left\langle u_{m}^{\prime }u_{n}^{\prime }\right\rangle =\left\langle u_{0}^{\prime }u_{n-m}^{\prime }\right\rangle \). Therefore, we can see that the covariance matrix becomes diagonal under the homogeneous assumption and each component takes the form

$$\begin{aligned} r_{j}=J\sum _{n}\left\langle u_{0}^{\prime }u_{n}^{\prime }\right\rangle e^{2\pi \imath \frac{nj}{J}}. \end{aligned}$$

In the same way, we can calculate the third moments as

$$\begin{aligned} \left\langle Z_{i}Z_{j}Z_{k}\right\rangle= & {} \sum _{m,n,s}\left\langle v_{i}^{m*}u_{m}^{\prime }v_{j}^{n*}u_{n}^{\prime }v_{k}^{s*}u_{s}^{\prime }\right\rangle \\= & {} \sum _{m,n,s}\left\langle u_{m}^{\prime }u_{n}^{\prime }u_{s}^{\prime }\right\rangle v_{i}^{m*}v_{j}^{n*}v_{k}^{s*}\\= & {} \sum _{n,s}\sum _{l}\left\langle u_{0}^{\prime }u_{n}^{\prime }u_{s}^{\prime }\right\rangle v_{i}^{l*}v_{j}^{n+l*}v_{k}^{s+l*}\\= & {} \sum _{n,s}\left\langle u_{0}^{\prime }u_{n}^{\prime }u_{s}^{\prime }\right\rangle e^{-2\pi \imath \frac{nj+ks}{J}}\sum _{l}e^{2\pi \imath l\frac{i+j+k}{J}}\\= & {} J\sum _{n,s}\left\langle u_{0}^{\prime }u_{n}^{\prime }u_{s}^{\prime }\right\rangle e^{-2\pi \imath \frac{nj+ks}{J}}\delta _{i+j+k}. \end{aligned}$$

The same we use the homogeneous property \(\left\langle u_{m}^{\prime }u_{n}^{\prime }u_{s}^{\prime }\right\rangle =\left\langle u_{0}^{\prime }u_{n-m}^{\prime }u_{s-m}^{\prime }\right\rangle \). Therefore, the third-order moments can only be nonzeros when the wavenumber satisfies \(i+j+k=0\), and the nonzero terms can be represented as

$$\begin{aligned} \left\langle Z_{-j-k}Z_{j}Z_{k}\right\rangle =J\sum _{n,s}\left\langle u_{0}^{\prime }u_{n}^{\prime }u_{s}^{\prime }\right\rangle e^{-2\pi \imath \frac{nj+ks}{J}}. \end{aligned}$$

With all these homogeneous properties, both the linear and nonlinear flux terms \(L_{v}\) and \(Q_{F}\) become diagonal. The detailed calculation is shown as follows. First, each part of the linear interaction term \(L_{v}\) simplifies

$$\begin{aligned} \sum _{i,j}R_{ij}\mathbf {B}\left( \mathbf {v}_{i},\mathbf {v}_{j}\right)= & {} \sum _{i}r_{i}\left( t\right) \mathbf {B}\left( \mathbf {v}_{i},\mathbf {v}_{i}\right) =\frac{1}{\sqrt{J}}\sum _{i}r_{i}\left( e^{2\pi \imath \frac{2i}{J}}-e^{-2\pi \imath \frac{i}{J}}\right) \mathbf {v}_{0},\\ \mathbf {B}\left( \mathbf {v}_{i},\bar{\mathbf {u}}\right)= & {} \bar{u}\left( t\right) \mathbf {B}\left( \mathbf {v}_{i},\mathbf {v}_{0}\right) =\mathbf {0},\\ \mathbf {B}\left( \bar{\mathbf {u}},\mathbf {v}_{i}\right)= & {} \bar{u}\left( t\right) \mathbf {B}\left( \mathbf {v}_{0},\mathbf {v}_{i}\right) =\bar{u}\left( t\right) \left( e^{2\pi \imath \frac{i}{J}}-e^{-2\pi \imath \frac{2i}{J}}\right) \mathbf {v}_{i},\\ L_{v,ij}= & {} -d\left( t\right) \delta _{ij}+\bar{u}\left( t\right) \left( e^{2\pi \imath \frac{j}{J}}-e^{-2\pi \imath \frac{2j}{J}}\right) \delta _{ij}. \end{aligned}$$

Then the simplified form of the nonlinear flux term \(Q_{F,}\) can also be calculated

$$\begin{aligned} Q_{F,ij}= & {} \frac{1}{\sqrt{J}}\sum _{m,n}\left\langle Z_{m}Z_{n}^{*}Z_{j}\right\rangle e^{-2\pi \imath \frac{m}{J}}\left( e^{-2\pi \imath \frac{n}{J}}-e^{2\pi \imath \frac{2n}{J}}\right) \delta _{n-m,i}\\&+\frac{1}{\sqrt{J}}\sum _{m,n}\left\langle Z_{m}^{*}Z_{n}Z_{i}^{*}\right\rangle e^{2\pi \imath \frac{m}{J}}\left( e^{2\pi \imath \frac{n}{J}}-e^{-2\pi \imath \frac{2n}{J}}\right) \delta _{n-m,j}\\= & {} \frac{1}{\sqrt{J}}\sum _{m,n}\left\langle Z_{m}Z_{-n}Z_{j}\right\rangle e^{-2\pi \imath \frac{m}{J}}\left( e^{-2\pi \imath \frac{n}{J}}-e^{2\pi \imath \frac{2n}{J}}\right) \delta _{n-m,i}\delta _{m-n+j}\\&+\frac{1}{\sqrt{J}}\sum _{m,n}\left\langle Z_{-m}Z_{n}Z_{-i}\right\rangle e^{2\pi \imath \frac{m}{J}}\left( e^{2\pi \imath \frac{n}{J}}-e^{-2\pi \imath \frac{2n}{J}}\right) \delta _{n-m,j}\delta _{-m+n-i}\\= & {} \frac{2}{\sqrt{J}}\sum _{m}\mathfrak {Re}\left\{ \left\langle Z_{m}Z_{-m-j}Z_{j}\right\rangle \left( e^{-2\pi \imath \frac{2m+j}{J}}-e^{2\pi \imath \frac{m+2j}{J}}\right) \right\} \delta _{ij}. \end{aligned}$$

Appendix 2: Numerical Strategies to Calculate the Kicked Response Operators

In Sect. 4, we use the kicked response theory to tune the imperfect model parameters in the training phase. Here we describe the details about calculating the kicked response operators for the mean and variance numerically. From the formula in (10), the response operators for the mean and variance can be achieved from the perturbation part of the probability density \(\delta \pi ^{\prime }\). And this density function is also used to measure the information distance between the truth and imperfect model result in the training phase. Below we describe the numerical procedure to get this distribution function \(\delta \pi ^{\prime }\) for the true system and the imperfect closure model separately.

  • Kicked response for the true model: For the true system, we want to achieve the most accurate possible estimation for the response operators both for comparison with the imperfect model results and for calculating the FDT linear prediction. Therefore, we use a Monte Carlo simulation with an ensemble size of 10, 000 particles to capture the response in density. The initial equilibrium ensemble is picked by sampling from a normal distribution with consistent equilibrium mean and variance of the true system. For the kicked response to the mean a constant perturbation with 10 percent of the equilibrium state, mean \(\delta \mathbf {u}=0.1\bar{\mathbf {u}}_{\infty }\) is added to each initial ensemble member (in fact, as observed in numerical experiments, this perturbation amplitude has little effect on the results in the response distribution as long as it is not too large); and the initial variance of the ensemble is kept unchanged. The response distribution \(\delta \pi ^{\prime }\) then is achieved by monitoring the decay of the ensemble particles back to equilibrium under unperturbed dynamics and uniformly perturbed initial value (and the length of the time window that we need to monitor depends on the mixing property of the turbulent system). See Bell (1980), Majda et al. (2005) for similar version of this alogrithm.

  • Kicked response for the imperfect model: For the imperfect model, we just need to run the closure equations to get the responses for the mean and variance. In the same way as the true model, the initial mean is taken from the equilibrium distribution and a perturbation with amplitude \(\delta \mathbf {u}=0.1\bar{\mathbf {u}}_{\infty }\) is added to the mean initial state. The initial value for the variance is taken the same as the equilibrium state value and kept unperturbed. Then using this initial mean and variance, the imperfect model with specific closure strategies is applied to monitor the decay of the mean and variance back to equilibrium.

One additional important point that requires attention is that even if the unperturbed equilibrium initial conditions are applied, the system will still deviate from the equilibrium state first and reapproach equilibrium again after some relaxation time. This is due to the insufficient characterization of the entire distribution of the true system with a Gaussian approximation (note that nonlinearities are also included in the imperfect closure methods). To eliminate this effect in computing the kicked response in both the true and imperfect models, we subtract the statistics computed using the unperturbed initial value from the statistics computed using the perturbed Gaussian initial condition to achieve more accurate characterization of the responses.

Appendix 3: Calculating Equilibrium Consistent Parameters for Statistical Closure Models

The improved statistical closure models discussed in Sect. 4 require equilibrium fidelity at stationary steady state as \(t\rightarrow \infty \) as a necessary condition. We proposed the climate consistent parameter values in (41) and (42) correspondingly for GC1 and GC2 models calculated through simple algebraic manipulations. In stationary steady state with uniform damping and forcing terms, each order of moments of the closure system (36) converges to the equilibrium state such that

$$\begin{aligned} \frac{{\text {d}}\bar{u}_{M,\infty }}{{\text {d}}t}=0,\quad \frac{{\text {d}}r_{k,\infty }^{M}}{{\text {d}}t}=0,\;\forall \, k. \end{aligned}$$

Applying (36) in the corresponding steady states, the equilibrium mean and variance \(\bar{u}_{M,\infty }\), \(r_{k,\infty }^{M}\) satisfies

$$\begin{aligned} 0= & {} -{\text {d}}\bar{u}_{M,\infty }+\frac{1}{J}\sum _{k=-{J}/{2}+1}^{{J}/{2}}r_{k,\infty }^{M}{\varGamma }_{k}+F,\\ 0= & {} 2\left[ -{\varGamma }_{k}\bar{u}_{M,\infty }-d\right] r_{k,\infty }^{M}+Q_{F,kk,\infty }^{M},\quad k=0,1,\dots ,{J}/{2}, \end{aligned}$$

with nonlinear interaction term in steady state \(Q_{F,kk,\infty }^{M}=-2d_{M,k}\left( R_{\infty }\right) r_{k,\infty }^{M}+\sigma _{M,k}^{2}\left( R_{\infty }\right) \). A similar procedure like that in Sect. 3 by summing up the modes in the second equations above and substituting the mean equation can be carried out so that

$$\begin{aligned} \sum _{k}Q_{F,kk,\infty }^{M}=2\bar{u}_{M,\infty }\sum _{k}{\varGamma }_{k}r_{k,\infty }^{M}+2d\mathrm {tr}R_{\infty }^{M}=2J\bar{u}_{M,\infty }\left( d\bar{u}_{M,\infty }-F\right) +2d\mathrm {tr}R_{\infty }^{M}. \end{aligned}$$

The right-hand side of the above relation between the equilibrium mean and total variance should equal to zero for climate consistent models. Furthermore, in order to enforce consistent equilibrium statistics with the truth along each direction, that is, \(\bar{u}_{M,\infty }=\bar{u}_{\infty }\), \(r_{k,\infty }^{M}=r_{k,\infty }\), one necessary condition requires that

$$\begin{aligned} 2\left[ {\varGamma }_{k}\bar{u}_{M,\infty }+d\right] r_{k,\infty }^{M}=-2d_{M,k}\left( R_{\infty }\right) r_{k,\infty }^{M}+\sigma _{M,k}^{2}\left( R_{\infty }\right) , \end{aligned}$$

for each mode \(k=0,\ldots ,J/2\) from the steady-state variance equations. Substituting the corresponding forms of GC1 and GC2 in (41) and (42) into the above equation, the model parameters satisfying equilibrium fidelity with the truth must follow the relations

$$\begin{aligned} -2d_{M,k}r_{k,\infty }^{M}+\sigma _{M}^{2}=2\left[ {\varGamma }_{k}\bar{u}_{M,\infty }+d\right] r_{k,\infty }^{M}, \end{aligned}$$

for GC1, and

$$\begin{aligned} -\varepsilon _{1,k}J\left( \mathrm {tr}R_{\infty }\right) ^{-1}r_{k,\infty }^{M}+\varepsilon _{M}=2\left[ {\varGamma }_{k}\bar{u}_{M,\infty }+d\right] r_{k,\infty }^{M}, \end{aligned}$$

for GC2. By solving the above equations, we find the consistent parameters \(d_{M,k}\) and \(\varepsilon _{1,k}\) for each mode in GC1 and GC2 respectively as shown in (41) and (42).

Rights and permissions

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Majda, A.J., Qi, D. Improving Prediction Skill of Imperfect Turbulent Models Through Statistical Response and Information Theory. J Nonlinear Sci 26, 233–285 (2016). https://doi.org/10.1007/s00332-015-9274-5

Download citation

  • Received:

  • Accepted:

  • Published:

  • Issue Date:

  • DOI: https://doi.org/10.1007/s00332-015-9274-5

Keywords

Mathematics Subject Classification

Navigation