refactor and progress on cost function coding

main.jl

@@ -5,6 +5,7 @@ using LightGraphs
 using NearestNeighbors
 using Statistics
 using BenchmarkTools
+using Plots
 
 
 struct contour
@@ -13,7 +14,7 @@ struct contour
 end
 
 
-struct problem
+struct contourdata
     contours::Vector{contour}
     G::SimpleDiGraph
     layers::Vector
@@ -22,6 +23,21 @@ struct problem
 end
 
 
+struct voxmap
+    seglen::Vector{Float64}
+    segoffset::Vector{Float64}
+    segcontours::Vector{Int}
+end
+
+
+struct voxdata
+    voxels::DataFrame
+    maps::Vector{voxmap}
+    below::Vector{Int}
+    width::Number
+end
+
+
 function vecvec_to_matrix(vecvec)
     # convert vector of vectors int a matrix
     dim1 = length(vecvec)
@@ -41,7 +57,7 @@ function contour(d::Dict)
 end
 
 
-function problem(cons::Vector{contour}, max_layers::Int, min_dist::Number)
+function contourdata(cons::Vector{contour}, max_layers::Int, min_dist::Number)
     G = LightGraphs.SimpleDiGraph(0)
 
     # separate contours into layers
@@ -82,14 +98,7 @@ function problem(cons::Vector{contour}, max_layers::Int, min_dist::Number)
         end
     end
 
-    return problem(cons, G, layers, Dict(), layer_heights[1])
-end
-
-
-struct voxmap
-    seglen::Vector{Float64}
-    segoffset::Vector{Float64}
-    segcontours::Vector{Int}
+    return contourdata(cons, G, layers, Dict(), layer_heights[1])
 end
 
 
@@ -136,11 +145,11 @@ interpolate(p1, p2, xi, axis) = p1 + (p2-p1)*(xi-p1[axis])/(p2[axis]-p1[axis])
 interpolate(p1, p2, x1, xi, x2) = p1 + (p2-p1)*(xi-x1)/(x2-x1)
 
 
-function voxmap(vox, vox_d, prob::problem)
+function voxmap(vox::Vector{Float64}, vox_d::Number, cdata::contourdata)
 
     # for one vox, get all contours which pass through it
     # only need to search contours in its layer
-    l = Int(round((vox[3] + prob.layer_height/2)/prob.layer_height))
+    l = Int(round((vox[3] + cdata.layer_height/2)/cdata.layer_height))
     voxx1 = vox[1] + vox_d/2
     voxx2 = vox[1] - vox_d/2
     voxy1 = vox[2] + vox_d/2
@@ -153,9 +162,9 @@ function voxmap(vox, vox_d, prob::problem)
     seglen_sofar = 0
     t_start = 0
 
-    for cid in prob.layers[l]
+    for cid in cdata.layers[l]
 
-        c = prob.contours[cid]
+        c = cdata.contours[cid]
 
         # check if contour passes thorough this vox
         for i in 2:size(c.pos)[1]
@@ -312,10 +321,21 @@ function voxmap(vox, vox_d, prob::problem)
 end
 
 
-function random_rollout(prob::problem)
+function voxdata(fname::String, cdata::contourdata)
+    voxels = DataFrames.DataFrame(CSV.File(fname))
+    w = dist(Vector(voxels[1, ["x","y","z"]]), Vector(voxels[2, ["x","y","z"]]))
+    vpos = [[v.x, v.y, v.z] for v in eachrow(voxels)]
+    voxms = [voxmap(v, w, cdata) for v in vpos]
+    below = indexin([v - [0,0,cdata.layer_height] for v in vpos], vpos)
+    replace!(below, nothing=>0)
+    return voxdata(voxels, voxms, below, w)
+end
+
+
+function random_rollout(cdata::contourdata)
     done_contours = Set{Int}()
-    avail_contours = Set(prob.layers[1])
-    todo_contours = Set(1:length(prob.contours))
+    avail_contours = Set(cdata.layers[1])
+    todo_contours = Set(1:length(cdata.contours))
     rollout = Vector{Int}()
 
     while length(avail_contours) > 0
@@ -331,13 +351,13 @@ function random_rollout(prob::problem)
         for i in todo_contours
             if i in avail_contours
                 continue
-            elseif length(inneighbors(prob.G, i)) == 0
+            elseif length(inneighbors(cdata.G, i)) == 0
                 push!(avail_contours, i)
                 continue
             end
 
             add = true
-            for j in inneighbors(prob.G, i)
+            for j in inneighbors(cdata.G, i)
                 if !(j in done_contours)
                     add = false
                     break
@@ -354,7 +374,7 @@ function random_rollout(prob::problem)
 end
 
 
