QuTiP 4— Simulating Quantum Sensors with QuTiP

Connecting phase, decoherence, and measurement in a simple two-level model

Infographic by NotebookLM

In the previous articles, we built up the basic ingredients of quantum dynamics in QuTiP. In Article 1, we introduced QuTiP as a practical Python toolbox for simulating quantum systems. In Article 2, we used a two-level system to study Rabi oscillations, the simplest picture of coherent quantum control. In Article 3, we made the model more realistic by adding relaxation and dephasing, showing how open quantum systems lose coherence over time.

In this article, we use those same ingredients for a new purpose: quantum sensing. A quantum sensor is a quantum system whose state changes in response to an external physical quantity. That quantity might be a magnetic field, an electric field, acceleration, gravity, temperature, time, or frequency. The central idea is simple: prepare a quantum state, let it interact with the quantity we want to estimate, and then measure how the state has changed.

Here, we will not try to model a complete experimental device. Instead, we will build the smallest useful picture: a spin-1/2 quantum system that accumulates phase in response to a weak external field. Then we will add dephasing and see how noise reduces the sensor signal.

A companion Google Colab notebook accompanies this article and contains every code example below in runnable form. You can open it here: https://colab.research.google.com/drive/1FQBD6w_h-uZbjwTVy2kjHtsdzxBtQ_BJ. Running the cells in order will reproduce each figure discussed in the text, so you are welcome to follow along, modify the parameters, and experiment with the model as you read.

1. The basic idea of a quantum sensor

A useful way to think about a quantum sensor is as a short sequence of steps: prepare a quantum state, expose it to the quantity we want to measure, allow the state to accumulate phase, measure the final state, and finally infer the unknown quantity from the measurement result.

For a two-level quantum system, the state can be visualised on the Bloch sphere. In Article 2, we saw that coherent dynamics move the state around the Bloch sphere, and in Article 3 we saw that decoherence shrinks or distorts that motion. For sensing, both effects matter. The signal comes from coherent phase accumulation, while the limitation often comes from noise, especially dephasing. A good sensor needs the state to remain coherent long enough for the unknown signal to leave a detectable imprint.

In this article, the unknown signal will be represented by a frequency shift. Physically, you can think of this as a spin experiencing a weak magnetic field. Mathematically, the Hamiltonian is H = (ω/2) σz, where σz is the Pauli-z operator and ω is the angular frequency caused by the external field. A larger field produces a faster phase accumulation. We will use arbitrary units so that the code stays focused on the concept rather than on the calibration of a particular sensor platform.

2. From Rabi oscillations to phase accumulation

In Article 2, we used a Hamiltonian proportional to σx to drive transitions between two states, which gave us Rabi oscillations. For sensing, we use a different but closely related idea: we prepare the system in a superposition and then let the relative phase between the two basis states evolve.

We return to the notation |0⟩ and |1⟩ here because this model does not depend on decay from an excited state to a ground state. The important point is the balanced superposition:

|+x⟩ = (|0⟩ + |1⟩)/√2.

This state points along the x-axis of the Bloch sphere. Under the Hamiltonian H = (ω/2) σz, it rotates around the z-axis, and the larger ω is, the faster the rotation. If ω depends on an external field, then observing the final state gives us information about that field.

Before plotting anything, we can already make a simple prediction. If ω = 1, then ⟨σx⟩ = cos(ωt) has period 2π/ω = 2π ≈ 6.3, so a time window from t = 0 to t = 10 should show about one and a half oscillations. That is the first sanity check for the simulation.

3. Installing and importing the libraries

As before, this article assumes you are working in Google Colab or a local Python environment.

!pip -q install qutip plotly

Now import the packages.

import numpy as np
import plotly.graph_objects as go

from qutip import basis, sigmax, sigmay, sigmaz, mesolve

We will continue to use Plotly for interactive visualisation.

4. Define the two-level sensor

First, we define the Pauli operators and the initial state.

# Pauli operators
sx = sigmax()
sy = sigmay()
sz = sigmaz()

# Basis states
zero = basis(2, 0)
one = basis(2, 1)

