{ "cells": [ { "cell_type": "code", "execution_count": 1, "id": "f52c5daf-9123-4b36-b83c-d5d2bcae5d91", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "10000" ] }, "execution_count": 1, "metadata": {}, "output_type": "execute_result" } ], "source": [ "using SpiDy\n", "using NPZ\n", "using DataFrames\n", "using CSV\n", "using ProgressMeter\n", "using Random\n", "using Statistics\n", "using LinearAlgebra\n", "using Plots\n", "\n", "########################\n", "########################\n", "\n", "Δt = .2\n", "N = 600\n", "tspan = (0., N*Δt)\n", "saveat = (0:1:N)*Δt\n", "\n", "α = 10. # 0.16\n", "ω0 = 7. # 1.4\n", "Γ = 5. # 0.5\n", "\n", "\n", "α = 0.16\n", "ω0 = 1.4\n", "Γ = 0.5\n", "\n", "\n", "J = LorentzianSD(α, ω0, Γ) # coloring the noise\n", "matrix = AnisoCoupling([-sin(π/4) 0. 0. # coupling to the environment\n", " 0. 0. 0.\n", " cos(π/4) 0. 0.]);\n", "\n", "matrix = IsoCoupling(1)\n", "\n", "T = .01\n", "noise = ClassicalNoise(T);\n", "\n", "nspin = 1 # number of spins\n", "\n", "s0 = [0, 1, 1] ./ sqrt(2)\n", "J0 = 1.\n", "JH = Nchain(nspin, J0)\n", "\n", "nruns = 10000" ] }, { "cell_type": "code", "execution_count": 2, "id": "6e3cdf02-c243-4622-a52a-a6d7f987c68d", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Starting...\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "\u001b[32mProgress: 100%|█████████████████████████████████████████| Time: 0:01:39\u001b[39m2:11\u001b[39m\n" ] } ], "source": [ "println(\"Starting...\")\n", "progress = Progress(nruns);\n", "\n", "sols = zeros(nruns, length(saveat), 3*nspin)\n", "\n", "Threads.@threads for i in 1:nruns\n", " bfields = [bfield(N, Δt, J, noise),\n", " bfield(N, Δt, J, noise),\n", " bfield(N, Δt, J, noise)];\n", " sol = diffeqsolver(s0, tspan, J, bfields, matrix; JH=JH, saveat=saveat);\n", " sols[i, :, :] = transpose(sol[:, :])\n", " next!(progress)\n", "end\n", "\n", "solavg = mean(sols, dims=1)[1, :, :];\n", "solstd = std(sols, dims=1)[1, :, :];" ] }, { "cell_type": "code", "execution_count": 3, "id": "d16b97dd-b9f5-415c-b388-2b1f3344be45", "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "plot(saveat, solavg[:, 1], ribbon=solstd[:, 1])" ] }, { "cell_type": "code", "execution_count": 4, "id": "a06addc4-47b2-4467-8148-3e68f3eae702", "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "plot(saveat, solavg[:, 2], ribbon=solstd[:, 2])" ] }, { "cell_type": "code", "execution_count": 5, "id": "bc1686fd-d86d-42b0-9d87-bcb911526d00", "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "plot(saveat, solavg[:, 3], ribbon=solstd[:, 3])" ] }, { "cell_type": "code", "execution_count": 6, "id": "59ac92d9-6a78-402b-91ae-3f5fd22cc6ff", "metadata": {}, "outputs": [], "source": [ "projected = zeros(nruns, length(saveat), 2)\n", "\n", "normalized_avg = [normalize(vec(row)) for row in eachrow(solavg)]\n", "\n", "Threads.@threads for i in 1:length(saveat)\n", " n = normalized_avg[i, :][1]\n", "\n", " u = normalize(cross(n, [0,0,1]))\n", " v = cross(u, n)\n", "\n", " for j in 1:nruns\n", " b = sols[j, i, :]\n", " proj = dot(b, n) * n\n", " b_ort = b - proj\n", "\n", " projected[j, i, 1] = dot(u,b)\n", " projected[j, i, 2] = dot(v,b)\n", " end\n", "end" ] }, { "cell_type": "code", "execution_count": 7, "id": "03467eca-3f0c-4bd4-b3f7-5649b68780bc", "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 7, "metadata": {}, "output_type": "execute_result" } ], "source": [ "histogram(projected[:, 280, 1])" ] }, { "cell_type": "code", "execution_count": 8, "id": "9c6484d3-3396-47bb-b063-98e8d013d020", "metadata": {}, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "histogram(projected[:, 280, 2])" ] }, { "cell_type": "code", "execution_count": 9, "id": "f6df949c-15b2-447f-a264-8ef48ce1725e", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "nth_moment (generic function with 1 method)" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "function nth_moment(data, N)\n", " return mean(data .^ N)\n", "end" ] }, { "cell_type": "code", "execution_count": 10, "id": "e4ad51fc-0dc7-4106-99e4-2b29e54b1930", "metadata": {}, "outputs": [], "source": [ "using HypothesisTests" ] }, { "cell_type": "code", "execution_count": 11, "id": "e5b9609e-f467-4bd7-9b51-9a4c1f175738", "metadata": {}, "outputs": [], "source": [ "JB = zeros(N)\n", "for i in 10:N\n", " JB[i] = HypothesisTests.JarqueBeraTest(projected[:, i, 1]).JB\n", "end" ] }, { "cell_type": "code", "execution_count": null, "id": "ebf486ee-5e68-45ba-b98d-b25e08d766b7", "metadata": { "scrolled": true }, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 21, "metadata": {}, "output_type": "execute_result" } ], "source": [ "plot(saveat[1:N], JB)" ] }, { "cell_type": "code", "execution_count": 33, "id": "2ba81fb8-bfed-40f2-8d33-4b1447a162b3", "metadata": { "scrolled": true }, "outputs": [ { "ename": "BoundsError", "evalue": "BoundsError: attempt to access 600-element StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64} at index [1:601]", "output_type": "error", "traceback": [ "BoundsError: attempt to access 600-element StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64} at index [1:601]", "", "Stacktrace:", " [1] throw_boundserror(A::StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, I::Tuple{UnitRange{Int64}})", " @ Base ./abstractarray.jl:691", " [2] checkbounds", " @ ./abstractarray.jl:656 [inlined]", " [3] getindex(r::StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, s::UnitRange{Int64})", " @ Base ./twiceprecision.jl:484", " [4] gr_draw_segments(series::Plots.Series, x::StepRangeLen{Float64, Base.TwicePrecision{Float64}, Base.TwicePrecision{Float64}, Int64}, y::Vector{Float64}, z::Nothing, fillrange::Tuple{Vector{Float64}, Vector{Float64}}, clims::Tuple{Float64, Float64})", " @ Plots ~/.julia/packages/Plots/ju9dp/src/backends/gr.jl:1825", " [5] gr_add_series(sp::Plots.Subplot{Plots.GRBackend}, series::Plots.Series)", " @ Plots ~/.julia/packages/Plots/ju9dp/src/backends/gr.jl:1742", " [6] gr_display(sp::Plots.Subplot{Plots.GRBackend}, w::Measures.AbsoluteLength, h::Measures.AbsoluteLength, vp_canvas::Plots.GRViewport{Float64})", " @ Plots ~/.julia/packages/Plots/ju9dp/src/backends/gr.jl:964", " [7] (::Plots.var\"#537#538\"{Int64, Int64, Plots.GRViewport{Float64}})(sp::Plots.Subplot{Plots.GRBackend})", " @ Plots ~/.julia/packages/Plots/ju9dp/src/backends/gr.jl:688", " [8] foreach(f::Plots.var\"#537#538\"{Int64, Int64, Plots.GRViewport{Float64}}, itr::Vector{Plots.Subplot})", " @ Base ./abstractarray.jl:2712", " [9] gr_display(plt::Plots.Plot{Plots.GRBackend}, dpi_factor::Int64)", " @ Plots ~/.julia/packages/Plots/ju9dp/src/backends/gr.jl:688", " [10] #582", " @ ~/.julia/packages/Plots/ju9dp/src/backends/gr.jl:2062 [inlined]", " [11] withenv(::Plots.var\"#582#583\"{Plots.Plot{Plots.GRBackend}, Int64}, ::Pair{String, String}, ::Vararg{Pair{String, String}})", " @ Base ./env.jl:172", " [12] _show(io::IOBuffer, #unused#::MIME{Symbol(\"image/svg+xml\")}, plt::Plots.Plot{Plots.GRBackend})", " @ Plots ~/.julia/packages/Plots/ju9dp/src/backends/gr.jl:2057", " [13] #invokelatest#2", " @ ./essentials.jl:716 [inlined]", " [14] invokelatest", " @ ./essentials.jl:714 [inlined]", " [15] show(io::IOBuffer, m::MIME{Symbol(\"image/svg+xml\")}, plt::Plots.Plot{Plots.GRBackend})", " @ Plots ~/.julia/packages/Plots/ju9dp/src/output.jl:232", " [16] sprint(::Function, ::MIME{Symbol(\"image/svg+xml\")}, ::Vararg{Any}; context::Nothing, sizehint::Int64)", " @ Base ./strings/io.jl:114", " [17] sprint", " @ ./strings/io.jl:108 [inlined]", " [18] _ijulia_display_dict(plt::Plots.Plot{Plots.GRBackend})", " @ Plots.IJuliaExt ~/.julia/packages/Plots/ju9dp/ext/IJuliaExt.jl:54", " [19] display_dict(plt::Plots.Plot{Plots.GRBackend})", " @ Plots.IJuliaExt ~/.julia/packages/Plots/ju9dp/ext/IJuliaExt.jl:70", " [20] #invokelatest#2", " @ ./essentials.jl:716 [inlined]", " [21] invokelatest", " @ ./essentials.jl:714 [inlined]", " [22] execute_request(socket::ZMQ.Socket, msg::IJulia.Msg)", " @ IJulia ~/.julia/packages/IJulia/Vo51o/src/execute_request.jl:112", " [23] #invokelatest#2", " @ ./essentials.jl:716 [inlined]", " [24] invokelatest", " @ ./essentials.jl:714 [inlined]", " [25] eventloop(socket::ZMQ.Socket)", " @ IJulia ~/.julia/packages/IJulia/Vo51o/src/eventloop.jl:8", " [26] (::IJulia.var\"#15#18\")()", " @ IJulia ./task.jl:429" ] } ], "source": [ "plot!(saveat, solavg[:, 3])" ] }, { "cell_type": "code", "execution_count": null, "id": "6ac5d053-d730-4118-8a73-f4af439a4058", "metadata": { "scrolled": true }, "outputs": [ { "data": { "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "execution_count": 13, "metadata": {}, "output_type": "execute_result" } ], "source": [ "histogram(projected[:, 200, 1])" ] }, { "cell_type": "code", "execution_count": 14, "id": "1f7d2316-131d-4faf-941e-eb4ca032727c", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "Jarque-Bera normality test\n", "--------------------------\n", "Population details:\n", " parameter of interest: skewness and kurtosis\n", " value under h_0: \"0 and 3\"\n", " point estimate: \"-0.17966985914759054 and 3.1431723658541975\"\n", "\n", "Test summary:\n", " outcome with 95% confidence: reject h_0\n", " one-sided p-value: <1e-13\n", "\n", "Details:\n", " number of observations: 10000\n", " JB statistic: 62.3431\n" ] }, "execution_count": 14, "metadata": {}, "output_type": "execute_result" } ], "source": [ "HypothesisTests.JarqueBeraTest(projected[:, 200, 1])" ] }, { "cell_type": "code", "execution_count": null, "id": "6c63a170-6eff-4c8b-8b80-542eb705666f", "metadata": {}, "outputs": [], "source": [] } ], "metadata": { "kernelspec": { "display_name": "Julia 1.7.3", "language": "julia", "name": "julia-1.7" }, "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.7.3" } }, "nbformat": 4, "nbformat_minor": 5 }