Receiving the Error: "MethodError: no method matching eps(::Int64)" When Using NonlinearProblem Solver

New to Julia
, , , ,
1

Hello, I am attempting to solve a non-linear optimization problem using the NonlinearProblem solver with the TrustRegion method. I kept getting the error saying there was no method for Float64(::ForwardDiff Dual etc). I figured that to resolve this I would make my primitive definitions type Real instead of type Float64. However I am now getting the error

ERROR: MethodError: no method matching eps(::Int64)
The function `eps` exists, but no method is defined for this combination of argument types.

I am getting this error even though I do not explicitly use the eps() function anywhere in my code. I assume that it is a part of the NonlinearProblem() solver but I do not know how it is used within the solver, or how to correct this issue.
I would really appreciate some help as I am still quite new to Julia. Thank you!

The relevant code is the following

using LinearAlgebra, Statistics, Distributions,Plots,FixedEffectModels,DataFrames, 
IJulia, Optimization,Optim,XLSX,Parameters,Printf,NonlinearSolve,RecursiveArrayTools, ForwardDiff


@with_kw struct primitives
    γ::Real = .0016093 ##  Parameter on Increasing Component of Costs 
    τ ::Real =.25 ## Tax Rate
    α::Real = -0.0097 ##Price coefficient 
    β::AbstractArray{Real,1} = [.011,0] ## Linear Characteristic Coefficients
    number_markets::Int64 = 16
    number_products::Int64 = 12
    pop::Int64 = 15000 ##Number of agents in the given market. 
    capacity_threshold::Real = .85
end 

mutable struct results
    prices::AbstractArray{Real,2}
    # prices::ArrayPartition{AbstractArray{Real,2}}
    # new_prices::ArrayPartition(AbstractArray{Real,2})
    new_prices::AbstractArray{Real,2} 
    avg_rooms::AbstractArray{Real,2} 
    shares::AbstractArray{Float64,2} 
    some_meals::AbstractArray{Float64,2}  ##Only Need One Because its the same throughout markets
    δ::AbstractArray{Float64,2}    
    δ_noprice::AbstractArray{Float64,2} 
    markups::AbstractArray{Real,2} 
    mc_linear::AbstractArray{Float64,2} 
    mc_increasing::AbstractArray{Float64,2} 
    mc_total::AbstractArray{Float64,2} 
    deviation::AbstractArray{Float64,2} 
    choice_probability::AbstractArray{Float64,3} 
    rc_price_draw::AbstractArray{Float64,2} 
    rc_constant_draw::AbstractArray{Float64,2} 
    rc_meals_draw::AbstractArray{Float64,2} 
    rc_rooms_draw::AbstractArray{Float64,2} 
    quantity_demanded::AbstractArray{Float64,2} 
    demand_prime::AbstractArray{Float64,2} 
    quantity_supplied::AbstractArray{Float64,2} 
    utility_numerators::AbstractArray{Float64,3} 
    consumer_surplus::AbstractArray{Float64,1} 
    weights::AbstractArray{Float64,2} 
    choice_probability_prime::AbstractArray{Float64,3}
    market_size::AbstractArray{Float64,2}
    count_establishments::AbstractArray{Float64,2}
    total_product_capacity::AbstractArray{Float64,2}
    supply_residuals::AbstractArray{Float64,2}
    accom_foc::AbstractArray{Real,2}
    capacity_constraint::AbstractArray{Float64,2}
    γ::Real  ##  Parameter on Increasing Component of Costs 
    τ ::Real ## Tax Rate
    α::Real ##Price coefficient 
    β::AbstractArray{Real,1}
    iter_index::Int64
end 

