3. Transport 4: Transshipment¶

ISE 754, Fall 2024

Package Used: No new packages used.

Examples¶

Ex: FTL Location¶

Where shSuld a DC be located to minimize transportation costs, given:

  1. FTLs containing a mix of products A and B shipped P2P from DC to customers in Winston-Salem, Durham, and Wilmington
  2. Each customer receives 20, 30, and 50% of total demand
  3. 100 tons/yr of A shipped FTL P2P to DC from a supplier in Asheville
  4. 380 tons/yr of B shipped FTL P2P to DC from Statesville
  5. Each carton of A weighs 30 lb, and occupies 10 ft3
  6. Each carton of B weighs 120 lb, and occupies 4 ft3
  7. Revenue per loaded truck mile is \$2
  8. Each truck’s cubic and weight capacity is 2,750 ft3 and 25 tons, respectively

image.png

image.png

In [1]:
using DataFrames, Optim

maxpayld(sh, tr) = min(tr.Kwt, sh.s*tr.Kcu/2000)
d1(x1, x2) = length(x1) == length(x2) ? sum(abs.(x1 .- x2)) : 
    error("Inputs not same length.")
D1(X₁, X₂) = [d1(i, j) for i in eachrow(X₁), j in eachrow(X₂)]

tr = (r = 2, Kwt = 25, Kcu = 2750)
uwt = [30, 120]
ucu = [10, 4]
shS = DataFrame(f=[100, 380], s=uwt./ucu)
pct = [20, 30, 50]/100
sh = vcat(shS, DataFrame(f=sum(shS.f)*pct, s=sum(shS.f)./sum(shS.f./shS.s)))
sh.qmax = [maxpayld(sh, tr) for sh in eachrow(sh)]
sh
Out[1]:
5×3 DataFrame
Rowfsqmax
Float64Float64Float64
1100.03.04.125
2380.030.025.0
396.010.434814.3478
4144.010.434814.3478
5240.010.434814.3478

FTL Min TC: Locate DC to min transport cost assuming FTL used for all transport

In [2]:
using Statistics

n = sh.f./sh.qmax                                                # TL/yr
w = n*tr.r                                                       # $/TL-mi
P = [50 150 190 270 420]'
xᵒ = optimize(x -> sum(w .* d1.(x, P)), [mean(P)]).minimizer[1]
Out[2]:
150.00000000228758
In [3]:
hcat(sh, DataFrame(n=n, w=w))
Out[3]:
5×5 DataFrame
Rowfsqmaxnw
Float64Float64Float64Float64Float64
1100.03.04.12524.242448.4848
2380.030.025.015.230.4
396.010.434814.34786.6909113.3818
4144.010.434814.347810.036420.0727
5240.010.434814.347816.727333.4545

(Extension) FTL Min TC w Freq Constr: Include monthly outbound frequency constraint: Outbound shipments must occur at least once each month

In [4]:
idxin, idxout = 1:2, 3:5
tmax = 1/12                                                  # yr/TL
nmin = 1/tmax                                                # TL/yr
n′ = [n[idxin]; max.(n[idxout], nmin)]                       # TL/yr
w′ = n′*tr.r                                                 # $/TL-mi
hcat(sh, DataFrame(n′=n′, w′=w′))
Out[4]:
5×5 DataFrame
Rowfsqmaxn′w′
Float64Float64Float64Float64Float64
1100.03.04.12524.242448.4848
2380.030.025.015.230.4
396.010.434814.347812.024.0
4144.010.434814.347812.024.0
5240.010.434814.347816.727333.4545
In [5]:
xᵒ′ = optimize(x -> sum(w′ .* d1.(x, P)), [mean(P)]).minimizer[1]
Out[5]:
189.9999999907548
In [6]:
# Increase in number of shipments:
sum(n), sum(n′), sum(n′) - sum(n)
Out[6]:
(72.89696969696969, 80.16969696969696, 7.272727272727266)

