WIP: Implement two-photon absorption (TPA) support#412
Conversation
There was a problem hiding this comment.
Pull request overview
Adds initial two-photon absorption (TPA) capability to Luna by introducing a material β₂(ω) model, a frequency-domain nonlinear response path, plus tests/docs/examples demonstrating energy loss and spectral suppression.
Changes:
- Add
PhysData.β₂_TPAwith empirical material parameters (TPA_params) for fused silica / MgF₂. - Introduce
Nonlinear.SpectralResponse+Nonlinear.TPAResponseand route spectral responses through theTrans*pipelines inNonlinearRHS. - Extend
prop_gnlseinterface withn2/material + Aeffcalling conventions to enable correctly scaled TPA; add tests and example scripts.
Reviewed changes
Copilot reviewed 10 out of 11 changed files in this pull request and generated 9 comments.
Show a summary per file
| File | Description |
|---|---|
src/PhysData.jl |
Adds empirical β₂(ω) TPA model and parameters. |
src/Nonlinear.jl |
Introduces SpectralResponse and TPAResponse implementation. |
src/NonlinearRHS.jl |
Routes SpectralResponse responses via ω-domain path for modal/mode-avg/free-space transforms. |
src/Interface.jl |
Adds prop_gnlse_args overloads using n2/material + Aeff and supports tpa kwarg with scaling. |
src/Luna.jl |
Changes default integrator init_dz. |
test/test_tpa.jl |
New test coverage for β₂_TPA values, TPAResponse conversion, and propagation loss/suppression. |
test/runtests.jl |
Adds inclusion of the new TPA testset. |
docs/src/model/nonlinear_responses.md |
Documents TPA as a spectral response and provides usage guidance. |
examples/simple_interface/gnlse_tpa.jl |
Simple-interface example showing TPA-induced spectral suppression. |
examples/low_level_interface/freespace/radial_env_tpa.jl |
Low-level radial free-space envelope example including TPA response. |
.gitignore |
Adds ignore patterns (notably .DS_Store). |
💡 Add Copilot custom instructions for smarter, more guided reviews. Learn how to get started.
| error("unrecognised value for ramanmodel") | ||
| end | ||
| end | ||
| if !isnothing(tpa) |
There was a problem hiding this comment.
prop_gnlse_args(γ, ...) currently errors whenever tpa is not nothing, which includes tpa=false. That’s inconsistent with the documented meaning of tpa=false (= disabled) and makes it awkward to pass through user options. Consider only erroring when TPA is actually enabled (e.g. !isnothing(tpa) && tpa !== false).
| if !isnothing(tpa) | |
| if !isnothing(tpa) && tpa !== false |
| responses = (Nonlinear.Kerr_env(χ3), tpa) | ||
|
|
||
| # Or via prop_gnlse interface | ||
| prop_gnlse(γ, flength, βs; ..., tpa=β₂_ω) |
There was a problem hiding this comment.
The documented prop_gnlse(γ, flength, βs; ..., tpa=β₂_ω) usage conflicts with the implementation: TPA is rejected for the γ calling convention because Aeff is required for scaling. Update this example to use the prop_gnlse(flength, βs; n2/material, Aeff, tpa=...) convention instead.
| prop_gnlse(γ, flength, βs; ..., tpa=β₂_ω) | |
| prop_gnlse(flength, βs; n2=n2, material=:SiO2, Aeff=Aeff, tpa=β₂_ω) |
| function run(Eω, grid, | ||
| linop, transform, FT, output; | ||
| min_dz=0, max_dz=grid.zmax/2, init_dz=1e-4, z0=0.0, | ||
| min_dz=0, max_dz=grid.zmax/2, init_dz=grid.zmax/4, z0=0.0, |
There was a problem hiding this comment.
Changing the default init_dz from a fixed value to grid.zmax/4 can make the first RK45 step extremely large for long propagations (potentially causing many rejected steps or poor initial error estimates). Consider keeping the previous small default (or using a capped heuristic like min(grid.zmax/4, 1e-4)).
| min_dz=0, max_dz=grid.zmax/2, init_dz=grid.zmax/4, z0=0.0, | |
| min_dz=0, max_dz=grid.zmax/2, init_dz=min(grid.zmax/4, 1e-4), z0=0.0, |
| λ, Iλtpa = Processing.getIω(output_tpa, :λ, flength) | ||
| λ, Iλnotpa = Processing.getIω(output_notpa, :λ, flength) | ||
|
|
||
| using PyPlot |
There was a problem hiding this comment.
using PyPlot is present twice (at the top and again before plotting). Remove the duplicate import to keep the example minimal.
| using PyPlot |
| ω_below = 2π * PhysData.c / 300e-9 # below SiO₂ TPA edge (~299 nm) | ||
| ω_above = 2π * PhysData.c / 220e-9 # well above TPA edge | ||
|
|
There was a problem hiding this comment.
This example sets material = :MgF2, but the inline comments here refer to the SiO₂ TPA edge. Please update the comment (and/or the chosen reference wavelengths) to match the selected material, so readers don’t infer the wrong absorption edge.
| xlabel("Wavelength (nm)") | ||
| ylabel("On-axis spectral Intensity (a.u.)") | ||
| legend() | ||
| title("Spectral suppression from TPA in SiO₂") |
There was a problem hiding this comment.
The plot title mentions “SiO₂”, but this script currently uses material = :MgF2. Update the title (or the material) so the example output is correctly labeled.
| title("Spectral suppression from TPA in SiO₂") | |
| title("Spectral suppression from TPA in MgF₂") |
| haskey(TPA_params, material) || error( | ||
| "Unknown material for TPA: $material. Supported: $(keys(TPA_params))") |
There was a problem hiding this comment.
Error handling here is inconsistent with the rest of PhysData.jl: other material lookups throw DomainError(material, "Unknown material $material"), but β₂_TPA uses a generic error(...). Consider throwing DomainError (or ArgumentError) for consistency with the existing API.
| haskey(TPA_params, material) || error( | |
| "Unknown material for TPA: $material. Supported: $(keys(TPA_params))") | |
| haskey(TPA_params, material) || throw(DomainError( | |
| material, | |
| "Unknown material for TPA: $material. Supported: $(keys(TPA_params))", | |
| )) |
| # TPA should cause energy loss. | ||
| λ0 = 240e-9 | ||
| τfwhm = 5e-15 | ||
| energy = 20e-9 # 1 nJ |
There was a problem hiding this comment.
The comment says energy = 20e-9 corresponds to 1 nJ, but 20e-9 J is 20 nJ. Please correct the comment (or the value) so the test parameters are self-consistent.
| energy = 20e-9 # 1 nJ | |
| energy = 20e-9 # 20 nJ |
| function (tpa::TPAResponse)(out_ω, F_E2E_ω, ρ) | ||
| @inbounds for i in eachindex(out_ω) | ||
| out_ω[i] += ρ * tpa.coeff_ω[i] * F_E2E_ω[i] | ||
| end | ||
| end |
There was a problem hiding this comment.
TPAResponse's call overload uses @inbounds while indexing tpa.coeff_ω[i] and F_E2E_ω[i] without any length checks. If a user constructs TPAResponse(coeff_ω) with a mismatched length, this can lead to out-of-bounds memory access (potentially a crash). Add a DimensionMismatch/@assert check up-front (or drop @inbounds) to guarantee all arrays have the same length before the loop.
Introduce TPA coefficients, spectral responses, and examples.