A JPA against WRspice at the signal frequency
This recipe compares the JPA's small-signal reflection gain from harmonic balance, a tangent calculation about a transient pump trajectory, and WRspice. It requires JosephsonCircuits, Plots, and XicTools_jll (or a local WRspice installation). WRspice comes with the XicTools_jll package on x86_64 Linux; elsewhere install it and use JosephsonCircuits.wrspice_cmd() for the executable.
1 2
o------[cc]--------o--------+
| | |
[p1] [jj] [cj]
| | |
o------------------o--------+
0The pump alone is solved once, for 200 ns at eighty steps per pump period with the 1 - sech(t/10 ns) rise of the WRspice input. The signals are then directions of the tangent about that recorded trajectory: a unit current at each signal frequency, its cosine and its sine, with the same rise, all on one pass. The reflection at a signal frequency comes from the last pump period exactly as JosephsonCircuits.wrspice_calcS_paramp reads it: the cosine and sine responses combined as cos + i sin carry the signal at its own frequency and the idler at the negative of its own, so demodulated at the signal frequency the idler is periodic in the pump and averages to nothing over one period, and a Norton current I into the 50 ohm port gives V = 25 I (1 + S11).
using JosephsonCircuits, XicTools_jll, Plots
circuit = Circuit([(:p1, 1, 0, Port(1; Z0 = 50.0)), (:cc, 1, 2, Capacitor(100e-15)),
(:jj, 2, 0, JosephsonJunction(1000e-12)), (:cj, 2, 0, Capacitor(1000e-15))])
fp, Ip = 4.75001e9, 0.00565e-6
fs = (4.5:0.01:5.0)*1e9
# harmonic balance
jpa = hbsolve(2pi*fs, (2pi*fp,), [(mode = (1,), port = 1, current = Ip)], (8,), (16,), circuit)
@assert jpa.nonlinear.solverinfo.converged
hbgain = 10*log10.(abs2.(jpa.linearized.S(outputmode = (0,), outputport = 1,
inputmode = (0,), inputport = 1, freqindex = :)))
# the pump alone in time, eighty steps per pump period, WRspice's rise
rise(t) = 1 - 2/(exp(t/10e-9) + exp(-t/10e-9))
drive(Ip, fp) = t -> 2Ip*rise(t)*cospi(2fp*t)
problem = transientproblem(circuit; sources = [TransientSource(1, drive(Ip, fp))])
steps = 80
dt = 1/(steps*fp)
pump = transientsolve(problem, (0.0, 76000dt); dt, method = GaussLegendre(), record = :phases)
# the signals as directions of the tangent: unit cosine and sine at each frequency
currents = zeros(1, length(pump.times), 2length(fs))
for (k, f) in enumerate(fs)
currents[1, :, 2k - 1] .= rise.(pump.times) .* cospi.(2f .* pump.times)
currents[1, :, 2k] .= rise.(pump.times) .* sinpi.(2f .* pump.times)
end
signal = transienttangent(pump, currents)
# the reflection over the last pump period
reflection(v, f, t) = 2*sum(v .* cispi.(-2f .* t))/length(t)/50 - 1
lastperiod = length(pump.times) - steps + 1:length(pump.times)
S11 = [reflection(signal.voltage[1, lastperiod, 2k - 1] .+ im .* signal.voltage[1, lastperiod, 2k], fs[k],
pump.times[lastperiod]) for k in eachindex(fs)]
tdgain = 10*log10.(abs2.(S11))
# WRspice: the pump alone, then the pump with a small sine and a small cosine at each frequency
netlist = JosephsonCircuits.exportnetlist(circuit)
input = JosephsonCircuits.wrspice_input_paramp(netlist.netlist, 2pi*fs, 2pi*fp, 2Ip, (0, 1), (0, 1);
stepsperperiod = steps)
output = JosephsonCircuits.spice_run(input, XicTools_jll.wrspice())
wrgain = 10*log10.(abs2.(JosephsonCircuits.wrspice_calcS_paramp(output, 2pi*fs, netlist.Nnodes;
stepsperperiod = steps).S11))
plot(fs/1e9, hbgain; label = "harmonic balance", xlabel = "Frequency (GHz)", ylabel = "Gain (dB)")
plot!(fs/1e9, tdgain; label = "JosephsonCircuits.jl transient", seriestype = :scatter)
plot!(fs/1e9, wrgain; label = "WRspice", seriestype = :scatter, marker = :x)Once compiled, the pump solve takes about 0.2 s and the tangent along its 102 directions about 1 s. One WRspice run of this input takes about 1.8 s, and spice_run runs one WRspice process per logical processor, so the 103 runs take about 103 × 1.8/8 ≈ 25 s on eight. The code gives these gains at seven of the signal frequencies:
| Signal (GHz) | Harmonic balance (dB) | Transient (dB) | WRspice (dB) |
|---|---|---|---|
| 4.50 | 0.0027 | 0.0027 | 0.0027 |
| 4.70 | 0.7108 | 0.7102 | 0.7120 |
| 4.73 | 4.1802 | 4.1781 | 4.1874 |
| 4.75 | 13.3023 | 13.1726 | 13.2484 |
| 4.77 | 4.1857 | 4.1843 | 4.1930 |
| 4.80 | 0.7115 | 0.7115 | 0.7130 |
| 4.90 | 0.0193 | 0.0194 | 0.0193 |
The two time domain results are within 0.03 dB of each other and of harmonic balance across the band, and at the gain peak within 0.08 dB of each other and 0.13 dB of harmonic balance, where the gain is most sensitive to the operating point. Check both settling time and time step before interpreting a peak-gain discrepancy.
The port response resolved in time at one signal frequency, 4.76 GHz, is the same demodulation slid along the record, one pump period at a time, through the rise of the pump. From WRspice it is the difference between the runs with and without the signal, which the tangent gives directly:
k = findfirst(==(4.76e9), fs)
Is = 1e-13
vw = ((output[2k + 1].values["V"][1, :] .- output[1].values["V"][1, :]) .+
im .* (output[2k].values["V"][1, :] .- output[1].values["V"][1, :])) ./ Is
vt = signal.voltage[1, :, 2k - 1] .+ im .* signal.voltage[1, :, 2k]
js = steps:steps:length(pump.times)
envelope(v) = [10*log10(abs2(reflection(v[j - steps + 1:j], fs[k], pump.times[j - steps + 1:j]))) for j in js]
plot(pump.times[js]*1e9, envelope(vt); label = "JosephsonCircuits.jl transient",
xlabel = "Time (ns)", ylabel = "Signal gain (dB)")
plot!(pump.times[js]*1e9, envelope(vw); label = "WRspice", linestyle = :dash)| Time (ns) | Transient (dB) | WRspice (dB) |
|---|---|---|
| 4.2 | -1.093 | -1.093 |
| 12.6 | -4.651 | -4.650 |
| 21.1 | -3.255 | -3.255 |
| 42.1 | 0.868 | 0.872 |
| 63.2 | 5.319 | 5.333 |
| 84.2 | 7.734 | 7.757 |
| 109.5 | 8.454 | 8.483 |
| 160.0 | 8.153 | 8.176 |
| 197.9 | 8.181 | 8.205 |
The signal is first attenuated; its gain then rises, overshoots, and settles. The two simulations differ by less than 0.03 dB at the listed times.