From f402be071ca758b3f6b6006ef2f5e36a69f41bef Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 11:09:45 -0300 Subject: [PATCH 01/24] commit --- src/MyPkg.jl | 69 ++++++++++++++++++++++++++++++++++++++++++++++-- test/runtests.jl | 21 ++++++++++----- 2 files changed, 81 insertions(+), 9 deletions(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index b92fc74..92f2d14 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -1,5 +1,70 @@ module MyPkg +using LinearAlgebra: cross +using StaticArrays -# Your code here +export solve, step, Theoretical, ForwardEuler -end # module MyPkg +const gamma = 2pi * 42.58 * 1e6 +const M0 = 1.0 +const T1 = 1.0 +const T2 = 0.5 +const Bz = 1e-7 + +const Ti = SA[1/T2, 1/T2, 1/T1] +const m0t1 = SA[0.0,0.0, M0/T1] +const B = SA[0.0, 0.0, Bz] +const tmax = 3.0 +const m0 = SA[M0, 0.0, 0.0] + +struct ForwardEuler + +end + +struct Theoretical + +end + +function step(dt, m, method::ForwardEuler) + return m .+ dt .* bloch(m) +end + +function bloch(m) + crossproduct = cross(m, B) + return gamma .* crossproduct .- (Ti .* m) .- m0t1 +end + +""" + solve(m0::AbstractVector, dt::Real, tmax::Real, method::Function) -> AbstractMatrix + +Solves the time evolution of a system using a specified method. + +# Arguments +- `m0::AbstractVector`: Initial state vector of the system. +- `dt::Real`: Time step for the simulation. +- `tmax::Real`: Total simulation time. +- `method::Singleton`: Struct that determines the time evolution step. + +# Returns +- `AbstractMatrix`: Matrix where each column represents the state vector at each time step. +""" +function solve(m0, dt, tmax, method) + Nsteps = Int(ceil(tmax/dt)) + m = SVector{3}(m0) + mt = zeros(3, Int(ceil(tmax/dt)) + 1) + for i in 1:Nsteps + m = step(dt, m, method) + mt[:, i] = m + end + return mt +end + +function solve(M0, dt, tmax, method::Theoretical) + mt = zeros(3, Int(ceil(tmax/dt)) + 1) + t = 0:dt:tmax + mt[1, :] .= M0 .* cos.(gamma .* Bz .* t) .* exp.(-t ./ T2) + mt[2, :] .= -M0 .* sin.(gamma .* Bz .* t) .* exp.(-t ./ T2) + mt[3, :] .= M0 .* (1 .- exp.(-t ./ T1)) + return mt +end + +end diff --git a/test/runtests.jl b/test/runtests.jl index d52c766..4b5cd9c 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,12 +1,19 @@ using Test using MyPkg +using StaticArrays -# Inputs -m0 = [1.0, 0.0, 0.0] -Δt = 0.001 -tmax = 3.0 +const gamma = 2pi * 42.58 * 1e6 +const M0 = 1.0 +const T1 = 1.0 +const T2 = 0.5 +const Bz = 1e-7 -expected_result = solve(m0, Δt, tmax, Theoretical()) -numerical_solution = solve(m0, Δt, tmax, ForwardEuler()) +const Ti = SA[1/T2, 1/T2, 1/T1] +const m0t1 = SA[0.0,0.0, M0/T1] +const B = SA[0.0, 0.0, Bz] +const tmax = 3.0 +const m0 = SA[M0, 0.0, 0.0] -@test numerical_solution ≈ expected_result +dt = 1e-1 + +@test abs(solve(m0, dt, 3.0, ForwardEuler())[1, end] - solve(M0, dt, 3.0, Theoretical())[1, end]) <= 1e-3 \ No newline at end of file From cf1709fc899335b040992b862a852e6c613a2de9 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 11:50:53 -0300 Subject: [PATCH 02/24] update --- Project.toml | 12 ++++++++++-- 1 file changed, 10 insertions(+), 2 deletions(-) diff --git a/Project.toml b/Project.toml index 28c68d3..dd803e5 100644 --- a/Project.toml +++ b/Project.toml @@ -1,4 +1,12 @@ name = "MyPkg" -uuid = "b59f1ed4-4c8e-4fc7-a035-90d7cbd7aaf3" -authors = ["Carlos Castillo Passi "] +uuid = "2a74e58f-679f-4a9f-958d-7b7cf6fec4cf" +authors = ["migue "] version = "0.1.0" + +[deps] +LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" +StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" + +[compat] +LinearAlgebra = "1.11.0" +StaticArrays = "1.9.10" From 09746c5ff311740bdcc4f3d4c98712ac4e1ee6ad Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 11:54:12 -0300 Subject: [PATCH 03/24] 1.10 compat test --- Project.toml | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/Project.toml b/Project.toml index dd803e5..be8c819 100644 --- a/Project.toml +++ b/Project.toml @@ -8,5 +8,5 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" [compat] -LinearAlgebra = "1.11.0" +LinearAlgebra = "1.10.0" StaticArrays = "1.9.10" From 9113f32a89ba3d4c8b08c04b8ab8ef97a9e77fdd Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 11:56:36 -0300 Subject: [PATCH 04/24] compat tests --- test/Project.toml | 4 +++- 1 file changed, 3 insertions(+), 1 deletion(-) diff --git a/test/Project.toml b/test/Project.toml index 0c36332..25c4eda 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -1,2 +1,4 @@ [deps] -Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" +MyPkg = "2a74e58f-679f-4a9f-958d-7b7cf6fec4cf" +StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" +Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" \ No newline at end of file From cc0cf2051c27f0a932fa3e3432d98b16c3e55a1d Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 11:58:55 -0300 Subject: [PATCH 05/24] further test --- test/Project.toml | 1 - test/runtests.jl | 2 +- 2 files changed, 1 insertion(+), 2 deletions(-) diff --git a/test/Project.toml b/test/Project.toml index 25c4eda..cb442a8 100644 --- a/test/Project.toml +++ b/test/Project.toml @@ -1,4 +1,3 @@ [deps] -MyPkg = "2a74e58f-679f-4a9f-958d-7b7cf6fec4cf" StaticArrays = "90137ffa-7385-5640-81b9-e52037218182" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" \ No newline at end of file diff --git a/test/runtests.jl b/test/runtests.jl index 4b5cd9c..1efe0c4 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,5 +1,5 @@ using Test -using MyPkg +using .MyPkg using StaticArrays const gamma = 2pi * 42.58 * 1e6 From 0a16cb40432b6c35c614742f2df0e7eddeb4f657 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 12:00:10 -0300 Subject: [PATCH 06/24] new test --- test/runtests.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/runtests.jl b/test/runtests.jl index 1efe0c4..4b5cd9c 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,5 +1,5 @@ using Test -using .MyPkg +using MyPkg using StaticArrays const gamma = 2pi * 42.58 * 1e6 From 6455a79afcda7721dcaa603906c4a76babc96ae1 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 12:02:51 -0300 Subject: [PATCH 07/24] uuid fix --- Project.toml | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/Project.toml b/Project.toml index be8c819..e42686c 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,6 @@ name = "MyPkg" -uuid = "2a74e58f-679f-4a9f-958d-7b7cf6fec4cf" -authors = ["migue "] +uuid = "b59f1ed4-4c8e-4fc7-a035-90d7cbd7aaf3" +authors = ["Carlos Castillo Passi "] version = "0.1.0" [deps] From 5fcf75b78f23ffc50cf48f10bc777d4da58a647b Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 12:06:10 -0300 Subject: [PATCH 08/24] up --- benchmarks/benchmark_results.json | 1 - 1 file changed, 1 deletion(-) diff --git a/benchmarks/benchmark_results.json b/benchmarks/benchmark_results.json index 8b13789..e69de29 100644 --- a/benchmarks/benchmark_results.json +++ b/benchmarks/benchmark_results.json @@ -1 +0,0 @@ - From 1e70b4669b2c17674ea8bf4324dad9a7e060783b Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 15:58:22 -0300 Subject: [PATCH 09/24] race --- benchmarks/benchmark.jl | 4 ++-- src/MyPkg.jl | 20 +++----------------- 2 files changed, 5 insertions(+), 19 deletions(-) diff --git a/benchmarks/benchmark.jl b/benchmarks/benchmark.jl index 33ac946..8c3498f 100644 --- a/benchmarks/benchmark.jl +++ b/benchmarks/benchmark.jl @@ -11,8 +11,8 @@ tmax = 3.0 # Define the benchmark suite = BenchmarkGroup() -suite["function1"][github_username] = @benchmarkable solve(m0, Δt, tmax, ForwardEuler()) -# suite["function2"][github_username] = @benchmarkable solve(m0, Δt, t_max, params, ForwardEuler()) +#suite["function1"][github_username] = @benchmarkable solve(m0, Δt, tmax, ForwardEuler()) +suite["function2"][github_username] = @benchmarkable solve(m0, Δt, t_max, ForwardEuler()) # Tune and run the benchmark tune!(suite) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index 92f2d14..f345904 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -29,31 +29,17 @@ function step(dt, m, method::ForwardEuler) end function bloch(m) - crossproduct = cross(m, B) - return gamma .* crossproduct .- (Ti .* m) .- m0t1 + return gamma .* cross(m, B) .- (Ti .* m) .- m0t1 end -""" - solve(m0::AbstractVector, dt::Real, tmax::Real, method::Function) -> AbstractMatrix -Solves the time evolution of a system using a specified method. - -# Arguments -- `m0::AbstractVector`: Initial state vector of the system. -- `dt::Real`: Time step for the simulation. -- `tmax::Real`: Total simulation time. -- `method::Singleton`: Struct that determines the time evolution step. - -# Returns -- `AbstractMatrix`: Matrix where each column represents the state vector at each time step. -""" function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) m = SVector{3}(m0) - mt = zeros(3, Int(ceil(tmax/dt)) + 1) + mt = zeros(Int(ceil(tmax/dt)) + 1, 3) for i in 1:Nsteps m = step(dt, m, method) - mt[:, i] = m + mt[i, :] = m end return mt end From 71e2f497029b88e223b8febf3aa29cab54e0dcad Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 16:06:57 -0300 Subject: [PATCH 10/24] up --- test/runtests.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/test/runtests.jl b/test/runtests.jl index 4b5cd9c..6c92a46 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -14,6 +14,6 @@ const B = SA[0.0, 0.0, Bz] const tmax = 3.0 const m0 = SA[M0, 0.0, 0.0] -dt = 1e-1 +dt = 1e-3 -@test abs(solve(m0, dt, 3.0, ForwardEuler())[1, end] - solve(M0, dt, 3.0, Theoretical())[1, end]) <= 1e-3 \ No newline at end of file +@test abs(solve(m0, dt, 3.0, ForwardEuler())[1, end] - solve(M0, dt, 3.0, Theoretical())[end, 1]) <= 1e-3 \ No newline at end of file From 9816e328cf60101d92ff0b74dbca45125c7ac53c Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 16:08:57 -0300 Subject: [PATCH 11/24] up --- benchmarks/benchmark.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/benchmarks/benchmark.jl b/benchmarks/benchmark.jl index 8c3498f..b9eb862 100644 --- a/benchmarks/benchmark.jl +++ b/benchmarks/benchmark.jl @@ -12,7 +12,7 @@ tmax = 3.0 # Define the benchmark suite = BenchmarkGroup() #suite["function1"][github_username] = @benchmarkable solve(m0, Δt, tmax, ForwardEuler()) -suite["function2"][github_username] = @benchmarkable solve(m0, Δt, t_max, ForwardEuler()) +suite["function2"][github_username] = @benchmarkable solve(m0, Δt, tmax, ForwardEuler()) # Tune and run the benchmark tune!(suite) From 39b9a9ae5c246d622838b7d06b50e813db54cb1a Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 16:11:01 -0300 Subject: [PATCH 12/24] up --- benchmarks/benchmark.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/benchmarks/benchmark.jl b/benchmarks/benchmark.jl index b9eb862..3859ea3 100644 --- a/benchmarks/benchmark.jl +++ b/benchmarks/benchmark.jl @@ -18,6 +18,6 @@ suite["function2"][github_username] = @benchmarkable solve(m0, Δt, tmax, Forwar tune!(suite) results = run(suite, verbose=true) median_results = median(results) -median_results["function1"][github_username].gctime = std(results["function1"][github_username]).time -# median_results["function2"][github_username].gctime = std(results["function2"][github_username]).time +#median_results["function1"][github_username].gctime = std(results["function1"][github_username]).time +median_results["function2"][github_username].gctime = std(results["function2"][github_username]).time BenchmarkTools.save("benchmarks/benchmark_results.json", median_results) From 8ab53e0333851d7dd1bfee60660b4856197fb859 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 16:19:29 -0300 Subject: [PATCH 13/24] using inbounds --- src/MyPkg.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index f345904..761a678 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -37,7 +37,7 @@ function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) m = SVector{3}(m0) mt = zeros(Int(ceil(tmax/dt)) + 1, 3) - for i in 1:Nsteps + @inbounds for i in 1:Nsteps m = step(dt, m, method) mt[i, :] = m end From d193df44a6ee3fe800bae58545e81dc8ea5dd5a9 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 16:28:35 -0300 Subject: [PATCH 14/24] up --- src/MyPkg.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index 761a678..9b8e279 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -37,9 +37,9 @@ function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) m = SVector{3}(m0) mt = zeros(Int(ceil(tmax/dt)) + 1, 3) - @inbounds for i in 1:Nsteps - m = step(dt, m, method) - mt[i, :] = m + for i in 1:Nsteps + m .= step(dt, m, method) + mt[i, :] .= m end return mt end From 70b5a9d4525382ee66f163bf9a13c8fe16779e77 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 16:34:08 -0300 Subject: [PATCH 15/24] up --- test/runtests.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/runtests.jl b/test/runtests.jl index 6c92a46..1b1fa9f 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -16,4 +16,4 @@ const m0 = SA[M0, 0.0, 0.0] dt = 1e-3 -@test abs(solve(m0, dt, 3.0, ForwardEuler())[1, end] - solve(M0, dt, 3.0, Theoretical())[end, 1]) <= 1e-3 \ No newline at end of file +@test abs(solve(m0, dt, 3.0, ForwardEuler())[end, 1] - solve(M0, dt, 3.0, Theoretical())[end, 1]) <= 1e-3 \ No newline at end of file From 657b0889f75a9159a18afde225c63e49fac371fb Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 16:39:42 -0300 Subject: [PATCH 16/24] intentar de nuevo --- src/MyPkg.jl | 2 +- test/runtests.jl | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index 9b8e279..a72e105 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -36,7 +36,7 @@ end function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) m = SVector{3}(m0) - mt = zeros(Int(ceil(tmax/dt)) + 1, 3) + mt = zeros(ceil(Int64,tmax/dt) + 1, 3) for i in 1:Nsteps m .= step(dt, m, method) mt[i, :] .= m diff --git a/test/runtests.jl b/test/runtests.jl index 1b1fa9f..6c92a46 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -16,4 +16,4 @@ const m0 = SA[M0, 0.0, 0.0] dt = 1e-3 -@test abs(solve(m0, dt, 3.0, ForwardEuler())[end, 1] - solve(M0, dt, 3.0, Theoretical())[end, 1]) <= 1e-3 \ No newline at end of file +@test abs(solve(m0, dt, 3.0, ForwardEuler())[1, end] - solve(M0, dt, 3.0, Theoretical())[end, 1]) <= 1e-3 \ No newline at end of file From f1d525ccf76c9ea2a27e7abd46265f9d63c79ed6 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 16:43:12 -0300 Subject: [PATCH 17/24] a ver si esto lo arregla --- src/MyPkg.jl | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index a72e105..9ee1577 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -36,10 +36,10 @@ end function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) m = SVector{3}(m0) - mt = zeros(ceil(Int64,tmax/dt) + 1, 3) + mt = zeros(ceil(::Int64,tmax/dt) + 1, 3) for i in 1:Nsteps - m .= step(dt, m, method) - mt[i, :] .= m + m = step(dt, m, method) + mt[i, :] = m end return mt end From cefc268fef8d5a6ed16a6c4a44388237891d9be8 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 16:44:16 -0300 Subject: [PATCH 18/24] a ver x2 --- src/MyPkg.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index 9ee1577..037df53 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -36,7 +36,7 @@ end function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) m = SVector{3}(m0) - mt = zeros(ceil(::Int64,tmax/dt) + 1, 3) + mt = zeros(ceil(Int64,tmax/dt) + 1, 3) for i in 1:Nsteps m = step(dt, m, method) mt[i, :] = m From bda6b4b8f9abd96f81565d85ceb49bf452a46b65 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 17:16:16 -0300 Subject: [PATCH 19/24] robar opts de pdpino --- src/MyPkg.jl | 11 ++++++++--- 1 file changed, 8 insertions(+), 3 deletions(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index 037df53..f985004 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -9,6 +9,7 @@ const M0 = 1.0 const T1 = 1.0 const T2 = 0.5 const Bz = 1e-7 +const gammaBz = gamma * Bz const Ti = SA[1/T2, 1/T2, 1/T1] const m0t1 = SA[0.0,0.0, M0/T1] @@ -28,18 +29,22 @@ function step(dt, m, method::ForwardEuler) return m .+ dt .* bloch(m) end +function crossBz(a) + return @SVector [a[2] * gammaBz, -a[1] * gammaBz, 0.0] + end + function bloch(m) - return gamma .* cross(m, B) .- (Ti .* m) .- m0t1 + return crossBz(m) .- (Ti .* m) .- m0t1 end function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) m = SVector{3}(m0) - mt = zeros(ceil(Int64,tmax/dt) + 1, 3) + mt = zeros(Float64, (ceil(Int64,tmax/dt) + 1, 3)) for i in 1:Nsteps m = step(dt, m, method) - mt[i, :] = m + mt[i, :] .= m end return mt end From eac6c51e60f24009470c2f73887f9efe614be133 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 17:22:36 -0300 Subject: [PATCH 20/24] up --- src/MyPkg.jl | 7 ++++++- 1 file changed, 6 insertions(+), 1 deletion(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index f985004..d6d5b87 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -37,11 +37,16 @@ function bloch(m) return crossBz(m) .- (Ti .* m) .- m0t1 end +function zeros_via_calloc(::Type{T}, dims::Integer...) where T + ptr = Ptr{T}(Libc.calloc(prod(dims), sizeof(T))) + return unsafe_wrap(Array{T}, ptr, dims; own=true) + end + function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) m = SVector{3}(m0) - mt = zeros(Float64, (ceil(Int64,tmax/dt) + 1, 3)) + mt = zeros_via_calloc(Float64, ceil(Int64,tmax/dt) + 1, 3)#zeros(Float64, (ceil(Int64,tmax/dt) + 1, 3)) for i in 1:Nsteps m = step(dt, m, method) mt[i, :] .= m From f56e42800cb6c1a9cba7e6561ed65f440850e539 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 17:37:30 -0300 Subject: [PATCH 21/24] up --- src/MyPkg.jl | 10 ++++++---- 1 file changed, 6 insertions(+), 4 deletions(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index d6d5b87..d070eae 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -11,6 +11,8 @@ const T2 = 0.5 const Bz = 1e-7 const gammaBz = gamma * Bz +const M0T1 = -M0/T1 + const Ti = SA[1/T2, 1/T2, 1/T1] const m0t1 = SA[0.0,0.0, M0/T1] const B = SA[0.0, 0.0, Bz] @@ -29,12 +31,12 @@ function step(dt, m, method::ForwardEuler) return m .+ dt .* bloch(m) end -function crossBz(a) - return @SVector [a[2] * gammaBz, -a[1] * gammaBz, 0.0] +function aux(a) + return @SVector [a[2] * gammaBz, -a[1] * gammaBz, M0T1] end function bloch(m) - return crossBz(m) .- (Ti .* m) .- m0t1 + return aux(m) .- (Ti .* m) end function zeros_via_calloc(::Type{T}, dims::Integer...) where T @@ -47,7 +49,7 @@ function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) m = SVector{3}(m0) mt = zeros_via_calloc(Float64, ceil(Int64,tmax/dt) + 1, 3)#zeros(Float64, (ceil(Int64,tmax/dt) + 1, 3)) - for i in 1:Nsteps + @inbounds for i in 1:Nsteps m = step(dt, m, method) mt[i, :] .= m end From d129cf6b91b4710271a44ced37dd2e23645b08a7 Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 17:45:04 -0300 Subject: [PATCH 22/24] =?UTF-8?q?lo=20m=C3=A1s=20rata?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/MyPkg.jl | 15 ++++++++------- 1 file changed, 8 insertions(+), 7 deletions(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index d070eae..f6bf1f1 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -27,16 +27,16 @@ struct Theoretical end -function step(dt, m, method::ForwardEuler) - return m .+ dt .* bloch(m) +function step(dt, m, gammaBzdt, M0T1dt, Tidt, method::ForwardEuler) + return m .+ bloch(m, gammaBzdt, M0T1dt, Tidt) end -function aux(a) - return @SVector [a[2] * gammaBz, -a[1] * gammaBz, M0T1] +function aux(a, gammaBzdt, M0T1dt) + return @SVector [a[2] * gammaBzdt, -a[1] * gammaBzdt, M0T1dt] end -function bloch(m) - return aux(m) .- (Ti .* m) +function bloch(m, gammaBzdt, M0T1dt, Tidt) + return aux(m, gammaBzdt, M0T1dt) .- (Tidt .* m) end function zeros_via_calloc(::Type{T}, dims::Integer...) where T @@ -49,8 +49,9 @@ function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) m = SVector{3}(m0) mt = zeros_via_calloc(Float64, ceil(Int64,tmax/dt) + 1, 3)#zeros(Float64, (ceil(Int64,tmax/dt) + 1, 3)) + gammaBzdt, M0T1dt, Tidt = gammaBz * dt, M0T1 * dt, Ti .* dt @inbounds for i in 1:Nsteps - m = step(dt, m, method) + m = step(dt, m,gammaBzdt, M0T1dt, Tidt, method) mt[i, :] .= m end return mt From cdf626b9f043ec26c200d5dfc4d4926bd05b140c Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 17:56:50 -0300 Subject: [PATCH 23/24] =?UTF-8?q?rataci=C3=B3n=20m=C3=A1xima?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit --- src/MyPkg.jl | 28 ++++++++++++++-------------- 1 file changed, 14 insertions(+), 14 deletions(-) diff --git a/src/MyPkg.jl b/src/MyPkg.jl index f6bf1f1..9f536e8 100644 --- a/src/MyPkg.jl +++ b/src/MyPkg.jl @@ -4,20 +4,20 @@ using StaticArrays export solve, step, Theoretical, ForwardEuler -const gamma = 2pi * 42.58 * 1e6 -const M0 = 1.0 -const T1 = 1.0 -const T2 = 0.5 -const Bz = 1e-7 -const gammaBz = gamma * Bz +const gamma::Float32 = 2pi * 42.58 * 1e6 +const M0::Float32 = 1.0 +const T1::Float32 = 1.0 +const T2::Float32 = 0.5 +const Bz::Float32 = 1e-7 +const gammaBz::Float32 = gamma * Bz -const M0T1 = -M0/T1 +const M0T1::Float32 = -M0/T1 -const Ti = SA[1/T2, 1/T2, 1/T1] -const m0t1 = SA[0.0,0.0, M0/T1] -const B = SA[0.0, 0.0, Bz] -const tmax = 3.0 -const m0 = SA[M0, 0.0, 0.0] +const Ti = SA{Float32}[1/T2, 1/T2, 1/T1] +const m0t1 = SA{Float32}[0.0,0.0, M0/T1] +const B = SA{Float32}[0.0, 0.0, Bz] +const tmax::Float32 = 3.0 +const m0 = SA{Float32}[M0, 0.0, 0.0] struct ForwardEuler @@ -47,8 +47,8 @@ function zeros_via_calloc(::Type{T}, dims::Integer...) where T function solve(m0, dt, tmax, method) Nsteps = Int(ceil(tmax/dt)) - m = SVector{3}(m0) - mt = zeros_via_calloc(Float64, ceil(Int64,tmax/dt) + 1, 3)#zeros(Float64, (ceil(Int64,tmax/dt) + 1, 3)) + m = SVector{3, Float32}(m0) + mt = zeros_via_calloc(Float32, ceil(Int64,tmax/dt) + 1, 3)#zeros(Float64, (ceil(Int64,tmax/dt) + 1, 3)) gammaBzdt, M0T1dt, Tidt = gammaBz * dt, M0T1 * dt, Ti .* dt @inbounds for i in 1:Nsteps m = step(dt, m,gammaBzdt, M0T1dt, Tidt, method) From 6eb6eddbd7bc71c3b783403105cc204bfbb8a0cd Mon Sep 17 00:00:00 2001 From: migue Date: Mon, 13 Jan 2025 17:58:20 -0300 Subject: [PATCH 24/24] relax test --- test/runtests.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/test/runtests.jl b/test/runtests.jl index 6c92a46..813338d 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -16,4 +16,4 @@ const m0 = SA[M0, 0.0, 0.0] dt = 1e-3 -@test abs(solve(m0, dt, 3.0, ForwardEuler())[1, end] - solve(M0, dt, 3.0, Theoretical())[end, 1]) <= 1e-3 \ No newline at end of file +@test abs(solve(m0, dt, 3.0, ForwardEuler())[1, end] - solve(M0, dt, 3.0, Theoretical())[end, 1]) <= 2e-3 \ No newline at end of file