histogram calculated differently

experiments.jl

@@ -4,11 +4,11 @@ tds_decay = [tempdecay(215,25, 0.02*i) for i in 1:10]
 push!(tds_decay, tempdecay(215,25, 0.4))
 push!(tds_decay, tempdecay(215,25, 0.6))
 decay_results = Dict[]
-for td in tds_decay[11:12]
+for td in tds_decay #[1:12]
     println(td)
     obj="obj" * string(rand(1:1000000))
     println(obj)
-    n_local_searches = 100
+    n_local_searches = 10
     max_iterations = 250
     k = 150
     # Load data
@@ -33,6 +33,7 @@ for td in tds_decay[11:12]
     println("Doing local search " * string(n_local_searches) * " times")
     #Threads.@threads 
     for i in 1:n_local_searches
+        println(i)
         rl = random_rollout(cdata)
         random_cost = cost_func(rl) # change to random
         update_result(results, obj, rl, random_cost, :random)
@@ -64,8 +65,8 @@ for speed in print_speeds
 
     obj="obj" * string(rand(1:1000000))
     println(obj)
-    n_local_searches = 100
-    max_iterations = 250
+    n_local_searches = 10
+    max_iterations = 50
     k = 150
     # Load data
     f_name = Symbol("cost_" * obj)
@@ -114,14 +115,16 @@ Plots.ylabel!("Cost")
 
 # different models
 
-models = ["M1", "M2", "M3", "M4"]
+models = ["M1", "M2", "M3"]
+speeds = [1,2,4]
 model_results = Dict[]
-for obj in models[3:end]
+for i in 1:length(models)#obj in models[3:end]
     # load data
+    obj = models[i]
+    speed = speeds[i]
     println("Loading data...")
     obj_add =  add * obj
-    contours = clean_contour.(contour.(JSON.parse(open(obj_add * "contours.json"))))
-    contours = [c for c in contours if c!=nothing]
+    contours = clean_contours(contour.(JSON.parse(open(obj_add * "contours.json"))))
     cdata = contourdata(contours, 20, 5) # contour data
     @time vd = voxdata(obj_add * "_voxels.csv", cdata)
     stress_multiplier!(vd.voxels,10)
@@ -129,7 +132,10 @@ for obj in models[3:end]
     n_local_searches = 100
     max_iterations = 250
     k = 150
-    # Load data
+
+    for c in cdata.contours
+        c.time *= speed
+    end
     f_name = Symbol("cost_" * obj)
     construct_cost(cdata, vd, ABS, td, f_name)
     cost_func = getfield(Main, f_name)
@@ -159,6 +165,9 @@ for obj in models[3:end]
     end
     println("best LS:" * string(results[obj]["cost_local"]))
     push!(model_results, results[obj])
+    for c in cdata.contours
+        c.time /= speed
+    end
 end
 ys_default = [r["cost_default"] for r in model_results]
 ys_local   = [r["cost_local"] for r in model_results]

main.jl

@@ -27,7 +27,7 @@ n_local_searches = 100
 max_iterations = 250
 k = 150
 
-obj = "M1"
+obj = "M3"
 
 # Load data
 f_name = Symbol("cost_" * obj)

utils.jl

@@ -777,16 +777,14 @@ function construct_cost_hist(cdata::contourdata, vd::voxdata, mat::material, td:
     contour_times = [cdata.contours[c].time[end] for c in 1:length(cdata.contours)]
 
     # considered voxels
-    println("n total voxels:")
-    println(length(vd.maps))
     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 = (1:length(vd.below))[(vd.below.!=0) .& not_empty_voxels .& (vd.above.!=0)]
     valid_voxels = valid_voxels[not_empty_voxels[vd.below[valid_voxels]]]
-    println("n non-empty voxels:")
-    println(sum(not_empty_voxels))
+    valid_voxels = valid_voxels[not_empty_voxels[vd.above[valid_voxels]]]
 
     considered_voxels = valid_voxels
     relbelows = vd.below[valid_voxels]
+    relaboves = vd.above[valid_voxels]
     rel_voxels = vd.voxels[considered_voxels,:]
     println("n valid voxels:")
     println(size(rel_voxels)[1])
@@ -819,18 +817,22 @@ function construct_cost_hist(cdata::contourdata, vd::voxdata, mat::material, td:
             voxtimes = sum($vox_seglen .* timestart[$vox_contour_id], dims=2) .+ $vox_c
 
             # voxel temps
-            Δt = voxtimes[$considered_voxels] - voxtimes[$relbelows]
-            ΔT = $(mat.T_cutoff.-td.ambient) .- $(td.extrusion-td.ambient).*ℯ.^(-$td.decay_rate.*Δt)
-            replace!(x-> x<0 ? 0 : x, ΔT)
+            Δt_below = voxtimes[$considered_voxels] - voxtimes[$relbelows]
+            ΔT_b = $(mat.T_cutoff.-td.ambient) .- $(td.extrusion-td.ambient).*ℯ.^(-$td.decay_rate.*Δt_below)
+            replace!(x-> x<0 ? 0 : x, ΔT_b)
+            Δt_above = voxtimes[$considered_voxels] - voxtimes[$relaboves]
+            ΔT_a = $(mat.T_cutoff.-td.ambient) .- $(td.extrusion-td.ambient).*ℯ.^(-$td.decay_rate.*-Δt_above)
+            replace!(x-> x<0 ? 0 : x, ΔT_a)
+            #println(ΔT_b - ΔT_a)
 
             # voxel stresses
             rel_v= $rel_voxels
-            σ11 = rel_v.Sx + $(mat.E*mat.α)*ΔT
-            σ22 = rel_v.Sy + $(mat.E*mat.α)*ΔT
+            σ11 = rel_v.Sx + $(mat.E*mat.α)*(ΔT_b-ΔT_a)
+            σ22 = rel_v.Sy + $(mat.E*mat.α)*(ΔT_b-ΔT_a)
             σ33 = rel_v.Sz
             σ12 = rel_v.Txy
-            σ23 = rel_v.Tyz + $((cdata.layer_height/vd.width)*mat.E*mat.α)*ΔT
-            σ31 = rel_v.Txy + $((cdata.layer_height/vd.width)*mat.E*mat.α)*ΔT
+            σ23 = rel_v.Tyz + $((cdata.layer_height/vd.width)*mat.E*mat.α)*(ΔT_b-ΔT_a)
+            σ31 = rel_v.Txy + $((cdata.layer_height/vd.width)*mat.E*mat.α)*(ΔT_b-ΔT_a)
             return ($vox_area_scaling .* ($F * (σ11 - σ22).^2 +
                 $G * ((σ33 - σ11).^2 + (σ33 - σ22).^2) +
                 $(2 * L) * (σ12).^2 +