Kubo Oscillator
GeometricProblems.KuboOscillator — Module
Kubo Oscillator
The Kubo oscillator is a unit-frequency harmonic oscillator, $\dot{q}_1 = q_2$, $\dot{q}_2 = -q_1$, driven by multiplicative (Stratonovich) noise. It is a standard test problem for stochastic geometric integrators. The module provides it as a stochastic differential equation (sdeproblem), a partitioned SDE (psdeproblem), a split partitioned SDE (spsdeproblem), and the underlying deterministic ODE (odeproblem). The *problem builders create a single problem from one initial condition; the *ensemble builders (sdeensemble, psdeensemble, spsdeensemble) build an ensemble from several initial conditions. The noise is a one-dimensional WienerProcess.
Because the diffusion is proportional to the drift, the solution is the deterministic one evaluated at the random time $\theta(t) = t + \nu W(t)$, so the energy is conserved exactly along every sample path. That is what makes this problem useful for checking a stochastic geometric integrator: any energy drift is the scheme's, not the problem's. exact_solution gives the closed form.
The damped_* builders add linear damping, following Kraus & Tyranowski, Variational integrators for stochastic dissipative Hamiltonian systems. There the system is $H = (p^2+q^2)/2$, $h = \nu (p^2+q^2)/2$ with forcing $F(q,p) = -\gamma p$ and $f(q,p) = -\nu \gamma p$. The energy then decays rather than being conserved, and exact_mean_energy gives $E(H)$ in closed form. Unlike the undamped problem, the damped one has genuinely non-zero $f_2$ and $G_2$, so it exercises the split half of the SPSDE.
System parameters: ν — the noise intensity; γ — the damping coefficient (0 for the undamped problems). Note the paper writes these as $\beta$ and $\nu$ respectively; the names here follow this module's existing convention, in which ν has always been the noise intensity.
The deterministic drift alone (no noise) is a plain harmonic oscillator, tracing a circle in phase space:
using GeometricProblems.KuboOscillator
using GeometricIntegrators: integrate, ImplicitMidpoint
using CairoMakie
sol = integrate(odeproblem(; timespan = (0.0, 2π)), ImplicitMidpoint())
fig = Figure()
ax = Axis(fig[1, 1]; xlabel = "q₁", ylabel = "q₂", aspect = DataAspect(), title = "Kubo oscillator (deterministic drift)")
lines!(ax, sol.q[:, 1], sol.q[:, 2])
fig
Library functions
GeometricProblems.KuboOscillator.damped_psdeproblem — Function
Damped Kubo oscillator as a partitioned SDE, with the damping folded into f and G.
Pairs with damped_spsdeproblem, which splits the same dynamics; the two must give the same trajectory on a common sample path.
GeometricProblems.KuboOscillator.damped_spsdeproblem — Function
Damped Kubo oscillator as a split partitioned SDE, with the Hamiltonian part in f1/G1 and the damping in f2/G2.
GeometricProblems.KuboOscillator.exact_mean_energy — Method
exact_mean_energy(t, q₀, p₀, t₀, params)Expected value of the Hamiltonian $H = (p^2 + q^2)/2$ of the damped Kubo oscillator at time t, in closed form (Kraus & Tyranowski §4.1), written in the elapsed time $s = t - t_0$:
\[E(H) = a \, e^{-\frac{\gamma (2 - \nu^2 \gamma)}{2} s} + e^{-((2 - \gamma^2)\nu^2 + \gamma) s} \Big[ b \cos \big( 2 (1 - \nu^2 \gamma) \omega s \big) + c \sin \big( 2 (1 - \nu^2 \gamma) \omega s \big) \Big],\]
with $\omega = \tfrac{1}{2}\sqrt{4 - \gamma^2}$ and
\[a = \frac{2 (p_0^2 + q_0^2 + \gamma p_0 q_0)}{4 - \gamma^2}, \qquad b = -\frac{\gamma^2 (p_0^2 + q_0^2) + 4 \gamma p_0 q_0}{2 (4 - \gamma^2)}, \qquad c = \frac{\gamma (q_0^2 - p_0^2)}{2 \sqrt{4 - \gamma^2}} .\]
This averages over the noise, so unlike exact_solution it needs no sample path. As there, $|\gamma| < 2$. Without damping it collapses to the constant $(p_0^2 + q_0^2)/2$, which is the exact conservation the undamped problems exhibit pathwise.
GeometricProblems.KuboOscillator.exact_solution — Method
exact_solution_q(t, W, q₀, p₀, t₀, params)
exact_solution_p(t, W, q₀, p₀, t₀, params)
exact_solution(t, W, q₀, p₀, t₀, params)
exact_solution(t, W, x₀, t₀, params)Exact solution of the (possibly damped) Kubo oscillator at time t, along the sample path whose Wiener increment $W(t) - W(t_0)$ is W.
Because the diffusion is proportional to the drift, the solution is the deterministic damped oscillator evaluated at the random time $\theta = (t - t_0) + \nu (W(t) - W(t_0))$:
\[\begin{aligned} q(t) &= e^{-\gamma \theta / 2} \left[ q_0 \cos \omega \theta + \tfrac{1}{\omega} \left( p_0 + \tfrac{\gamma}{2} q_0 \right) \sin \omega \theta \right], \\ p(t) &= e^{-\gamma \theta / 2} \left[ p_0 \cos \omega \theta - \tfrac{1}{\omega} \left( q_0 + \tfrac{\gamma}{2} p_0 \right) \sin \omega \theta \right], \end{aligned}\]
with $\omega = \tfrac{1}{2} \sqrt{4 - \gamma^2}$. params needs a noise intensity ν; a damping coefficient γ is taken as zero when absent, which reduces the above to a rotation by $\theta$. This is the underdamped solution, so $|\gamma| < 2$; anything else is a DomainError.
q₀ and p₀ are numbers or one-element vectors, and the result matches: a tuple of numbers or a tuple of vectors. The form taking a single state vector x₀ serves the sdeproblem formulation, where $q = (q_1, q_2)$ plays the role of $(q, p)$, and returns [q, p].
GeometricProblems.KuboOscillator.exact_solution_p — Method
exact_solution_p(t, W, q₀, p₀, t₀, params)Momentum component of exact_solution.
GeometricProblems.KuboOscillator.exact_solution_q — Method
exact_solution_q(t, W, q₀, p₀, t₀, params)Position component of exact_solution.