General Example

Contents

General Example#

This is a generic example demonstrating basic usage patterns of the latex_parser library. It serves as a minimal end-to-end demo.

Source#

  1#!/usr/bin/env python
  2# coding: utf-8
  3# flake8: noqa
  4
  5import sys
  6from pathlib import Path
  7
  8import matplotlib.pyplot as plt
  9import numpy as np
 10from qutip import (
 11    about,
 12    basis,
 13    destroy,
 14    mesolve,
 15    qeye,
 16    tensor,
 17)
 18
 19ROOT = Path(__file__).resolve().parents[1]
 20if str(ROOT) not in sys.path:
 21    sys.path.insert(0, str(ROOT))
 22
 23from latex_parser.latex_api import compile_model as compile_latex_model
 24
 25
 26def run_vacuum_rabi():
 27    """
 28    Reproduce the QuTiP vacuum Rabi oscillations tutorial, but build
 29    H and c_ops from your LaTeX DSL instead of manual Qobj expressions.
 30    """
 31    # === Problem parameters (same physics as the QuTiP example) ===
 32    N = 15  # number of cavity Fock states
 33    omega_c = 1.0 * 2 * np.pi  # \omega_c: cavity frequency
 34    omega_a = 1.0 * 2 * np.pi  # \omega_a: atom frequency
 35    g = 0.05 * 2 * np.pi  # g: coupling strength
 36    kappa = 0.005  # \kappa: cavity dissipation rate
 37    gamma = 0.05  # \gamma: atom dissipation rate
 38    n_th_a = 0.0  # thermal photon number (environment)
 39    use_rwa = True
 40
 41    tlist = np.linspace(0, 40, 100)
 42
 43    # === Build model at T = 0 using the simple API ===
 44    if use_rwa:
 45        H_latex = r"""
 46            \omega_c \hat{n}_{1}
 47            + \frac{\omega_a}{2} \sigma_{z,1}
 48            + g \left( a_{1}^{\dagger} \sigma_{-,1}
 49                     + a_{1} \sigma_{+,1} \right)
 50        """
 51    else:
 52        H_latex = r"""
 53            \omega_c \hat{n}_{1}
 54            + \frac{\omega_a}{2} \sigma_{z,1}
 55            + g \left( a_{1}^{\dagger} + a_{1} \right)
 56                \left( \sigma_{-,1} + \sigma_{+,1} \right)
 57        """
 58
 59    c_ops_latex = [
 60        r"\sqrt{\kappa (1 + n_{th})} \, a_{1}",
 61        r"\sqrt{\kappa n_{th}} \, a_{1}^{\dagger}",
 62        r"\sqrt{\gamma} \, \sigma_{-,1}",
 63    ]
 64
 65    params = {
 66        "omega_c": omega_c,
 67        "omega_a": omega_a,
 68        "g": g,
 69        "kappa": kappa,
 70        "gamma": gamma,
 71        "n_th": n_th_a,
 72    }
 73
 74    compiled_0 = compile_latex_model(
 75        H_latex,
 76        params,
 77        c_ops_latex=c_ops_latex,
 78        qubits=1,
 79        bosons=[N],
 80        t_name="t",
 81    )
 82
 83    H_0 = compiled_0.H  # Qobj (static)
 84    c_ops_0 = compiled_0.c_ops
 85    args_0 = compiled_0.args  # currently just the params dict
 86
 87    # === Initial state ===
 88    # Hilbert ordering: [qubit, boson] -> 2 ⊗ N
 89    # "Atom excited, cavity ground" -> |e> ⊗ |0>.
 90    # QuTiP tutorial uses tensor(basis(N, 0), basis(2, 0)) with different ordering;
 91    # we swap to match [qubit, boson].
 92    psi0 = tensor(basis(2, 0), basis(N, 0))
 93
 94    # === Operators for expectation values (cavity & atom excitation) ===
 95    # Cavity annihilation in our ordering: I_qubit ⊗ a
 96    a = tensor(qeye(2), destroy(N))
 97    n_cav = a.dag() * a
 98
 99    # Atomic excitation projector (σ_+ σ_- in our ordering)
100    sm_raising = tensor(destroy(2).dag(), qeye(N))  # σ_+
101    atom_excited_proj = sm_raising.dag() * sm_raising  # σ_- σ_+
102
103    e_ops = [n_cav, atom_excited_proj]
104
105    # === Evolve at T = 0 ===
106    output_0 = mesolve(H_0, psi0, tlist, c_ops_0, e_ops, args=args_0)
107
108    # === Now finite temperature: n_th_a > 0 ===
109    n_th_a = 2.0
110
111    params_T = dict(params)
112    params_T["n_th"] = n_th_a
113    compiled_T = compile_latex_model(
114        H_latex,
115        params_T,
116        c_ops_latex=c_ops_latex,
117        qubits=1,
118        bosons=[N],
119        t_name="t",
120    )
121
122    H_T = compiled_T.H
123    c_ops_T = compiled_T.c_ops
124    args_T = compiled_T.args
125
126    output_T = mesolve(H_T, psi0, tlist, c_ops_T, e_ops, args=args_T)
127
128    # === Plot T = 0 ===
129    fig, ax = plt.subplots(figsize=(8, 5))
130    ax.plot(tlist, output_0.expect[0], label="Cavity")
131    ax.plot(tlist, output_0.expect[1], label="Atom excited state")
132    ax.legend()
133    ax.set_xlabel("Time")
134    ax.set_ylabel("Occupation probability")
135    ax.set_title("Vacuum Rabi oscillations at T = 0")
136    plt.tight_layout()
137    plt.show()
138
139    # === Plot T > 0 ===
140    fig, ax = plt.subplots(figsize=(8, 5))
141    ax.plot(tlist, output_T.expect[0], label="Cavity")
142    ax.plot(tlist, output_T.expect[1], label="Atom excited state")
143    ax.legend()
144    ax.set_xlabel("Time")
145    ax.set_ylabel("Occupation probability")
146    ax.set_title(f"Vacuum Rabi oscillations at T = {n_th_a}")
147    plt.tight_layout()
148    plt.show()
149
150    # === Optional: analytic test in RWA, no collapse operators ===
151    if use_rwa:
152        output_no_cops = mesolve(H_0, psi0, tlist, [], e_ops, args=args_0)
153        freq = 0.25 * np.sqrt(g**2 * (N + 1))
154        atom_pop = np.array(output_no_cops.expect[1])
155        analytic = (np.cos(tlist * freq)) ** 2
156        assert np.allclose(atom_pop, analytic, atol=1e-3)
157
158    print()
159    print("QuTiP / environment info:")
160    about()
161
162
163if __name__ == "__main__":
164    run_vacuum_rabi()

Run#

python examples/example.py

Notes#

  • This example is a good starting point for understanding the overall workflow.

Contents