Impact of Freq Constr: In addition to the change of location, the total transport costs increase due to the increase in the number of shipments. Could also calculate the increase in cost, but this would require having to determine distances, which is tedious if working by hand.

Ex: TLC Location with Transshipment¶

TL Min TLC: Continuing with Ex: FTL Location, locate DC to minimize TLC assuming all shipments are TL and using either (a) uncoordinated inventory or (b) perfect cross-docking for transshipment at the DC.

(a) Uncoordinated inventory (UC):

In [7]:
uval = [300, 450]                                     # unit value ($)
shS.v .= uval./(uwt/2000)
αS = 0                                                # batch prod
shS.α .= αS + 0.5                                     # uncoord DC
shC = DataFrame(f=sum(shS.f)*pct,
                s=sum(shS.f)./sum(shS.f./shS.s), 
                v=sum(shS.f.*shS.v)/sum(shS.f))
shC.α .= 0.5                                          # constant consumption
sh = vcat(shS, shC)
sh.h .= 0.3
sh
Out[7]:
5×5 DataFrame
Rowfsvαh
Float64Float64Float64Float64Float64
1100.03.020000.00.50.3
2380.030.07500.00.50.3
396.010.434810104.20.50.3
4144.010.434810104.20.50.3
5240.010.434810104.20.50.3
In [8]:
# Add 'd' to sh based on location `x`
shd(x) = hcat(sh, DataFrame(d=[d1.(x, P[idxin]); d1.(x, P[idxout])]))
shd([0.])
Out[8]:
5×6 DataFrame
Rowfsvαhd
Float64Float64Float64Float64Float64Float64
1100.03.020000.00.50.350.0
2380.030.07500.00.50.3150.0
396.010.434810104.20.50.3190.0
4144.010.434810104.20.50.3270.0
5240.010.434810104.20.50.3420.0

Re-define functions from Tran 3, where tranchgTL is 0 if the distance of a shipment is 0 to allow an EF to be located at an EF and thus incur no transport charges. This assumes all EFs are actual facilities and not aggregate demand points.

In [9]:
maxpayld(sh, tr) = min(tr.Kwt, sh.s*tr.Kcu/2000)

# TL transport charge: charge is 0 if distance ≈ 0
tranchgTL(q, sh, tr) = isapprox(sh.d, 0., atol=1e-4) ? 0. :
    max(ceil(q/maxpayld(sh, tr))*sh.d*tr.r, 45tr.r/2)

rateLTL(q,s,d,ppi) = ppi*(s^2/8 + 14)/((q^(1/7) * d^(15/29) - 7/2) * (s^2 + 2*s + 14))
tranchgLTL(q, sh, ppi) = max(rateLTL(q, sh.s, sh.d, ppi) * q * sh.d, 
    (ppi/104.2)*(45 + sh.d^(28/19)/1625))

totlogcost(q, c, sh) = sh.f*c/q + sh.α*sh.v*sh.h*q

# Determine independent shipment size that minimizes TLC
function minTLC(sh, tr = nothing, ppi = nothing)
    qᵒ, TLCᵒ, isLTL = nothing, Inf, false
    
    if tr !== nothing    # Check TL option
        qTL = min(sqrt((sh.f * max(tr.r*sh.d, 45tr.r/2))/(sh.α*sh.v*sh.h)),
            maxpayld(sh, tr))
        TLCtl = totlogcost(qTL, tranchgTL(qTL, sh, tr), sh)
        qᵒ, TLCᵒ = qTL, TLCtl
    end

    if ppi !== nothing   # Check LTL option
        qLTL = optimize(q -> totlogcost(q, tranchgLTL(q,sh,ppi), sh),
            150/2000, min(5,650sh.s/2000)).minimizer
        TLCltl = totlogcost(qLTL, tranchgLTL(qLTL, sh, ppi), sh)
        if TLCltl < TLCᵒ
            qᵒ, TLCᵒ, isLTL = qLTL, TLCltl, true
        end
    end

    return (qᵒ = qᵒ, TLCᵒ = TLCᵒ, isLTL = isLTL)
end;
In [10]:
TLCh(x) = sum(sh -> minTLC(sh, tr).TLCᵒ, eachrow(shd(x)))
TLCh([0.])
Out[10]:
105147.58927589661
In [11]:
using CairoMakie
dcf() = display(current_figure())

x = 0:.1:450
lines(x, [TLCh([i]) for i in x])
dcf();
No description has been provided for this image

Continuous optimization has difficulty finding xᵒ = 150 (Statesville).

In [12]:
xˣ = optimize(TLCh, [mean(P)]).minimizer[1]
TLCh(xˣ), xˣ
Out[12]:
(70920.5813093357, 172.4999999998286)
In [13]:
xᵒ = 150
TLCh(xᵒ), xᵒ
Out[13]:
(69206.62608617966, 150)
In [14]:
sh = shd(xᵒ)
sh.qmax = [maxpayld(sh, tr) for sh in eachrow(sh)]
transform!(sh, AsTable(:) => ByRow(row -> minTLC(row, tr)) => AsTable)
Out[14]:
5×10 DataFrame
RowfsvαhdqmaxqᵒTLCᵒisLTL
Float64Float64Float64Float64Float64Int64Float64Float64Float64Bool
1100.03.020000.00.50.31004.1252.5819915491.9false
2380.030.07500.00.50.3025.03.898724386.06false
396.010.434810104.20.50.34014.34782.251056823.49false
4144.010.434810104.20.50.312014.34784.7751914474.8false
5240.010.434810104.20.50.327014.34789.2471228030.3false