function Initialize()
    prim= primitives()
    prices = zeros(prim.number_markets,prim.number_products)
    new_prices = zeros(prim.number_markets,prim.number_products)
    avg_rooms = zeros(prim.number_markets,prim.number_products)
    shares = zeros(prim.number_markets,prim.number_products)
    some_meals = zeros(prim.number_markets, prim.number_products)
    δ = zeros(prim.number_markets,prim.number_products)
    δ_noprice = zeros(prim.number_markets,prim.number_products)
    markups = zeros(prim.number_markets,prim.number_products)
    mc_linear = zeros(prim.number_markets,prim.number_products)
    mc_increasing = zeros(prim.number_markets,prim.number_products)
    mc_total=  zeros(prim.number_markets,prim.number_products)
    deviation = zeros(prim.number_markets,prim.number_products)
    choice_probability = zeros(prim.number_markets,prim.number_products,prim.pop)
    rc_price_draw = zeros(prim.number_markets,prim.pop)
    rc_constant_draw = zeros(prim.number_markets,prim.pop)
    rc_meals_draw = zeros(prim.number_markets,prim.pop)
    rc_rooms_draw = zeros(prim.number_markets,prim.pop)
    quantity_demanded = zeros(prim.number_markets,prim.number_products)
    demand_prime = zeros(prim.number_markets,prim.number_products)
    quantity_supplied = zeros(prim.number_markets,prim.number_products)
    utility_numerators = zeros(prim.number_markets,prim.number_products,prim.pop)
    consumer_surplus = zeros(prim.number_markets)
    weights = zeros(prim.number_markets,prim.pop)
    choice_probability_prime = zeros(prim.number_markets,prim.number_products,prim.pop)
    market_size = zeros(prim.number_markets,prim.number_products)
    count_establishments = zeros(prim.number_markets,prim.number_products)
    total_product_capacity = zeros(prim.number_markets,prim.number_products)
    supply_residuals = zeros(prim.number_markets,prim.number_products)
    accom_foc = zeros(prim.number_markets,prim.number_products)
    capacity_constraint = zeros(prim.number_markets,prim.number_products)
    iter_index = 1


    γ = .0016093 ##  Parameter on Increasing Component of Costs 
    τ =.70 ## Tax Rate
    α = -0.0097 ##Price coefficient 
    β = [0.011,0]  





    res = results(prices,new_prices,avg_rooms,shares,some_meals,δ,δ_noprice,markups,
    mc_linear,mc_increasing,mc_total,deviation,choice_probability,rc_price_draw,
    rc_constant_draw,rc_meals_draw,rc_rooms_draw,quantity_demanded,demand_prime,
    quantity_supplied,utility_numerators,consumer_surplus,weights,choice_probability_prime,
    market_size,count_establishments,total_product_capacity,supply_residuals,
    accom_foc,capacity_constraint,γ,τ,α,β,iter_index)
   prim,res

end

prim, res = Initialize()

function compute_equilibrium(prices,res) 
    #  @unpack α,β,τ,γ = prim 
    # @unpack demand_prime = res 
    # out_test = []
    # out_test = test_random_coeff(prim,res,df_deltas_o,market_ids_vector,df_agents,loop_index)

    res.iter_index =2

    loop_index = res.iter_index

    
    res.markups[loop_index,:] = (res.shares[loop_index,:]./(1+ res.τ)).*(1 ./ (res.count_establishments[loop_index,:]).*res.demand_prime[loop_index,:])
    #pre_deviation = (((res.shares[loop_index,:].*res.market_size[loop_index,:]) .-  (res.γ.*res.total_product_capacity[loop_index,:])) ./ res.count_establishments[loop_index,:])
     res.deviation[loop_index,:] = max.((((res.shares[loop_index,:].*res.market_size[loop_index,:]) .-  (res.γ.*res.total_product_capacity[loop_index,:])) ./ res.count_establishments[loop_index,:]),0)
    #  deviation = max.(pre_deviation,0)
     println(res.deviation[loop_index,:])
    # println(prices)
    # print(res.deviation[loop_index,:])
    println("So far everything is running!! :D. The deviation has now been calculated.")
    # pre_accom_foc = (prices ./ (1+res.τ)) .+ res.markups[loop_index,:] .- (res.supply_residuals[loop_index,:] .+ res.mc_linear[loop_index,:] .+ (res.γ.*res.deviation[loop_index,:]))
    # println("Pre Accom FOC dims are: ",size(pre_accom_foc))
    println(size(prices))
    # xyz = (prices ./ (1+res.τ)) .+ markups[loop_index,:] .- (res.supply_residuals[loop_index,:] .+ res.mc_linear[loop_index,:] .+ (res.γ.*res.deviation[loop_index,:]))
    # print(xyz)
    println(prices)
    
    res.accom_foc[loop_index,:]  = (prices ./ (1+res.τ)) .+ res.markups[loop_index,:] .- (res.supply_residuals[loop_index,:] .+ res.mc_linear[loop_index,:] .+ (res.γ.*res.deviation[loop_index,:]))
   println(res.accom_foc)
    println("The Accommodation FOC has been calculated and Added to the Results Struct")

    # pre_capacity_constraint = res.shares[loop_index,:] #.*res.market_size[loop_index,:]) #.-  (res.total_product_capacity[loop_index,:])
    # println(size(pre_capacity_constraint))
    res.capacity_constraint[loop_index,:] =  (res.shares[loop_index,:] .*res.market_size[loop_index,:]) .-  (res.total_product_capacity[loop_index,:])
    println("The Capacity Constraint has been calculated and Added to the Results Struct!!")

    println("This Portion of Finding an Equilibrium is Completed!!")
    error_vector = Float64.(res.accom_foc[loop_index,:] .+ res.capacity_constraint[loop_index,])
    println("I have now obtained the Error Vector!!")
    


    return error_vector
