diff --git a/.gitignore b/.gitignore new file mode 100644 index 0000000..2251642 --- /dev/null +++ b/.gitignore @@ -0,0 +1 @@ +Manifest.toml \ No newline at end of file diff --git a/Examples/qo_singlephoton-DESKTOP-O66V958.jl b/Examples/qo_singlephoton-DESKTOP-O66V958.jl new file mode 100644 index 0000000..be07d87 --- /dev/null +++ b/Examples/qo_singlephoton-DESKTOP-O66V958.jl @@ -0,0 +1,85 @@ +using QuantumOptics +using WaveguideQED +using LinearAlgebra +using PyPlot +pygui(false) +pygui(true) +include("singlecavity_simulations.jl") +#Parameter structure imported from singlecavity_simulations (bit overkill for now, but can prove usefull) +param = parameters() + +#Set detuning: +param.δ = 0 +#Set non linearity: +param.x3 = 0 + +param.times = 0:0.1:10 + +dt = param.times[2] - param.times[1] + +#Create operators for two photons interacting with cavity +bc = FockBasis(1) +bw = WaveguideBasis(1,param.times) +a = destroy(bc) +ad = create(bc); +n = ad*a ⊗ identityoperator(bw) +#$w†a and a†w efficient implementation$ +wda = emission(bc,bw) +adw = absorption(bc,bw) +H = param.δ*n + im*sqrt(param.γ/dt)*(adw-wda) + param.x3/4*(n*n+n) + + +#Define input onephoton state shape +ξfun(t,σ,t0) = complex(sqrt(2/σ)* (log(2)/pi)^(1/4)*exp(-2*log(2)*(t-t0)^2/σ^2)) +ξvec = ξfun.(param.times,param.σ,param.t0) + +ψ_cw = onephoton(bw,ξvec) +psi = fockstate(bc,0) ⊗ ψ_cw + +#Solve +@btime ψ = waveguide_evolution(param.times, psi, H) + +#REFERENCE SOLUTION +sol1 = solve_differentialeq(param,ξfun) +ref_sol = ξfun.(sol1.t,param.σ,param.t0)-sqrt(param.γ)*sol1 + +#Plot single photon waveguide state +ψ_single = OnePhotonView(ψ)/sqrt(dt) + +fig,ax = subplots(1,1,figsize=(9,4.5)) +ax.plot(param.times,abs.(ξvec).^2,"g-",label="Input pulse") +ax.plot(param.times,abs.(ψ_single).^2,"ro",label="δ = 0",fillstyle="none") +ax.plot(sol1.t,abs.(ref_sol).^2,"r-") +ax.set_xlabel("time [1/γ]") +ax.set_ylabel(L"$\xi_{out}^{(1)}$") +ax.legend() +plt.tight_layout() + +using DifferentialEquations +using LinearAlgebra +using PyPlot +function dpsi!(dpsi,psi,p,t) + y,d,nsteps,dt = p + timeindex = round(Int,t/dt,RoundDown) + 1 + dpsi .= 0 + dpsi[2+nsteps] = -sqrt(y/dt)*psi[1+timeindex] + dpsi[1+timeindex] = sqrt(y/dt)*psi[2+nsteps] +end + +#Define parameters +y,d,dt = 1,0,0.1 +times = 0:dt:10 +N = length(times) +p = (y,d,N,dt) + +#Define input gaussian state with width s = 1 and arrivial time t0=5 +xi(t,s,t0) = complex(sqrt(2/s)* (log(2)/pi)^(1/4)*exp(-2*log(2)*(t-t0)^2/s^2)) +psi = zeros(ComplexF64,2*(N+1)) +psi[2:N+1] .= sqrt(dt)*xi.(times,1,5) +prob = ODEProblem(dpsi!, psi, (times[1], times[end]+dt), p) +sol = solve(prob, OrdinaryDiffEq.DP5();reltol = 1.0e-8,abstol = 1.0e-10); +psi_out = sol[end][2:N+1] + +#Plot the output +fig,ax = subplots(1,1,figsize=(9,4.5)) +ax.plot(times,abs.(psi_out).^2,"r-",label="Input pulse") \ No newline at end of file diff --git a/Examples/qo_singlephoton.jl b/Examples/qo_singlephoton.jl index caa3df0..16b65a2 100644 --- a/Examples/qo_singlephoton.jl +++ b/Examples/qo_singlephoton.jl @@ -9,7 +9,7 @@ include("singlecavity_simulations.jl") param = parameters() #Set detuning: -param.δ = 0 +param.δ = 2 #Set non linearity: param.x3 = 0 @@ -53,4 +53,5 @@ ax.plot(sol1.t,abs.(ref_sol).^2,"r-") ax.set_xlabel("time [1/γ]") ax.set_ylabel(L"$\xi_{out}^{(1)}$") ax.legend() -plt.tight_layout() \ No newline at end of file +plt.tight_layout() + diff --git a/Examples/qo_twophoton.jl b/Examples/qo_twophoton.jl index c931d47..7ecaa7d 100644 --- a/Examples/qo_twophoton.jl +++ b/Examples/qo_twophoton.jl @@ -41,6 +41,7 @@ psi = fockstate(bc,0) ⊗ ψ_cw ψ_double = TwoPhotonView(ψ) + fig,ax = subplots(1,1,figsize=(4.5,4.5)) xgrid = repeat(param.times',length(param.times),1) ygrid = repeat(param.times,1,length(param.times)) diff --git a/Manifest.toml b/Manifest.toml index ff41b8a..b8c57ff 100644 --- a/Manifest.toml +++ b/Manifest.toml @@ -1,20 +1,28 @@ # This file is machine-generated - editing it directly is not advised -julia_version = "1.8.5" +julia_version = "1.9.0" manifest_format = "2.0" -project_hash = "17108abd943556615807072c82e1cef39046804b" +project_hash = "4dde56137bf92cdd4b71288bb369d001941db81e" [[deps.AbstractFFTs]] -deps = ["ChainRulesCore", "LinearAlgebra"] +deps = ["LinearAlgebra"] git-tree-sha1 = "16b6dbc4cf7caee4e1e75c49485ec67b667098a0" uuid = "621f4979-c628-5d54-868e-fcf4e3e8185c" version = "1.3.1" +weakdeps = ["ChainRulesCore"] + + [deps.AbstractFFTs.extensions] + AbstractFFTsChainRulesCoreExt = "ChainRulesCore" [[deps.Adapt]] deps = ["LinearAlgebra", "Requires"] git-tree-sha1 = "cc37d689f599e8df4f464b2fa3870ff7db7492ef" uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" version = "3.6.1" +weakdeps = ["StaticArrays"] + + [deps.Adapt.extensions] + AdaptStaticArraysExt = "StaticArrays" [[deps.ArgTools]] uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" @@ -104,12 +112,6 @@ git-tree-sha1 = "c6d890a52d2c4d55d326439580c3b8d0875a77d9" uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" version = "1.15.7" -[[deps.ChangesOfVariables]] -deps = ["ChainRulesCore", "LinearAlgebra", "Test"] -git-tree-sha1 = "485193efd2176b88e6622a39a246f8c5b600e74e" -uuid = "9e997f8a-9a97-42d5-a9f1-ce6bfc15e2c0" -version = "0.1.6" - [[deps.CloseOpenIntervals]] deps = ["ArrayInterface", "Static"] git-tree-sha1 = "d61300b9895f129f4bd684b2aff97cf319b6c493" @@ -128,15 +130,19 @@ uuid = "bbf7d656-a473-5ed7-a52c-81e309532950" version = "0.3.0" [[deps.Compat]] -deps = ["Dates", "LinearAlgebra", "UUIDs"] +deps = ["UUIDs"] git-tree-sha1 = "7a60c856b9fa189eb34f5f8a6f6b5529b7942957" uuid = "34da2185-b29b-5c13-b0c7-acf172513d20" version = "4.6.1" +weakdeps = ["Dates", "LinearAlgebra"] + + [deps.Compat.extensions] + CompatLinearAlgebraExt = "LinearAlgebra" [[deps.CompilerSupportLibraries_jll]] deps = ["Artifacts", "Libdl"] uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae" -version = "1.0.1+0" +version = "1.0.2+0" [[deps.ConstructionBase]] deps = ["LinearAlgebra"] @@ -144,6 +150,14 @@ git-tree-sha1 = "89a9db8d28102b094992472d333674bd1a83ce2a" uuid = "187b0558-2788-49d3-abe0-74a17ed4e7c9" version = "1.5.1" + [deps.ConstructionBase.extensions] + IntervalSetsExt = "IntervalSets" + StaticArraysExt = "StaticArrays" + + [deps.ConstructionBase.weakdeps] + IntervalSets = "8197267c-284f-5f27-9208-e0e47529a953" + StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" + [[deps.CpuId]] deps = ["Markdown"] git-tree-sha1 = "fcbb72b032692610bfbdb15018ac16a36cf2e406" @@ -170,18 +184,34 @@ version = "1.0.0" deps = ["Printf"] uuid = "ade2ca70-3891-5945-98fb-dc099432e06a" -[[deps.DensityInterface]] -deps = ["InverseFunctions", "Test"] -git-tree-sha1 = "80c3e8639e3353e5d2912fb3a1916b8455e2494b" -uuid = "b429d917-457f-4dbc-8f4c-0cc954292b1d" -version = "0.4.0" - [[deps.DiffEqBase]] -deps = ["ArrayInterfaceCore", "ChainRulesCore", "DataStructures", "Distributions", "DocStringExtensions", "EnumX", "FastBroadcast", "ForwardDiff", "FunctionWrappers", "FunctionWrappersWrappers", "LinearAlgebra", "Logging", "Markdown", "MuladdMacro", "Parameters", "PreallocationTools", "Printf", "RecursiveArrayTools", "Reexport", "Requires", "SciMLBase", "Setfield", "SparseArrays", "Static", "StaticArraysCore", "Statistics", "Tricks", "ZygoteRules"] +deps = ["ArrayInterfaceCore", "ChainRulesCore", "DataStructures", "DocStringExtensions", "EnumX", "FastBroadcast", "ForwardDiff", "FunctionWrappers", "FunctionWrappersWrappers", "LinearAlgebra", "Logging", "Markdown", "MuladdMacro", "Parameters", "PreallocationTools", "Printf", "RecursiveArrayTools", "Reexport", "Requires", "SciMLBase", "Setfield", "SparseArrays", "Static", "StaticArraysCore", "Statistics", "Tricks", "ZygoteRules"] git-tree-sha1 = "81470904b958f3f24173fa013d4a54198842cb8d" uuid = "2b5f629d-d688-5b77-993f-72d75c75574e" version = "6.119.0" + [deps.DiffEqBase.extensions] + DiffEqBaseDistributionsExt = "Distributions" + DiffEqBaseGeneralizedGeneratedExt = "GeneralizedGenerated" + DiffEqBaseMPIExt = "MPI" + DiffEqBaseMeasurementsExt = "Measurements" + DiffEqBaseMonteCarloMeasurementsExt = "MonteCarloMeasurements" + DiffEqBaseReverseDiffExt = "ReverseDiff" + DiffEqBaseTrackerExt = "Tracker" + DiffEqBaseUnitfulExt = "Unitful" + DiffEqBaseZygoteExt = "Zygote" + + [deps.DiffEqBase.weakdeps] + Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" + GeneralizedGenerated = "6b9d7cbe-bcb9-11e9-073f-15a7a543e2eb" + MPI = "da04e1cc-30fd-572f-bb4f-1f8673147195" + Measurements = "eff96d63-e80a-5855-80a2-b1b0885c5ab7" + MonteCarloMeasurements = "0987c9cc-fe09-11e8-30f0-b96dd679fdca" + ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" + Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" + Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" + Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" + [[deps.DiffEqCallbacks]] deps = ["DataStructures", "DiffEqBase", "ForwardDiff", "LinearAlgebra", "Markdown", "NLsolve", "Parameters", "RecipesBase", "RecursiveArrayTools", "SciMLBase", "StaticArraysCore"] git-tree-sha1 = "b497f63a13fe37e03ed7ac72d71b72aad17b46c4" @@ -194,6 +224,12 @@ git-tree-sha1 = "2c4ed3eedb87579bfe9f20ecc2440de06b9f3b89" uuid = "77a26b50-5914-5dd7-bc55-306e6241c503" version = "5.16.0" + [deps.DiffEqNoiseProcess.extensions] + DiffEqNoiseProcessReverseDiffExt = "ReverseDiff" + + [deps.DiffEqNoiseProcess.weakdeps] + ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" + [[deps.DiffResults]] deps = ["StaticArraysCore"] git-tree-sha1 = "782dd5f4561f5d267313f23853baaaa4c52ea621" @@ -217,11 +253,19 @@ deps = ["Random", "Serialization", "Sockets"] uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" [[deps.Distributions]] -deps = ["ChainRulesCore", "DensityInterface", "FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SparseArrays", "SpecialFunctions", "Statistics", "StatsBase", "StatsFuns", "Test"] +deps = ["FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SparseArrays", "SpecialFunctions", "Statistics", "StatsBase", "StatsFuns", "Test"] git-tree-sha1 = "da9e1a9058f8d3eec3a8c9fe4faacfb89180066b" uuid = "31c24e10-a181-5473-b8eb-7969acd0382f" version = "0.25.86" + [deps.Distributions.extensions] + DistributionsChainRulesCoreExt = "ChainRulesCore" + DistributionsDensityInterfaceExt = "DensityInterface" + + [deps.Distributions.weakdeps] + ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" + DensityInterface = "b429d917-457f-4dbc-8f4c-0cc954292b1d" + [[deps.DocStringExtensions]] deps = ["LibGit2"] git-tree-sha1 = "2fb1e02f2b635d0845df5d7c167fec4dd739b00d" @@ -300,10 +344,14 @@ uuid = "6a86dc24-6348-571c-b903-95158fe2bd41" version = "2.17.0" [[deps.ForwardDiff]] -deps = ["CommonSubexpressions", "DiffResults", "DiffRules", "LinearAlgebra", "LogExpFunctions", "NaNMath", "Preferences", "Printf", "Random", "SpecialFunctions", "StaticArrays"] +deps = ["CommonSubexpressions", "DiffResults", "DiffRules", "LinearAlgebra", "LogExpFunctions", "NaNMath", "Preferences", "Printf", "Random", "SpecialFunctions"] git-tree-sha1 = "00e252f4d706b3d55a8863432e742bf5717b498d" uuid = "f6369f11-7733-5829-9624-2563aa707210" version = "0.10.35" +weakdeps = ["StaticArrays"] + + [deps.ForwardDiff.extensions] + ForwardDiffStaticArraysExt = "StaticArrays" [[deps.FunctionWrappers]] git-tree-sha1 = "d62485945ce5ae9c0c48f124a84998d755bae00e" @@ -380,12 +428,6 @@ version = "2018.0.3+2" deps = ["Markdown"] uuid = "b77e0a4c-d291-57a0-90e8-8db25a27a240" -[[deps.InverseFunctions]] -deps = ["Test"] -git-tree-sha1 = "49510dfcb407e572524ba94aeae2fced1f3feb0f" -uuid = "3587e190-3f89-42d0-90ee-14403ec27112" -version = "0.1.8" - [[deps.IrrationalConstants]] git-tree-sha1 = "630b497eafcc20001bba38a4651b327dcfc491d2" uuid = "92d709cd-6900-40b7-9082-c6be49f344b6" @@ -488,7 +530,7 @@ uuid = "d3d80556-e9d4-5f37-9878-2ab0fcc64255" version = "7.2.0" [[deps.LinearAlgebra]] -deps = ["Libdl", "libblastrampoline_jll"] +deps = ["Libdl", "OpenBLAS_jll", "libblastrampoline_jll"] uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" [[deps.LinearMaps]] @@ -503,20 +545,41 @@ git-tree-sha1 = "d1fce810e9a4213607f0182cf25ffd6ce13e19b6" uuid = "7ed4a6bd-45f5-4d41-b270-4a48e9bafcae" version = "1.37.0" + [deps.LinearSolve.extensions] + LinearSolveHYPRE = "HYPRE" + + [deps.LinearSolve.weakdeps] + HYPRE = "b5ffcf37-a2bd-41ab-a3da-4bd9bc8ad771" + [[deps.LogExpFunctions]] -deps = ["ChainRulesCore", "ChangesOfVariables", "DocStringExtensions", "InverseFunctions", "IrrationalConstants", "LinearAlgebra"] +deps = ["DocStringExtensions", "IrrationalConstants", "LinearAlgebra"] git-tree-sha1 = "0a1b7c2863e44523180fdb3146534e265a91870b" uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688" version = "0.3.23" + [deps.LogExpFunctions.extensions] + LogExpFunctionsChainRulesCoreExt = "ChainRulesCore" + LogExpFunctionsChangesOfVariablesExt = "ChangesOfVariables" + LogExpFunctionsInverseFunctionsExt = "InverseFunctions" + + [deps.LogExpFunctions.weakdeps] + ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" + ChangesOfVariables = "9e997f8a-9a97-42d5-a9f1-ce6bfc15e2c0" + InverseFunctions = "3587e190-3f89-42d0-90ee-14403ec27112" + [[deps.Logging]] uuid = "56ddb016-857b-54e1-b83d-db4d58db5568" [[deps.LoopVectorization]] -deps = ["ArrayInterface", "ArrayInterfaceCore", "ArrayInterfaceOffsetArrays", "ArrayInterfaceStaticArrays", "CPUSummary", "ChainRulesCore", "CloseOpenIntervals", "DocStringExtensions", "ForwardDiff", "HostCPUFeatures", "IfElse", "LayoutPointers", "LinearAlgebra", "OffsetArrays", "PolyesterWeave", "SIMDTypes", "SLEEFPirates", "SnoopPrecompile", "SpecialFunctions", "Static", "ThreadingUtilities", "UnPack", "VectorizationBase"] +deps = ["ArrayInterface", "ArrayInterfaceCore", "ArrayInterfaceOffsetArrays", "ArrayInterfaceStaticArrays", "CPUSummary", "CloseOpenIntervals", "DocStringExtensions", "HostCPUFeatures", "IfElse", "LayoutPointers", "LinearAlgebra", "OffsetArrays", "PolyesterWeave", "SIMDTypes", "SLEEFPirates", "SnoopPrecompile", "Static", "ThreadingUtilities", "UnPack", "VectorizationBase"] git-tree-sha1 = "9696a80c21a56b937e3fd89e972f8db5db3186e2" uuid = "bdcacae8-1622-11e9-2a5c-532679323890" version = "0.12.150" +weakdeps = ["ChainRulesCore", "ForwardDiff", "SpecialFunctions"] + + [deps.LoopVectorization.extensions] + ForwardDiffExt = ["ChainRulesCore", "ForwardDiff"] + SpecialFunctionsExt = "SpecialFunctions" [[deps.MKL_jll]] deps = ["Artifacts", "IntelOpenMP_jll", "JLLWrappers", "LazyArtifacts", "Libdl", "Pkg"] @@ -542,7 +605,7 @@ uuid = "d6f4376e-aef5-505a-96c1-9c027394607a" [[deps.MbedTLS_jll]] deps = ["Artifacts", "Libdl"] uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1" -version = "2.28.0+0" +version = "2.28.2+0" [[deps.Missings]] deps = ["DataAPI"] @@ -555,7 +618,7 @@ uuid = "a63ad114-7e13-5084-954f-fe012c677804" [[deps.MozillaCACerts_jll]] uuid = "14a3606d-f60d-562e-9121-12d972cd8159" -version = "2022.2.1" +version = "2022.10.11" [[deps.MuladdMacro]] git-tree-sha1 = "cac9cc5499c25554cba55cd3c30543cff5ca4fab" @@ -599,7 +662,7 @@ version = "1.12.9" [[deps.OpenBLAS_jll]] deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"] uuid = "4536629a-c528-5b80-bd46-f80d51c5b363" -version = "0.3.20+0" +version = "0.3.21+4" [[deps.OpenLibm_jll]] deps = ["Artifacts", "Libdl"] @@ -642,9 +705,9 @@ uuid = "d96e819e-fc66-5662-9728-84c9c7592b0a" version = "0.12.3" [[deps.Pkg]] -deps = ["Artifacts", "Dates", "Downloads", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"] +deps = ["Artifacts", "Dates", "Downloads", "FileWatching", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"] uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" -version = "1.8.0" +version = "1.9.0" [[deps.PoissonRandom]] deps = ["Random"] @@ -676,6 +739,18 @@ git-tree-sha1 = "2c7658dd593e3adc118b00429e1048829f1abb8c" uuid = "d236fae5-4411-538c-8e31-a6e3d9e00b46" version = "0.4.11" + [deps.PreallocationTools.extensions] + PreallocationToolsReverseDiffExt = "ReverseDiff" + + [deps.PreallocationTools.weakdeps] + ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" + +[[deps.PrecompileTools]] +deps = ["Preferences"] +git-tree-sha1 = "259e206946c293698122f63e2b513a7c99a244e8" +uuid = "aea7be01-6a6a-4083-8856-8a6e6704d82a" +version = "1.1.1" + [[deps.Preferences]] deps = ["TOML"] git-tree-sha1 = "47e5f437cc0e7ef2ce8406ce1e7e24d44915f88d" @@ -705,15 +780,15 @@ version = "0.1.0" [[deps.QuantumOptics]] deps = ["Arpack", "DiffEqBase", "DiffEqCallbacks", "FFTW", "ForwardDiff", "IterativeSolvers", "LinearAlgebra", "LinearMaps", "OrdinaryDiffEq", "QuantumOpticsBase", "Random", "RecursiveArrayTools", "Reexport", "SparseArrays", "StochasticDiffEq", "WignerSymbols"] -git-tree-sha1 = "81d7cf56eddb8b7c6b65c9da638d10bf2f10d7c2" +git-tree-sha1 = "ff110fa403b1cd68b97527e42fb1333ab0b67969" uuid = "6e0679c1-51ea-5a7c-ac74-d61b76210b0c" -version = "1.0.9" +version = "1.0.10" [[deps.QuantumOpticsBase]] -deps = ["Adapt", "FFTW", "LRUCache", "LinearAlgebra", "QuantumInterface", "Random", "SparseArrays", "Strided", "UnsafeArrays"] -git-tree-sha1 = "fd7c7704d09c2ab619b04634fade1e33bfce8602" +deps = ["Adapt", "FFTW", "FillArrays", "LRUCache", "LinearAlgebra", "QuantumInterface", "Random", "SparseArrays", "Strided", "UnsafeArrays"] +git-tree-sha1 = "0d916919aa637558cc8e146dfaffa60130219db7" uuid = "4f57444f-1401-5e15-980d-4471b28d5678" -version = "0.3.9" +version = "0.4.0" [[deps.REPL]] deps = ["InteractiveUtils", "Markdown", "Sockets", "Unicode"] @@ -752,6 +827,12 @@ git-tree-sha1 = "54e055256bbd41fd10566880bc4baa5316bca6fe" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" version = "2.37.0" + [deps.RecursiveArrayTools.extensions] + RecursiveArrayToolsTrackerExt = "Tracker" + + [deps.RecursiveArrayTools.weakdeps] + Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" + [[deps.RecursiveFactorization]] deps = ["LinearAlgebra", "LoopVectorization", "Polyester", "SnoopPrecompile", "StrideArraysCore", "TriangularSolve"] git-tree-sha1 = "9088515ad915c99026beb5436d0a09cd8c18163e" @@ -845,6 +926,12 @@ git-tree-sha1 = "54c78ac3cc0343a16785adabe5bbf4063c737967" uuid = "727e6d20-b764-4bd8-a329-72de5adea6c7" version = "0.1.14" + [deps.SimpleNonlinearSolve.extensions] + SimpleBatchedNonlinearSolveExt = "NNlib" + + [deps.SimpleNonlinearSolve.weakdeps] + NNlib = "872c559c-99b0-510c-b3b7-b6c96a88d5cd" + [[deps.SimpleTraits]] deps = ["InteractiveUtils", "MacroTools"] git-tree-sha1 = "5d7e3f4e11935503d3ecaf7186eac40602e7d231" @@ -867,7 +954,7 @@ uuid = "a2af1166-a08f-5f64-846c-94a0d3cef48c" version = "1.1.0" [[deps.SparseArrays]] -deps = ["LinearAlgebra", "Random"] +deps = ["Libdl", "LinearAlgebra", "Random", "Serialization", "SuiteSparse_jll"] uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" [[deps.SparseDiffTools]] @@ -883,10 +970,14 @@ uuid = "e56a9233-b9d6-4f03-8d0f-1825330902ac" version = "0.3.9" [[deps.SpecialFunctions]] -deps = ["ChainRulesCore", "IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"] +deps = ["IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"] git-tree-sha1 = "ef28127915f4229c971eb43f3fc075dd3fe91880" uuid = "276daf66-3868-5448-9aa4-cd146d93841b" version = "2.2.0" +weakdeps = ["ChainRulesCore"] + + [deps.SpecialFunctions.extensions] + SpecialFunctionsChainRulesCoreExt = "ChainRulesCore" [[deps.Static]] deps = ["IfElse"] @@ -908,6 +999,7 @@ version = "1.4.0" [[deps.Statistics]] deps = ["LinearAlgebra", "SparseArrays"] uuid = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" +version = "1.9.0" [[deps.StatsAPI]] deps = ["LinearAlgebra"] @@ -922,11 +1014,19 @@ uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" version = "0.33.21" [[deps.StatsFuns]] -deps = ["ChainRulesCore", "HypergeometricFunctions", "InverseFunctions", "IrrationalConstants", "LogExpFunctions", "Reexport", "Rmath", "SpecialFunctions"] +deps = ["HypergeometricFunctions", "IrrationalConstants", "LogExpFunctions", "Reexport", "Rmath", "SpecialFunctions"] git-tree-sha1 = "f625d686d5a88bcd2b15cd81f18f98186fdc0c9a" uuid = "4c63d2b9-4356-54db-8cca-17b64c39e42c" version = "1.3.0" + [deps.StatsFuns.extensions] + StatsFunsChainRulesCoreExt = "ChainRulesCore" + StatsFunsInverseFunctionsExt = "InverseFunctions" + + [deps.StatsFuns.weakdeps] + ChainRulesCore = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" + InverseFunctions = "3587e190-3f89-42d0-90ee-14403ec27112" + [[deps.StochasticDiffEq]] deps = ["Adapt", "ArrayInterface", "DataStructures", "DiffEqBase", "DiffEqNoiseProcess", "DocStringExtensions", "FillArrays", "FiniteDiff", "ForwardDiff", "JumpProcesses", "LevyArea", "LinearAlgebra", "Logging", "MuladdMacro", "NLsolve", "OrdinaryDiffEq", "Random", "RandomNumbers", "RecursiveArrayTools", "Reexport", "SciMLBase", "SparseArrays", "SparseDiffTools", "StaticArrays", "UnPack"] git-tree-sha1 = "c6b4b802d4d830e0e958f5f2098d8dea0a935f4b" @@ -940,10 +1040,16 @@ uuid = "7792a7ef-975c-4747-a70f-980b88e8d1da" version = "0.4.7" [[deps.Strided]] -deps = ["LinearAlgebra", "TupleTools"] -git-tree-sha1 = "a7a664c91104329c88222aa20264e1a05b6ad138" +deps = ["LinearAlgebra", "StridedViews", "TupleTools"] +git-tree-sha1 = "b32eadf6ac726a790567fdc872b63117712e16a8" uuid = "5e0ebb24-38b0-5f93-81fe-25c709ecae67" -version = "1.2.3" +version = "2.0.1" + +[[deps.StridedViews]] +deps = ["LinearAlgebra"] +git-tree-sha1 = "59cc024139c20d1ed8400c419c6fe608637d583d" +uuid = "4db3bf67-4bd7-4b4e-b153-31dc3fb37143" +version = "0.1.2" [[deps.SuiteSparse]] deps = ["Libdl", "LinearAlgebra", "Serialization", "SparseArrays"] @@ -952,7 +1058,7 @@ uuid = "4607b0f0-06f3-5cda-b6b1-a6196a1729e9" [[deps.SuiteSparse_jll]] deps = ["Artifacts", "Libdl", "Pkg", "libblastrampoline_jll"] uuid = "bea87d4a-7f5b-5778-9afe-8cc45184846c" -version = "5.10.1+0" +version = "5.10.1+6" [[deps.SymbolicIndexingInterface]] deps = ["DocStringExtensions"] @@ -963,7 +1069,7 @@ version = "0.2.2" [[deps.TOML]] deps = ["Dates"] uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76" -version = "1.0.0" +version = "1.0.3" [[deps.TableTraits]] deps = ["IteratorInterfaceExtensions"] @@ -980,7 +1086,7 @@ version = "1.10.1" [[deps.Tar]] deps = ["ArgTools", "SHA"] uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e" -version = "1.10.1" +version = "1.10.0" [[deps.Test]] deps = ["InteractiveUtils", "Logging", "Random", "Serialization"] @@ -1052,7 +1158,7 @@ version = "2.0.0" [[deps.Zlib_jll]] deps = ["Libdl"] uuid = "83775a58-1f1d-513f-b197-d71354ab007a" -version = "1.2.12+3" +version = "1.2.13+0" [[deps.ZygoteRules]] deps = ["ChainRulesCore", "MacroTools"] @@ -1061,9 +1167,9 @@ uuid = "700de1a5-db45-46bc-99cf-38207098b444" version = "0.2.3" [[deps.libblastrampoline_jll]] -deps = ["Artifacts", "Libdl", "OpenBLAS_jll"] +deps = ["Artifacts", "Libdl"] uuid = "8e850b90-86db-534c-a0d3-1478176c7d93" -version = "5.1.1+0" +version = "5.7.0+0" [[deps.nghttp2_jll]] deps = ["Artifacts", "Libdl"] diff --git a/Project.toml b/Project.toml index 9ca6809..997db60 100644 --- a/Project.toml +++ b/Project.toml @@ -6,16 +6,18 @@ version = "0.1.0" [deps] LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" Parameters = "d96e819e-fc66-5662-9728-84c9c7592b0a" +PrecompileTools = "aea7be01-6a6a-4083-8856-8a6e6704d82a" QuantumOptics = "6e0679c1-51ea-5a7c-ac74-d61b76210b0c" -SnoopPrecompile = "66db9d55-30c0-4569-8b51-7e840670fc0c" +QuantumOpticsBase = "4f57444f-1401-5e15-980d-4471b28d5678" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" Strided = "5e0ebb24-38b0-5f93-81fe-25c709ecae67" UnsafeArrays = "c4a57d5a-5b31-53a6-b365-19f8c011fbd6" [compat] Parameters = "0.12" +PrecompileTools = "1" QuantumOptics = "1" -SnoopPrecompile = "1" -Strided = "1.2" +QuantumOpticsBase = "0.4" +Strided = "1,2" UnsafeArrays = "1" -julia = "1.8" \ No newline at end of file +julia = "1.9" diff --git a/README.md b/README.md index e4e25f3..70fcdd5 100644 --- a/README.md +++ b/README.md @@ -1,12 +1,12 @@ # WaveguideQED.jl - + -A julia package for simulating quantum states of photon wavepackets using a discrete-time formalism [Phys. Rev. A 101, 042322](https://journals.aps.org/pra/abstract/10.1103/PhysRevA.101.042322). Package works as an extension to [QuantumOptics.jl](https://qojulia.org/) where basises and operators from WaveguideQED.jl can be used together with operators and basises from QuantumOpics.jl. +A Julia package for simulating quantum states of photon wavepackets using a discrete-time formalism [Phys. Rev. A 101, 042322](https://journals.aps.org/pra/abstract/10.1103/PhysRevA.101.042322). The package works as an extension to [QuantumOptics.jl](https://qojulia.org/) where bases and operators from WaveguideQED.jl can be used together with operators and bases from QuantumOpics.jl. ### Example of usage: -Define a waveguide basis, containing a two photon wavepacket for a time interval 0 to 20 with 0.2 timesteps: +Define a waveguide basis, containing a two-photon wavepacket for a time interval from 0 to 20 with timesteps of 0.2: ```julia @@ -32,19 +32,19 @@ wda = a ⊗ wd adw = ad ⊗ w ``` -Finally, we can define an initial twophoton gaussian wavepacket state with view_twophoton and zero photons in the cavity, an Hamiltonian, and simulate the evolution: +Finally, we can define an initial two-photon Gaussian wavepacket state with view_twophoton and zero photons in the cavity, a Hamiltonian, and simulate the evolution: ```julia ξfun(t1,t2,σ1,σ2,t0) = sqrt(2/σ1)* (log(2)/pi)^(1/4)*exp(-2*log(2)*(t1-t0)^2/σ1^2)*sqrt(2/σ2)* (log(2)/pi)^(1/4)*exp(-2*log(2)*(t2-t0)^2/σ2^2) -ψ_cw = twophoton(bw,ξfun,times,1,1,5) +ψ_cw = twophoton(bw,ξfun,1,1,5) psi = fockstate(bc,0) ⊗ ψ_cw dt = times[2] - times[1] H = im*sqrt(1/dt)*(adw-wda) ψ = waveguide_evolution(times, psi, H) ``` -Plotting the twophoton state is also simple: +Plotting the two-photon state is also simple: ```julia diff --git a/docs/.gitignore b/docs/.gitignore index a303fff..aa030f0 100644 --- a/docs/.gitignore +++ b/docs/.gitignore @@ -1,2 +1,3 @@ build/ site/ +Manifest.toml \ No newline at end of file diff --git a/docs/Manifest.toml b/docs/Manifest.toml index 0f3d902..27e2bc4 100644 --- a/docs/Manifest.toml +++ b/docs/Manifest.toml @@ -1,8 +1,8 @@ # This file is machine-generated - editing it directly is not advised -julia_version = "1.8.5" +julia_version = "1.9.0" manifest_format = "2.0" -project_hash = "7dff462257cbbec9dfee017dc62a921b33e39405" +project_hash = "d6a5993505edb88a10fbee1b0a4d47aa45a4e462" [[deps.ANSIColoredPrinters]] git-tree-sha1 = "574baf8110975760d391c710b6341da1afa48d8c" @@ -10,16 +10,24 @@ uuid = "a4c015fc-c6ff-483c-b24f-f7ea428134e9" version = "0.0.1" [[deps.AbstractFFTs]] -deps = ["ChainRulesCore", "LinearAlgebra"] +deps = ["LinearAlgebra"] git-tree-sha1 = "16b6dbc4cf7caee4e1e75c49485ec67b667098a0" uuid = "621f4979-c628-5d54-868e-fcf4e3e8185c" version = "1.3.1" +weakdeps = ["ChainRulesCore"] + + [deps.AbstractFFTs.extensions] + AbstractFFTsChainRulesCoreExt = "ChainRulesCore" [[deps.Adapt]] deps = ["LinearAlgebra", "Requires"] git-tree-sha1 = "cc37d689f599e8df4f464b2fa3870ff7db7492ef" uuid = "79e6a3ab-5dfb-504d-930d-738a2a938a0e" version = "3.6.1" +weakdeps = ["StaticArrays"] + + [deps.Adapt.extensions] + AdaptStaticArraysExt = "StaticArrays" [[deps.ArgTools]] uuid = "0dad84c5-d112-42e6-8d28-ef12dabb789f" @@ -49,6 +57,22 @@ git-tree-sha1 = "38911c7737e123b28182d89027f4216cfc8a9da7" uuid = "4fba245c-0d91-5ea0-9b3e-6abc04ee57a9" version = "7.4.3" + [deps.ArrayInterface.extensions] + ArrayInterfaceBandedMatricesExt = "BandedMatrices" + ArrayInterfaceBlockBandedMatricesExt = "BlockBandedMatrices" + ArrayInterfaceCUDAExt = "CUDA" + ArrayInterfaceGPUArraysCoreExt = "GPUArraysCore" + ArrayInterfaceStaticArraysCoreExt = "StaticArraysCore" + ArrayInterfaceTrackerExt = "Tracker" + + [deps.ArrayInterface.weakdeps] + BandedMatrices = "aae01518-5342-5314-be14-df237901396f" + BlockBandedMatrices = "ffab5731-97b5-5995-9138-79e8c1846df0" + CUDA = "052768ef-5323-5732-b1bb-66c8b64840ba" + GPUArraysCore = "46192b85-c4d5-4398-a991-12ede77f4527" + StaticArraysCore = "1e83bf80-4336-4d27-bf5d-d5a4f845583c" + Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" + [[deps.ArrayInterfaceCore]] deps = ["LinearAlgebra", "SnoopPrecompile", "SparseArrays", "SuiteSparse"] git-tree-sha1 = "e5f08b5689b1aad068e01751889f2f615c7db36d" @@ -108,10 +132,14 @@ uuid = "d360d2e6-b24c-11e9-a2a3-2a2ae2dbcce4" version = "1.15.7" [[deps.ChangesOfVariables]] -deps = ["ChainRulesCore", "LinearAlgebra", "Test"] +deps = ["LinearAlgebra", "Test"] git-tree-sha1 = "485193efd2176b88e6622a39a246f8c5b600e74e" uuid = "9e997f8a-9a97-42d5-a9f1-ce6bfc15e2c0" version = "0.1.6" +weakdeps = ["ChainRulesCore"] + + [deps.ChangesOfVariables.extensions] + ChangesOfVariablesChainRulesCoreExt = "ChainRulesCore" [[deps.CloseOpenIntervals]] deps = ["Static", "StaticArrayInterface"] @@ -119,6 +147,12 @@ git-tree-sha1 = "70232f82ffaab9dc52585e0dd043b5e0c6b714f1" uuid = "fb6a15b2-703c-40df-9091-08a04967cfa9" version = "0.1.12" +[[deps.CodeTracking]] +deps = ["InteractiveUtils", "UUIDs"] +git-tree-sha1 = "d730914ef30a06732bdd9f763f6cc32e92ffbff1" +uuid = "da1fd8a2-8d9e-5ec2-8556-3022fb5608a2" +version = "1.3.1" + [[deps.CodecZlib]] deps = ["TranscodingStreams", "Zlib_jll"] git-tree-sha1 = "ded953804d019afa9a3f98981d99b33e3db7b6da" @@ -157,7 +191,7 @@ version = "4.5.0" [[deps.CompilerSupportLibraries_jll]] deps = ["Artifacts", "Libdl"] uuid = "e66e0078-7015-5450-92f7-15fbd957f2ae" -version = "1.0.1+0" +version = "1.0.2+0" [[deps.Conda]] deps = ["Downloads", "JSON", "VersionParsing"] @@ -171,6 +205,14 @@ git-tree-sha1 = "89a9db8d28102b094992472d333674bd1a83ce2a" uuid = "187b0558-2788-49d3-abe0-74a17ed4e7c9" version = "1.5.1" + [deps.ConstructionBase.extensions] + IntervalSetsExt = "IntervalSets" + StaticArraysExt = "StaticArrays" + + [deps.ConstructionBase.weakdeps] + IntervalSets = "8197267c-284f-5f27-9208-e0e47529a953" + StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" + [[deps.CpuId]] deps = ["Markdown"] git-tree-sha1 = "fcbb72b032692610bfbdb15018ac16a36cf2e406" @@ -204,11 +246,33 @@ uuid = "b429d917-457f-4dbc-8f4c-0cc954292b1d" version = "0.4.0" [[deps.DiffEqBase]] -deps = ["ArrayInterface", "ChainRulesCore", "DataStructures", "Distributions", "DocStringExtensions", "EnumX", "FastBroadcast", "ForwardDiff", "FunctionWrappers", "FunctionWrappersWrappers", "LinearAlgebra", "Logging", "Markdown", "MuladdMacro", "Parameters", "PreallocationTools", "Printf", "RecursiveArrayTools", "Reexport", "Requires", "SciMLBase", "Setfield", "SparseArrays", "Static", "StaticArraysCore", "Statistics", "Tricks", "TruncatedStacktraces", "ZygoteRules"] +deps = ["ArrayInterface", "ChainRulesCore", "DataStructures", "DocStringExtensions", "EnumX", "FastBroadcast", "ForwardDiff", "FunctionWrappers", "FunctionWrappersWrappers", "LinearAlgebra", "Logging", "Markdown", "MuladdMacro", "Parameters", "PreallocationTools", "Printf", "RecursiveArrayTools", "Reexport", "Requires", "SciMLBase", "Setfield", "SparseArrays", "Static", "StaticArraysCore", "Statistics", "Tricks", "TruncatedStacktraces", "ZygoteRules"] git-tree-sha1 = "117b2d02e737aeefd58cd4a4803abecadd37c8cc" uuid = "2b5f629d-d688-5b77-993f-72d75c75574e" version = "6.122.2" + [deps.DiffEqBase.extensions] + DiffEqBaseDistributionsExt = "Distributions" + DiffEqBaseGeneralizedGeneratedExt = "GeneralizedGenerated" + DiffEqBaseMPIExt = "MPI" + DiffEqBaseMeasurementsExt = "Measurements" + DiffEqBaseMonteCarloMeasurementsExt = "MonteCarloMeasurements" + DiffEqBaseReverseDiffExt = "ReverseDiff" + DiffEqBaseTrackerExt = "Tracker" + DiffEqBaseUnitfulExt = "Unitful" + DiffEqBaseZygoteExt = "Zygote" + + [deps.DiffEqBase.weakdeps] + Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f" + GeneralizedGenerated = "6b9d7cbe-bcb9-11e9-073f-15a7a543e2eb" + MPI = "da04e1cc-30fd-572f-bb4f-1f8673147195" + Measurements = "eff96d63-e80a-5855-80a2-b1b0885c5ab7" + MonteCarloMeasurements = "0987c9cc-fe09-11e8-30f0-b96dd679fdca" + ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" + Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" + Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d" + Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" + [[deps.DiffEqCallbacks]] deps = ["DataStructures", "DiffEqBase", "ForwardDiff", "LinearAlgebra", "Markdown", "NLsolve", "Parameters", "RecipesBase", "RecursiveArrayTools", "SciMLBase", "StaticArraysCore"] git-tree-sha1 = "63b6be7b396ad395825f3cc48c56b53bfaf7e69d" @@ -221,6 +285,12 @@ git-tree-sha1 = "2c4ed3eedb87579bfe9f20ecc2440de06b9f3b89" uuid = "77a26b50-5914-5dd7-bc55-306e6241c503" version = "5.16.0" + [deps.DiffEqNoiseProcess.extensions] + DiffEqNoiseProcessReverseDiffExt = "ReverseDiff" + + [deps.DiffEqNoiseProcess.weakdeps] + ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" + [[deps.DiffResults]] deps = ["StaticArraysCore"] git-tree-sha1 = "782dd5f4561f5d267313f23853baaaa4c52ea621" @@ -244,10 +314,15 @@ deps = ["Random", "Serialization", "Sockets"] uuid = "8ba89e20-285c-5b6f-9357-94700520ee1b" [[deps.Distributions]] -deps = ["ChainRulesCore", "DensityInterface", "FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SparseArrays", "SpecialFunctions", "Statistics", "StatsBase", "StatsFuns", "Test"] +deps = ["FillArrays", "LinearAlgebra", "PDMats", "Printf", "QuadGK", "Random", "SparseArrays", "SpecialFunctions", "Statistics", "StatsBase", "StatsFuns", "Test"] git-tree-sha1 = "da9e1a9058f8d3eec3a8c9fe4faacfb89180066b" uuid = "31c24e10-a181-5473-b8eb-7969acd0382f" version = "0.25.86" +weakdeps = ["ChainRulesCore", "DensityInterface"] + + [deps.Distributions.extensions] + DistributionsChainRulesCoreExt = "ChainRulesCore" + DistributionsDensityInterfaceExt = "DensityInterface" [[deps.DocStringExtensions]] deps = ["LibGit2"] @@ -351,10 +426,14 @@ uuid = "53c48c17-4a7d-5ca2-90c5-79b7896eea93" version = "0.8.4" [[deps.ForwardDiff]] -deps = ["CommonSubexpressions", "DiffResults", "DiffRules", "LinearAlgebra", "LogExpFunctions", "NaNMath", "Preferences", "Printf", "Random", "SpecialFunctions", "StaticArrays"] +deps = ["CommonSubexpressions", "DiffResults", "DiffRules", "LinearAlgebra", "LogExpFunctions", "NaNMath", "Preferences", "Printf", "Random", "SpecialFunctions"] git-tree-sha1 = "00e252f4d706b3d55a8863432e742bf5717b498d" uuid = "f6369f11-7733-5829-9624-2563aa707210" version = "0.10.35" +weakdeps = ["StaticArrays"] + + [deps.ForwardDiff.extensions] + ForwardDiffStaticArraysExt = "StaticArrays" [[deps.FunctionWrappers]] git-tree-sha1 = "d62485945ce5ae9c0c48f124a84998d755bae00e" @@ -488,6 +567,12 @@ git-tree-sha1 = "8d928db71efdc942f10e751564e6bbea1e600dfe" uuid = "7d188eb4-7ad8-530c-ae41-71a32a6d4692" version = "1.0.1" +[[deps.JuliaInterpreter]] +deps = ["CodeTracking", "InteractiveUtils", "Random", "UUIDs"] +git-tree-sha1 = "6a125e6a4cb391e0b9adbd1afa9e771c2179f8ef" +uuid = "aa1ae85d-cabe-5617-a682-6adf51b2e16a" +version = "0.9.23" + [[deps.JumpProcesses]] deps = ["ArrayInterface", "DataStructures", "DiffEqBase", "DocStringExtensions", "FunctionWrappers", "Graphs", "LinearAlgebra", "Markdown", "PoissonRandom", "Random", "RandomNumbers", "RecursiveArrayTools", "Reexport", "SciMLBase", "StaticArrays", "TreeViews", "UnPack"] git-tree-sha1 = "740c685ba3d7f218663436b2152041563c19db6e" @@ -579,7 +664,7 @@ uuid = "d3d80556-e9d4-5f37-9878-2ab0fcc64255" version = "7.2.0" [[deps.LinearAlgebra]] -deps = ["Libdl", "libblastrampoline_jll"] +deps = ["Libdl", "OpenBLAS_jll", "libblastrampoline_jll"] uuid = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" [[deps.LinearMaps]] @@ -594,11 +679,23 @@ git-tree-sha1 = "1d3e720d603557d697fedc036bd1af43fe7b3474" uuid = "7ed4a6bd-45f5-4d41-b270-4a48e9bafcae" version = "1.41.1" + [deps.LinearSolve.extensions] + LinearSolveHYPRE = "HYPRE" + + [deps.LinearSolve.weakdeps] + HYPRE = "b5ffcf37-a2bd-41ab-a3da-4bd9bc8ad771" + [[deps.LogExpFunctions]] -deps = ["ChainRulesCore", "ChangesOfVariables", "DocStringExtensions", "InverseFunctions", "IrrationalConstants", "LinearAlgebra"] +deps = ["DocStringExtensions", "IrrationalConstants", "LinearAlgebra"] git-tree-sha1 = "0a1b7c2863e44523180fdb3146534e265a91870b" uuid = "2ab3a3ac-af41-5b50-aa03-7779005ae688" version = "0.3.23" +weakdeps = ["ChainRulesCore", "ChangesOfVariables", "InverseFunctions"] + + [deps.LogExpFunctions.extensions] + LogExpFunctionsChainRulesCoreExt = "ChainRulesCore" + LogExpFunctionsChangesOfVariablesExt = "ChangesOfVariables" + LogExpFunctionsInverseFunctionsExt = "InverseFunctions" [[deps.Logging]] uuid = "56ddb016-857b-54e1-b83d-db4d58db5568" @@ -610,10 +707,21 @@ uuid = "e6f89c97-d47a-5376-807f-9c37f3926c36" version = "1.0.0" [[deps.LoopVectorization]] -deps = ["ArrayInterface", "ArrayInterfaceCore", "CPUSummary", "ChainRulesCore", "CloseOpenIntervals", "DocStringExtensions", "ForwardDiff", "HostCPUFeatures", "IfElse", "LayoutPointers", "LinearAlgebra", "OffsetArrays", "PolyesterWeave", "SIMDTypes", "SLEEFPirates", "SnoopPrecompile", "SpecialFunctions", "Static", "StaticArrayInterface", "ThreadingUtilities", "UnPack", "VectorizationBase"] +deps = ["ArrayInterface", "ArrayInterfaceCore", "CPUSummary", "CloseOpenIntervals", "DocStringExtensions", "HostCPUFeatures", "IfElse", "LayoutPointers", "LinearAlgebra", "OffsetArrays", "PolyesterWeave", "SIMDTypes", "SLEEFPirates", "SnoopPrecompile", "Static", "StaticArrayInterface", "ThreadingUtilities", "UnPack", "VectorizationBase"] git-tree-sha1 = "a282dbdbc2860134d6809acd951543ce359bcf15" uuid = "bdcacae8-1622-11e9-2a5c-532679323890" version = "0.12.155" +weakdeps = ["ChainRulesCore", "ForwardDiff", "SpecialFunctions"] + + [deps.LoopVectorization.extensions] + ForwardDiffExt = ["ChainRulesCore", "ForwardDiff"] + SpecialFunctionsExt = "SpecialFunctions" + +[[deps.LoweredCodeUtils]] +deps = ["JuliaInterpreter"] +git-tree-sha1 = "60168780555f3e663c536500aa790b6368adc02a" +uuid = "6f1432cf-f94c-5a45-995e-cdbf5db27b0b" +version = "2.3.0" [[deps.MKL_jll]] deps = ["Artifacts", "IntelOpenMP_jll", "JLLWrappers", "LazyArtifacts", "Libdl", "Pkg"] @@ -645,7 +753,7 @@ version = "1.1.7" [[deps.MbedTLS_jll]] deps = ["Artifacts", "Libdl"] uuid = "c8ffd9c3-330d-5841-b78e-0817d7145fa1" -version = "2.28.0+0" +version = "2.28.2+0" [[deps.Missings]] deps = ["DataAPI"] @@ -658,7 +766,7 @@ uuid = "a63ad114-7e13-5084-954f-fe012c677804" [[deps.MozillaCACerts_jll]] uuid = "14a3606d-f60d-562e-9121-12d972cd8159" -version = "2022.2.1" +version = "2022.10.11" [[deps.MuladdMacro]] git-tree-sha1 = "cac9cc5499c25554cba55cd3c30543cff5ca4fab" @@ -702,7 +810,7 @@ version = "1.12.9" [[deps.OpenBLAS_jll]] deps = ["Artifacts", "CompilerSupportLibraries_jll", "Libdl"] uuid = "4536629a-c528-5b80-bd46-f80d51c5b363" -version = "0.3.20+0" +version = "0.3.21+4" [[deps.OpenLibm_jll]] deps = ["Artifacts", "Libdl"] @@ -763,9 +871,9 @@ uuid = "69de0a69-1ddd-5017-9359-2bf0b02dc9f0" version = "2.5.3" [[deps.Pkg]] -deps = ["Artifacts", "Dates", "Downloads", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"] +deps = ["Artifacts", "Dates", "Downloads", "FileWatching", "LibGit2", "Libdl", "Logging", "Markdown", "Printf", "REPL", "Random", "SHA", "Serialization", "TOML", "Tar", "UUIDs", "p7zip_jll"] uuid = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" -version = "1.8.0" +version = "1.9.0" [[deps.PoissonRandom]] deps = ["Random"] @@ -797,6 +905,18 @@ git-tree-sha1 = "f739b1b3cc7b9949af3b35089931f2b58c289163" uuid = "d236fae5-4411-538c-8e31-a6e3d9e00b46" version = "0.4.12" + [deps.PreallocationTools.extensions] + PreallocationToolsReverseDiffExt = "ReverseDiff" + + [deps.PreallocationTools.weakdeps] + ReverseDiff = "37e2e3b7-166d-5795-8a7a-e32c996b4267" + +[[deps.PrecompileTools]] +deps = ["Preferences"] +git-tree-sha1 = "259e206946c293698122f63e2b513a7c99a244e8" +uuid = "aea7be01-6a6a-4083-8856-8a6e6704d82a" +version = "1.1.1" + [[deps.Preferences]] deps = ["TOML"] git-tree-sha1 = "47e5f437cc0e7ef2ce8406ce1e7e24d44915f88d" @@ -885,6 +1005,16 @@ git-tree-sha1 = "140cddd2c457e4ebb0cdc7c2fd14a7fbfbdf206e" uuid = "731186ca-8d62-57ce-b412-fbd966d074cd" version = "2.38.3" + [deps.RecursiveArrayTools.extensions] + RecursiveArrayToolsMeasurementsExt = "Measurements" + RecursiveArrayToolsTrackerExt = "Tracker" + RecursiveArrayToolsZygoteExt = "Zygote" + + [deps.RecursiveArrayTools.weakdeps] + Measurements = "eff96d63-e80a-5855-80a2-b1b0885c5ab7" + Tracker = "9f7883ad-71c0-57eb-9f7f-b5c9e6d3789c" + Zygote = "e88e6eb3-aa80-5325-afca-941959d7151f" + [[deps.RecursiveFactorization]] deps = ["LinearAlgebra", "LoopVectorization", "Polyester", "SnoopPrecompile", "StrideArraysCore", "TriangularSolve"] git-tree-sha1 = "9088515ad915c99026beb5436d0a09cd8c18163e" @@ -908,6 +1038,12 @@ git-tree-sha1 = "256eeeec186fa7f26f2801732774ccf277f05db9" uuid = "ae5879a3-cd67-5da8-be7f-38c6eb64a37b" version = "1.1.1" +[[deps.Revise]] +deps = ["CodeTracking", "Distributed", "FileWatching", "JuliaInterpreter", "LibGit2", "LoweredCodeUtils", "OrderedCollections", "Pkg", "REPL", "Requires", "UUIDs", "Unicode"] +git-tree-sha1 = "feafdc70b2e6684314e188d95fe66d116de834a7" +uuid = "295af30f-e4ad-537b-8983-00126c2a3abe" +version = "3.5.2" + [[deps.Rmath]] deps = ["Random", "Rmath_jll"] git-tree-sha1 = "f65dcb5fa46aee0cf9ed6274ccbd597adc49aa7b" @@ -983,6 +1119,12 @@ git-tree-sha1 = "54c78ac3cc0343a16785adabe5bbf4063c737967" uuid = "727e6d20-b764-4bd8-a329-72de5adea6c7" version = "0.1.14" + [deps.SimpleNonlinearSolve.extensions] + SimpleBatchedNonlinearSolveExt = "NNlib" + + [deps.SimpleNonlinearSolve.weakdeps] + NNlib = "872c559c-99b0-510c-b3b7-b6c96a88d5cd" + [[deps.SimpleTraits]] deps = ["InteractiveUtils", "MacroTools"] git-tree-sha1 = "5d7e3f4e11935503d3ecaf7186eac40602e7d231" @@ -1010,7 +1152,7 @@ uuid = "a2af1166-a08f-5f64-846c-94a0d3cef48c" version = "1.1.0" [[deps.SparseArrays]] -deps = ["LinearAlgebra", "Random"] +deps = ["Libdl", "LinearAlgebra", "Random", "Serialization", "SuiteSparse_jll"] uuid = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" [[deps.SparseDiffTools]] @@ -1026,10 +1168,14 @@ uuid = "e56a9233-b9d6-4f03-8d0f-1825330902ac" version = "0.3.9" [[deps.SpecialFunctions]] -deps = ["ChainRulesCore", "IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"] +deps = ["IrrationalConstants", "LogExpFunctions", "OpenLibm_jll", "OpenSpecFun_jll"] git-tree-sha1 = "ef28127915f4229c971eb43f3fc075dd3fe91880" uuid = "276daf66-3868-5448-9aa4-cd146d93841b" version = "2.2.0" +weakdeps = ["ChainRulesCore"] + + [deps.SpecialFunctions.extensions] + SpecialFunctionsChainRulesCoreExt = "ChainRulesCore" [[deps.StableRNGs]] deps = ["Random", "Test"] @@ -1048,6 +1194,11 @@ deps = ["ArrayInterface", "Compat", "IfElse", "LinearAlgebra", "Requires", "Snoo git-tree-sha1 = "fd5f417fd7e103c121b0a0b4a6902f03991111f4" uuid = "0d7ed370-da01-4f52-bd93-41d350b8b718" version = "1.3.0" +weakdeps = ["OffsetArrays", "StaticArrays"] + + [deps.StaticArrayInterface.extensions] + StaticArrayInterfaceOffsetArraysExt = "OffsetArrays" + StaticArrayInterfaceStaticArraysExt = "StaticArrays" [[deps.StaticArrays]] deps = ["LinearAlgebra", "Random", "StaticArraysCore", "Statistics"] @@ -1063,6 +1214,7 @@ version = "1.4.0" [[deps.Statistics]] deps = ["LinearAlgebra", "SparseArrays"] uuid = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" +version = "1.9.0" [[deps.StatsAPI]] deps = ["LinearAlgebra"] @@ -1077,10 +1229,15 @@ uuid = "2913bbd2-ae8a-5f71-8c99-4fb6c76f3a91" version = "0.33.21" [[deps.StatsFuns]] -deps = ["ChainRulesCore", "HypergeometricFunctions", "InverseFunctions", "IrrationalConstants", "LogExpFunctions", "Reexport", "Rmath", "SpecialFunctions"] +deps = ["HypergeometricFunctions", "IrrationalConstants", "LogExpFunctions", "Reexport", "Rmath", "SpecialFunctions"] git-tree-sha1 = "f625d686d5a88bcd2b15cd81f18f98186fdc0c9a" uuid = "4c63d2b9-4356-54db-8cca-17b64c39e42c" version = "1.3.0" +weakdeps = ["ChainRulesCore", "InverseFunctions"] + + [deps.StatsFuns.extensions] + StatsFunsChainRulesCoreExt = "ChainRulesCore" + StatsFunsInverseFunctionsExt = "InverseFunctions" [[deps.StochasticDiffEq]] deps = ["Adapt", "ArrayInterface", "DataStructures", "DiffEqBase", "DiffEqNoiseProcess", "DocStringExtensions", "FillArrays", "FiniteDiff", "ForwardDiff", "JumpProcesses", "LevyArea", "LinearAlgebra", "Logging", "MuladdMacro", "NLsolve", "OrdinaryDiffEq", "Random", "RandomNumbers", "RecursiveArrayTools", "Reexport", "SciMLBase", "SparseArrays", "SparseDiffTools", "StaticArrays", "UnPack"] @@ -1113,7 +1270,7 @@ uuid = "4607b0f0-06f3-5cda-b6b1-a6196a1729e9" [[deps.SuiteSparse_jll]] deps = ["Artifacts", "Libdl", "Pkg", "libblastrampoline_jll"] uuid = "bea87d4a-7f5b-5778-9afe-8cc45184846c" -version = "5.10.1+0" +version = "5.10.1+6" [[deps.SymbolicIndexingInterface]] deps = ["DocStringExtensions"] @@ -1124,7 +1281,7 @@ version = "0.2.2" [[deps.TOML]] deps = ["Dates"] uuid = "fa267f1f-6049-4f14-aa54-33bafae1ed76" -version = "1.0.0" +version = "1.0.3" [[deps.TableTraits]] deps = ["IteratorInterfaceExtensions"] @@ -1141,7 +1298,7 @@ version = "1.10.1" [[deps.Tar]] deps = ["ArgTools", "SHA"] uuid = "a4e569a6-e804-4fa4-b0f3-eef7a1d5b13e" -version = "1.10.1" +version = "1.10.0" [[deps.Test]] deps = ["InteractiveUtils", "Logging", "Random", "Serialization"] @@ -1227,7 +1384,7 @@ uuid = "19fa3120-7c27-5ec5-8db8-b0b0aa330d6f" version = "0.2.0" [[deps.WaveguideQED]] -deps = ["LinearAlgebra", "Parameters", "QuantumOptics", "SnoopPrecompile", "SparseArrays", "Strided", "UnsafeArrays"] +deps = ["LinearAlgebra", "Parameters", "PrecompileTools", "QuantumOptics", "QuantumOpticsBase", "SparseArrays", "Strided", "UnsafeArrays"] path = ".." uuid = "c4555495-0e8d-488d-8e6a-965ecefe9055" version = "0.1.0" @@ -1247,7 +1404,7 @@ version = "0.4.8" [[deps.Zlib_jll]] deps = ["Libdl"] uuid = "83775a58-1f1d-513f-b197-d71354ab007a" -version = "1.2.12+3" +version = "1.2.13+0" [[deps.ZygoteRules]] deps = ["ChainRulesCore", "MacroTools"] @@ -1256,9 +1413,9 @@ uuid = "700de1a5-db45-46bc-99cf-38207098b444" version = "0.2.3" [[deps.libblastrampoline_jll]] -deps = ["Artifacts", "Libdl", "OpenBLAS_jll"] +deps = ["Artifacts", "Libdl"] uuid = "8e850b90-86db-534c-a0d3-1478176c7d93" -version = "5.1.1+0" +version = "5.7.0+0" [[deps.nghttp2_jll]] deps = ["Artifacts", "Libdl"] diff --git a/docs/Project.toml b/docs/Project.toml index c150958..727a775 100644 --- a/docs/Project.toml +++ b/docs/Project.toml @@ -3,5 +3,6 @@ Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" DocumenterCitations = "daee34ce-89f3-4625-b898-19384cb65244" PyPlot = "d330b81b-6aea-500a-939a-2ce795aea3ee" QuantumOptics = "6e0679c1-51ea-5a7c-ac74-d61b76210b0c" +Revise = "295af30f-e4ad-537b-8983-00126c2a3abe" StableRNGs = "860ef19b-820b-49d6-a774-d7a799459cd3" WaveguideQED = "c4555495-0e8d-488d-8e6a-965ecefe9055" diff --git a/docs/make.jl b/docs/make.jl index 1937b03..d37a38c 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -24,7 +24,7 @@ pages = [ "multiplewaveguides.md", "Examples" => ["Scattering on two level system" => "example_lodahl.md"], "API" => "API.md", -"References" => "references.md", +"References and suggested readings" => "references.md", ] ) diff --git a/docs/src/detection.md b/docs/src/detection.md index a13a692..440e912 100644 --- a/docs/src/detection.md +++ b/docs/src/detection.md @@ -26,7 +26,7 @@ $$\begin{align*} &= \frac{1}{2} \left ( W_c^\dagger(\xi_a) W_c^\dagger(\xi_b) \ket{0}_c - W_d^\dagger(\xi_a) W_d^\dagger(\xi_b) \ket{0}_d + \int_{t_0}^{t_{end}} \mathrm{d}t \int_{t_0}^{t_{end}} \mathrm{d}t' \left [ \xi_a(t)\xi_b(t') - \xi_a(t')\xi_b(t) \right] \ket{1}_c\ket{1}_d \right) \end{align*}$$ -where we introduced $$W_{c/d}^\dagger(\xi_a) W_{c/d}^\dagger(\xi_b) \ket{0}_{c/d} = int_{t_0}^{t_{end}} \mathrm{d}t \int_{t_0}^{t_{end}} \mathrm{d}t' \xi_a(t)\xi_b(t') w_{c/d}^\dagger(t) w_{c/d}^\dagger(t') \ket{0}_{c/d}$$. $$W_{c/d}^\dagger(\xi_a) W_{c/d}^\dagger(\xi_b) \ket{0}_{c/d}$$ thus corresponds to both photons going into the same direction. It is also evident that if $$\xi_a(t)\xi_b(t') - \xi_a(t')\xi_b(t) \right = 0$$ then we will have no photons in waveguide c and d simultaneously. This condition is exactly fulfilled if the photon in waveguide a is indistinguishable from the photon in waveguide b. This also means that if the photons ARE distinguishable, we will start to see photon occurring in waveguide c and d simultaneously. All this and more can be simulated in the code, and in the next section, we walk through how to set the above example up in the code. +where we introduced $$W_{c/d}^\dagger(\xi_a) W_{c/d}^\dagger(\xi_b) \ket{0}_{c/d} = int_{t_0}^{t_{end}} \mathrm{d}t \int_{t_0}^{t_{end}} \mathrm{d}t' \xi_a(t)\xi_b(t') w_{c/d}^\dagger(t) w_{c/d}^\dagger(t') \ket{0}_{c/d}$$. $$W_{c/d}^\dagger(\xi_a) W_{c/d}^\dagger(\xi_b) \ket{0}_{c/d}$$ thus corresponds to both photons going into the same direction. It is also evident that if $$\xi_a(t)\xi_b(t') - \xi_a(t')\xi_b(t) = 0$$ then we will have no photons in waveguide c and d simultaneously. This condition is exactly fulfilled if the photon in waveguide a is indistinguishable from the photon in waveguide b. This also means that if the photons ARE distinguishable, we will start to see photons occurring in waveguides c and d simultaneously. All this and more can be simulated in the code, and in the next section, we walk through how to set the above example up in the code. ## Beamsplitter and detection in WaveguideQED.jl @@ -38,8 +38,8 @@ using QuantumOptics #hide times = 0:0.1:20 bw = WaveguideBasis(1,times) ξfun(t,σ,t0) = complex(sqrt(2/σ)* (log(2)/pi)^(1/4)*exp(-2*log(2)*(t-t0)^2/σ^2)) -waveguide_a = onephoton(bw,ξfun,times,1,10) -waveguide_b = onephoton(bw,ξfun,times,1,10) +waveguide_a = onephoton(bw,ξfun,1,10) +waveguide_b = onephoton(bw,ξfun,1,10) wa = destroy(bw) wb = destroy(bw) nothing #hide @@ -60,7 +60,7 @@ Dminus = Detector(wa/sqrt(2),-wb/sqrt(2)) nothing #hide ``` -The [`Detector`](@ref) applies the first operator (`wa/sqrt(2)`) to the first `Ket` in LazyTensorKet (`waveguide_a`) and the second operator (L`$\pm$ wb/sqrt(2)`) to the second `Ket` in `LazyTensorKet` (waveguide_b). The probability of detecting a photon in the detectors can then be calculated by: +The [`Detector`](@ref) applies the first operator (`wa/sqrt(2)`) to the first `Ket` in LazyTensorKet (`waveguide_a`) and the second operator ($$\pm $$ `wb/sqrt(2)`) to the second `Ket` in `LazyTensorKet` (waveguide_b). The probability of detecting a photon in the detectors can then be calculated by: ```@repl detection p_plus = Dplus * ψ_total @@ -84,8 +84,8 @@ p_minus_plus = Dminus * Dplus * ψ_total As expected, the resulting probabilities are zero. If we instead displace the photons in time so that one is centered around $t = 5$ and another around $t = 15$ we get: ```@repl detection -waveguide_a = onephoton(bw,ξfun,times,1,5); -waveguide_b = onephoton(bw,ξfun,times,1,15); +waveguide_a = onephoton(bw,ξfun,1,5); +waveguide_b = onephoton(bw,ξfun,1,15); ψ_total = LazyTensorKet(waveguide_a,waveguide_b); p_plus_plus = Dplus * Dplus * ψ_total p_minus_minus = Dminus * Dminus * ψ_total diff --git a/docs/src/example_lodahl.md b/docs/src/example_lodahl.md index 18da5d6..4436681 100644 --- a/docs/src/example_lodahl.md +++ b/docs/src/example_lodahl.md @@ -58,8 +58,8 @@ We can now study how single or two-photon states scatter on the atom. We define ξ₂(t1,t2,σ1,σ2,t0) = ξ₁(t1,σ1,t0) * ξ₁(t2,σ2,t0) w = 1 t0 = 5 -ψ1 = onephoton(bw,1,ξ₁,times,w,t0) ⊗ fockstate(be,0) -ψ2 = twophoton(bw,1,ξ₂,times,w,w,t0) ⊗ fockstate(be,0) +ψ1 = onephoton(bw,1,ξ₁,w,t0) ⊗ fockstate(be,0) +ψ2 = twophoton(bw,1,ξ₂,w,w,t0) ⊗ fockstate(be,0) ψScat1 = waveguide_evolution(times,ψ1,H) ψScat2 = waveguide_evolution(times,ψ2,H) nothing #hide @@ -93,7 +93,7 @@ end nothing #hide ``` -Finally, we can plot the scattered wavefunctions, and we note that this matches Fig. 3 in Ref. ^[1]: +Finally, we can plot the scattered wavefunctions, and we note that this matches Fig. 3 in Ref.[^1]: ```@example lodahl using PyPlot; #hide diff --git a/docs/src/index.md b/docs/src/index.md index 532bf36..bc55ec2 100644 --- a/docs/src/index.md +++ b/docs/src/index.md @@ -2,10 +2,10 @@ WaveguideQED.jl is a package for simulating continuous fockstates in waveguides. It expands on [`QuantumOptics.jl`](https://qojulia.org/) by adding a custom basis, operators, and routines for doing detection. -### Dev docs +## Dev docs Added functionalities: -* [`WaveguideBasis`](@ref) for representing the waveguide space and the related generator functions: [`zerophoton`](@ref), [`onephoton`](@ref), and [`twophoton`](@ref). Also see [`OnePhotonView`](@ref), [`TwoPhotonView`](@ref), and [`plot_twophoton!`](@ref) for viewing the waveguide data for plotting. Note that [`WaveguideBasis`](@ref) can contain multiple waveguides. -* [`WaveguideOperator`](@ref) which are specialized operators allowing efficient annihilation and creation operators at each time-bin in the waveguide. They are created by giving a basis to [`WaveguideQED.destroy`](@ref) and [`WaveguideQED.create`](@ref) +* [`WaveguideBasis`](@ref) for representing the waveguide Hilbert space and the related functions for generating states in this Hilbert space: [`zerophoton`](@ref), [`onephoton`](@ref), and [`twophoton`](@ref). Also see [`OnePhotonView`](@ref), [`TwoPhotonView`](@ref), and [`plot_twophoton!`](@ref) for viewing the waveguide states and plotting them. Note that [`WaveguideBasis`](@ref) can contain multiple waveguides. +* [`WaveguideOperator`](@ref) are specialized operators allowing efficient annihilation and creation operators at each time-bin in the waveguide. They are created by giving a basis to [`WaveguideQED.destroy`](@ref) and [`WaveguideQED.create`](@ref) * Since the interaction between the waveguide time-bin mode $k$ and cavity/emitter is given as: $a^\dagger w_k - a w_k^\dagger$ there are specially optimized functions for doing these operations called [`CavityWaveguideOperator`](@ref) which are created using a fockbasis and a waveguide basis and the functions [`emission`](@ref) and [`absorption`](@ref). * [`Detector`](@ref), [`LazyTensorKet`](@ref), and [`LazySumKet`](@ref), together with [`detect_single_click`](@ref) and [`detect_double_click`](@ref) allow one to do a beamsplitter interference and subsequent detection of photons coming from two waveguides. diff --git a/docs/src/multiplewaveguides.md b/docs/src/multiplewaveguides.md index aa6020c..9e7b9fc 100644 --- a/docs/src/multiplewaveguides.md +++ b/docs/src/multiplewaveguides.md @@ -34,17 +34,17 @@ Similarly, initializing one or two-photon states in the first or second waveguid ```@example multiple ξ(t,σ,t0) = sqrt(2/σ)* (log(2)/pi)^(1/4)*exp(-2*log(2)*(t-t0)^2/σ^2) ξ2(t1,t2,σ,t0) = ξ(t1,σ,t0)*ξ(t2,σ,t0) -ψ_single_1 = onephoton(bw,1,ξ,times,2,5) -ψ_double_1 = twophoton(bw,1,ξ2,times,2,5) -ψ_single_2 = onephoton(bw,2,ξ,times,2,5) -ψ_double_2 = twophoton(bw,2,ξ2,times,2,5) +ψ_single_1 = onephoton(bw,1,ξ,2,5) +ψ_double_1 = twophoton(bw,1,ξ2,2,5) +ψ_single_2 = onephoton(bw,2,ξ,2,5) +ψ_double_2 = twophoton(bw,2,ξ2,2,5) nothing #hide ``` If we want to describe a simultaneous excitation in both waveguides (states like $\ket{1_i}_\mathrm{1}\ket{1_j }_\mathrm{2}$ where the subscript $$\ket{1_i}_\mathrm{i}$ means waveguide i) we specify both indices of the waveguides: ```@example multiple -ψ_single_1_and_2 = twophoton(bw,[1,2],ξ2,times,2,5) +ψ_single_1_and_2 = twophoton(bw,[1,2],ξ2,2,5) nothing #hide ``` diff --git a/docs/src/references.bib b/docs/src/references.bib index 9a10ad9..e03a997 100644 --- a/docs/src/references.bib +++ b/docs/src/references.bib @@ -27,6 +27,147 @@ @article{Heuck2020Photon-photonCavities arxivId = {1905.02134} } +@article{Heuck2020, + abstract = {We show that relatively simple integrated photonic circuits have the potential to realize a high fidelity deterministic controlled-phase gate between photonic qubits using bulk optical nonlinearities. The gate is enabled by converting travelling continuous-mode photons into stationary cavity modes using strong classical control fields that dynamically change the effective cavity-waveguide coupling rate. This architecture succeeds because it reduces the wave packet distortions that otherwise accompany the action of optical nonlinearities [J. Shapiro, Phys. Rev. A 73, 062305 (2006)PLRAAN1050-294710.1103/PhysRevA.73.062305; J. Gea-Banacloche, Phys. Rev. A 81, 043823 (2010)PLRAAN1050-294710.1103/PhysRevA.81.043823]. We show that high-fidelity gates can be achieved with self-phase modulation in χ(3) materials as well as second-harmonic generation in χ(2) materials. The gate fidelity asymptotically approaches unity with increasing storage time for an incident photon wave packet with fixed duration. We also show that dynamically coupled cavities enable a trade-off between errors due to loss and wave packet distortion. Our proposed architecture represents a new approach to practical implementation of quantum gates that is roomerature compatible and only relies on components that have been individually demonstrated.}, + author = {Mikkel Heuck and Kurt Jacobs and Dirk R. Englund}, + doi = {10.1103/PhysRevLett.124.160501}, + issn = {10797114}, + issue = {16}, + journal = {Physical Review Letters}, + month = {4}, + pmid = {32383940}, + publisher = {American Physical Society}, + title = {Controlled-Phase Gate Using Dynamically Coupled Cavities and Optical Nonlinearities}, + volume = {124}, + year = {2020}, +} + +@article{Krastanov2022, + abstract = {
We propose an architecture for achieving high-fidelity deterministic quantum logic gates on dual-rail encoded photonic qubits by letting photons interact with a two-level emitter (TLE) inside an optical cavity. The photon wave packets that define the qubit are preserved after the interaction due to a quantum control process that actively loads and unloads the photons from the cavity and dynamically alters their effective coupling to the TLE. The controls rely on nonlinear wave mixing between cavity modes enhanced by strong externally modulated electromagnetic fields or on AC Stark shifts of the TLE transition energy. We numerically investigate the effect of imperfections in terms of loss and dephasing of the TLE as well as control field miscalibration. Our results suggest that III-V quantum dots in GaAs membranes is a promising platform for photonic quantum information processing.
}, + author = {Stefan Krastanov and Kurt Jacobs and Gerald Gilbert and Dirk R. Englund and Mikkel Heuck}, + doi = {10.1038/s41534-022-00604-5}, + issn = {2056-6387}, + issue = {1}, + journal = {npj Quantum Information}, + month = {9}, + pages = {103}, + title = {Controlled-phase gate by dynamic coupling of photons to a two-level emitter}, + volume = {8}, + url = {https://www.nature.com/articles/s41534-022-00604-5}, + year = {2022}, +} + +@article{Fischer2018, + abstract = {In this work, we discuss connections between different theoretical physics communities and their works, all related to systems that act as sources of particles such as photons, phonons, or electrons. Our interest is to understand how a low-dimensional quantum system driven by coherent fields, e.g., a two-level system, a Jaynes-Cummings system, or a photon pair source driven by a laser pulse, emits photons into a waveguide. Of particular relevance to solid-state sources is that we provide a way to include dissipation into the formalism for temporal-mode quantum optics. We will discuss the connections between temporal-mode quantum optics, scattering matrices, quantum stochastic calculus, continuous matrix product states and operators, and very traditional quantum optical concepts such as the Mandel photon counting formula and the Lindblad form of the quantum-optical master equation. We close with an example of how our formalism relates to single-photon sources with dephasing.}, + author = {Kevin A. Fischer and Rahul Trivedi and Daniil Lukin}, + doi = {10.1103/PhysRevA.98.023853}, + issn = {24699934}, + issue = {2}, + journal = {Physical Review A}, + month = {8}, + publisher = {American Physical Society}, + title = {Particle emission from open quantum systems}, + volume = {98}, + year = {2018}, +} + + +@article{Ciccarello2018, + abstract = {Quantum collision models (CMs) provide advantageous case studies for investigating major issues in open quantum systems theory, and especially quantum non-Markovianity. After reviewing their general definition and distinctive features, we illustrate the emergence of a CM in a familiar quantum optics scenario. This task is carried out by highlighting the close connection between the well-known input-output formalism and CMs. Within this quantum optics framework, usual assumptions in the CMs’ literature - such as considering a bath of noninteracting yet initially correlated ancillas - have a clear physical origin.}, + author = {Francesco Ciccarello}, + doi = {10.1515/qmetro-2017-0007}, + issue = {1}, + journal = {Quantum Measurements and Quantum Metrology}, + month = {1}, + publisher = {Portico}, + title = {Collision models in quantum optics}, + volume = {4}, + year = {2018}, +} + +@article{HughesWQEDMPP2021, + abstract = {We present two different methods for modeling non-Markovian quantum light-matter interactions in waveguide QED systems, using matrix product states (MPSs) and a space-discretized waveguide (SDW) model. After describing the general theory and implementation of both approaches, we compare and contrast these methods directly on three topical problems of interest in waveguide-QED, including (i) a two-level system (TLS) coupled to an infinite (one-dimensional) waveguide, (ii) a TLS coupled to a terminated waveguide with a time-delayed coherent feedback, and (iii) two spatially separated TLSs coupled within an infinite waveguide. Both approaches are shown to efficiently describe multiphoton nonlinear dynamics in highly non-Markovian regimes, and we highlight the advantages and disadvantages of these methods for modeling waveguide QED interactions, including their implementation in python, computational run times, and ease of conceptual understanding. We explore both vacuum dynamics as well as regimes of strong optical pumping, where a weak excitation approximation cannot be applied. The MPS approach scales better when modeling multiphoton dynamics and long delay times and explicitly includes non-Markovian memory effects. In contrast, the SDW model accounts for non-Markovian effects through space discretization and solves Markovian equations of motion, yet rigorously includes the effects of retardation. The SDW model, based on an extension of recent collisional pictures in quantum optics, is solved through quantum trajectory techniques and can more easily add in additional dissipation processes, including off-chip decay and TLS pure dephasing. The impact of these processes is shown directly on feedback-induced population trapping and TLS entanglement between spatially separated TLSs.}, + author = {Sofia Arranz Regidor and Gavin Crowder and Howard Carmichael and Stephen Hughes}, + doi = {10.1103/PhysRevResearch.3.023030}, + issn = {2643-1564}, + issue = {2}, + journal = {Physical Review Research}, + month = {4}, + pages = {023030}, + publisher = {American Physical Society}, + title = {Modeling quantum light-matter interactions in waveguide QED with retardation, nonlinear interactions, and a time-delayed feedback: Matrix product states versus a space-discretized waveguide model}, + volume = {3}, + url = {https://link.aps.org/doi/10.1103/PhysRevResearch.3.023030}, + year = {2021}, +} + + +@article{Kiilerich2020, + abstract = {This paper presents a general master equation formalism for the interaction between traveling pulses of quantum radiation and localized quantum systems. Traveling fields populate a continuum of free-space radiation modes and the Jaynes-Cummings model, valid for a discrete eigenmode of a cavity, does not apply. We develop a complete input-output theory to describe the driving of quantum systems by arbitrary incident pulses of radiation and the quantum state of the field emitted into any desired outgoing temporal mode. Our theory is applicable to the transformation and interaction of pulses of radiation by their coupling to a wide class of material quantum systems. We discuss the most essential differences between quantum interactions with pulses and with discrete radiative eigenmodes and present examples relevant to quantum information protocols with optical, microwave, and acoustic waves.}, + author = {Alexander Holm Kiilerich and Klaus Mølmer}, + doi = {10.1103/PhysRevA.102.023717}, + issn = {2469-9926}, + issue = {2}, + journal = {Physical Review A}, + month = {8}, + pages = {023717}, + publisher = {American Physical Society}, + title = {Quantum interactions with pulses of radiation}, + volume = {102}, + url = {https://link.aps.org/doi/10.1103/PhysRevA.102.023717}, + year = {2020}, +} + +@article{Kiilerich2019, + abstract = {We present a formalism that accounts for the interaction of a local quantum system, such as an atom or a cavity, with traveling pulses of quantized radiation. We assume Markovian coupling of the stationary system to the input and output fields and nondispersive asymptotic propagation of the pulses before and after the interaction. This permits derivation of a master equation where the input and output pulses are treated as single oscillator modes that both couple to the local system in a cascaded manner. As examples of our theory, we analyze reflection by an empty cavity with phase noise, stimulated atomic emission by a quantum light pulse, and formation of a Schrödinger-cat state by the dispersive interaction of a coherent pulse and a single atom in a cavity.}, + author = {Alexander Holm Kiilerich and Klaus Mølmer}, + doi = {10.1103/PhysRevLett.123.123604}, + issn = {0031-9007}, + issue = {12}, + journal = {Physical Review Letters}, + month = {9}, + pages = {123604}, + publisher = {American Physical Society}, + title = {Input-Output Theory with Quantum Pulses}, + volume = {123}, + url = {https://link.aps.org/doi/10.1103/PhysRevLett.123.123604}, + year = {2019}, +} + + +@article{Yang2022, + abstract = {Sorting quantum fields into different modes according to their Fock-space quantum numbers is a highly desirable quantum operation. In this Letter, we show that a pair of two-level emitters, chirally coupled to a waveguide, may scatter single- and two-photon components of an input pulse into orthogonal temporal modes with a fidelity ≳0.9997. We develop a general theory to characterize and optimize this process and reveal that such a high fidelity is enabled by an interesting two-photon scattering dynamics: while the first emitter gives rise to a complex multimode field, the second emitter recombines the field amplitudes, and the net two-photon scattering induces a self-time reversal of the input pulse mode. The presented scheme can be employed to construct logic elements for propagating photons, such as a deterministic nonlinear-sign gate with a fidelity ≳0.9995.}, + author = {Fan Yang and Mads M. Lund and Thomas Pohl and Peter Lodahl and Klaus Mølmer}, + doi = {10.1103/PhysRevLett.128.213603}, + issn = {10797114}, + issue = {21}, + journal = {Physical Review Letters}, + month = {5}, + pmid = {35687472}, + publisher = {American Physical Society}, + title = {Deterministic Photon Sorting in Waveguide QED Systems}, + volume = {128}, + year = {2022}, +} + +@article{Christiansen2023, + abstract = {The interaction of a propagating pulse of quantum radiation with a localized quantum system can be described by a cascaded master equation with a distinct initially populated input and a finally populated output field mode [Kiilerich and Mølmer, Phys. Rev. Lett. 123, 123604 (2019)0031-900710.1103/PhysRevLett.123.123604]. By transformation to an appropriate interaction picture, we break the cascaded nature of the master equation and recover an effective time-dependent interaction with a lossless single mode and a supplementary lossy mode. The former closely represents the traveling pulse, while the latter constitutes a non-Markovian component in the exchange of quanta between the scatterer and the quantized field. The transformed master equation offers important insights into the system dynamics, and it permits numerically efficient solutions.}, + author = {Victor Rueskov Christiansen and Alexander Holm Kiilerich and Klaus Mølmer}, + doi = {10.1103/PhysRevA.107.013706}, + issn = {2469-9926}, + issue = {1}, + journal = {Physical Review A}, + month = {1}, + pages = {013706}, + publisher = {American Physical Society}, + title = {Interactions of quantum systems with pulses of quantized radiation: From a cascaded master equation to a traveling mode perspective}, + volume = {107}, + url = {https://link.aps.org/doi/10.1103/PhysRevA.107.013706}, + year = {2023}, +} + + + @book{Gerry2004, author = {Christopher Gerry and Peter Knight}, doi = {10.1017/CBO9780511791239}, diff --git a/docs/src/references.md b/docs/src/references.md index 78f29bc..f1afe10 100644 --- a/docs/src/references.md +++ b/docs/src/references.md @@ -1,4 +1,18 @@ -# References +# Suggested readings + +### Theory and background +The theory of time-bin continuous fockstates is introduced in [Heuck2020Photon-photonCavities](@cite) where it is used to derive the equations of motion for Waveguide QED systems. The numerical method in this library is heavily based on this approach, where instead of deriving the equations, the systems are solved by applying the operators themselves. The time-bin method is also used in [Heuck2020](@cite) and [Krastanov2022](@cite). + +### Similar approaches +The following is a list of approaches that are trying to solve problems that can also be treated with this library. +* [HughesWQEDMPP2021](@cite) considers feedback in waveguide systems and uses a space-discretized waveguide picture with Monte Carlo trajectories +* [Fischer2018](@cite) relates many approaches to solving WaveguideQED problems with each other and also introduces a framework that aims to deal with similar problems through master equations, photon counting and tracing out the waveguide. +* The SLH formalism introduced in [Kiilerich2019](@cite) and [Kiilerich2019](@cite) uses cascaded cavities to simulate quantum pulses. Further work also includes: [Yang2022](@cite) [Christiansen2023](@cite) + +### Papers where we reproduce results from +* The theoretical results in [DynamicalPhotonLodahl2022](@cite) are reproduced in [Scattering on a two-level system](https://qojulia.github.io/WaveguideQED.jl/dev/example_lodahl/) + +# References ```@bibliography ``` diff --git a/docs/src/theoreticalbackground.md b/docs/src/theoreticalbackground.md index 4a4d0b2..9f0578b 100644 --- a/docs/src/theoreticalbackground.md +++ b/docs/src/theoreticalbackground.md @@ -2,7 +2,7 @@ In this section, we go over the necessary theory to work with continuous fockstates in the **WaveguideQED.jl** -## continuous Fock States +## Continuous Fock States The single photon continuous fock state can be defined as: @@ -10,7 +10,7 @@ $$\begin{equation*} \ket{\psi} = W^\dagger(\xi) \ket{0} = \int_{t_0}^{t_{end}} \mathrm{d}t \ \xi(t) w^\dagger(t) \ket{\emptyset} \end{equation*}$$ -here $W^\dagger(\xi)$ creates a photon with the wavefunction $\xi(t)$. $w^\dagger(t)$ is the creation operator for a photon at time $t$, and it obeys the commutation relation: $\left[w(t),w(t')\right ] = \delta(t-t')$. The probability of observing a photon at time $t$ is given by: $\bra{0} w(t) \ket{\psi} = |\xi(t)|^2$. The wavefunction $\xi(t)$ thus describes the temporal distribution of the photon. +here $W^\dagger(\xi)$ creates a photon with the wavefunction $\xi(t)$. $w^\dagger(t)$ is the creation operator for a photon at time $t$, and it obeys the commutation relation: $\left[w(t),w(t')\right ] = \delta(t-t')$. The probability of observing a photon at time $t$ is given by: $\bra{\psi} w^\dagger(t) w(t) \ket{\psi} = |\xi^{(1)}(t)|^2$. The interpretation of the wavefunction $\xi^{(1)}(t)$. The wavefunction $\xi(t)$ thus describes the temporal distribution of the photon. The heart of the photon time-binning is discretizing the continuous fock state into time-bins of width $\Delta t$. The interaction with the emitter/cavity is then assumed to span only one time-bin at a time, corresponding to a spectrally flat interaction between the waveguide and emitter/cavity. We thus discretize the annihilation and creation operators by taking[^1]: @@ -18,14 +18,14 @@ $$\begin{equation*} w(t_k) = w(k \Delta t) \rightarrow \frac{w_k}{\sqrt{\Delta t}} \ \ \ \text{with} \ \left[ w_j, w_k^\dagger \right ] = \delta_{jk} \end{equation*}$$ -where $w_k$ is the descritized operator and the factor of $1/\sqrt{\Delta t}$ assures the commutator relation in the limit of $\Delta t \rightarrow 0$. We denote the action of the discretized creation operator as: $w_k\dagger \ket{\emptyset} = \ket{1_k}$ meaning a single photon in time-bin $k$. This means that the single photon continuous fock state becomes: +where $w_k$ is the descritized operator and the factor of $1/\sqrt{\Delta t}$ assures the commutator relation in the limit of $\Delta t \rightarrow 0$. We denote the action of the discretized creation operator as: $w_k^\dagger \ket{\emptyset} = \ket{1_k}$ meaning a single photon in time-bin $k$. This means that the single photon continuous fock state becomes: $$\begin{equation*} \ket{\psi} = \int_{t_0}^{t_{end}} \mathrm{d}t \ \xi(t) w^\dagger(t) \ket{\emptyset} \rightarrow \sum_{k=1}^N \sqrt{\Delta t} \xi(t_k) w_k^\dagger \ket{\emptyset} \end{equation*}$$ -In `WaveguideQED.jl, the time-bins above are represented as elements in arrays corresponding to each time-bin: +In `WaveguideQED.jl`, the time-bins above are represented as elements in arrays corresponding to each time-bin: ![Alt text](./illustrations/onephoton_array.png) @@ -38,12 +38,12 @@ bw = WaveguideBasis(1,times) nothing #hide ``` -Notice that the input for WaveguideBasis is `1` and `times`. `1` denotes the maximum excitation number of fockstates (currently can only be 1 or 2), and `times` the is the time interval over which the continuous fock state is defined. To define the continuous fockstate, we need to define a wavefunction $\xi$. In the following, we define a Gaussian wavefunction located around $t=5$ with a width of $\sigma = 1$: +Notice that the input for WaveguideBasis is `1` and `times`. `1` denotes the maximum excitation number of fockstates (currently can only be 1 or 2), and `times` the is the time interval over which the continuous fockstate is defined. To define the continuous fockstate, we need to define a wavefunction $\xi$. In the following, we define a Gaussian wavefunction located around $t=5$ with a width of $\sigma = 1$: ```@example theory ξ(t,σ,t0) = sqrt(2/σ)* (log(2)/pi)^(1/4)*exp(-2*log(2)*(t-t0)^2/σ^2) σ,t0 = 1,5 -ψ = onephoton(bw,ξ,times,σ,t0) +ψ = onephoton(bw,ξ,σ,t0) nothing #hide ``` @@ -99,7 +99,7 @@ This is also seen if we plot the creation operator acting on the vacuum: ψ = wd*zerophoton(bw) viewed_state = OnePhotonView(ψ) fig,ax = subplots(1,1,figsize=(9,4.5)) -ax.plot(times,viewed_state,"r-"); +ax.plot(times,real.(viewed_state),"r-"); ax.set_xlabel("Time [a.u]") ax.set_ylabel(L"$\xi(t)$") plt.tight_layout() @@ -145,15 +145,15 @@ We can define a two-photon basis and corresponding operator by: ```@example theory bw = WaveguideBasis(2,times) +w = destroy(bw) wd = create(bw) -println(pwd()) nothing #hide ``` The creation operator can then be visualized by acting on [`onephoton`](@ref) filled with ones. This is seen in the following. Note that the state is visualized as a contour plot mirrored around the diagonal. ```@example theory set_waveguidetimeindex!(wd,50) -psi_plot = wd*onephoton(bw,x->1,times) +psi_plot = wd*onephoton(bw,x->1) fig,ax = subplots(1,1,figsize=(4.5,4.5)) plot_twophoton!(ax,psi_plot,times) plt.savefig("twophoton_created.svg") #hide @@ -168,7 +168,7 @@ If we want to create a two-photon Gaussian state, we instead do: ξ(t,σ,t0) = sqrt(2/σ)* (log(2)/pi)^(1/4)*exp(-2*log(2)*(t-t0)^2/σ^2) ξ2(t1,t2,σ,t0) = ξ(t1,σ,t0)*ξ(t2,σ,t0) σ,t0 = 1,5 -ψ = twophoton(bw,ξ2,times,σ,t0) +ψ = twophoton(bw,ξ2,σ,t0) nothing #hide ``` diff --git a/docs/src/tutorial.md b/docs/src/tutorial.md index a5932cc..ff465ad 100644 --- a/docs/src/tutorial.md +++ b/docs/src/tutorial.md @@ -18,7 +18,7 @@ bc = FockBasis(1) nothing #hide ``` -Next, we want to create the Hamiltonian for the system. The interaction between the waveguide and cavity is at timestep k given by[^1] $$H_k = i \hbar \sqrt{\gamma / \Delta t}( a^\dagger w_k - a w_k^\dagger)$$, where $$a$$ ($$a^\dagger$$) is the cavity annihilation (creation) operator, $$w_k$$($$w_k^\dagger$$) is the waveguide annihilation (creation) operator, $$\gamma$$ is the leakage rate of the cavity, and `\Delta t = times[2]-times[1]` is the width of the time-bin. `WaveguideQED.jl` follows the same syntax as [`QuantumOptics.jl`](https://qojulia.org/), and operators are defined from a basis. Operators of different Hilbert spaces are then combined using ⊗ (``\otimes``) or `tensor`: +Next, we want to create the Hamiltonian for the system. The interaction between the waveguide and cavity is at timestep k given by[^1] $$H_k = i \hbar \sqrt{\gamma / \Delta t}( a^\dagger w_k - a w_k^\dagger)$$, where $$a$$ ($$a^\dagger$$) is the cavity annihilation (creation) operator, $$w_k$$($$w_k^\dagger$$) is the waveguide annihilation (creation) operator, $$\gamma$$ is the leakage rate of the cavity, and $$\Delta t = \mathrm{times[2]}-\mathrm{times[1]}$$ is the width of the time-bin. `WaveguideQED.jl` follows the same syntax as [`QuantumOptics.jl`](https://qojulia.org/), and operators are defined from a basis. Operators of different Hilbert spaces are then combined using ⊗ (`\otimes`) or `tensor`: ```@example tutorial a = destroy(bc) @@ -36,7 +36,7 @@ With this, we can now simulate the scattering of a single photon with a Gaussian ```@example tutorial ξ(t,σ,t0) = sqrt(2/σ)* (log(2)/pi)^(1/4)*exp(-2*log(2)*(t-t0)^2/σ^2) σ,t0 = 1,5 -ψ_waveguide = onephoton(bw,ξ,times,σ,t0) +ψ_waveguide = onephoton(bw,ξ,σ,t0) nothing #hide ``` @@ -77,7 +77,7 @@ nothing #hide ``` ![alt text](scat_onephoton.svg) -We see that the wavefunction has changed after the interaction with the cavity. More specifically, we see how the pulse gets absorbed into the cavity leading and a corresponding phase change of the wave. This phase change also leads to destructive interference between the photon being emitted from the cavity and the reflection of the incoming photon. This leads to a dip in the photon wavefunction after the interaction. +We see that the wavefunction has changed after the interaction with the cavity. More specifically, we see how the pulse gets absorbed into the cavity leading to a phase change in the wavefunction. This phase change also leads to destructive interference between the photon being emitted from the cavity and the reflection of the incoming photon. This leads to a dip in the photon wavefunction after the interaction. ## Expectation values diff --git a/src/CavityWaveguideOperator.jl b/src/CavityWaveguideOperator.jl index a81968f..368f68b 100644 --- a/src/CavityWaveguideOperator.jl +++ b/src/CavityWaveguideOperator.jl @@ -3,6 +3,7 @@ Abstract type used for operators on acting on a combined [`WaveguideBasis`](@ref """ abstract type CavityWaveguideOperator{BL,BR} <: AbstractOperator{BL,BR} end +#TO DO: CLEAN UP #Indexing structure used to loop over state. Needs cleanup and can possible be removed with the addition of eachslice(a,dims=(1,3,...)) in julia 1.9 struct WaveguideIndexing{N} ndims::Int @@ -17,7 +18,6 @@ struct WaveguideIndexing{N} new{length(iter)}(length(strides),k_idx,end_idx,strides,idx_vec1,idx_vec2,range_idx,iter) end end - function WaveguideIndexing(b::Basis,loc) dims = b.shape alldims = [1:length(dims)...] @@ -109,14 +109,18 @@ end function QuantumOpticsBase.:+(a::CavityWaveguideOperator,b::CavityWaveguideOperator) @assert a.basis_l == b.basis_l @assert a.basis_r == b.basis_r - LazySum(a,b) + LazySum(a) +LazySum(b) end function QuantumOpticsBase.:-(a::CavityWaveguideOperator,b::CavityWaveguideOperator) @assert a.basis_l == b.basis_l @assert a.basis_r == b.basis_r - LazySum([1,-1],[a,b]) + LazySum(a) - LazySum(b) +end +function QuantumOpticsBase.:-(a::CavityWaveguideOperator) + out = copy(a) + out.factor = -a.factor + return out end - """ absorption(b1::WaveguideBasis{T},b2::FockBasis) where T @@ -225,10 +229,10 @@ end Return identityoperator(a.basis_l). """ -function QuantumOptics.identityoperator(a::CavityWaveguideOperator) +function QuantumOpticsBase.identityoperator(a::CavityWaveguideOperator) identityoperator(a.basis_l) end -function QuantumOptics.identityoperator(::Type{T}, b1::Basis, b2::Basis) where {T<:CavityWaveguideOperator} +function QuantumOpticsBase.identityoperator(::Type{T}, b1::Basis, b2::Basis) where {T<:CavityWaveguideOperator} @assert b1==b2 identityoperator(b1) end @@ -339,7 +343,7 @@ function is_destroy(data,basis::CompositeBasis) for k = 1:N ind .= 0 ind[k] = 1 - if isequal(data,tensor([i == 0 ? identityoperator(basis.bases[j]) : destroy(basis.bases[j]) for (j,i) in enumerate(ind)]...).data) + if isequal(data,tensor([ (i==0 || !isa(basis.bases[j],FockBasis)) ? identityoperator(basis.bases[j]) : destroy(basis.bases[j]) for (j,i) in enumerate(ind)]...).data) return k end end @@ -351,7 +355,7 @@ function is_create(data,basis::CompositeBasis) for k = 1:N ind .= 0 ind[k] = 1 - if isequal(data,tensor([i == 0 ? identityoperator(basis.bases[j]) : create(basis.bases[j]) for (j,i) in enumerate(ind)]...).data) + if isequal(data,tensor([(i==0 || !isa(basis.bases[j],FockBasis)) ? identityoperator(basis.bases[j]) : create(basis.bases[j]) for (j,i) in enumerate(ind)]...).data) return k end end @@ -409,7 +413,8 @@ function tensor(a::Operator,b::T) where {T<:WaveguideOperator} end - +#TO DO: CLEAN UP +#Currently indexing is very complicated using the CavityIndexing structure and could possible be done smoother and faster. """ mul!(result::Ket{B1}, a::CavityWaveguideEmission, b::Ket{B2}, alpha, beta) where {B1<:Basis,B2<:Basis} mul!(result::Ket{B1}, a::CavityWaveguideAbsorption, b::Ket{B2}, alpha, beta) where {B1<:Basis,B2<:Basis} @@ -419,38 +424,41 @@ Fast in-place multiplication of operators/state vectors. Updates `result` as `re function mul!(result::Ket{B1}, a::CavityWaveguideEmission, b::Ket{B2}, alpha, beta) where {B1<:Basis,B2<:Basis} dims = basis(result).shape i,j = a.loc[1],a.loc[2] + beta1 = beta + if iszero(beta) + fill!(result.data,beta) + beta1 = one(beta) + end if length(dims) == 2 - loop_destroy_ax!(result.data,a,b.data,alpha,beta,a.indexing) - a.indexing.idx_vec1[j] = dims[j] - li1 = linear_index(a.indexing.ndims,a.indexing.idx_vec1,a.indexing.strides) - rmul!(view(result.data,li1:a.indexing.strides[a.loc[1]]:li1+(a.indexing.end_idx-1)*a.indexing.strides[a.loc[1]]),beta) + loop_destroy_ax!(result.data,a,b.data,alpha,beta1,a.indexing) elseif length(dims) == 3 - loop_destroy_third_axis!(result.data,a,b.data,alpha,beta,a.indexing,1) + loop_destroy_third_axis!(result.data,a,b.data,alpha,beta1,a.indexing,1) a.indexing.idx_vec1[j] = dims[j] - loop_rmul_axis!(result.data,a,b.data,alpha,beta,a.indexing,1) + loop_rmul_axis!(result.data,a,b.data,alpha,beta1,a.indexing,1) else - iterate_over_iter!(result.data,a,b.data,alpha,beta,a.indexing,1,loop_destroy_third_axis!) + iterate_over_iter!(result.data,a,b.data,alpha,beta1,a.indexing,1,loop_destroy_third_axis!) a.indexing.idx_vec1[j] = dims[j] - iterate_over_iter!(result.data,a,b.data,alpha,beta,a.indexing,1,loop_rmul_axis!) + iterate_over_iter!(result.data,a,b.data,alpha,beta1,a.indexing,1,loop_rmul_axis!) end return result end function mul!(result::Ket{B1}, a::CavityWaveguideAbsorption, b::Ket{B2}, alpha, beta) where {B1<:Basis,B2<:Basis} dims = basis(result).shape i,j = a.loc[1],a.loc[2] + beta1 = beta + if iszero(beta) + fill!(result.data,beta) + beta1 = one(beta) + end if length(dims) == 2 - loop_create_ax!(result.data,a,b.data,alpha,beta,a.indexing) - a.indexing.idx_vec1[j] = 1 - li1 = linear_index(a.indexing.ndims,a.indexing.idx_vec1,a.indexing.strides) - rmul!(view(result.data,li1:a.indexing.strides[a.loc[1]]:li1+(a.indexing.end_idx-1)*a.indexing.strides[a.loc[1]]),beta) + loop_create_ax!(result.data,a,b.data,alpha,beta1,a.indexing) elseif length(dims) == 3 - loop_create_third_axis!(result.data,a,b.data,alpha,beta,a.indexing,1) + loop_create_third_axis!(result.data,a,b.data,alpha,beta1,a.indexing,1) a.indexing.idx_vec1[j] = 1 - loop_rmul_axis!(result.data,a,b.data,alpha,beta,a.indexing,1) else - iterate_over_iter!(result.data,a,b.data,alpha,beta,a.indexing,1,loop_create_third_axis!) + iterate_over_iter!(result.data,a,b.data,alpha,beta1,a.indexing,1,loop_create_third_axis!) a.indexing.idx_vec1[j] = 1 - iterate_over_iter!(result.data,a,b.data,alpha,beta,a.indexing,1,loop_rmul_axis!) + iterate_over_iter!(result.data,a,b.data,alpha,beta1,a.indexing,1,loop_rmul_axis!) end return result end diff --git a/src/WaveguideInteraction.jl b/src/WaveguideInteraction.jl index 700399b..9c95a59 100644 --- a/src/WaveguideInteraction.jl +++ b/src/WaveguideInteraction.jl @@ -1,3 +1,38 @@ +#= +#Indexing structure used to loop over state. Needs cleanup and can possible be removed with the addition of eachslice(a,dims=(1,3,...)) in julia 1.9 +struct WaveguideStride + looplength::Int + arraylength::Int + s1::Int + s2::Int +end + +function WaveguideStride(b::Basis,loc) + dims = b.shape + alldims = [1:length(dims)...] + i = loc[1] + exclude_dims = [i] + otherdims = setdiff(alldims, exclude_dims) + looplength = prod(dims[otherdims]) + arraylength = dims[i] + + s1 = 1 + s2 = 1 + for k in 1:i-1 + s1 *= dims[k] + end + for k in 1:i + s2 *= dims[k] + end + return WaveguideStride(looplength,arraylength,s1,s2) +end + st::WaveguideStride + function WaveguideInteraction(basis_l::B1,basis_r::B2,factor::ComplexF64,op1::AbstractOperator,op2::AbstractOperator,loc) where{B1,B2} + new{B1,B2}(basis_l,basis_r,factor,op1,op2,loc,WaveguideStride(basis_l,loc)) + end + +=# + mutable struct WaveguideInteraction{B1,B2} <: AbstractOperator{B1,B2} basis_l::B1 basis_r::B2 @@ -30,12 +65,12 @@ Base.:*(b::WaveguideInteraction,a::Number)=*(a,b) function QuantumOpticsBase.:+(a::WaveguideInteraction,b::WaveguideInteraction) @assert a.basis_l == b.basis_l @assert a.basis_r == b.basis_r - LazySum(a,b) + LazySum(a) + LazySum(b) end function QuantumOpticsBase.:-(a::WaveguideInteraction,b::WaveguideInteraction) @assert a.basis_l == b.basis_l @assert a.basis_r == b.basis_r - LazySum([1,-1],[a,b]) + LazySum(a) - LazySum(b) end function set_waveguidetimeindex!(op::WaveguideInteraction,timeindex) @@ -67,7 +102,6 @@ end function mul!(result::Ket{B1}, a::WaveguideInteraction{B1,B2}, b::Ket{B2}, alpha, beta) where {B1<:Basis,B2<:Basis} dims = basis(result).shape - alldims = [1:length(dims)...] i = a.loc[1] exclude_dims = [i] @@ -89,9 +123,22 @@ function mul!(result::Ket{B1}, a::WaveguideInteraction{B1,B2}, b::Ket{B2}, alpha end_idx = dims[i] waveguide_interaction_mul!(view(result.data,li1:stride:li1+(end_idx-1)*stride),a.op1,a.op2,view(b.data,li1:stride:li1+(end_idx-1)*stride),alpha,beta) else - iterate_over_interaction!(result.data,b.data,a,a.factor*alpha,beta,iter,idx_vec1,range_idx,1,stride,dims[i],dims) + iterate_over_interaction!(result.data,b.data,a,alpha,beta,iter,idx_vec1,range_idx,1,stride,dims[i],dims) end + #= + #println(a.st.looplength) + #println(a.st.arraylength) + #println(a.st.s1) + #println(a.st.s2) + for i in 1:a.st.looplength + idx = 1+(i-1)*a.st.s2 + println(idx) + #println(idx+(a.st.arraylength-1)*a.st.s2) + #println(length(view(result.data,idx:a.st.s1:idx+(a.st.arraylength)*a.st.s1-1))) + waveguide_interaction_mul!(view(result.data,idx:a.st.s1:idx+(a.st.arraylength-1)*a.st.s1),a.op1,a.op2,view(b.data,idx:a.st.s1:idx+(a.st.arraylength-1)*a.st.s1),alpha*a.factor,beta) + end + =# return result end @@ -101,7 +148,7 @@ function iterate_over_interaction!(result,b,a,alpha,beta,iter::Tuple,idx_vec1,ra @inbounds idx_vec1[range_idx[idx]] = i li1 = linear_index(dims,idx_vec1) @uviews result b begin - waveguide_interaction_mul!(view(result,li1:stride:li1+(end_idx-1)*stride),a.op1,a.op2,view(b,li1:stride:li1+(end_idx-1)*stride),alpha,beta) + waveguide_interaction_mul!(view(result,li1:stride:li1+(end_idx-1)*stride),a.op1,a.op2,view(b,li1:stride:li1+(end_idx-1)*stride),alpha*a.factor,beta) end end else @@ -116,10 +163,12 @@ function waveguide_interaction_mul!(result,a::WaveguideCreate{B1,B2,1,idx1},bop: if !isone(beta) rmul!(result,beta) end - timeindex = a.timeindex + timeindex_a = a.timeindex + timeindex_b = bop.timeindex nsteps = a.basis_l.nsteps - - @inbounds result[1+(idx1-1)*nsteps+timeindex] += alpha*a.factor*bop.factor*b[1+(idx2-1)*nsteps+timeindex] + if timeindex_a>0 && timeindex_b>0 + @inbounds result[1+(idx1-1)*nsteps+timeindex_a] += alpha*a.factor*bop.factor*b[1+(idx2-1)*nsteps+timeindex_b] + end result end diff --git a/src/WaveguideOperator.jl b/src/WaveguideOperator.jl index 34832ef..a1713a6 100644 --- a/src/WaveguideOperator.jl +++ b/src/WaveguideOperator.jl @@ -153,10 +153,10 @@ end Return identityoperator(a.basis_l). """ -function QuantumOptics.identityoperator(a::WaveguideOperator) +function QuantumOpticsBase.identityoperator(a::WaveguideOperator) identityoperator(a.basis_l) end -function QuantumOptics.identityoperator(::Type{T}, b1::Basis, b2::Basis) where {T<:WaveguideOperator} +function QuantumOpticsBase.identityoperator(::Type{T}, b1::Basis, b2::Basis) where {T<:WaveguideOperator} @assert b1==b2 identityoperator(b1) end @@ -174,13 +174,28 @@ function mul!(result::Ket{B1}, a::LazyTensor{B1,B2,F,I,T}, b::Ket{B2}, alpha, be b_data = Base.ReshapedArray(b.data, QuantumOpticsBase._comp_size(basis(b)), ()) result_data = Base.ReshapedArray(result.data, QuantumOpticsBase._comp_size(basis(result)), ()) - tp_ops = (tuple(( (isa(op,DataOperator) ? op.data : op) for op in a.operators)...), a.indices) - iso_ops = QuantumOpticsBase._explicit_isometries(a.indices, a.basis_l, a.basis_r) + tp_ops = QuantumOpticsBase._tpops_tuple(a) + iso_ops = QuantumOpticsBase._explicit_isometries(eltype(a), a.indices, a.basis_l, a.basis_r) QuantumOpticsBase._tp_sum_matmul!(result_data, tp_ops, iso_ops, b_data, alpha * a.factor, beta) result end +function QuantumOpticsBase._tpops_tuple(a::LazyTensor{B1,B2,F,I,T}; shift=0, op_transform=identity) where {B1<:Basis,B2<:Basis, F,I,T<:Tuple{Vararg{AbstractOperator}}} + length(a.operators) == 0 == length(a.indices) && return () + op_pairs = tuple((((isa(op,DataOperator) ? op_transform(op.data) : op_transform(op)), i + shift) for (op, i) in zip(a.operators, a.indices))...) + + # Filter out identities: + # This induces a non-trivial cost only if _is_square_eye is not inferrable. + # This happens if we have Eyes that are not SquareEyes. + # This can happen if the user constructs LazyTensor operators including + # explicit identityoperator(b,b). + filtered = filter(p->!QuantumOpticsBase._is_square_eye(p[1]), op_pairs) + return filtered +end +QuantumOpticsBase._is_square_eye(a::AbstractOperator) = false + + #Called from _tp_sum_matmul! #Makes sure operator works on correct part of tensor. function QuantumOpticsBase.:_tp_matmul!(result, a::WaveguideOperator, loc::Integer, b, α::Number, β::Number) @@ -283,40 +298,49 @@ end #Destroy 1 waveguide photon function waveguide_mul!(result,a::WaveguideDestroy{B,B,1,idx},b,alpha,beta) where {B,idx} - if !isone(beta) - rmul!(result,beta) + if iszero(beta) + fill!(result, beta) + elseif !isone(beta) + rmul!(result, beta) end - add_zerophoton_onephoton!(result,b,alpha*a.factor,a.timeindex+(idx-1)*a.basis_l.nsteps) + @inbounds result[1] += alpha*a.factor*b[a.timeindex+(idx-1)*a.basis_l.nsteps+1] return end #Destroy 2 waveguide photon function waveguide_mul!(result,a::WaveguideDestroy{B,B,2,1},b,alpha,beta) where {B<:SingleWaveguideBasis} - if !isone(beta) - rmul!(result,beta) + if iszero(beta) + fill!(result, beta) + elseif !isone(beta) + rmul!(result, beta) end timeindex = a.timeindex nsteps = a.basis_l.nsteps - add_zerophoton_onephoton!(result,b,alpha*a.factor,timeindex) - twophotonview = TwoPhotonTimestepView(b,timeindex,nsteps,nsteps+1) - axpy!(alpha*a.factor,twophotonview,view(result,2:1:nsteps+1)) + @inbounds result[1] += alpha*a.factor*b[timeindex+1] + #twophotonview = TwoPhotonTimestepView(b,timeindex,nsteps,nsteps+1) + #axpy!(alpha*a.factor,twophotonview,view(result,2:1:nsteps+1)) + twophoton_destroy!(view(result,2:1:nsteps+1),b,alpha*a.factor,timeindex,nsteps,nsteps+1) return end #Destroy 2 waveguide photon function waveguide_mul!(result,a::WaveguideDestroy{B,B,2,idx},b,alpha,beta) where {B<:MultipleWaveguideBasis,idx} - if !isone(beta) - rmul!(result,beta) + if iszero(beta) + fill!(result, beta) + elseif !isone(beta) + rmul!(result, beta) end timeindex = a.timeindex nsteps = a.basis_l.nsteps Nw = get_number_of_waveguides(a.basis_l) - add_zerophoton_onephoton!(result,b,alpha*a.factor,timeindex+(idx-1)*nsteps) - twophotonview = TwoPhotonTimestepView(b,timeindex,nsteps,Nw*nsteps+1+(idx-1)*(nsteps*(nsteps+1))÷2) - axpy!(alpha*a.factor,twophotonview,view(result,2+(idx-1)*nsteps:1:idx*nsteps+1)) - @simd for k in filter(x -> x != idx, 1:Nw) + @inbounds result[1] += alpha*a.factor*b[timeindex+(idx-1)*a.basis_l.nsteps+1] + #twophotonview = TwoPhotonTimestepView(b,timeindex,nsteps,Nw*nsteps+1+(idx-1)*(nsteps*(nsteps+1))÷2) + #axpy!(alpha*a.factor,twophotonview,view(result,2+(idx-1)*nsteps:1:idx*nsteps+1)) + twophoton_destroy!(view(result,2+(idx-1)*nsteps:1:idx*nsteps+1),b,alpha*a.factor,timeindex,nsteps,Nw*nsteps+1+(idx-1)*(nsteps*(nsteps+1))÷2) + for k in filter(x -> x != idx, 1:Nw) i,j = min(k,idx),max(k,idx) index = (i-1)*Nw + j - (i*(i+1))÷2 - twophotonview = TwoWaveguideTimestepView(b,timeindex,nsteps,1+Nw*nsteps+Nw*(nsteps*(nsteps+1))÷2+(index-1)*nsteps^2,i==idx) - axpy!(alpha*a.factor,twophotonview,view(result,2+(k-1)*nsteps:1:k*nsteps+1)) + #twophotonview = TwoWaveguideTimestepView(b,timeindex,nsteps,1+Nw*nsteps+Nw*(nsteps*(nsteps+1))÷2+(index-1)*nsteps^2,i==idx) + #axpy!(alpha*a.factor,twophotonview,view(result,2+(k-1)*nsteps:1:k*nsteps+1)) + twowaveguide_destroy!(view(result,2+(k-1)*nsteps:1:k*nsteps+1),b,alpha*a.factor,timeindex,nsteps,1+Nw*nsteps+Nw*(nsteps*(nsteps+1))÷2+(index-1)*nsteps^2,i==idx) end return end @@ -324,53 +348,105 @@ end #Create 1 waveguide photon function waveguide_mul!(result,a::WaveguideCreate{B,B,1,idx},b,alpha,beta) where {B,idx} - if !isone(beta) - rmul!(result,beta) + if iszero(beta) + fill!(result, beta) + elseif !isone(beta) + rmul!(result, beta) end - add_onephoton_zerophoton!(result,b,alpha*a.factor,a.timeindex+(idx-1)*a.basis_l.nsteps) + @inbounds result[a.timeindex+(idx-1)*a.basis_l.nsteps+1] += alpha*a.factor*b[1] return end #Create 2 waveguide photon function waveguide_mul!(result,a::WaveguideCreate{B,B,2,1},b,alpha,beta) where {B<:SingleWaveguideBasis} - if !isone(beta) - rmul!(result,beta) + if iszero(beta) + fill!(result, beta) + elseif !isone(beta) + rmul!(result, beta) end timeindex = a.timeindex nsteps = a.basis_l.nsteps - add_onephoton_zerophoton!(result,b,alpha*a.factor,timeindex) - twophotonview = TwoPhotonTimestepView(result,timeindex,nsteps,nsteps+1) - axpy!(alpha*a.factor,view(b,2:1:nsteps+1),twophotonview) + @inbounds result[timeindex+1] += alpha*a.factor*b[1] + #twophotonview = TwoPhotonTimestepView(result,timeindex,nsteps,nsteps+1) + #axpy!(alpha*a.factor,view(b,2:1:nsteps+1),twophotonview) + twophoton_create!(result,view(b,2:1:nsteps+1),alpha*a.factor,timeindex,nsteps,nsteps+1) return end #Create 2 waveguide photon function waveguide_mul!(result,a::WaveguideCreate{B,B,2,idx},b,alpha,beta) where {B<:MultipleWaveguideBasis,idx} - if !isone(beta) - rmul!(result,beta) + if iszero(beta) + fill!(result, beta) + elseif !isone(beta) + rmul!(result, beta) end timeindex = a.timeindex nsteps = a.basis_l.nsteps Nw = get_number_of_waveguides(a.basis_l) - add_onephoton_zerophoton!(result,b,alpha*a.factor,timeindex+(idx-1)*nsteps) - twophotonview = TwoPhotonTimestepView(result,timeindex,nsteps,Nw*nsteps+1+(idx-1)*(nsteps*(nsteps+1))÷2) - axpy!(alpha*a.factor,view(b,2+(idx-1)*nsteps:1:idx*nsteps+1),twophotonview) + @inbounds result[timeindex+(idx-1)*a.basis_l.nsteps+1] += alpha*a.factor*b[1] + #twophotonview = TwoPhotonTimestepView(result,timeindex,nsteps,Nw*nsteps+1+(idx-1)*(nsteps*(nsteps+1))÷2) + #axpy!(alpha*a.factor,view(b,2+(idx-1)*nsteps:1:idx*nsteps+1),twophotonview) + twophoton_create!(result,view(b,2+(idx-1)*nsteps:1:idx*nsteps+1),alpha*a.factor,timeindex,nsteps,Nw*nsteps+1+(idx-1)*(nsteps*(nsteps+1))÷2) @simd for k in filter(x -> x != idx, 1:Nw) i,j = min(k,idx),max(k,idx) index = (i-1)*Nw + j - (i*(i+1))÷2 - twophotonview = TwoWaveguideTimestepView(result,timeindex,nsteps,1+Nw*nsteps+Nw*(nsteps*(nsteps+1))÷2+(index-1)*nsteps^2,i==idx) - axpy!(alpha*a.factor,view(b,2+(k-1)*nsteps:1:k*nsteps+1),twophotonview) + #twophotonview = TwoWaveguideTimestepView(result,timeindex,nsteps,1+Nw*nsteps+Nw*(nsteps*(nsteps+1))÷2+(index-1)*nsteps^2,i==idx) + #axpy!(alpha*a.factor,view(b,2+(k-1)*nsteps:1:k*nsteps+1),twophotonview) + twowaveguide_create!(result,view(b,2+(k-1)*nsteps:1:k*nsteps+1),alpha*a.factor,timeindex,nsteps,1+Nw*nsteps+Nw*(nsteps*(nsteps+1))÷2+(index-1)*nsteps^2,i==idx) end return end +function twophoton_create!(result,b,alpha,timeindex,nsteps,offset) + @simd for i in 1:timeindex-1 + @inbounds result[offset + twophoton_index(i,nsteps,timeindex)] += alpha*b[i] + end + @simd for i in timeindex+1:lastindex(b) + @inbounds result[offset + twophoton_index(timeindex,nsteps,i)] += alpha*b[i] + end + @inbounds result[offset + twophoton_index(timeindex,nsteps,timeindex)] += sqrt(2)*alpha*b[timeindex] +end +function twophoton_destroy!(result,b,alpha,timeindex,nsteps,offset) + @simd for i in 1:timeindex-1 + @inbounds result[i] += alpha*b[offset + twophoton_index(i,nsteps,timeindex)] + end + @simd for i in timeindex+1:lastindex(result) + @inbounds result[i] += alpha*b[offset + twophoton_index(timeindex,nsteps,i)] + end + @inbounds result[timeindex] += sqrt(2)*alpha*b[offset + twophoton_index(timeindex,nsteps,timeindex)] +end + +function twowaveguide_create!(result,b,alpha,timeindex,nsteps,offset,order) + if order + @simd for i in eachindex(b) + @inbounds result[offset + (i-1)*nsteps + timeindex] += alpha*b[i] + end + else + @simd for i in eachindex(b) + @inbounds result[offset + (timeindex-1)*nsteps + i] += alpha*b[i] + end + end +end +function twowaveguide_destroy!(result,b,alpha,timeindex,nsteps,offset,order) + if order + @simd for i in eachindex(result) + @inbounds result[i] += alpha*b[offset + (i-1)*nsteps + timeindex] + end + else + @simd for i in eachindex(result) + @inbounds result[i] += alpha*b[offset + (timeindex-1)*nsteps + i] + end + end +end + + -function add_zerophoton_onephoton!(a,b,alpha,timeindex::Int) - a[1] += alpha*b[timeindex+1] +@inline function add_zerophoton_onephoton!(a,b,alpha,timeindex::Int) + @inbounds a[1] += alpha*b[timeindex+1] end -function add_onephoton_zerophoton!(a,b,alpha,timeindex::Int) - a[1+timeindex] += alpha*b[1] +@inline function add_onephoton_zerophoton!(a,b,alpha,timeindex::Int) + @inbounds a[1+timeindex] += alpha*b[1] end function add_twophoton_onephoton!(a,b,alpha) diff --git a/src/WaveguideQED.jl b/src/WaveguideQED.jl index 7da28b4..2c43beb 100644 --- a/src/WaveguideQED.jl +++ b/src/WaveguideQED.jl @@ -4,10 +4,10 @@ using QuantumOptics using Strided using UnsafeArrays import LinearAlgebra: axpy!, dot, mul!, rmul! -import QuantumOptics: create, dagger, destroy, expect, identityoperator, tensor +import QuantumOpticsBase: create, dagger, destroy, expect, identityoperator, tensor export TwoPhotonTimestepView,TwoWaveguideTimestepView,OnePhotonView,TwoPhotonView,TwoWaveguideView, - WaveguideBasis,zerophoton,onephoton,twophoton,view_waveguide,get_waveguidetimeindex,set_waveguidetimeindex!,get_nsteps,get_waveguide_location,get_waveguide_basis,get_number_of_waveguides,get_waveguide_operators, + WaveguideBasis,zerophoton,onephoton,twophoton,view_waveguide,get_waveguidetimeindex,set_waveguidetimeindex!,get_dt,get_nsteps,get_waveguide_location,get_waveguide_basis,get_number_of_waveguides,get_waveguide_operators, WaveguideOperator,WaveguideDestroy,WaveguideCreate, CavityWaveguideAbsorption,CavityWaveguideEmission,emission,absorption, WaveguideInteraction, @@ -22,7 +22,6 @@ include("WaveguideOperator.jl") include("CavityWaveguideOperator.jl") include("WaveguideInteraction.jl") include("solver.jl") -include("should_upstream.jl") include("detection.jl") include("plotting.jl") include("precompile.jl") diff --git a/src/basis.jl b/src/basis.jl index 348816a..1ce6ea3 100644 --- a/src/basis.jl +++ b/src/basis.jl @@ -5,24 +5,25 @@ Basis for time binned Waveguide where `Np` is the number of photons in the waveguide and `Nw` the number of waveguides (default is 1). Currently number of photons is restricted to either 1 or 2. Times is timeinterval over which the photon state should be binned. """ -mutable struct WaveguideBasis{Np,Nw} <: QuantumOptics.Basis +mutable struct WaveguideBasis{Np,Nw} <: QuantumOpticsBase.Basis shape::Vector{Int} N::Int offset::Int nsteps::Int + dt::Float64 function WaveguideBasis(Np::Int,Nw::Int,times) dim = 0 N = length(times) if Np == 1 - dim = 1 + N*Nw + dim = N*Nw elseif Np == 2 #Number of combinations of waveguides ([1,2],[1,3],[2,3] and so on) combinations = Nw*(Nw-1)/2 - dim = 1 + Nw*N + Nw*(N*(N+1))÷2 + N^2*combinations + dim = Nw*N + Nw*(N*(N+1))÷2 + N^2*combinations else error("Currently no more than two simultanoues photons are allowed") end - new{Np,Nw}([dim+1], dim, 0,N) + new{Np,Nw}([dim+1], dim, 0,N,times[2]-times[1]) end end function WaveguideBasis(Np::Int,times) @@ -32,9 +33,9 @@ end Base.:(==)(b1::WaveguideBasis,b2::WaveguideBasis) = (b1.N==b2.N && b1.offset==b2.offset && b1.nsteps==b2.nsteps) #Type unions to dispatch on. -const SingleWaveguideBasis{Np} = Union{CompositeBasis{Vector{Int64}, T},WaveguideBasis{Np,1}} where {T<:Tuple{Vararg{Union{FockBasis,WaveguideBasis{Np,1}}}},Np} +const SingleWaveguideBasis{Np} = Union{CompositeBasis{Vector{Int64}, T},WaveguideBasis{Np,1}} where {T<:Tuple{Vararg{Union{NLevelBasis,FockBasis,WaveguideBasis{Np,1}}}},Np} const SingleWaveguideKet = Ket{T, Vector{ComplexF64}} where T<:SingleWaveguideBasis -const MultipleWaveguideBasis{Np,Nw} = Union{CompositeBasis{Vector{Int64},T},WaveguideBasis{Np,Nw},WaveguideBasis{Np,Nw}} where {T<:Tuple{Vararg{Union{FockBasis,WaveguideBasis{Np,Nw}}}},Np,Nw} +const MultipleWaveguideBasis{Np,Nw} = Union{CompositeBasis{Vector{Int64},T},WaveguideBasis{Np,Nw},WaveguideBasis{Np,Nw}} where {T<:Tuple{Vararg{Union{NLevelBasis,FockBasis,WaveguideBasis{Np,Nw}}}},Np,Nw} const MultipleWaveguideKet = Ket{T, Vector{ComplexF64}} where T<:MultipleWaveguideBasis @@ -50,7 +51,7 @@ using QuantumOptics; times=0:1:10; bw = WaveguideBasis(1,times); bc = FockBasis(2); -ψ_waveguide = Ket(bw,ones(length(times)+1)); +ψ_waveguide = Ket(bw,ones(length(times))); ψ_total = ψ_waveguide ⊗ fockstate(bc,0) ⊗ fockstate(bc,0); ψ_view = view_waveguide(ψ_total); ψ_view_index = view_waveguide(ψ_total,[:,1,1]); @@ -95,7 +96,7 @@ See [`onephoton`](@ref) on how to create onephoton wavepackets and [`view_wavegu using QuantumOptics; times = 0:1:10; bw = WaveguideBasis(1,times); -ψ1 = onephoton(bw,x->1,times,norm=false) +ψ1 = onephoton(bw,x->1,norm=false); OnePhotonView(ψ1) == ones(length(times)) ``` @@ -107,7 +108,7 @@ OnePhotonView(ψCavity,[3,:]) == ones(length(times)) ```@example onephotonview bw = WaveguideBasis(1,3,times); -ψ2 = onephoton(bw,2,x->1,times,norm=false) +ψ2 = onephoton(bw,2,x->1,norm=false) OnePhotonView(ψ,2) == ones(length(times)) ``` @@ -184,14 +185,15 @@ OnePhotonView(ψ::T,index::I) where {T<:MultipleWaveguideKet,I<:Union{Vector{An """ -Create a onephoton wavepacket in waveguide of the form ``W^\\dagger(\\xi) |0 \\rangle = \\int_{t_0}^{t_{end}} dt \\xi(t) w_{\\mathrm{i}}^\\dagger(t) |\\emptyset \\rangle`` where ``i`` is the index of the waveguide and return it as a `Ket`. +Create a onephoton wavepacket in waveguide of the form ``W^\\dagger(\\xi) |0 \\rangle = \\int_{t_0}^{t_{end}} dt \\xi(t) w_{\\mathrm{i}}^\\dagger(t) |\\emptyset \\rangle = \\sum_{k} `` where ``i`` is the index of the waveguide and return it as a `Ket`. See also [`WaveguideBasis`](@ref) and [`OnePhotonView`](@ref) for how to view the state. # Examples ```@example onewaveguide -times = 1:1:10; +dt = 1 +times = 1:dt:10; bw = WaveguideBasis(1,times); -ψ = onephoton(bw,x->1,times,norm=false); +ψ = onephoton(bw,x->dt,norm=false); OnePhotonView(ψ) == ones(length(times)) ``` ```@example onewaveguide @@ -201,25 +203,24 @@ OnePhotonView(ψ) == vec ``` ```@example onewaveguide bw = WaveguideBasis(1,3,times); -ψ = onephoton(bw,2,x->1,times); +ψ = onephoton(bw,2,x->dt); OnePhotonView(ψ,2) == ones(length(times)) ``` """ function onephoton end """ - onephoton(b::WaveguideBasis{T,1},ξ::Function,times,args...,norm=True) where {T} + onephoton(b::WaveguideBasis{T,1},ξ::Function,args...,norm=True) where {T} * Since `b` only contains a single waveguide, the index of the waveguide is not needed. -* `ξ` should be broadcastable as ξ.(times,args...). -* `times`: A vector or range of times where the wavefunction is evaluated. +* `ξ` should be broadcastable as ξ.(times,args...), where `times = 0:b.dt:(b.nsteps-1)*b.dt` (the times used to define the waveguidebasis) * `args...`: additional arguments to be passed to ξ if it is a function. * `norm::Bool=true`: If true normalize the resulting wavepacket. """ -function onephoton(b::WaveguideBasis{T,1},ξ::Function,times,args...; norm=true) where {T} +function onephoton(b::WaveguideBasis{T,1},ξ::Function,args...;norm=true) where {T} state = Ket(b) view = OnePhotonView(state) - view .= ξ.(times,args...) + view .= ξ.(0:b.dt:(b.nsteps-1)*b.dt,args...)*sqrt(b.dt) if norm normalize!(state) end @@ -230,13 +231,13 @@ end onephoton(b::WaveguideBasis{T,1},ξ::AbstractArray;norm=true) where {T} * Since `b` only contains a single waveguide, the index of the waveguide is not needed. -* `ξ` should be `AbstractArray` with length(ξ) == length(b.times) +* `ξ` should be `AbstractArray` with length(ξ) == b.nsteps * `norm::Bool=true`: If true normalize the resulting wavepacket. """ function onephoton(b::WaveguideBasis{T,1},ξ::AbstractArray;norm=true) where {T} state = Ket(b) view = OnePhotonView(state) - view .= ξ + view .= ξ*sqrt(b.dt) if norm normalize!(state) end @@ -244,19 +245,18 @@ function onephoton(b::WaveguideBasis{T,1},ξ::AbstractArray;norm=true) where {T} end """ - onephoton(b::WaveguideBasis{T,Nw},i::Int,ξ::Function,times,args...; norm=true) where {T,Nw} + onephoton(b::WaveguideBasis{T,Nw},i::Int,ξ::Function,args...; norm=true) where {T,Nw} * Since `b` contains `Nw` waveguides, the index of the waveguide needed. * `i` is the index of the waveguide at which the onephoton wavepacket is created in -* `ξ` should be broadcastable as ξ.(times,args...). -* `times`: A vector or range of times where the wavefunction is evaluated. +* `ξ` should be broadcastable as ξ.(times,args...), where `times = 0:b.dt:(b.nsteps-1)*b.dt` (the times used to define the waveguidebasis) * `args...`: additional arguments to be passed to ξ if it is a function. * `norm::Bool=true`: If true normalize the resulting wavepacket. """ -function onephoton(b::WaveguideBasis{T,Nw},i::Int,ξ::Function,times,args...; norm=true) where {T,Nw} +function onephoton(b::WaveguideBasis{T,Nw},i::Int,ξ::Function,args...; norm=true) where {T,Nw} state = Ket(b) view = OnePhotonView(state,i) - view .= ξ.(times,args...) + view .= ξ.(0:b.dt:(b.nsteps-1)*b.dt,args...)*sqrt(b.dt) if norm normalize!(state) end @@ -268,19 +268,19 @@ end * Since `b` contains a `Nw` waveguides, the index of the waveguide needed. * `i` is the index of the waveguide at which the onephoton wavepacket is created in -* `ξ` should be `AbstractArray` with length(ξ) == length(b.times) +* `ξ` should be `AbstractArray` with length(ξ) == b.nsteps * `norm::Bool=true`: If true normalize the resulting wavepacket. """ function onephoton(b::WaveguideBasis{T,Nw},i::Int,ξ::AbstractArray;norm=true) where {T,Nw} state = Ket(b) view = OnePhotonView(state,i) - view .= ξ + view .= ξ*sqrt(b.dt) if norm normalize!(state) end return state end -onephoton(b::WaveguideBasis{T,Nw},ξ::Function,times,args...;norm=true) where {T,Nw} = throw(ArgumentError("WaveguideBasis contains multiple waveguides. Please provide the index of the waveguide in which the excitation should be created.")) +onephoton(b::WaveguideBasis{T,Nw},ξ::Function,args...;norm=true) where {T,Nw} = throw(ArgumentError("WaveguideBasis contains multiple waveguides. Please provide the index of the waveguide in which the excitation should be created.")) onephoton(b::WaveguideBasis{T,Nw},ξ::AbstractArray;norm=true) where {T,Nw} = throw(ArgumentError("WaveguideBasis contains multiple waveguides. Please provide the index of the waveguide in which the excitation should be created.")) @@ -298,7 +298,7 @@ Basic viewing: using LinearAlgebra; #hide times = 1:1:10; bw = WaveguideBasis(2,times); -ψ = twophoton(bw,(t1,t2)->1,times,norm=false); +ψ = twophoton(bw,(t1,t2)->1,norm=false); ψview = TwoPhotonView(ψ); ψview == ones((length(times),length(times))) ``` @@ -315,7 +315,7 @@ bc = FockBasis(2); Viewing twophoton state in waveguide 2 with multiple waveguides ```@example twophotview bw = WaveguideBasis(2,3,times) -ψ = twophoton(bw,2,(t1,t2)->1,times,norm=false); +ψ = twophoton(bw,2,(t1,t2)->1,norm=false); ψview = TwoPhotonView(ψ,2); ψview == ones((length(times),length(times))) ``` @@ -330,7 +330,7 @@ Viewing twophoton state in waveguide 2 with multiple waveguides combined with ot Viewing twophotons across waveguide 1 and 2 ```@example twophotview bw = WaveguideBasis(2,3,times) -ψ = twophoton(bw,1,2,(t1,t2)->1,times,norm=false); +ψ = twophoton(bw,1,2,(t1,t2)->1,norm=false); ψview = TwoPhotonView(ψ,1,2); ψview == ones((length(times),length(times))) ``` @@ -358,15 +358,14 @@ function Base.getindex(x::TwoPhotonView,i::Int,j::Int) 1/sqrt(2)*x.data[x.offset + twophoton_index(j,x.nsteps,i)] end end -function Base.setindex!(x::TwoPhotonView,Left,i::Int,j::Int) +function Base.setindex!(x::TwoPhotonView,input,i::Int,j::Int) if i