(b) Perfect cross-docking (XD): To simplify, we will use the optimal UC location, which may not be optimal for XD. As a result, we can use the final sh from UC since it has a column d with the distance from the optimal UC location to each EF. sh is copied to a new data frame shXD so that the α column can reflect cross-docking (can't just add another column to sh because totlogcost is looking for the α column).

In [15]:
shXD = sh[:, Not([:qᵒ, :TLCᵒ, :isLTL])]
shXD[idxin, :α] .= αS   # Outbound α is unchanged
shXD
Out[15]:
5×7 DataFrame
Rowfsvαhdqmax
Float64Float64Float64Float64Float64Int64Float64
1100.03.020000.00.00.31004.125
2380.030.07500.00.00.3025.0
396.010.434810104.20.50.34014.3478
4144.010.434810104.20.50.312014.3478
5240.010.434810104.20.50.327014.3478
In [16]:
TLC_XDh(t) = sum(sh -> totlogcost(sh.f*t[1], tranchgTL(sh.f*t[1], sh, tr), sh), 
    eachrow(shXD))
t₀ = maximum(shXD.qmax./shXD.f)
TLC_XDh(t₀), t₀, 365.25t₀
Out[16]:
(128668.52865612647, 0.14945652173913043, 54.58899456521739)
In [17]:
tᵒ = optimize(TLC_XDh, [t₀]).minimizer[1]    # 1-D Nelder-Mead
TLC_XDh(tᵒ), tᵒ, 365.25tᵒ
Out[17]:
(-Inf, -4.2075327623855226e302, -1.5368013414613122e305)
In [18]:
tᵒ = optimize(TLC_XDh, 0, t₀).minimizer[1]   # 1-D bounded search (univariate opt)
TLC_XDh(tᵒ), tᵒ, 365.25tᵒ
Out[18]:
(55539.17536298139, 0.03817125453381726, 13.942050718476755)
In [19]:
shXD.qXD = [sh.f*tᵒ for sh in eachrow(shXD)]
shXD.TLC_XD = [totlogcost(sh.f*tᵒ, tranchgTL(sh.f*tᵒ, sh, tr), sh) 
    for sh in eachrow(shXD)];

Determining Optimal Location: Determining the optimal location would turn the problem into a 2-D optimization problem even though the locations are in 1-D, where both the location x and the shipment interval t that minimizes total logistics cost would need to be determined.

Allocated FTL: Determine TLC assuming allocated FTL charge applied to qXD from optimal XD solution.

In [20]:
shXD.TLC_AFTL = [sh.f*tr.r*sh.d/sh.qmax + sh.α*sh.v*sh.h*sh.qXD for sh in eachrow(shXD)]
sum(shXD.TLC_AFTL)
Out[20]:
44594.79979456418

Report results:

In [21]:
hcat(sh[:, Not([:s, :v, :α, :h, :d, :isLTL])],
    shXD[:, Not([:f, :s, :v, :α, :h, :d, :qmax])])
Out[21]:
5×7 DataFrame
RowfqmaxqᵒTLCᵒqXDTLC_XDTLC_AFTL
Float64Float64Float64Float64Float64Float64Float64
1100.04.1252.5819915491.93.817135239.544848.48
2380.025.03.898724386.0614.50510.00.0
396.014.34782.251056823.493.664447649.746089.19
4144.014.34784.7751914474.85.4966614618.310739.6
5240.014.34789.2471228030.39.161128031.622917.5
In [22]:
using Printf

function prt(value, label, units="", digits=2)
    fmtval = split(Printf.format(Printf.Format("%."*string(digits)*"f"),value),".")
    fmtval[1] = reverse(join(Iterators.partition(reverse(fmtval[1]), 3), ","))
    println("$label\t= $(join(fmtval, ".")) $units")
    return value
end

prt( sum(sh.TLCᵒ), "Uncooridinated   ", "(\$/yr)")
prt( sum(shXD.TLC_XD), "Perfect Cross-Docking", "(\$/yr)") 
prt( sum(shXD.TLC_AFTL), "Allocated FTL   ", "(\$/yr)");

Uncooridinated   	= 69,206.63 ($/yr)
Perfect Cross-Docking	= 55,539.18 ($/yr)
Allocated FTL   	= 44,594.80 ($/yr)

Ex: Direct vs. Transshipment¶

In [23]:
using Logjam.DataTools, DataFrames

function dgc(xy₁, xy₂; unit=:mi)
    length(xy₁) == length(xy₂) == 2 || error("Inputs must have length 2.")
    unit in [:mi, :km] || error("Unit must be :mi or :km")

    Δx, Δy = xy₂[1] - xy₁[1], xy₂[2] - xy₁[2]
    a = sind(Δy / 2)^2 + cosd(xy₁[2]) * cosd(xy₂[2]) * sind(Δx / 2)^2
    2 * asin(min(sqrt(a), 1.0)) * (unit == :mi ? 3958.75 : 6371.00)
end

Dgc(X₁, X₂) = [dgc(i, j) for i in eachrow(X₁), j in eachrow(X₂)]

function name2lonlat(name, st, df)
    idx = findfirst(r -> startswith(r[:NAME], name) && r.ST == st, eachrow(df))
    if idx === nothing
        error("'$name', '$st' not found in $df")
    end
    return collect(df[idx, [:LON, :LAT]])
end

# Determine supplier (S) and customer (C) locations
Scity = ["Rochester", "Orlando", "San Antonio"]
Sst = ["NY", "FL", "TX"]
SXY = vcat([name2lonlat(i, j, usplace())' for (i,j) in zip(Scity, Sst)]...)
Ccity = ["Portland", "Tucson", "Denver", "Milwaukee"]
Cst = ["OR", "AZ", "CO", "WI"]
CXY = vcat([name2lonlat(i, j, usplace())' for (i,j) in zip(Ccity, Cst)]...)
Out[23]:
4×2 Matrix{Float64}:
 -122.65    45.537
 -110.871   32.153
 -104.881   39.7619
  -87.9667  43.0633
In [24]:
# Data: 3 Different Products Supplied to 4 Customers

ppiTL = 131.4  # Jan 2018 (P)
ppiLTL = 179.4 # Jan 2018 (P)
tr = (r = 2ppiTL/102.7, Kwt = 25, Kcu = 2750)
fS = [500, 75, 300]            # demand (ton/yr)
pct = [20, 35, 25, 20]/100     # customer demand percentage
sS = [32, 3, 12]               # density (lb/ft^3)
v = [50_000, 25_000, 10_000]   # product cost (in $/ton)
αS = 0                         # batch production at suppliers
αC = 0.5;                      # constant consumption at customers

1. Direct Shipments

In [25]:
F = fS[:] * pct[:]'   # (3 x 1) * (1 x 4) = (3 x 4) matrix
Out[25]:
3×4 Matrix{Float64}:
 100.0  175.0   125.0   100.0
  15.0   26.25   18.75   15.0
  60.0  105.0    75.0    60.0
In [26]:
S = repeat(sS[:], 1, size(CXY, 1))
Out[26]:
3×4 Matrix{Int64}:
 32  32  32  32
  3   3   3   3
 12  12  12  12
In [27]:
V = repeat(v[:], 1, size(CXY, 1))
D = Dgc(SXY, CXY) * 1.2
Out[27]:
3×4 Matrix{Float64}:
 2641.83  2346.59   1709.51   626.061
 3043.26  2136.29   1867.03  1295.17
 2057.81   906.461   957.05  1327.78
In [28]:
sh = DataFrame(f=F[:], s=S[:], α=αS+αC, v=V[:], h=0.3, d=D[:])
sh.qmax = [maxpayld(sh, tr) for sh in eachrow(sh)]
transform!(sh, AsTable(:) => ByRow(row -> minTLC(row, tr, ppiLTL)) => AsTable)
sh.t = 365.25sh.qᵒ./sh.f
sh[:, Not(:h)]
Out[28]:
12×10 DataFrame
RowfsαvdqmaxqᵒTLCᵒisLTLt
Float64Int64Float64Int64Float64Float64Float64Float64BoolFloat64
1100.0320.5500002641.8325.01.999051.14292e5true7.30154
215.030.5250003043.264.1254.12543786.7false100.444
360.0120.5100002057.8116.514.513143539.3false88.3485
4175.0320.5500002346.5925.03.102251.77379e5true6.47485
526.2530.5250002136.294.1254.12550256.0false57.3964
6105.0120.510000906.46116.512.742438227.2false44.3253
7125.0320.5500001709.5125.02.072471.17048e5true6.05576
818.7530.5250001867.034.1254.12537184.9false80.355
975.0120.510000957.0516.511.065733197.2false53.8901
10100.0320.550000626.06125.01.238866432.0true4.52473
1115.030.5250001295.174.1253.64127307.5false88.6583
1260.0120.5100001327.7816.511.657934973.7false70.9675
In [29]:
using Statistics
out = DataFrame(Method="Direct Shipments", TLC=round(Int, sum(sh.TLCᵒ)), 
    t̄=365.25mean(sh.qᵒ./sh.f), LTL=sum(sh.isLTL))
Out[29]:
1×4 DataFrame
RowMethodTLCt̄LTL
StringInt64Float64Int64
1Direct Shipments78362350.72854

2. Uncoordinated: Use a DC located at Memphis, TN, with no coordination

In [30]:
DCcity = "Memphis"
DCXY = name2lonlat(DCcity, "TN", usplace())'

shS = DataFrame(f=fS, s=sS, α=αS+0.5, v=v)

shC = DataFrame(f=sum(shS.f)*pct,
                s=sum(shS.f)./sum(shS.f./shS.s),
                α=αC,
                v=sum(shS.f.*shS.v)/sum(shS.f))

sh = vcat(shS, shC)
sh.h .= 0.3
sh.d = Dgc(DCXY, [SXY; CXY])[:] * 1.2
sh.qmax = [maxpayld(sh, tr) for sh in eachrow(sh)]
transform!(sh, AsTable(:) => ByRow(row -> minTLC(row, tr)) => AsTable)
sh.t = 365.25sh.qᵒ./sh.f
sh[:, Not(:h)]
Out[30]:
7×10 DataFrame
RowfsαvdqmaxqᵒTLCᵒisLTLt
Float64Float64Float64Float64Float64Float64Float64Float64BoolFloat64
1500.032.00.550000.01036.0125.013.29421.99414e5false9.71144
275.03.00.525000.0827.4594.1254.12553966.8false20.0888
3300.012.00.510000.0760.31416.516.560124.1false20.0888
4175.013.33330.534142.92220.2318.333313.93321.42715e5false29.0805
5306.2513.33330.534142.91460.9618.333314.95161.53148e5false17.8321
6218.7513.33330.534142.91053.3218.333310.72961.09902e5false17.9154
7175.013.33330.534142.9671.8918.33337.6647878509.3false15.9975
In [31]:
out = vcat(out, DataFrame(Method="Uncoord Inv at DC", TLC=round(Int, sum(sh.TLCᵒ)), 
    t̄=365.25mean(sh.qᵒ./sh.f), LTL=sum(sh.isLTL)))
Out[31]:
2×4 DataFrame
RowMethodTLCt̄LTL
StringInt64Float64Int64
1Direct Shipments78362350.72854
2Uncoord Inv at DC79777918.67350

3. Cross-Docking: Use a DC located at Memphis, TN, with perfect cross-docking

In [32]:
shXD = sh[:, Not([:qᵒ, :TLCᵒ, :isLTL, :t])]
shXD[idxin, :α] .= αS   # Outbound α is unchanged
TLC_XDh(t) = sum(sh -> totlogcost(sh.f*t[1], tranchgTL(sh.f*t[1], sh, tr), sh), 
    eachrow(shXD))
t₀ = maximum(shXD.qmax./shXD.f)
tᵒ = optimize(TLC_XDh, 0, t₀).minimizer[1]   # 1-D bounded search (univariate opt)
shXD.qXD = [sh.f*tᵒ for sh in eachrow(shXD)]
shXD.TLC_XD = [totlogcost(sh.f*tᵒ, tranchgTL(sh.f*tᵒ, sh, tr), sh) 
    for sh in eachrow(shXD)]
shXD.t = 365.25shXD.qXD./shXD.f
shXD[:, Not(:h)]
Out[32]:
7×9 DataFrame
RowfsαvdqmaxqXDTLC_XDt
Float64Float64Float64Float64Float64Float64Float64Float64Float64
1500.032.00.050000.01036.0125.025.053021.018.2625
275.03.00.025000.0827.4594.1253.7542347.818.2625
3300.012.00.510000.0760.31416.515.061411.518.2625
4175.013.33330.534142.92220.2318.33338.751.5844e518.2625
5306.2513.33330.534142.91460.9618.333315.31251.53191e518.2625
6218.7513.33330.534142.91053.3218.333310.93751.09922e518.2625
7175.013.33330.534142.9671.8918.33338.7579198.618.2625
In [33]:
out = vcat(out, DataFrame(Method="Cross-Docking at DC", 
        TLC=round(Int, sum(shXD.TLC_XD)), t̄=mean(shXD.t), LTL=0))
Out[33]:
3×4 DataFrame
RowMethodTLCt̄LTL
StringInt64Float64Int64
1Direct Shipments78362350.72854
2Uncoord Inv at DC79777918.67350
3Cross-Docking at DC65753218.26250
In [34]:
using GeoMakie, Logjam.MapTools

XY = vcat(SXY, DCXY, CXY)
fig, ax = makemap(eachcol(XY)..., xexpand=0.25, yexpand=0.25)

h = []
push!(h, scatter!(ax, eachcol(SXY)..., 
        marker=:circle, color=:red, label="Suppliers"))
push!(h, scatter!(ax, eachcol(DCXY)..., 
        marker=:utriangle, markersize=14, label="DC"))
push!(h, scatter!(ax, eachcol(CXY)..., marker=:rect, label="Customers"))

text!(ax, eachcol(vcat(SXY, DCXY, CXY))..., 
    text=vcat(Scity, DCcity, Ccity); aligntext(eachcol(vcat(SXY, DCXY, CXY))...)...)

Legend(fig[1,2], h, [x.label for x in h])
fig
Out[34]:
No description has been provided for this image