add greedy rollout, identify issue with dependency graph

main.jl

@@ -30,46 +30,70 @@ add = "/Users/jayant/phd/tempaware/models/"
 
 ### MAIN LOOP
 results = JSON.parse(open(add * "results.json"))
-n_local_searches = 5
+n_local_searches = 10
 n_models = 10
-max_iterations = 20
-k = 50
+max_iterations = 50
+k = 150
 for i in 1:n_models
     s = time()
     obj = "M" * string(i)
-    if obj in keys(results) || obj == "M4" || obj == "M9"
+    # if obj in keys(results) || obj == "M4" || obj == "M9"
+    if obj == "M4" || obj == "M9"
         continue # skip is already done
     end
     println("Starting " * obj)
 
     # Load data
-    println("Loading data...")
-    obj_add =  add * obj
-    contours = clean_contour.(contour.(JSON.parse(open(obj_add * "contours.json"))))
-    cdata = contourdata(contours, 40, 5) # contour data
-    vd = voxdata(obj_add * "_voxels.csv", cdata)
-    construct_cost(cdata, vd, ABS, td, :cost_f)
-    construct_best_neighbor(cdata, cost_f, k)
+    f_name = Symbol("cost_" * obj)
+    try
+        cost_func = getfield(Main, f_name)
+    catch
+        println("Loading data...")
+        obj_add =  add * obj
+        contours = clean_contour.(contour.(JSON.parse(open(obj_add * "contours.json"))))
+        cdata = contourdata(contours, 40, 5) # contour data
+        println(string(length(cdata.contours)) * " contours")
+        vd = voxdata(obj_add * "_voxels.csv", cdata)
+        println("Constructing Cost function...")
+        construct_cost(cdata, vd, ABS, td, f_name)
+        cost_func = getfield(Main, f_name)
+        println(cost_func)
+        construct_best_neighbor(cdata, cost_func, k)
+    end
+    cost_func = getfield(Main, f_name)
 
     # default rollout cost
-    rl = Vector(1:length(cdata.contours))
-    c = cost_f(rl)
-    update_result(results, obj, rl, c, :default)
+    rl_d = Vector(1:length(cdata.contours))
+    c = cost_func(rl_d)
+    update_result(results, obj, rl_d, c, :default)
+
+    # greedy rollout cost
+    rl_g = greedy_rollout(cdata)
+    c = cost_func(rl_g)
+    update_result(results, obj, rl_g, c, :greedy)
 
     # local search
-    println("Doing local search " * string(n_local_searches) * " times")
-    opt_costs = zeros(n_local_searches)
-    Threads.@threads for i in 1:n_local_searches
-        rl = random_rollout(cdata)
-        random_cost = cost_f(rl)
-        update_result(results, obj, rl, random_cost, :random)
-        local_cost = local_search!(rl, max_iterations)
-        update_result(results, obj, rl, local_cost, :local)
-    end
+    # println("Doing local search " * string(n_local_searches) * " times")
+    # opt_costs = zeros(n_local_searches)
+    # Threads.@threads for i in 1:n_local_searches
+    #     rl = random_rollout(cdata)
+    #     random_cost = cost_f(rl)
+    #     update_result(results, obj, rl, random_cost, :random)
+    #     local_cost = local_search!(rl, max_iterations)
+    #     update_result(results, obj, rl, local_cost, :local)
+    # end
     println("Completed in " * string(time()-s) * " seconds")
 end
 save_result(results, add * "results.json")
 
+### Load individual models for the plots
+obj = "M1"
+obj_add = add * obj
+contours = clean_contour.(contour.(JSON.parse(open(obj_add * "contours.json"))))
+cdata = contourdata(contours, 40, 5) # contour data
+vd = voxdata(obj_add * "_voxels.csv", cdata)
+
+
 ### Visualisations
 # distribution of stresses
 data_stress = .√( vd.voxels.Sx.^2 + vd.voxels.Sy.^2 + vd.voxels.Sz.^2)
@@ -79,16 +103,21 @@ histogram(data) # vast majority of voxels near 0 stress - can ignore
 histogram!(sort(data, rev=true)[1:4000]) # vast majority of voxels near 0 stress - can ignore
 
 # voxmap sanity check
-plot(vd, 46, cdata)
+plot(vd, 1500, cdata)
 
 # rollout order 3d
-plot(Vector(1:length(rl)), cdata)
-plot(rl, cdata)
+rl_d = Vector(1:length(cdata.contours))
+rl_r = random_rollout(cdata)
+rl_g = greedy_rollout(cdata)
+rl_l = Vector{Int}(results[obj]["best_rollout"])
+
+plot(rl_d, cdata)
+plot(rl_r, cdata)
+plot(rl_g, cdata)
+plot(rl_l, cdata)
 
 # rollout order animated
-rl = random_rollout(cdata)
-rl = Vector(1:length(rl))
-@time plot_animate(rl, cdata)
+@time plot_animate(rl_g, cdata, rate=5)
 
 # plot results
 plot(results)