end
price_vector0 = df_deltas_o[df_deltas_o.market_ids.==current_market,"prices"]
price_vector3::AbstractArray{Real,1} = price_vector0

prob = NonlinearProblem(compute_equilibrium, price_vector3, res)
sol = solve(prob, TrustRegion())

The stack trace is the following.

ERROR: MethodError: no method matching eps(::Int64)
The function `eps` exists, but no method is defined for this combination of argument types.

Closest candidates are:
  eps(::Type{Dates.DateTime})
   @ Dates C:\Users\AppData\Local\Programs\Julia-1.11.1\share\julia\stdlib\v1.11\Dates\src\types.jl:452
  eps(::Type{Dates.Date})
   @ Dates C:\Users\AppData\Local\Programs\Julia-1.11.1\share\julia\stdlib\v1.11\Dates\src\types.jl:453
  eps(::Type{Dates.Time})
   @ Dates C:\Users\mcket\AppData\Local\Programs\Julia-1.11.1\share\julia\stdlib\v1.11\Dates\src\types.jl:454
  
The 
Stacktrace:
  [1] get_tolerance(::Nothing, ::Type{Real})
    @ NonlinearSolveBase C:\Users\\.julia\packages\NonlinearSolveBase\sduQ7\src\common_defaults.jl:39
  [2] get_tolerance(::Vector{Real}, η::Nothing, ::Type{Real})
    @ NonlinearSolveBase C:\Users\\.julia\packages\NonlinearSolveBase\sduQ7\src\common_defaults.jl:43
  [3] init_termination_cache(prob::NonlinearProblem{…}, abstol::Nothing, reltol::Nothing, du::Vector{…}, u::Vector{…}, tc::AbsNormSafeBestTerminationMode{…}, ::Val{…})
    @ NonlinearSolveBase C:\Users\.julia\packages\NonlinearSolveBase\sduQ7\src\termination_conditions.jl:291
  [4] init_termination_cache(prob::NonlinearProblem{…}, abstol::Nothing, reltol::Nothing, du::Vector{…}, u::Vector{…}, ::Nothing, callee::Val{…})
    @ NonlinearSolveBase C:\Users\\.julia\packages\NonlinearSolveBase\sduQ7\src\termination_conditions.jl:283
  [5] __init(::NonlinearProblem{…}, ::GeneralizedFirstOrderAlgorithm{…}; stats::SciMLBase.NLStats, alias_u0::Bool, maxiters::Int64, abstol::Nothing, reltol::Nothing, maxtime::Nothing, termination_condition::Nothing, internalnorm::Function, linsolve_kwargs::@NamedTuple{}, kwargs::@Kwargs{})
    @ NonlinearSolveFirstOrder C:\Users\\.julia\packages\NonlinearSolveFirstOrder\XQs9h\src\solve.jl:154
  [6] __init
    @ C:\Users\.julia\packages\NonlinearSolveFirstOrder\XQs9h\src\solve.jl:119 [inlined]      
  [7] __solve(::NonlinearProblem{…}, ::GeneralizedFirstOrderAlgorithm{…}; kwargs::@Kwargs{})        
    @ NonlinearSolveBase C:\Users\\.julia\packages\NonlinearSolveBase\sduQ7\src\solve.jl:5     
  [8] __solve
    @ C:\Users\\.julia\packages\NonlinearSolveBase\sduQ7\src\solve.jl:1 [inlined]
  [9] #solve_call#44
    @ C:\Users\\.julia\packages\DiffEqBase\HW4ge\src\solve.jl:632 [inlined]
 [10] solve_call
    @ C:\Users\\.julia\packages\DiffEqBase\HW4ge\src\solve.jl:589 [inlined]
 [11] #solve_up#53
    @ C:\Users\t\.julia\packages\DiffEqBase\HW4ge\src\solve.jl:1105 [inlined]
 [12] solve_up
    @ C:\Users\\.julia\packages\DiffEqBase\HW4ge\src\solve.jl:1099 [inlined]
 [13] #solve#52
    @ C:\Users\\.julia\packages\DiffEqBase\HW4ge\src\solve.jl:1093 [inlined]
 [14] solve(prob::NonlinearProblem{…}, args::GeneralizedFirstOrderAlgorithm{…})
    @ DiffEqBase C:\Users\\.julia\packages\DiffEqBase\HW4ge\src\solve.jl:1083
 [15] top-level scope
    @ c:\Users\Documents\n\Demand_Tasks\Creating_BLP_Dataset\Code\firmbehavior_pemodel_v2.jl:376