-function valid_swap(rollout::Vector{Int}, i::Int, j::Int, prob::problem)
+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
     # IT ONLY CHECKS DEPENDENCIES TO A DEPTH OF 1
@@ -368,13 +388,13 @@ function valid_swap(rollout::Vector{Int}, i::Int, j::Int, prob::problem)
 
     c1 = rollout[i]
     c2 = rollout[j]
-    c2_dependson = inneighbors(prob.G, c2)
+    c2_dependson = inneighbors(cdata.G, c2)
 
     if c1 in c2_dependson
         return false
     end
 
-    c1_dependents = outneighbors(prob.G, c1)
+    c1_dependents = outneighbors(cdata.G, c1)
     c_between = rollout[i+1:j-1]
 
     for c in c_between
@@ -387,12 +407,12 @@ function valid_swap(rollout::Vector{Int}, i::Int, j::Int, prob::problem)
 end
 
 
-function check_validity(rollout::Vector{Int}, prob::problem)
+function check_validity(rollout::Vector{Int}, cdata::contourdata)
     # make sure a given rollout is valid
     done_contours = Set{Int}()
 
     for c in rollout
-        c_dependson = inneighbors(prob.G, c)
+        c_dependson = inneighbors(cdata.G, c)
 
         if !issubset(c_dependson, done_contours)
             return false
@@ -425,26 +445,42 @@ function test_voxmap()
     contour5 = contour(pos5,time5)
 
     contours = [contour1, contour2, contour3, contour4, contour5]
-    prob = problem(contours, 1, 1)
-    vm = voxmap(vox, vox_d, prob)
+    cdata = contourdata(contours, 1, 1)
+    vm = voxmap(vox, vox_d, cdata)
     return vm
 end
 
 
-function rollout2time(rollout::Vector{Int}, prob)
+function rollout2time(rollout::Vector{Int}, cdata::contourdata)
     # start time of each contour, assuming no travel time
-    return cumsum([prob.contours[c].time[end] for c in rollout])
+    return cumsum([cdata.contours[c].time[end] for c in rollout])
 end
 
+# Temperature function
+T0 = 215 # extrusion temp
+Tc = 25 # room temp
+Tcutoff = 100 # temperature above which strain isn't happening
+k = 0.005 #value unknown  - assumes cooling to 50C from 215 in 8 minutes
+Temp(t::Number) = Tc + (T0-Tc)*ℯ^(-k*t)
+# see shape of temp function
+#x = 1:100
+#y = Temp.(x)
+#plot(x,y)
 
-function calc_cost(rollout::Vector{Int}, prob::problem, voxms::Vector{voxmap})
+function calc_cost(rollout::Vector{Int}, cdata::contourdata, vd::voxdata)
 
     # go from rollout to timestart for each contour
-    timestart = rollout2time(rollout, prob)
+    timestart = rollout2time(rollout, cdata)
 
     # calculate time at each voxel
+    voxtimes = [sum(v.seglen.*(v.segoffset + timestart[v.segcontours]))/sum(v.seglen) for v in vd.maps]
+    println("Empty voxels: ", sum(isnan.(voxtimes)))
+    #replace!(voxtimes, NaN=>0)
 
     # calculate temp difference from voxel below it
+    considered_voxels = (1:length(vd.below))[(below.!=0) .& (.!isnan.(voxtimes))] # cannot calculate cost is no voxel underneath
+    considered_voxels = considered_voxels[.!isnan.(voxtimes[below[considered_voxels]])]
+    Δt = voxtimes[considered_voxels] - voxtimes[below[considered_voxels]]
 
     # calculate stresses at each voxel
 
@@ -454,7 +490,7 @@ function calc_cost(rollout::Vector{Int}, prob::problem, voxms::Vector{voxmap})
 end
 
 
-function construct_cost(prob::problem)
+function construct_cost(cdata::contourdata)
     a = 1
     b = :($a + 5) # using $ 'interpolates' literal expression into the 'quoted' expression
     c = :(5*5 + a)
@@ -471,9 +507,8 @@ function construct_cost(prob::problem)
 end
 
 
-voxels = DataFrames.DataFrame(CSV.File("tensile-1-1.csv"))
 contours = contour.(JSON.parse(open("tensilecontours.json")))
-dt = problem(contours, 5, 5)
-@benchmark random_rollout(dt)
-rl = random_rollout(dt)
-voxms = [voxmap([v.x, v.y, v.z], 1, dt) for v in eachrow(voxels)]
\ No newline at end of file
+cdata = contourdata(contours, 5, 5) # contour data
+vd = voxdata("tensile-1-1.csv", cdata)
+@benchmark random_rollout(cdata)
+rl = random_rollout(cdata)
\ No newline at end of file