# Initial state: |+x> = (|0> + |1>) / sqrt(2)
psi_plus_x = (zero + one).unit()

# Time points
tlist = np.linspace(0, 10, 500)

The initial state is a balanced superposition. This matters because phase accumulation is visible through the relative phase between |0⟩ and |1⟩. If the system started in |0⟩ or |1⟩, the Hamiltonian H = (ω/2) σz would only add a phase to that basis state, and that would not create the same oscillating measurement signal in the x-basis. In short, a useful phase sensor needs coherence.

5. A single simulation function

Now we define one simulation function and reuse it throughout the article. This avoids having separate pieces of code that quietly implement the same physics in slightly different ways.

def simulate_sensor(time_points, initial_state, omega, gamma_phi=0.0, observables=None):
    """
    Simulate a two-level quantum sensor.
    Parameters
    ----------
    time_points : array-like
        Times at which the state is evaluated.
    initial_state : qutip.Qobj
        Initial quantum state.
    omega : float
        Signal-dependent angular frequency.
    gamma_phi : float
        Pure dephasing rate.
    observables : list of qutip.Qobj
        Operators whose expectation values are returned.
    Returns
    -------
    result : qutip.solver.Result
        QuTiP solver result containing expectation values.
    """
    if observables is None:
        observables = [sx]
    H = 0.5 * omega * sz
    c_ops = []
    if gamma_phi > 0:
        c_ops = [np.sqrt(gamma_phi / 2) * sz]
    result = mesolve(
        H,
        initial_state,
        time_points,
        c_ops=c_ops,
        e_ops=observables
    )
    return result

The Hamiltonian contains the signal, the collapse operators contain the noise, and the observables define what we read out from the simulation. This is the modelling pattern we will reuse for the rest of the article.

6. Simulate ideal phase accumulation

Let us compare three different signal strengths.

omegas = [0.6, 1.0, 1.4]

fig = go.Figure()
for omega in omegas:
    result = simulate_sensor(
        time_points=tlist,
        initial_state=psi_plus_x,
        omega=omega,
        gamma_phi=0.0,
        observables=[sx]
    )
    sx_expect = result.expect[0]
    fig.add_trace(
        go.Scatter(
            x=tlist,
            y=sx_expect,
            mode="lines",
            name=f"omega = {omega}"
        )
    )

fig.update_layout(
    title="Ideal quantum sensor: phase accumulation for different signal strengths",
    xaxis_title="Time",
    yaxis_title="⟨σx⟩",
    template="plotly_white"
)

fig.show()

The expectation value ⟨σx⟩ oscillates because the state rotates around the z-axis of the Bloch sphere, and a larger value of ω produces faster oscillations. For ω = 1, the curve completes about one and a half oscillations across the plotted range. That matches the period estimate we made before running the code: 2π ≈ 6.3, so the interval from t = 0 to t = 10 contains roughly 10 / 6.3 ≈ 1.6 cycles.

This is the simplest version of quantum sensing. The external field changes the phase accumulation rate, and the final quantum state carries information about that field.

7. Turn the quantum state into a measurement signal

In a real experiment, we do not directly see the quantum state; we measure it. For this simple model, imagine measuring in the x-basis. The probability of finding the system in the |+x⟩ state is related to ⟨σx⟩ by

P(+x) = (1 + ⟨σx⟩) / 2.

Let us calculate this probability over time for one value of ω.

omega = 1.0

result = simulate_sensor(
    time_points=tlist,
    initial_state=psi_plus_x,
    omega=omega,
    gamma_phi=0.0,
    observables=[sx]
)

sx_expect = result.expect[0]

p_plus_x = 0.5 * (1 + sx_expect)

fig = go.Figure()

fig.add_trace(
    go.Scatter(
        x=tlist,
        y=p_plus_x,
        mode="lines",
        name="P(+x)"
    )
)

fig.update_layout(
    title="Measurement probability in the x-basis",
    xaxis_title="Time",
    yaxis_title="Probability",
    template="plotly_white"
)

fig.show()

The probability oscillates between high and low values. This oscillating response is often called a fringe. The position of the fringe depends on the accumulated phase, and the accumulated phase depends on the signal. In a real sensing experiment, we would repeat the measurement many times, estimate this probability statistically, and infer the unknown field from it.

8. Add dephasing

The ideal model is useful, but it is incomplete. In Article 3, we saw that dephasing destroys the phase coherence between two basis states. For a phase-based sensor this is especially important, because the signal itself is stored in the relative phase.

Get Ken Moriwaki’s stories in your inbox

Join Medium for free to get updates from this writer.

Remember me for faster sign in

We can add pure dephasing using a collapse operator proportional to σz:

Cφ = √(γφ / 2) σz.

With this convention, the coherence decays exactly at rate γφ in this pure-dephasing model. The ⟨σx⟩ oscillation is therefore multiplied by the envelope exp(−γφt).

Before plotting, we can predict the general behaviour: the oscillation frequency will still be determined by ω, but the amplitude will decrease exponentially with time. Larger values of γφ should therefore produce a faster loss of contrast.

omega = 1.0
gamma_values = [0.0, 0.05, 0.15, 0.30]

fig = go.Figure()

for gamma_phi in gamma_values:
    result = simulate_sensor(
        time_points=tlist,
        initial_state=psi_plus_x,
        omega=omega,
        gamma_phi=gamma_phi,
        observables=[sx]
    )
    sx_expect = result.expect[0]
    fig.add_trace(
        go.Scatter(
            x=tlist,
            y=sx_expect,
            mode="lines",
            name=f"gamma_phi = {gamma_phi}"
        )
    )

fig.update_layout(
    title="Effect of dephasing on the sensor signal",
    xaxis_title="Time",
    yaxis_title="⟨σx⟩",
    template="plotly_white",
    legend=dict(
        orientation="h",
        yanchor="bottom", y=-0.3,
        xanchor="center", x=0.5,
        title_text="",
    ),
    margin=dict(b=110),
)

fig.show()

The oscillation frequency is still determined by ω, but the contrast decays as dephasing increases. This is a central trade-off in phase-based sensing: a longer sensing time allows more phase to accumulate, but it also gives noise more time to destroy the coherence that carries the signal.

9. Visualise the shrinking coherence

The previous plot shows the signal as a function of time, but there is another useful way to see the same physics. Under H = (ω/2) σz, the Bloch vector rotates in the x-y plane. With pure dephasing, it still rotates, but the radius of that rotation shrinks, so the state traces an inward spiral.

This picture is useful because it shows what dephasing is doing geometrically. The signal is not just “getting smaller” in an abstract sense; the coherence represented by the transverse Bloch-vector components is being reduced.

omega = 1.0
gamma_values = [0.0, 0.15, 0.30]

fig = go.Figure()

for gamma_phi in gamma_values:
    result = simulate_sensor(
        time_points=tlist,
        initial_state=psi_plus_x,
        omega=omega,
        gamma_phi=gamma_phi,
        observables=[sx, sy]
    )
    sx_expect = result.expect[0]
    sy_expect = result.expect[1]
    fig.add_trace(
        go.Scatter(
            x=sx_expect,
            y=sy_expect,
            mode="lines",
            name=f"gamma_phi = {gamma_phi}"
        )
    )

fig.update_layout(
    title="Bloch-vector spiral in the x-y plane",
    xaxis_title="⟨σx⟩",
    yaxis_title="⟨σy⟩",
    template="plotly_white"
)

fig.update_yaxes(
    scaleanchor="x",
    scaleratio=1
)

fig.show()

The ideal curve stays on the outside of the circle, while with dephasing the trajectory spirals inward toward the centre. This is one of the clearest pictures of the sensing problem: the external field is encoded in the phase angle, but dephasing reduces the radius of the vector that carries that phase information.

10. Compare ideal and noisy sensor response curves

Now we can look directly at the readout probability after a fixed sensing time. Suppose we prepare the state, let it evolve for an interrogation time T, and then measure in the x-basis. We repeat this calculation for different possible values of ω, which gives us a sensor response curve.

For the response curves below, we use an interrogation time of T = 4. Before running the scan, we can predict the reduction in contrast caused by dephasing. With γφ = 0.15,