Some type information was truncated. Use `show(err)` to see complete types.
2

I am not an expert of any of this but I think it is the way that you use types.

ForwardDiff works by passing a special Dual type to your function instead of a Float. When this type is run through the function it outputs both the value and the derivative. Having to strict type definitions and using explicit type casting will break this.

error_vector = Float64.(res.accom_foc[loop_index,:] .+ res.capacity_constraint[loop_index,])

Is an example of a line that probably won’t work with ForwardDiff due to the explicit conversion.

I also think that integers can cause some issues if you want to do a smooth optimization. You could try to replace some int values with float. E.g. replace 0 with 0.0 in this line.

res.deviation[loop_index,:] = max.((((res.shares[loop_index,:].*res.market_size[loop_index,:]) .-  (res.γ.*res.total_product_capacity[loop_index,:])) ./ res.count_establishments[loop_index,:]),0.0)

Finally, I don’t think this is the best way to use struct in Julia. I think there is a simpler and more clear solution to this problem without structs. The structs here are just used as a one time container of data.

1 Like
3

I agree with lupemba. There is almost always surprises with a Vector{Real}. I don’t see exactly what is the problem in your code, but note that a Vector{Real} can contain both Float64 and Int, as well as other things:

julia> v = Real[2.34, 0, pi]
3-element Vector{Real}:
 2.34
 0
 π = 3.1415926535897...

julia> typeof.(v)
3-element Vector{DataType}:
 Float64
 Int64
 Irrational{:π}

I guess somehow an Int has sneaked into one of your vectors.

I realize that you use Real to accommodate the Dual type of ForwardDiff, but it is also generally a bad idea to have structs with abstract data types like AbstractVector/AbstractArray. The reason is mainly that it is slow.

If you need temporary storage inside your objective function which is compatible with ForwardDiff, there are tools for it in Home · PreallocationTools.jl.

1 Like
4

Thanks so much for this response. I will go ahead and start looking into some alternatives to explicit type conversion. In general have you found that it is best to just not specify types like how I have specified mine.

I will also look into some better ways to use structs. This is the only way I have learned so far so it is useful to know that this is not the best format to use. Do you think that reading through the introductory documentation would be my best bet for learning about alternatives to using structs as I have done?

5

Thanks for this reply. That is very useful to know about Vector{Real}. I was not aware that it was known to be tricky. I also did not know that AbstractVector types generally have trouble with structs. I will definitely check out pre-allocation tools a little later today, that may be exactly what I need!

In general is it not wise to specify data types in structs? I am still trying to figure out the best way to incorporate them into my code.

6

You could look at the examples in Optimization.jl documentation since you are optimizing something. Solving the Rosenbrock Problem in >10 Ways · Optimization.jl

Here are some resource on dual numbers if you want to know more Forward-Mode Automatic Differentiation (AD) via High Dimensional Algebras - MIT Parallel Computing and Scientific Machine Learning (SciML)

The package documentation is often a good resource

1 Like
7

On the contrary. The fields of structs which are used in the hard working part of your program should have a specified datatype. But it should be concrete datatypes, not abstract ones. On the other hand, the arguments to functions need not have a datatype. This is typically only necessary if you do multiple dispatch (e.g. write separate versions of the function for Float64 and Float32).

That is, your @kwdef struct is fine, it’s not used in the hard working part of your program. @kwdef struct are typically used for global parameters in some initialization.

But the results struct can be problematic. It’s used all the time in your compute_equilibrium function and should have concrete types like Vector{Float64}, Array{Float64, 3}, Float64 etc.

The reason is that when the compile_equilibrium function is compiled, the compiler only knows the types of the actual arguments, and the declared types in the struct. The compiler does not know what kind of array is actually stored in an AbstractArray{Float64, 2} field, so it can’t create efficient code. It might be a Matrix{Float64}, or a SparseMatrix{Float64} or even some kind of matrix type which hasn’t been created yet. So it must insert code to figure out the actual type at run time, every time you access the field, which is wasting a mountain of cpu-cycles.

As you have seen, using concrete types creates problems for ForwardDiff, the differentiation package. That’s one of the reasons the PreallocationTools package was created.

Apart from all that, do you need all the intermediate results, i.e. for every loop_index, or do you just need to store thing temporarily inside the compute_equilibrium? And fetch them for the final loop index? In such a case there are packages for quick temporary allocations (mainly Bumper.jl), and you could get rid of the pre allocations altogether. If needed, you could recompute the intermediate arrays once you have obtained the solution.

2 Likes
8

Thank you for these recommendations! They helped a lot and I was able to update my code and resolve the issues.

1 Like
9

This explanation was incredible helpful, thank you. I was able to resolve the issues I was dealing with.