Skip to content

RFC: Additions to Compressor.jl #9

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Open
jtravs wants to merge 5 commits into
base: master
Choose a base branch
from
Open

RFC: Additions to Compressor.jl #9

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
80dade0
include gradient factor
jtravs Feb 11, 2023
eb96980
add max pressure
jtravs Feb 11, 2023
fd34825
add maximum flow
jtravs Feb 11, 2023
ddadf87
fix max flow
jtravs Feb 11, 2023
e6831d6
add polarisation
jtravs Feb 12, 2023
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
125 changes: 90 additions & 35 deletions src/Compressor.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -4,31 +4,52 @@ import HiSol.Solitons: T0P0
import HiSol.HCF: Aeff0, Leff, transmission, intensity_modeavg, α
import HiSol.Focusing: max_flength, diverged_beam
import HiSol.Data: n2_gas
import Luna.PhysData: pressure
import HiSol.GasFlow: tube_PVflow_choked, bar, slpm
import Luna.PhysData: pressure, density
import Luna.Tools: τfw_to_τ0
import Logging: @debug, @info
import Printf: @sprintf
import Roots: find_zero
import PyPlot: plt

function fix_pressure(pr, a, flength, gas, max_pressure, max_flow, ngradient)
pr = min(pr, max_pressure)
flow = 0.0
(flength <= 0) && return 0.0, 0.0
if ngradient > 0
flow = ngradient*slpm(tube_PVflow_choked(a, flength/ngradient, gas, pr*bar))
if flow > max_flow
pr = find_zero(pr) do prg
flow = ngradient*slpm(tube_PVflow_choked(a, flength/ngradient,
gas, prg*bar))
max_flow - flow
end
end
end
flow = ngradient*slpm(tube_PVflow_choked(a, flength/ngradient, gas, pr*bar))
pr, flow
end

function optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
thickness=1e-3, material=:SiO2, Bmax=0.2,
LIDT=2000, S_fluence=5, S_sf=1.5, S_ion=10,
entrance_window=true, exit_window=true)
entrance_window=true, exit_window=true, ngradient=0,
max_pressure=60.0, max_flow=1.0, polarisation=:linear)
factor = τfwhm_in/τfwhm_out

ρcrit = critical_density(gas, λ0, τfwhm_in, energy)
ρ = ρcrit/S_sf
pr = pressure(gas, ρ)
n2 = n2_gas(gas, pr)
prm = pressure(gas, ρ)
@debug(@sprintf("ρcrit: %.4e m⁻³", ρcrit))
@debug(@sprintf("ρ: %.4e m⁻³", ρ))
@debug(@sprintf("pressure: %.3f bar", pr))
@debug(@sprintf("n₂: %.4e m²/W", n2))
@debug(@sprintf("crit pressure: %.3f bar", pr))
@debug(@sprintf("crit n₂: %.4e m²/W", n2))

_, P0 = T0P0(τfwhm_in, energy)

φm = nonlinear_phase(factor)
φm *= (ngradient > 0 ? 3/2 : 1.0)
φm *= (polarisation == :circular ? 3/2 : 1.0)