exp(−γφT) = exp(−0.15 × 4) = exp(−0.6) ≈ 0.55.

The dephased response should therefore retain about 55% of the ideal fringe contrast.

We can reuse the same simulation function and simply return the final readout probability.

def readout_probability(omega, gamma_phi, interrogation_time, initial_state):
    """
    Return P(+x) after a fixed interrogation time.
    """
    local_tlist = np.linspace(0, interrogation_time, 300)

    result = simulate_sensor(
        time_points=local_tlist,
        initial_state=initial_state,
        omega=omega,
        gamma_phi=gamma_phi,
        observables=[sx]
    )
    sx_final = result.expect[0][-1]
    p_plus_x = 0.5 * (1 + sx_final)
    return p_plus_x

Now scan over possible signal strengths.

interrogation_time = 4.0
omega_scan = np.linspace(0, 2.5, 120)

p_ideal = []
p_dephased = []

for omega in omega_scan:
    p_ideal.append(
        readout_probability(
            omega=omega,
            gamma_phi=0.0,
            interrogation_time=interrogation_time,
            initial_state=psi_plus_x
        )
    )
    p_dephased.append(
        readout_probability(
            omega=omega,
            gamma_phi=0.15,
            interrogation_time=interrogation_time,
            initial_state=psi_plus_x
        )
    )

fig = go.Figure()

fig.add_trace(
    go.Scatter(
        x=omega_scan,
        y=p_ideal,
        mode="lines",
        name="Ideal"
    )
)

fig.add_trace(
    go.Scatter(
        x=omega_scan,
        y=p_dephased,
        mode="lines",
        name="With dephasing"
    )
)

fig.update_layout(
    title=f"Sensor response after interrogation time T = {interrogation_time}",
    xaxis_title="Signal-dependent frequency omega",
    yaxis_title="P(+x)",
    template="plotly_white"
)

fig.show()

As predicted, the ideal curve has high contrast, while the dephased curve is flatter. For γφ = 0.15 and T = 4, the noisy fringe retains approximately 55% of the ideal contrast.

A flatter curve is less useful for sensing because different values of ω become harder to distinguish. This is where the open-system physics from Article 3 becomes directly relevant to quantum technology. Decoherence is not just a loss of mathematical purity; in a sensor, it reduces the information available in the measurement.

11. Reading sensitivity from the slope

A sensor is useful when a small change in the physical quantity produces a measurable change in the output. In this example, the output is P(+x) and the unknown parameter is ω. The steepest part of the response curve is usually the most sensitive operating region, because around that point a small change in ω produces a relatively large change in measurement probability.

We can see this numerically by estimating the slope of the response curve.

p_ideal = np.array(p_ideal)
p_dephased = np.array(p_dephased)

slope_ideal = np.gradient(p_ideal, omega_scan)
slope_dephased = np.gradient(p_dephased, omega_scan)

fig = go.Figure()

fig.add_trace(
    go.Scatter(
        x=omega_scan,
        y=slope_ideal,
        mode="lines",
        name="Ideal slope"
    )
)

fig.add_trace(
    go.Scatter(
        x=omega_scan,
        y=slope_dephased,
        mode="lines",
        name="Dephased slope"
    )
)

fig.update_layout(
    title="Sensitivity proxy: slope of the sensor response",
    xaxis_title="Signal-dependent frequency omega",
    yaxis_title="dP(+x) / domega",
    template="plotly_white"
)

fig.show()

This slope is not yet a full precision analysis. A proper quantum metrology treatment would also consider projection noise, technical noise, repeated measurements, Fisher information, and experimental constraints. But as a first simulation it captures the main intuition: the sensor response is most informative where the measurement probability changes rapidly with the parameter we want to estimate, and dephasing reduces that slope by reducing the fringe contrast.

12. What this simple model leaves out

This article uses a highly simplified spin-1/2 sensor. It is useful for learning, but it is not a complete device model. A real quantum sensor may require a calibrated relationship between the physical field and ω, realistic control pulses, finite preparation and measurement fidelity, relaxation as well as dephasing, technical noise, ensemble averaging, sensor readout noise, environmental drift, and a statistical estimation procedure.

