From f4438434662852c661d79e7af39d1274f89f9b87 Mon Sep 17 00:00:00 2001 From: MRIDUL JAIN <105979087+Spinachboul@users.noreply.github.com> Date: Sun, 24 Mar 2024 23:54:06 +0530 Subject: [PATCH 1/7] Add multi-threading --- src/Radials.jl | 34 +++++++++++----------------------- 1 file changed, 11 insertions(+), 23 deletions(-) diff --git a/src/Radials.jl b/src/Radials.jl index f84d67b7..ec1d8d77 100644 --- a/src/Radials.jl +++ b/src/Radials.jl @@ -1,5 +1,7 @@ using LinearAlgebra using ExtendableSparse +using Surrogates +using Base.Threads _copy(t::Tuple) = t _copy(t) = copy(t) @@ -196,33 +198,19 @@ end _ret_copy(v::Base.RefValue) = v[] _ret_copy(v) = copy(v) -function _approx_rbf(val, rad::RadialBasis) +function _approx_rbf_threaded(val, rad::RadialBasis) n = length(rad.x) - - # make sure @inbounds is safe - if n > size(rad.coeff, 2) - throw("Length of model's x vector exceeds number of calculated coefficients ($n != $(size(rad.coeff, 2))).") - end - approx = _make_approx(val, rad) - - if rad.phi === linearRadial().phi - for i in 1:n - tmp = zero(eltype(val)) - @inbounds @simd ivdep for j in 1:length(val) - tmp += ((val[j] - rad.x[i][j]) / rad.scale_factor)^2 - end - tmp = sqrt(tmp) - _add_tmp_to_approx!(approx, i, tmp, rad) - end - else - tmp = collect(val) - @inbounds for i in 1:n - tmp = (val .- rad.x[i]) ./ rad.scale_factor - _add_tmp_to_approx!(approx, i, tmp, rad; f = rad.phi) + + Threads.@threads for i in 1:n + tmp = zero(eltype(val)) + @inbounds @simd ivdep for j in 1:length(val) + tmp += ((val[j] - rad.x[i][j]) / rad.scale_factor)^2 end + tmp = sqrt(tmp) + _add_tmp_to_approx!(approx, i, tmp, rad) end - + return _ret_copy(approx) end From 84cd68c280807f9408aecba5872820c942dfcc36 Mon Sep 17 00:00:00 2001 From: MRIDUL JAIN <105979087+Spinachboul@users.noreply.github.com> Date: Mon, 25 Mar 2024 00:22:07 +0530 Subject: [PATCH 2/7] Update Radials.jl --- src/Radials.jl | 1 - 1 file changed, 1 deletion(-) diff --git a/src/Radials.jl b/src/Radials.jl index ec1d8d77..345073a0 100644 --- a/src/Radials.jl +++ b/src/Radials.jl @@ -1,6 +1,5 @@ using LinearAlgebra using ExtendableSparse -using Surrogates using Base.Threads _copy(t::Tuple) = t From 2ee48d3bd044774db5d2514b631a2a668db2398a Mon Sep 17 00:00:00 2001 From: MRIDUL JAIN <105979087+Spinachboul@users.noreply.github.com> Date: Mon, 25 Mar 2024 00:48:33 +0530 Subject: [PATCH 3/7] Update Radials.jl --- src/Radials.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/Radials.jl b/src/Radials.jl index 345073a0..c18ba4f2 100644 --- a/src/Radials.jl +++ b/src/Radials.jl @@ -197,7 +197,7 @@ end _ret_copy(v::Base.RefValue) = v[] _ret_copy(v) = copy(v) -function _approx_rbf_threaded(val, rad::RadialBasis) +function _approx_rbf(val, rad::RadialBasis) n = length(rad.x) approx = _make_approx(val, rad) From fec3ef51e6cb960417540071eb37bab547b5dc62 Mon Sep 17 00:00:00 2001 From: MRIDUL JAIN <105979087+Spinachboul@users.noreply.github.com> Date: Wed, 27 Mar 2024 01:45:42 +0530 Subject: [PATCH 4/7] Update Radials.jl --- src/Radials.jl | 54 +++++++++++++++++++++++++++++++++----------------- 1 file changed, 36 insertions(+), 18 deletions(-) diff --git a/src/Radials.jl b/src/Radials.jl index c18ba4f2..137beca9 100644 --- a/src/Radials.jl +++ b/src/Radials.jl @@ -2,6 +2,8 @@ using LinearAlgebra using ExtendableSparse using Base.Threads +abstract type AbstractSurrogate end + _copy(t::Tuple) = t _copy(t) = copy(t) @@ -44,8 +46,15 @@ of a polynomial term. References: https://en.wikipedia.org/wiki/Polyharmonic_spline """ -function RadialBasis(x, y, lb, ub; rad::RadialFunction = linearRadial(), - scale_factor::Real = 0.5, sparse = false) +function RadialBasis( + x, + y, + lb, + ub; + rad::RadialFunction = linearRadial(), + scale_factor::Real = 0.5, + sparse = false +) q = rad.q phi = rad.phi coeff = _calc_coeffs(x, y, lb, ub, phi, q, scale_factor, sparse) @@ -110,10 +119,9 @@ using ChainRulesCore: @non_differentiable function _make_combination(n, d, ix) exponents_combinations = [e - for e - in collect(Iterators.product(Iterators.repeated(0:n, - d)...))[:] - if sum(e) <= n] + for + e in collect(Iterators.product(Iterators.repeated( + 0:n, d)...))[:] if sum(e) <= n] return exponents_combinations[ix + 1] end @@ -146,9 +154,7 @@ function multivar_poly_basis(x, ix, d, n) if n == 0 return one(eltype(x)) else - prod(a^d - for (a, d) - in zip(x, _make_combination(n, d, ix))) + prod(a^d for (a, d) in zip(x, _make_combination(n, d, ix))) end end @@ -182,13 +188,17 @@ function _add_tmp_to_approx!(approx, i, tmp, rad::RadialBasis; f = identity) end end # specialise when only single output dimension -function _make_approx(val, - ::RadialBasis{F, Q, X, <:AbstractArray{<:Number}}) where {F, Q, X} +function _make_approx( + val, ::RadialBasis{F, Q, X, <:AbstractArray{<:Number}}) where {F, Q, X} return Ref(zero(eltype(val))) end -function _add_tmp_to_approx!(approx::Base.RefValue, i, tmp, +function _add_tmp_to_approx!( + approx::Base.RefValue, + i, + tmp, rad::RadialBasis{F, Q, X, <:AbstractArray{<:Number}}; - f = identity) where {F, Q, X} + f = identity +) where {F, Q, X} @inbounds @simd ivdep for j in 1:size(rad.coeff, 1) approx[] += rad.coeff[j, i] * f(tmp) end @@ -200,16 +210,16 @@ _ret_copy(v) = copy(v) function _approx_rbf(val, rad::RadialBasis) n = length(rad.x) approx = _make_approx(val, rad) - + Threads.@threads for i in 1:n tmp = zero(eltype(val)) - @inbounds @simd ivdep for j in 1:length(val) + @inbounds @simd ivdep for j in eachindex(val) tmp += ((val[j] - rad.x[i][j]) / rad.scale_factor)^2 end tmp = sqrt(tmp) _add_tmp_to_approx!(approx, i, tmp, rad) end - + return _ret_copy(approx) end @@ -231,7 +241,15 @@ function add_point!(rad::RadialBasis, new_x, new_y) append!(rad.x, new_x) append!(rad.y, new_y) end - rad.coeff = _calc_coeffs(rad.x, rad.y, rad.lb, rad.ub, rad.phi, rad.dim_poly, - rad.scale_factor, rad.sparse) + rad.coeff = _calc_coeffs( + rad.x, + rad.y, + rad.lb, + rad.ub, + rad.phi, + rad.dim_poly, + rad.scale_factor, + rad.sparse + ) nothing end From 070749cbce5ae990c23cef9eb496b6c17ba6758f Mon Sep 17 00:00:00 2001 From: MRIDUL JAIN <105979087+Spinachboul@users.noreply.github.com> Date: Wed, 27 Mar 2024 01:56:13 +0530 Subject: [PATCH 5/7] Update Radials.jl --- src/Radials.jl | 2 -- 1 file changed, 2 deletions(-) diff --git a/src/Radials.jl b/src/Radials.jl index 137beca9..1139d80a 100644 --- a/src/Radials.jl +++ b/src/Radials.jl @@ -2,8 +2,6 @@ using LinearAlgebra using ExtendableSparse using Base.Threads -abstract type AbstractSurrogate end - _copy(t::Tuple) = t _copy(t) = copy(t) From 01a8e675157bb8229b795ad0ec528255fbee158c Mon Sep 17 00:00:00 2001 From: MRIDUL JAIN <105979087+Spinachboul@users.noreply.github.com> Date: Thu, 28 Mar 2024 15:55:40 +0530 Subject: [PATCH 6/7] Update Radials.jl --- src/Radials.jl | 12 ++++++++++-- 1 file changed, 10 insertions(+), 2 deletions(-) diff --git a/src/Radials.jl b/src/Radials.jl index 1139d80a..a1f6bc3f 100644 --- a/src/Radials.jl +++ b/src/Radials.jl @@ -209,12 +209,20 @@ function _approx_rbf(val, rad::RadialBasis) n = length(rad.x) approx = _make_approx(val, rad) - Threads.@threads for i in 1:n + # Define a function to compute tmp for a single index i + function compute_tmp(i) tmp = zero(eltype(val)) @inbounds @simd ivdep for j in eachindex(val) tmp += ((val[j] - rad.x[i][j]) / rad.scale_factor)^2 end - tmp = sqrt(tmp) + return sqrt(tmp) + end + + # Use pmap to parallelize the computation of tmp + tmp_values = pmap(compute_tmp, 1:n) + + # Update the approx using the computed tmp values + for (i, tmp) in enumerate(tmp_values) _add_tmp_to_approx!(approx, i, tmp, rad) end From da3846cedbfa976240eeb53dddd20838fd5833bf Mon Sep 17 00:00:00 2001 From: MRIDUL JAIN <105979087+Spinachboul@users.noreply.github.com> Date: Thu, 28 Mar 2024 16:08:37 +0530 Subject: [PATCH 7/7] Update Radials.jl --- src/Radials.jl | 1 + 1 file changed, 1 insertion(+) diff --git a/src/Radials.jl b/src/Radials.jl index a1f6bc3f..39b40b5c 100644 --- a/src/Radials.jl +++ b/src/Radials.jl @@ -1,6 +1,7 @@ using LinearAlgebra using ExtendableSparse using Base.Threads +using Distributed _copy(t::Tuple) = t _copy(t) = copy(t)