Commit 3c05268b authored by Jayant Khatkar's avatar Jayant Khatkar
Browse files

multithreading main local search loop 2x speedup

parent de99aeb9
......@@ -524,68 +524,6 @@ function calc_cost(rollout::Vector{Int}, cdata::contourdata, vd::voxdata, mat::m
end
function construct_cost(cdata::contourdata, vd::voxdata, mat::material)
# store contour time lengths inside function
contour_times = [cdata.contours[c].time[end] for c in 1:length(cdata.contours)]
# considered voxels
not_empty_voxels = length.([m.seglen for m in vd.maps]) .>0
valid_voxels = (1:length(vd.below))[(vd.below.!=0) .& not_empty_voxels]
valid_voxels = valid_voxels[not_empty_voxels[vd.below[valid_voxels]]]
#stress_size = .√( vd.voxels.Sx.^2 + vd.voxels.Sy.^2 + vd.voxels.Sz.^2)
#top_voxes = Set(sortperm(stress_size)[1:Int(length(stress_size)/4)]) # top stressors only
#top_valid_voxels = sort(collect(Set(valid_voxels) ∩ top_voxes))
#rel_voxels = vd.voxels[top_valid_voxels,:]
#relmaps = sort(collect(Set(vd.below[top_valid_voxels]) ∪ Set(top_valid_voxels)))
#considered_voxels = indexin(top_valid_voxels, relmaps)
# relbelows = indexin(vd.below[top_valid_voxels], relmaps)
considered_voxels = valid_voxels
relbelows = vd.below[valid_voxels]
rel_voxels = vd.voxels[considered_voxels,:]
relmaps = vd.maps
# calculate stresses at each voxel
F = (2/mat.σ̄[1,1]^2 - 1/mat.σ̄[3,3]^2)/2
G = 1/(2*mat.σ̄[3,3]^2)
L = 1/(2*mat.σ̄[1,2]^2)
M = 1/(2*mat.σ̄[1,3]^2)
return quote
function cost_f(rl::Vector{Int})
# println("Cumulative Sum") # 100k
timestart = cumsum([$contour_times[c] for c in rl])
#println("Vox times") # 600
voxtimes = [v.seglen timestart[v.segcontours] + v.c for v in relmaps]
# println("time diff") # 28k
Δt = voxtimes[considered_voxels] - voxtimes[relbelows]
#println("Temp diff") # 6k
ΔT = Tcutoff .- (Tc .+ $(T0-Tc).*.^(-k.*Δt))
# println("Temp diff clean") # 82k
replace!(x-> x<0 ? 0 : x, ΔT)
# println("Stresses") # 6k
σ11 = rel_voxels.Sx + $(mat.E*mat.α)*ΔT
σ22 = rel_voxels.Sy + $(mat.E*mat.α)*ΔT
σ33 = rel_voxels.Sz
σ12 = rel_voxels.Txy
σ23 = rel_voxels.Tyz + $((cdata.layer_height/vd.width)*mat.E*mat.α)*ΔT
σ31 = rel_voxels.Txy + $((cdata.layer_height/vd.width)*mat.E*mat.α)*ΔT
cost = sum($F * (σ11 - σ22).^2 +
$G * ((σ33 - σ11).^2 + (σ33 - σ22).^2) +
$(2 * L) * (σ12).^2 +
$(2 * M) * (σ23 + σ31).^2)
return cost
end
end
end
function clean_contour(c::contour)
# remove first element of array if second element is the same
while c.pos[1,:] == c.pos[2,:]
......@@ -607,7 +545,7 @@ end
# ABS material properties
α_ABS = 10e-6#100 # 78 - 108 #10-6
α_ABS = 63e-6
E_ABS = 2.5e9 # 1.19e9 - 2.9e9
σ̄ = 6e7 # mega
σ̄_ABS = [
......@@ -618,41 +556,131 @@ E_ABS = 2.5e9 # 1.19e9 - 2.9e9
ABS = material(α_ABS, E_ABS, σ̄_ABS)
### LOAD IN DATA
obj = "/Users/jayant/phd/tempaware/" * "M1"
contours = clean_contour.(contour.(JSON.parse(open(obj * "contours.json"))))
cdata = contourdata(contours, 5, 5) # contour data
@time vd = voxdata(obj * "_voxels.csv", cdata)
rl = random_rollout(cdata)
calc_cost(rl, cdata, vd, ABS) # 500/second
@time calc_cost(rl, cdata, vd, ABS) # 500/second
#####################################
###### CONSTRUCT COST FUNCTION
mat = ABS
# store contour time lengths inside function
contour_times = [cdata.contours[c].time[end] for c in 1:length(cdata.contours)]
# considered voxels
not_empty_voxels = length.([m.seglen for m in vd.maps]) .>0
valid_voxels = (1:length(vd.below))[(vd.below.!=0) .& not_empty_voxels]
valid_voxels = valid_voxels[not_empty_voxels[vd.below[valid_voxels]]]
#stress_size = .√( vd.voxels.Sx.^2 + vd.voxels.Sy.^2 + vd.voxels.Sz.^2)
#top_voxes = Set(sortperm(stress_size)[1:Int(length(stress_size)/4)]) # top stressors only
#top_valid_voxels = sort(collect(Set(valid_voxels) ∩ top_voxes))
#rel_voxels = vd.voxels[top_valid_voxels,:]
#relmaps = sort(collect(Set(vd.below[top_valid_voxels]) ∪ Set(top_valid_voxels)))
#considered_voxels = indexin(top_valid_voxels, relmaps)
# relbelows = indexin(vd.below[top_valid_voxels], relmaps)
considered_voxels = valid_voxels
relbelows = vd.below[valid_voxels]
rel_voxels = vd.voxels[considered_voxels,:]
relmaps = vd.maps
# calculate stresses at each voxel
F = (2/mat.σ̄[1,1]^2 - 1/mat.σ̄[3,3]^2)/2
G = 1/(2*mat.σ̄[3,3]^2)
L = 1/(2*mat.σ̄[1,2]^2)
M = 1/(2*mat.σ̄[1,3]^2)
a = quote
function cost_f(rl::Vector{Int})
# println("Cumulative Sum") # 100k
timestart = cumsum([$contour_times[c] for c in rl])
#println("Vox times") # 600
voxtimes = [v.seglen timestart[v.segcontours] + v.c for v in relmaps]
# println("time diff") # 28k
Δt = voxtimes[considered_voxels] - voxtimes[relbelows]
#println("Temp diff") # 6k
ΔT = Tcutoff .- (Tc .+ $(T0-Tc).*.^(-k.*Δt))
# println("Temp diff clean") # 82k
replace!(x-> x<0 ? 0 : x, ΔT)
# println("Stresses") # 6k
σ11 = rel_voxels.Sx + $(mat.E*mat.α)*ΔT
σ22 = rel_voxels.Sy + $(mat.E*mat.α)*ΔT
σ33 = rel_voxels.Sz
σ12 = rel_voxels.Txy
σ23 = rel_voxels.Tyz + $((cdata.layer_height/vd.width)*mat.E*mat.α)*ΔT
σ31 = rel_voxels.Txy + $((cdata.layer_height/vd.width)*mat.E*mat.α)*ΔT
cost = sum($F * (σ11 - σ22).^2 +
$G * ((σ33 - σ11).^2 + (σ33 - σ22).^2) +
$(2 * L) * (σ12).^2 +
$(2 * M) * (σ23 + σ31).^2)
return cost
end
end
eval(a)
cost_f(rl)
# local search loop
o_rl = copy(rl)
n = length(rl)
k = 5
cost_val = Inf
o=0
@time for l in 1:100
println(l)
println(cost_val)
k = 5 # TODO k is hardcoded in function below (uses global variable)
function best_neighbor!(rl::Vector{Int}, current_cost::Number)
costs = Dict()
for i in 1:n-1
for j in i+1:min(i+k,n)
for i in 1:length(rl)-1
for j in i+1:min(i+k,length(rl)) # TODO length rl is constant
if valid_swap(rl, i, j, cdata)
swap!(rl, i, j)
costs[i,j] = cost_f(rl)
costs[i,j] = cost_f(rl) # TODO cost function global vairable being used here
swap!(rl, i, j)
o += 1
end
end
end
v, (i,j) = findmin(costs)
if v < cost_val
cost_val = v
if abs(v - current_cost) < current_cost/1e6
return 0
elseif v < current_cost
swap!(rl, i, j)
return v
end
return 0
end
function local_search!(rl::Vector{Int}, max_iter::Int)
cost_val = Inf
for l in 1:max_iter
c = best_neighbor!(rl, cost_val)
if c 0
cost_val = c
else
break
end
end
return cost_val
end
rl = random_rollout(cdata)
@time local_search!(rl, 5)
n_starts=10
my_costs = zeros(n_starts)
threads = zeros(n_starts)
@time Threads.@threads for i in 1:n_starts
my_costs[i] = local_search!(random_rollout(cdata), 5)
threads[i] = Threads.threadid()
end
@time for i in 1:n_starts
my_costs[i] = local_search!(random_rollout(cdata), 5)
threads[i] = Threads.threadid()
end
println(o)
### Visualise distribution of stresses
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment