.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "introduction/lessons/plot_02_soliton.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_introduction_lessons_plot_02_soliton.py: GNLSE propagation: an N = 3 soliton =================================== A pulse in anomalous dispersion can balance dispersive spreading against nonlinear self-phase modulation. A fundamental soliton keeps its shape in the ideal cubic nonlinear Schrödinger equation. A higher-order soliton instead changes shape periodically. Here an **N = 3 soliton** compresses, develops several peaks, and returns to its initial shape after one soliton period. The GNLSE interface lets us specify dispersion and nonlinear strength directly. For this lesson we retain only second-order dispersion and the instantaneous Kerr response: no loss, Raman response, higher-order dispersion, or self-steepening. These choices give an ideal recurrence problem. They are useful for understanding the balance of effects, but are not a complete model of a particular fibre. The preceding :doc:`plot_01_rdw` lesson retained the broader physics needed for ultraviolet generation. Removing those perturbations here makes it easier to connect pulse parameters to dynamics. The full equation and available extensions are in :doc:`/guide/waveguides/gnlse`. .. GENERATED FROM PYTHON SOURCE LINES 28-43 Choose the pulse and its length scales -------------------------------------- The input amplitude is :math:`A(t)=\sqrt{P_0}\,\mathrm{sech}(t/\tau_0)`. Its intensity FWHM is :math:`2\operatorname{arcosh}(\sqrt{2})\tau_0`, about :math:`1.763\tau_0`; ``width`` below is :math:`\tau_0`, not that FWHM. Define the dispersion length :math:`L_D=\tau_0^2/|\beta_2|` and nonlinear length :math:`L_{NL}=1/(\gamma P_0)`. Their ratio fixes the soliton order, :math:`N^2=L_D/L_{NL}`. Choosing N = 3 therefore sets the peak power rather than leaving it as an independent parameter. The soliton period is :math:`z_0=\pi L_D/2`. The example coefficients set useful scales; changing them while holding N fixed rescales distance and power in this ideal equation. The 4 ps time window keeps the pulse tails far from the boundaries, and the wavelength band accommodates the much broader spectrum at compression. See :doc:`/guide/setup/grids`. .. GENERATED FROM PYTHON SOURCE LINES 43-68 .. code-block:: Python import time import matplotlib.pyplot as plt import numpy as np from solprop.driver import prop_gnlse from solprop.fourier import inverse from solprop.grid import EnvelopeGrid started = time.perf_counter() BETA2 = -11.83e-27 # s^2/m, anomalous GAMMA = 0.11 # 1/(W m) TAU0 = 100e-15 # s, the natural width of sech(t / tau0) ORDER = 3 dispersion_length = TAU0**2 / abs(BETA2) period = np.pi / 2 * dispersion_length peak_power = ORDER**2 * abs(BETA2) / (GAMMA * TAU0**2) grid = EnvelopeGrid.from_limits(period, 835e-9, (500e-9, 2e-6), 4e-12) print(f"dispersion length {dispersion_length:.3f} m, soliton period {period:.3f} m") print(f"N = {ORDER}: peak power {peak_power:.1f} W") .. rst-class:: sphx-glr-script-out .. code-block:: none dispersion length 0.845 m, soliton period 1.328 m N = 3: peak power 96.8 W .. GENERATED FROM PYTHON SOURCE LINES 69-75 Run one complete period ------------------------ The controller is asked for a local tolerance of :math:`10^{-9}` and the run stores two hundred planes over the period, plus the one at :math:`z = 0`. Without an ``output`` argument the planes are kept in memory; passing a path would stream them to an HDF5 file instead, and that file could be resumed if the run were interrupted. .. GENERATED FROM PYTHON SOURCE LINES 75-95 .. code-block:: Python output = prop_gnlse( GAMMA, period, [0.0, 0.0, BETA2], lambda0=835e-9, grid=grid, width=TAU0, peak_power=peak_power, pulseshape="sech", shock=False, rtol=1e-9, n_save=201, progress="none", ) z = output["z"] attempts = output["step_history/counts"].sum() print(f"{output.n_saved} planes saved, {attempts} steps taken") .. rst-class:: sphx-glr-script-out .. code-block:: none 201 planes saved, 1128 steps taken .. GENERATED FROM PYTHON SOURCE LINES 96-102 Check the field, including phase --------------------------------- What the handler holds is the stored spectrum at every plane, in the optics convention, and the statistics evaluated there. The recurrence error is the relative distance between the final and the initial spectrum once the expected global phase is removed; in the optics convention that phase is :math:`+\pi/4`. .. GENERATED FROM PYTHON SOURCE LINES 102-111 .. code-block:: Python stored = output["field"][:, 0, 0] recurrence = np.linalg.norm(stored[-1] - np.exp(1j * np.pi / 4) * stored[0]) recurrence /= np.linalg.norm(stored[0]) peak_bin = int(np.argmax(np.abs(stored[0]))) phase = float(np.angle(stored[-1][peak_bin] / stored[0][peak_bin])) print(f"recurrence error after one period: {recurrence:.2e}") print(f"global phase at the peak bin: {phase:+.5f} rad (expected {np.pi / 4:+.5f})") .. rst-class:: sphx-glr-script-out .. code-block:: none recurrence error after one period: 2.65e-09 global phase at the peak bin: +0.78540 rad (expected +0.78540) .. GENERATED FROM PYTHON SOURCE LINES 112-118 Follow compression and recurrence ----------------------------------- The evolution over one period. The pulse compresses to a fraction of its width, throws out side lobes, and comes back; the map is symmetric about the half period. The power comes from the same spectra in the internal convention, which is what an inverse transform expects. .. GENERATED FROM PYTHON SOURCE LINES 118-150 .. code-block:: Python power = np.abs(inverse(output.fields("fft")[:, 0, 0], kind="envelope")) ** 2 window = np.abs(grid.t) < 400e-15 fig, (left, right) = plt.subplots( 1, 2, figsize=(8.0, 3.4), width_ratios=(1.6, 1.0), layout="constrained" ) mesh = left.pcolormesh( grid.t[window] * 1e15, z / period, power[:, window] / peak_power, shading="nearest", rasterized=True, ) left.set_xlabel("time (fs)") left.set_ylabel("$z / z_0$") left.set_title(f"N = {ORDER} soliton over one period") fig.colorbar(mesh, ax=left, label="$P / P_0$") right.plot(grid.t[window] * 1e15, power[0, window] / peak_power, lw=2.4, label="z = 0") right.plot( grid.t[window] * 1e15, power[-1, window] / peak_power, lw=1.2, ls="--", label="$z = z_0$", ) right.set_xlabel("time (fs)") right.set_ylabel("$P / P_0$") right.set_title("recurrence") legend = right.legend(fontsize="small", frameon=False) .. image-sg:: /introduction/lessons/images/sphx_glr_plot_02_soliton_001.png :alt: N = 3 soliton over one period, recurrence :srcset: /introduction/lessons/images/sphx_glr_plot_02_soliton_001.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 151-159 Interpret the recorded statistics ---------------------------------- The statistics are recorded at every plane, so the compression is a line rather than a second pass over the fields: the peak power rises about sixfold a quarter of the way in and again three quarters of the way, with the split pulse between them, and the duration follows it down. A single FWHM becomes less informative when the pulse has several peaks, so read it beside the temporal map. Energy should remain constant: measure the largest departure over all saved planes, not just the exit. .. GENERATED FROM PYTHON SOURCE LINES 159-173 .. code-block:: Python fig, axes = plt.subplots(1, 2, figsize=(8.0, 3.0), layout="constrained", sharex=True) axes[0].plot(z / period, output["stats/peak_power"][:, 0] / peak_power) axes[0].set_ylabel("$P_\\mathrm{peak} / P_0$") axes[1].plot(z / period, output["stats/fwhm_t"] * 1e15) axes[1].set_ylabel("duration FWHM (fs)") for axis in axes: axis.set_xlabel("$z / z_0$") axes[0].set_title("compression") axes[1].set_title("duration") energy = output["stats/energy"][:, 0, 0] print(f"energy drift over the period: {np.max(np.abs(energy / energy[0] - 1)):.2e}") .. image-sg:: /introduction/lessons/images/sphx_glr_plot_02_soliton_002.png :alt: compression, duration :srcset: /introduction/lessons/images/sphx_glr_plot_02_soliton_002.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-script-out .. code-block:: none energy drift over the period: 5.34e-10 .. GENERATED FROM PYTHON SOURCE LINES 174-188 Assumptions and next steps -------------------------- Recurrence tests the complete field rather than only its energy. Small energy drift is necessary but cannot by itself prove that the pulse shape is accurate. The propagation tolerance controls longitudinal integration, while time-window and spectral-band refinement test different numerical errors. See :doc:`/guide/numerics/index` for how to vary these separately. In real fibres, higher-order dispersion, Raman scattering, loss and self-steepening can break this recurrence. Compare the worked examples of :doc:`/auto_examples/processes/plot_02_soliton_self_compression` and :doc:`/auto_examples/processes/plot_05_raman_self_frequency_shift` before adding those effects to a model. Next, :doc:`plot_03_filament` leaves the guide for a beam that focuses and ionises in free space. .. GENERATED FROM PYTHON SOURCE LINES 188-190 .. code-block:: Python print(f"{time.perf_counter() - started:.1f} s; backend {output.metadata.backend.name}") .. rst-class:: sphx-glr-script-out .. code-block:: none 4.2 s; backend numba_cpu .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 4.172 seconds) .. _sphx_glr_download_introduction_lessons_plot_02_soliton.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_02_soliton.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_02_soliton.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_02_soliton.zip `