{
"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",
"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:06:29\u001b[39m3:53\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"
]
},
"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"
]
},
"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"
]
},
"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"
]
},
"execution_count": 7,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"histogram(projected[:, 500, 1])"
]
},
{
"cell_type": "code",
"execution_count": 8,
"id": "9c6484d3-3396-47bb-b063-98e8d013d020",
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n"
]
},
"execution_count": 8,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"histogram(projected[:, 500, 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": 12,
"id": "3ad3f882-caa1-4dfc-b83c-3b57c214ecea",
"metadata": {},
"outputs": [
{
"data": {
"image/svg+xml": [
"\n",
"\n"
]
},
"execution_count": 12,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"plot(saveat[10:N], JB[10:N])"
]
},
{
"cell_type": "code",
"execution_count": 13,
"id": "bb3c6b84-5b21-4d79-b5bc-3345a8f2d893",
"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.06640955564275836 and 2.878302889621456\"\n",
"\n",
"Test summary:\n",
" outcome with 95% confidence: reject h_0\n",
" one-sided p-value: 0.0012\n",
"\n",
"Details:\n",
" number of observations: 10000\n",
" JB statistic: 13.5213\n"
]
},
"execution_count": 13,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"HypothesisTests.JarqueBeraTest(projected[:, 400, 1])"
]
},
{
"cell_type": "code",
"execution_count": 14,
"id": "e331f6c7-fb4b-4f21-84d3-6a3f203a95b6",
"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.009755563978670194 and 3.0706123545537634\"\n",
"\n",
"Test summary:\n",
" outcome with 95% confidence: fail to reject h_0\n",
" one-sided p-value: 0.3269\n",
"\n",
"Details:\n",
" number of observations: 10000\n",
" JB statistic: 2.23616\n"
]
},
"execution_count": 14,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"HypothesisTests.JarqueBeraTest(projected[:, 2, 2])"
]
},
{
"cell_type": "code",
"execution_count": 15,
"id": "240660fb-d008-4dab-8bc1-aaf530885102",
"metadata": {},
"outputs": [
{
"data": {
"text/plain": [
"64"
]
},
"execution_count": 15,
"metadata": {},
"output_type": "execute_result"
}
],
"source": [
"using Base.Threads\n",
"nthreads()"
]
},
{
"cell_type": "code",
"execution_count": 16,
"id": "8150abab-cacc-4617-9d95-bebdc4c7cac8",
"metadata": {},
"outputs": [
{
"name": "stdout",
"output_type": "stream",
"text": [
" PID COMMAND VSZ RSS SIZE\n",
" 54147 julia 20209472 1568276 15911276\n"
]
}
],
"source": [
"run(`ps -p $(getpid()) -o pid,comm,vsize,rss,size`);"
]
},
{
"cell_type": "code",
"execution_count": null,
"id": "e21d4205-2f8e-4b39-b1ea-ba722ea3fb13",
"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
}