.. DO NOT EDIT.
.. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY.
.. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE:
.. "api/gallery/climate/plot_11_hasselmann_stochastic_climate.py"
.. LINE NUMBERS ARE GIVEN BELOW.

.. only:: html

    .. note::
        :class: sphx-glr-download-link-note

        :ref:`Go to the end <sphx_glr_download_api_gallery_climate_plot_11_hasselmann_stochastic_climate.py>`
        to download the full example code or to run this example in your browser via JupyterLite.

.. rst-class:: sphx-glr-example-title

.. _sphx_glr_api_gallery_climate_plot_11_hasselmann_stochastic_climate.py:


Hasselmann's stochastic climate (1976)
======================================

Klaus Hasselmann showed that slow climate variability needs no slow cause.
The ocean integrates fast, random weather forcing, just as a heavy particle
integrates molecular kicks in Brownian motion. The result is red noise: a
spectrum flat below the frequency ``lambda / (2 pi C)`` and falling as
frequency squared above it. Here a seeded white-noise forcing drives a 50 m
mixed layer, sampled every month for 4000 years.

.. GENERATED FROM PYTHON SOURCE LINES 12-44



.. image-sg:: /api/gallery/climate/images/sphx_glr_plot_11_hasselmann_stochastic_climate_001.png
   :alt: Weather noise reddened by ocean heat capacity
   :srcset: /api/gallery/climate/images/sphx_glr_plot_11_hasselmann_stochastic_climate_001.png
   :class: sphx-glr-single-img





.. code-block:: Python


    import matplotlib.pyplot as plt
    import numpy as np
    from scipy.signal import welch

    from earthsciencekit.climate import heat_capacity_column, red_noise_spectrum
    from earthsciencekit.units import convert

    restoring, capacity = 1.2, heat_capacity_column(1025, 3990, 50)
    forcing_spectrum = 4e7  # W^2 m^-4 Hz^-1, two-sided
    dt = convert(1 / 12, "year", "s")
    decay = np.exp(-restoring * dt / capacity)
    variance = forcing_spectrum / (2 * restoring * capacity)
    rng = np.random.default_rng(1976)
    kicks = rng.normal(0, np.sqrt(variance * (1 - decay**2)), 48000)
    temperature = np.zeros(kicks.size)
    for i in range(1, kicks.size):
        temperature[i] = decay * temperature[i - 1] + kicks[i]
    assert np.isclose(temperature.var(), variance, rtol=0.15)
    assert np.isclose(np.corrcoef(temperature[:-1], temperature[1:])[0, 1], decay, atol=0.01)

    frequency, density = welch(temperature, fs=1 / dt, nperseg=4096)
    per_year = convert(1, "year", "s")
    fig, ax = plt.subplots()
    ax.loglog(frequency[1:] * per_year, density[1:] / per_year, lw=0.8, label="simulated")
    theory = 2 * red_noise_spectrum(frequency[1:], forcing_spectrum, restoring, capacity)
    ax.loglog(frequency[1:] * per_year, theory / per_year, label="Hasselmann")
    ax.axvline(restoring / (2 * np.pi * capacity) * per_year, color="0.5", ls=":", label="corner frequency")
    ax.set(xlabel="frequency [cycles/year]", ylabel="one-sided spectrum [K² year]", title="Weather noise reddened by ocean heat capacity")
    ax.legend()

    plt.show()


.. rst-class:: sphx-glr-timing

   **Total running time of the script:** (0 minutes 0.634 seconds)


.. _sphx_glr_download_api_gallery_climate_plot_11_hasselmann_stochastic_climate.py:

.. only:: html

  .. container:: sphx-glr-footer sphx-glr-footer-example

    .. container:: lite-badge

      .. image:: images/jupyterlite_badge_logo.svg
        :target: ../../../lite/lab/index.html?path=api/gallery/climate/plot_11_hasselmann_stochastic_climate.ipynb
        :alt: Launch JupyterLite
        :width: 150 px

    .. container:: sphx-glr-download sphx-glr-download-jupyter

      :download:`Download Jupyter notebook: plot_11_hasselmann_stochastic_climate.ipynb <plot_11_hasselmann_stochastic_climate.ipynb>`

    .. container:: sphx-glr-download sphx-glr-download-python

      :download:`Download Python source code: plot_11_hasselmann_stochastic_climate.py <plot_11_hasselmann_stochastic_climate.py>`

    .. container:: sphx-glr-download sphx-glr-download-zip

      :download:`Download zipped: plot_11_hasselmann_stochastic_climate.zip <plot_11_hasselmann_stochastic_climate.zip>`


.. only:: html

 .. rst-class:: sphx-glr-signature

    `Gallery generated by Sphinx-Gallery <https://sphinx-gallery.github.io>`_