Different platforms would also require different modelling choices. An NV-centre magnetometer, trapped-ion sensor, superconducting flux qubit, atomic vapour sensor, and atomic clock all use quantum coherence, but they do not all share the same Hamiltonian, noise channels, or readout model. The point of the minimal model is not to replace those details. It is to make the core mechanism visible before the platform-specific physics is added.

13. Where QuTiP fits in

This example shows why QuTiP is useful for quantum sensing. With the same workflow, we can change the Hamiltonian, add noise channels, choose observables, and compare ideal and noisy dynamics. In this article, mesolve was enough because we modelled ensemble-averaged dynamics with a Hamiltonian and a pure-dephasing channel. More advanced sensing models can include time-dependent control pulses, multi-level structure, relaxation, and continuous-measurement effects. Those are natural extensions, but they are not needed for the first sensing model.

The important point is that QuTiP lets us move from a physical sensing idea to an executable model without changing frameworks. Once the state, dynamics, noise, and readout are specified, the same structure can be reused across many quantum-technology examples.

14. Summary

In this article, we connected the ideas from the previous QuTiP articles to quantum sensing. We started with a two-level system prepared in a balanced superposition, allowed it to accumulate phase under a Hamiltonian proportional to σz, and interpreted that phase evolution as the response to an external field.

We then added dephasing through a collapse operator and saw how noise reduces the contrast of the measurement signal. The time-domain plot showed the decaying oscillation, the x-y plot showed the shrinking Bloch-vector spiral, and the response-curve slope showed why dephasing reduces sensitivity.

The main lesson is that a quantum sensor is not mysterious at the simulation level. It is a quantum system whose state changes in response to a parameter, and QuTiP helps us model how that change appears in observables. In the next article, we can build on this by connecting the same ideas to atomic clocks, Ramsey interrogation, pulse design, and coherent control.

References and further reading

• Nation, P.D. (2026). QuTiP — Quantum Toolbox in Python. [online] Qutip.org. Available at: https://qutip.org [Accessed 14 Jun. 2026].
• ‌QuTiP developers and contributors (2026). QuTiP: Quantum Toolbox in Python — QuTiP 5.3 Documentation. [online] Available at: https://qutip.readthedocs.io [Accessed 14 Jun. 2026].
• Lambert, N., Giguère, E., Menczel, P., Li, B., Hopf, P., Suárez, G., Gali, M., Lishman, J., Rushiraj Gadhvi, Agarwal, R., Galicia, A., Shammah, N., Nation, P., Johansson, J.R., Ahmed, S., Cross, S., Pitchford, A. and Nori, F. (2025). QuTiP 5: The Quantum Toolbox in Python. Physics Reports, 1153, pp.1–62. doi:https://doi.org/10.1016/j.physrep.2025.10.001.
• Johansson, J.R., Nation, P.D. and Nori, F. (2012). QuTiP: An open-source Python framework for the dynamics of open quantum systems. Computer Physics Communications, 183(8), pp.1760–1772. doi:https://doi.org/10.1016/j.cpc.2012.02.021.
• Johansson, J.R., Nation, P.D. and Nori, F. (2013). QuTiP 2: A Python framework for the dynamics of open quantum systems. Computer Physics Communications, 184(4), pp.1234–1240. doi:https://doi.org/10.1016/j.cpc.2012.11.019.
• Degen, C. L., Reinhard, F. and Cappellaro, P. (2017). Quantum sensing. Reviews of Modern Physics, 89(3). doi:https://doi.org/10.1103/revmodphys.89.035002.
• Acín, A., Bloch, I., Buhrman, H., Calarco, T., Eichler, C., Eisert, J., Esteve, D., Gisin, N., Glaser, S.J., Jelezko, F., Kuhr, S., Lewenstein, M., Riedel, M.F., Schmidt, P.O., Thew, R., Wallraff, A., Walmsley, I. and Wilhelm, F.K. (2018). The quantum technologies roadmap: a European community view. New Journal of Physics, 20(8), p.080201. doi:https://doi.org/10.1088/1367-2630/aad1ea.