From 65a0afbd1dd6f51765459168f18152357ece7c21 Mon Sep 17 00:00:00 2001 From: oloft Date: Sat, 13 Jan 2018 15:20:21 +0100 Subject: [PATCH] Changed notation to \pi. --- examples/example2.jl | 17 ++--------------- src/solve.jl | 32 ++++++++++++++------------------ src/types/PiecewiseQuadratic.jl | 21 +++++++++++---------- src/types/QuadraticPolynomial.jl | 9 --------- test/piecewise_quadratics.jl | 10 +++++----- 5 files changed, 32 insertions(+), 57 deletions(-) diff --git a/examples/example2.jl b/examples/example2.jl index a058956..2dfe6c2 100644 --- a/examples/example2.jl +++ b/examples/example2.jl @@ -1,33 +1,20 @@ using IterTools using Plots - -# The problem Σerror^2 + card(I) -# Heuristic - -include(joinpath(Pkg.dir("DynamicApproximations"),"src","DynamicApproximations.jl")) -## -DA = DynamicApproximations +using DynamicApproximations data = readdlm(joinpath(Pkg.dir("DynamicApproximations"),"examples","data","snp500.txt")) data = data[1:300] N = length(data) - - @time ℓ = DA.compute_discrete_transition_costs(data); - -#@time I, Y, f = brute_force_search(ℓ, K-1); - cost_last = DA.QuadraticPolynomial(1.0, -2*data[end], data[end]^2) start_time = time() gc() -@time Λ = DA.pwq_dp_constrained(ℓ, cost_last, 10, 1.65); +@time Λ = pwq_dp_constrained(ℓ, cost_last, 10, 1.65); println("Time: ", time()-start_time) -#println(counter1, " ", counter2) - for k=3:10 println(k, " : ", DA.recover_optimal_index_set(Λ[k, 1])[3]) end diff --git a/src/solve.jl b/src/solve.jl index 31a650c..d862392 100644 --- a/src/solve.jl +++ b/src/solve.jl @@ -24,9 +24,9 @@ function add_quadratic!{T}(Λ::PiecewiseQuadratic{T}, ρ::QuadraticPolynomial{T} right_endpoint = get_right_endpoint(λ_curr) - Δa = ρ.a - λ_curr.p.a - Δb = ρ.b - λ_curr.p.b - Δc = ρ.c - λ_curr.p.c + Δa = ρ.a - λ_curr.π.a + Δb = ρ.b - λ_curr.π.b + Δc = ρ.c - λ_curr.π.c b2_minus_4ac = Δb^2 - 4*Δa*Δc @@ -142,11 +142,11 @@ end @inline function update_segment_new(λ_prev, λ_curr, ρ) ρ.has_been_used = true # TODO: Could be moved to the second clause - if λ_prev.p === ρ + if λ_prev.π === ρ λ_prev.next = λ_curr.next v1, v2 = λ_prev, λ_curr.next else - λ_curr.p = ρ + λ_curr.π = ρ v1, v2 = λ_curr, λ_curr.next end return v1, v2 @@ -162,7 +162,7 @@ end @inline function update_segment_new_old(λ_prev, λ_curr, ρ, root) # ... λ_prev | (new segment) | λ_curr ... - if λ_prev.p === ρ + if λ_prev.π === ρ λ_curr.left_endpoint = root else ρ.has_been_used = true @@ -185,7 +185,7 @@ end # The second new piece contains a copy of the polynomial in λ_curr ρ.has_been_used = true - λ_2 = PiecewiseQuadratic(λ_curr.p, root2, λ_curr.next) + λ_2 = PiecewiseQuadratic(λ_curr.π, root2, λ_curr.next) λ_1 = PiecewiseQuadratic(ρ, root1, λ_2) λ_curr.next = λ_1 return λ_2, λ_2.next @@ -268,11 +268,8 @@ function pwq_dp_constrained{T}(ℓ::AbstractTransitionCost{T}, V_N::QuadraticPol for ip=i+1:N-m+1 DEBUG && println("(m:$m, i:$i, ip:$ip)") for λ in Λ[m-1, ip] - p = λ.p - #counter1 += - - minimize_wrt_x2(ℓ[i,ip], p, ρ) + minimize_wrt_x2(ℓ[i,ip], λ.π, ρ) DEBUG && println("Obtained ρ = $ρ") @@ -288,7 +285,7 @@ function pwq_dp_constrained{T}(ℓ::AbstractTransitionCost{T}, V_N::QuadraticPol if ρ.has_been_used == true ρ.time_index = ip - ρ.ancestor = p + ρ.ancestor = λ.π ρ = QuadraticPolynomial{T}() ρ.has_been_used = false end @@ -336,10 +333,9 @@ function pwq_dp_regularized{T}(ℓ::AbstractTransitionCost{T}, V_N::QuadraticPol ζ_level_insertion = false for λ in Λ[ip] - p = λ.p #counter1 += 1 - minimize_wrt_x2(ℓ[i,ip], p, ρ) + minimize_wrt_x2(ℓ[i,ip], λ.π, ρ) ρ.c += ζ # add cost for break point @@ -347,7 +343,7 @@ function pwq_dp_regularized{T}(ℓ::AbstractTransitionCost{T}, V_N::QuadraticPol if ρ.has_been_used ρ.time_index = ip - ρ.ancestor = p + ρ.ancestor = λ.π ρ = QuadraticPolynomial{T}() ρ.has_been_used = false @@ -469,9 +465,9 @@ function poly_minus_constant_is_greater{T}(Λ::PiecewiseQuadratic{T}, ρ::Quadra left_endpoint = λ_curr.left_endpoint right_endpoint = get_right_endpoint(λ_curr) - Δa = ρ.a - λ_curr.p.a - Δb = ρ.b - λ_curr.p.b - Δc = (ρ.c - ζ) - λ_curr.p.c + Δa = ρ.a - λ_curr.π.a + Δb = ρ.b - λ_curr.π.b + Δc = (ρ.c - ζ) - λ_curr.π.c b2_minus_4ac = Δb^2 - 4*Δa*Δc diff --git a/src/types/PiecewiseQuadratic.jl b/src/types/PiecewiseQuadratic.jl index bdc2f31..8be6696 100644 --- a/src/types/PiecewiseQuadratic.jl +++ b/src/types/PiecewiseQuadratic.jl @@ -11,9 +11,10 @@ Note that the first element of the list is sometimes assumed to be a list head w The list/element `x` is the last element in the list when `x.next === x` or `x.next.left_endpoint==Inf`. """ mutable struct PiecewiseQuadratic{T} - p::QuadraticPolynomial{T} + π::QuadraticPolynomial{T} left_endpoint::Float64 next::PiecewiseQuadratic{T} + PiecewiseQuadratic{T}() where T = (x=new(); x.left_endpoint=Inf; x.next = x) PiecewiseQuadratic{T}(p, left_endpoint, next) where T = new(p, left_endpoint, next) end @@ -126,7 +127,7 @@ function show{T}(io::IO, Λ::PiecewiseQuadratic{T}) end print(io, "\t : ") - show(io, λ.p) + show(io, λ.π) print(io, "\n") end return @@ -135,7 +136,7 @@ end function (Λ::PiecewiseQuadratic)(x::Number) for λ in Λ if λ.left_endpoint <= x < get_right_endpoint(λ) - return λ.p(x) + return λ.π(x) end end return NaN @@ -145,7 +146,7 @@ function (Λ::PiecewiseQuadratic)(x::AbstractArray) y = zeros(x) for λ in Λ inds = λ.left_endpoint .<= x .< get_right_endpoint(λ) - y[inds] .= λ.p.(x[inds]) + y[inds] .= λ.π.(x[inds]) end return y end @@ -184,7 +185,7 @@ function get_vals(Λ::PiecewiseQuadratic) end y_grid = linspace(left_endpoint, right_endpoint) - vals = λ.p.(y_grid) + vals = λ.π.(y_grid) append!(x, y_grid) append!(y, vals) push!(x_all, y_grid) @@ -216,25 +217,25 @@ function find_minimum(Λ::PiecewiseQuadratic) # TODO: True? f_opt = Inf - p_opt = Λ.p # The segment containing the smallest polynimal + π_opt = similar(Λ.π) x_opt = NaN for λ in Λ - x, f = find_minimum(λ.p) + x, f = find_minimum(λ.π) if f < f_opt f_opt = f x_opt = x - p_opt = λ.p + π_opt = λ.π end end - return p_opt, x_opt, f_opt + return π_opt, x_opt, f_opt end function find_minimum_value(Λ::PiecewiseQuadratic) f_opt = Inf for λ in Λ - _, f = find_minimum(λ.p) + _, f = find_minimum(λ.π) if f < f_opt f_opt = f end diff --git a/src/types/QuadraticPolynomial.jl b/src/types/QuadraticPolynomial.jl index 13779aa..6369c77 100644 --- a/src/types/QuadraticPolynomial.jl +++ b/src/types/QuadraticPolynomial.jl @@ -37,14 +37,6 @@ Base.zero{T<:QuadraticPolynomial}(S::T) = zero(T) -(p1::QuadraticPolynomial, p2::QuadraticPolynomial) = QuadraticPolynomial(p1.a-p2.a, p1.b-p2.b, p1.c-p2.c) -# # x .= y .- z -# # x .= y .+ z -# function Base.broadcast!{T}(op::Function, x::QuadraticPolynomial{T}, y::QuadraticPolynomial{T}, z::QuadraticPolynomial{T}) -# x.a = op(y.a,z.a) -# x.b = op(y.b,z.b) -# x.c = op(y.c,z.c) -# end - ==(p1::QuadraticPolynomial, p2::QuadraticPolynomial) = (p1.a==p2.a) && (p1.b==p2.b) && (p1.c==p2.c) @inline function (p::QuadraticPolynomial)(x) @@ -52,7 +44,6 @@ Base.zero{T<:QuadraticPolynomial}(S::T) = zero(T) end - # Finds the minimum of a positive definite quadratic one variable polynomial # the find_minimum fcn returns (opt_x, opt_val) function find_minimum(p::QuadraticPolynomial) diff --git a/test/piecewise_quadratics.jl b/test/piecewise_quadratics.jl index 1797fe4..61d2cee 100644 --- a/test/piecewise_quadratics.jl +++ b/test/piecewise_quadratics.jl @@ -45,8 +45,8 @@ function verify_piecewise_quadratic(c_mat, show_plot=false) break end - TEST_PRINT_LEVEL > 0 && println("Continutiy diff: ", λ.p(x) - λ.next.p(x)) - @test λ.p(x) ≈ λ.next.p(x) + TEST_PRINT_LEVEL > 0 && println("Continutiy diff: ", λ.π(x) - λ.next.π(x)) + @test λ.π(x) ≈ λ.next.π(x) end # Check that the piecewise quadratic is smaller than all the quadratics @@ -100,11 +100,11 @@ add_quadratic!(pwq2, p3) add_quadratic!(pwq2, p1) @test length(pwq2) == 4 -@test pwq1[1].p === pwq2[1].p == p1 +@test pwq1[1].π === pwq2[1].π == p1 @test pwq1[2].left_endpoint == pwq2[2].left_endpoint == -1 -@test pwq1[2].p === pwq2[2].p == p2 +@test pwq1[2].π === pwq2[2].π == p2 @test pwq1[3].left_endpoint == pwq2[3].left_endpoint == 0 -@test pwq1[3].p === pwq2[3].p +@test pwq1[3].π === pwq2[3].π #---