# φm = γP0Leff
γLeff = φm/P0
Expand All @@ -53,6 +74,8 @@ function optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
aguess *= 0.99
continue
end
pr, flow = fix_pressure(prm, aguess, flength, gas, max_pressure, max_flow, ngradient)
n2 = n2_gas(gas, pr)
γLeff_this = n2*k0/(A0*aguess^2)*Leff(flength, aguess, λ0)
@debug(
@sprintf(
Expand All @@ -64,9 +87,11 @@ function optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
enough = γLeff_this > γLeff
if aguess < 10e-6
error("""Could not find core radius: not enough nonlinearity.
Decrease the compression factor or increase the maximum length.""")
Decrease the compression factor or increase the maximum length,
maximum pressure or maximum flow.""")
end
if P0/(A0*aguess^2) > Isupp/S_ion
I_this = P0/(A0*aguess^2) * (polarisation == :circular ? 2/3 : 1.0)
if I_this > Isupp/S_ion
error("""Could not find core radius: intensity too high.
Decrease the energy or increase the maximum length.""")
end
Expand All @@ -78,37 +103,40 @@ function optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
global flength = max_flength(a, λ0, energy, τfwhm_in, maxlength;
thickness, material, Bmax, LIDT, S_fluence,
entrance_window, exit_window)
global pr, flow = fix_pressure(prm, a, flength, gas, max_pressure, max_flow, ngradient)
n2 = n2_gas(gas, pr)
γLeff_this = n2*k0/(A0*a^2)*Leff(flength, a, λ0)
γLeff_this - γLeff
end
@debug("Optimised", 1e6*aopt)
t = transmission(flength, aopt, λ0)
@debug(@sprintf("Results: a = %.3f μm, %.2f m long, transmission %.4f",
1e6aopt, flength, t))
return aopt, flength, pr, t
return aopt, flength, pr, flow, t
end

function params_maxlength(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
thickness=1e-3, material=:SiO2, Bmax=0.2, S_sf=1.5,
entrance_window=true, exit_window=true, LIDT=2000, S_fluence=5)

ρcrit = critical_density(gas, λ0, τfwhm_in, energy)
ρ = ρcrit/S_sf
pr = pressure(gas, ρ)
n2 = n2_gas(gas, pr)

entrance_window=true, exit_window=true, LIDT=2000, S_fluence=5,
ngradient=0, max_pressure=60.0, max_flow=1.0, polarisation=:linear)
_, P0 = T0P0(τfwhm_in, energy)

broadfac_req = τfwhm_in/τfwhm_out
φnl_req = nonlinear_phase(broadfac_req)

φnl_req *= (ngradient > 0 ? 3/2 : 1.0)
φnl_req *= (polarisation == :circular ? 3/2 : 1.0)
k0 = 2π/λ0
A0 = Aeff0()

function params(a)
maxflength = max_flength(a, λ0, energy, τfwhm_in, maxlength;
thickness, material, Bmax, LIDT, S_fluence,
entrance_window, exit_window)
thickness, material, Bmax, LIDT, S_fluence,
entrance_window, exit_window)
ρcrit = critical_density(gas, λ0, τfwhm_in, energy)
ρ = ρcrit/S_sf
prm = pressure(gas, ρ)
pr, flow = fix_pressure(prm, a, maxflength, gas, max_pressure, max_flow, ngradient)
n2 = n2_gas(gas, pr)

γthis = n2*k0/(A0*a^2)
Leff_req = φnl_req/P0/γthis
maxLeff = Leff(maxflength, a, λ0)
Expand All @@ -121,26 +149,27 @@ function params_maxlength(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;

γLeff = γthis * Leff(flength, a, λ0)
t = transmission.(flength, a, λ0)
intensity = P0/(A0*a^2)
intensity = P0/(A0*a^2) * (polarisation == :circular ? 2/3 : 1.0)
φnl = P0*γLeff
window_distance = (maxlength-maxflength)/2
beamsize_w0 = diverged_beam(a, λ0, window_distance)
(;φnl=P0*γLeff, transmission=t, intensity, flength,
broadening_factor=broadening_factor(φnl), pressure=pr, n2=n2,
P0=P0, window_distance, beamsize_w0)
P0=P0, window_distance, beamsize_w0, flow)
end
end

function plot_optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
thickness=1e-3, material=:SiO2, Bmax=0.2, S_sf=1.5, S_ion=10,
entrance_window=true, exit_window=true, LIDT=2000, S_fluence=5,
amin=25e-6, amax=500e-6, Na=512)
amin=25e-6, amax=500e-6, Na=512,
ngradient=0, max_pressure=60.0, max_flow=1.0, polarisation=:linear)
factor = τfwhm_in/τfwhm_out

f = params_maxlength(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
thickness, material, Bmax, S_sf,
entrance_window, exit_window, LIDT, S_fluence)

entrance_window, exit_window, LIDT, S_fluence, ngradient,
max_pressure, max_flow, polarisation)

Isupp = barrier_suppression_intensity(gas)
A0 = Aeff0()
Expand All @@ -152,21 +181,23 @@ function plot_optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
t = getindex.(params, :transmission)
intensity = getindex.(params, :intensity)
P0 = params[1].P0
pr = params[1].pressure
pr = getindex.(params, :pressure)
flow = getindex.(params, :flow)

t[flength .== 0] .= NaN

try
global aopt, flopt, _, topt = optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
thickness, material, Bmax, S_sf, S_ion, LIDT, S_fluence, entrance_window, exit_window)
global intopt = P0/(A0*aopt^2)
global aopt, flopt, propt, flowopt, topt = optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
thickness, material, Bmax, S_sf, S_ion, LIDT, S_fluence, entrance_window, exit_window,
ngradient, max_pressure, max_flow, polarisation)
global intopt = P0/(A0*aopt^2) * (polarisation == :circular ? 2/3 : 1.0)
catch e
global aopt = missing
end

fig = plt.figure()
fig.set_size_inches(15, 4)
plt.subplot(1, 4, 1)
fig.set_size_inches(15, 7)
plt.subplot(2, 4, 1)
plt.plot(1e6a, broadfac)
if ~ismissing(aopt)
plt.plot(1e6aopt, factor, "o"; color="k", label=@sprintf("%.1f μm", 1e6aopt))
Expand All @@ -179,7 +210,7 @@ function plot_optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
plt.ylabel("Broadening factor")


plt.subplot(1, 4, 2)
plt.subplot(2, 4, 2)
plt.plot(1e6a, flength)
if ~ismissing(aopt)
plt.plot(1e6aopt, flopt, "o"; color="k", label=@sprintf("%.3f m", flopt))
Expand All @@ -190,7 +221,7 @@ function plot_optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
plt.xlabel("Core radius (μm)")
plt.ylabel("Fibre length (m)")

plt.subplot(1, 4, 3)
plt.subplot(2, 4, 3)
plt.plot(1e6a, 100*t)
if ~ismissing(aopt)
plt.plot(1e6aopt, 100*topt, "o"; color="k", label=@sprintf("%.1f %%", 100*topt))
Expand All @@ -201,7 +232,7 @@ function plot_optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
plt.xlabel("Core radius (μm)")
plt.ylabel("Transmission (%)")

plt.subplot(1, 4, 4)
plt.subplot(2, 4, 4)
plt.plot(1e6a, intensity*1e-4)
plt.axhline(Isupp*1e-4/S_ion; linestyle="--", color="r", label="Limit")
if ~ismissing(aopt)
Expand All @@ -214,9 +245,33 @@ function plot_optimise(τfwhm_in, τfwhm_out, gas, λ0, energy, maxlength;
plt.xlabel("Core radius (μm)")
plt.ylabel("Intensity (W/cm\$^{-2}\$)")

plt.subplot(2, 4, 5)
plt.plot(1e6a, flow)
plt.axhline(max_flow; linestyle="--", color="r", label="Limit")
if ~ismissing(aopt)
plt.plot(1e6aopt, flowopt, "o"; color="k", label=@sprintf("%.3f slpm", flowopt))
plt.legend()
end
plt.xlim(1e6.*extrema(a))
plt.ylim(ymin=0)
plt.xlabel("Core radius (μm)")
plt.ylabel("Gas flow (slpm)")

plt.subplot(2, 4, 6)
plt.plot(1e6a, pr)
plt.axhline(max_pressure; linestyle="--", color="r", label="Limit")
if ~ismissing(aopt)
plt.plot(1e6aopt, propt, "o"; color="k", label=@sprintf("%.3f bar", propt))
plt.legend()
end
plt.xlim(1e6.*extrema(a))
plt.ylim(ymin=0)
plt.xlabel("Core radius (μm)")
plt.ylabel("Pressure (bar)")

fig.tight_layout()
plt.subplots_adjust(top=0.9)
plt.suptitle(@sprintf("Input peak power: %.2f GW | Pressure: %.2f bar", P0*1e-9, pr))
plt.suptitle(@sprintf("Input peak power: %.2f GW", P0*1e-9))

return fig, f
end
Expand Down
2 changes: 1 addition & 1 deletion src/HiSol.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -5,8 +5,8 @@ include("HCF.jl")
include("Solitons.jl")
include("Limits.jl")
include("Focusing.jl")
include("Compressor.jl")
include("GasFlow.jl")
include("Compressor.jl")
include("Design.jl")

end