\ No newline at end of file

utils.jl

@@ -413,6 +413,61 @@ function random_rollout(cdata::contourdata)
 end
 
 
+function greedy_rollout(cdata::contourdata)
+    done_contours = Set{Int}()
+    avail_contours = Set(cdata.layers[1])
+    todo_contours = Set(1:length(cdata.contours))
+    rollout = Vector{Int}()
+    contour_order = zeros(length(cdata.contours))
+    dep_times = zeros(length(cdata.contours))
+
+    while length(avail_contours) > 0
+        # get average times of dependency completion
+        for c in avail_contours
+            deps = inneighbors(cdata.G, c)
+            if length(deps) == 0
+                dep_times[c] = 1
+                continue
+            end
+            dep_times[c] = mean(contour_order[deps])
+        end
+        temp_avail_list = collect(avail_contours)
+        _,i = findmax(dep_times[temp_avail_list])
+        c = temp_avail_list[i] # contour with deps printed most recently
+        push!(rollout, c)
+        contour_order[c] = length(rollout)
+
+        # remove selected contour from todo and avail, add to done
+        delete!(avail_contours, c)
+        delete!(todo_contours, c)
+        push!(done_contours, c)
+
+        # update available contours
+        for i in todo_contours
+            if i in avail_contours
+                continue
+            elseif length(inneighbors(cdata.G, i)) == 0
+                push!(avail_contours, i)
+                continue
+            end
+
+            add = true
+            for j in inneighbors(cdata.G, i)
+                if !(j in done_contours)
+                    add = false
+                    break
+                end
+            end
+
+            if add
+                push!(avail_contours, i)
+            end
+        end
+    end
+    return rollout
+end
+
+
 function valid_swap(rollout::Vector{Int}, i::Int, j::Int, cdata::contourdata)
     # would swapping indices i and j in rollout result in another valid rollout?
     # NOTE THIS FUNCTION DOESNT WORK
@@ -578,7 +633,7 @@ end
 
 function construct_best_neighbor(cdata::contourdata, cost_f::Function, k::Int=50)
     n_contours=length(cdata.contours)
-    bn_func = quote 
+    bn_func = quote
         function best_neighbor!(rl::Vector{Int}, current_cost::Number)
             costs = Dict()
             for i in 1:1:$n_contours-1
@@ -627,7 +682,7 @@ function plot(vd::voxdata, i::Int, cdata::contourdata)
     sq_y = [loc[2]+w, loc[2]-w, loc[2]-w, loc[2]+w, loc[2]+w]
 
     traces = Vector{GenericTrace}()
-    push!(traces, scatter(x=sq_x, y=sq_y, 
+    push!(traces, scatter(x=sq_x, y=sq_y,
         mode="lines",
         name="Voxel",
         line=attr(color="black", width=4)
@@ -651,7 +706,8 @@ function plot(rl::Vector{Int}, cdata::contourdata)
     cols = range(HSV(720,1,1), stop=HSV(-720,1,1), length=n)
     traces = Vector{GenericTrace}()
     layout = Layout(
-        yaxis=attr(scaleanchor="x", scaleratio=1)
+        scene_aspect_ratio="data",
+        showlegend=false
     )
 
     for i in 1:n
@@ -692,10 +748,12 @@ function plot(results::Dict)
     rkey = keys(results)
     default_costs = [results[k]["cost_default"] for k in rkey]
     random_costs = [results[k]["cost_random"] for k in rkey]
+    greedy_costs = [results[k]["cost_greedy"] for k in rkey]
     local_costs = [results[k]["cost_local"] for k in rkey]
     trace2 = scatter(x=rkey,y=random_costs./default_costs, name="Random")
     trace3 = scatter(x=rkey,y=local_costs./default_costs, name="Local Search")
-    traces = [trace2, trace3]
+    trace4 = scatter(x=rkey,y=greedy_costs./default_costs, name="Greedy")
+    traces = [trace2, trace3, trace4]
     plot(traces)
 end
 
@@ -717,14 +775,16 @@ function update_result(
             # TODO: also save parameters, save all costs
             )
     end
-
-    if type==:default && results[obj]["cost_default"]>cost
-        results[obj]["cost_default"] = cost
-    elseif type==:random && results[obj]["cost_random"]>cost
-        results[obj]["cost_random"] = cost
-    elseif type==:local && results[obj]["cost_local"]>cost
+    if type==:local && results[obj]["cost_local"]>cost
         results[obj]["cost_local"] = cost
         results[obj]["best_rollout"] = rl
+    else
+        if !("cost_" * string(type) in keys(results[obj]))
+            results[obj]["cost_" * string(type)] = Inf
+        end
+        if results[obj]["cost_" * string(type)]>cost
+            results[obj]["cost_" * string(type)] = cost
+        end
     end
 end