那些遇到过的问题

Julia 中的加速递归系统

amr*_*ods 3 recursion julia

我需要加速这个递归系统。任何人都可以帮忙吗？

我注释了所有变量并且只使用了 Float64-consistent 函数。我还尝试通过 Memoize.jl 使用记忆化（因为系统是递归的），但是我必须在优化算法中计算 U、V 和 Z 数千次，并且我的计算机很快就会耗尽内存。

const Tage = 60.0
const initT = 1975.0
const Ttime = 1990.0
const K = 6.0
const ? = 0.95
const s = 0.5

?0 = [
    0.0011
 -0.0045
 -0.0075
 -0.0013
  3.60   
  0.2587
 -0.0026  
 -0.0528
  0.0060
  1.3772
  0.2932
 -0.0030
  0.0021
  0.0044
 10.2593    
  0.7498
  1.1765
]./10000

function wagem(t::Float64, agem::Float64, edum::Float64, ?::Vector{Float64})
    ?[5] + ?[6]*agem + ?[7]*agem^2 + ?[8]*edum + ?[9]*edum^2
end

function wagef(t::Float64, agef::Float64, eduf::Float64, ?::Vector{Float64})
    ?[10] + ?[11]*agef + ?[12]*agef^2 + ?[13]*eduf + ?[14]*eduf^2
end

function ?(agem::Float64, edum::Float64, agef::Float64, eduf::Float64, ?::Vector{Float64})
    ?[1]*agem*agef + ?[2]*agem*eduf + ?[3]*edum*agef + ?[4]*edum*eduf
end

function z(t::Float64, agem::Float64, edum::Float64, agef::Float64, eduf::Float64, ?::Vector{Float64})
    ?(agem, edum, agef, eduf, ?) + wagem(t, agem, edum, ?) + wagef(t, agef, eduf, ?)
end

function dc(?::Vector{Float64})
    ret::Float64 = 0.0
    ret = ?[15]
    return ret
end

function Z(t::Float64, agem::Float64, edum::Float64, agef::Float64, eduf::Float64, ?::Vector{Float64})
    ret::Float64 = 0.0
    if agem >= Tage || agef >= Tage # terminal age conditions
        ret = 0.0
    elseif t >= Ttime # terminal time condition
        ret = 1.5
    else
        ret = log(
        exp(z(t, agem, edum, agef, eduf, ?) + ?*Z(t+1, agem+1, edum, agef+1, eduf, ?))
        + exp(-dc(?) + ?*(V(t+1, agem+1, edum, ?) + U(t+1, agef+1, eduf, ?)))
        )
    end
    return ret
end

function V(t::Float64, agem::Float64, edum::Float64, ?::Vector{Float64})
    ret::Float64 = 0.0
    if agem >= Tage
        ret = 0.0
    elseif t >= Ttime
        ret = 1.5
    else
        suma::Float64 = 0.0
        for agef in 16.0:Tage-1, eduf in 1.0:K
            suma += exp(s*z(t, agem, edum, agef, eduf, ?) + wagem(t, agem, edum, ?) + ?*(s*(Z(t+1, agem+1, edum, agef+1, eduf, ?) - V(t+1, agem+1, edum, ?) - U(t+1, agef+1, eduf, ?)) + V(t+1, agem+1, edum, ?)))
        end
        ret = log(
                exp(wagem(t, agem, edum, ?) + ?*(V(t+1, agem+1, edum, ?))) + suma
                )
    end
    return ret
end

function U(t::Float64, agef::Float64, eduf::Float64, ?::Vector{Float64})
    ret::Float64 = 0.0
    if agef >= Tage
        ret = 0.0
    elseif t >= Ttime
        ret = 1.5
    else
        suma::Float64 = 0.0
        for agem in 16.0:Tage-1, edum in 1.0:K
            suma += exp((1-s)*z(t, agem, edum, agef, eduf, ?) + wagef(t, agef, eduf, ?) + ?*((1-s)*(Z(t+1, agem+1, edum, agef+1, eduf, ?) - V(t+1, agem+1, edum, ?) - U(t+1, agef+1, eduf, ?)) + U(t+1, agef+1, eduf, ?)))
        end
        ret = log(
        exp(wagef(t, agef, eduf, ?) + ?*(U(t+1, agef+1, eduf, ?))) + suma
        )
    end
    return ret
end
Run Code Online (Sandbox Code Playgroud)

我希望能够非常快速地计算 U(1984.0, 17.0, 3.0, ?0) 或 V(1975.0, 16.0, 1.0, ?0) 之类的东西。

任何帮助，将不胜感激。

Rez*_*lan 6

我认为代码就像一个怪物！
之后你已经有 3 个递归函数U(), V(), Z()，U()调用V()反之亦然，然后V()调用Z()它本身调用U()......
尝试重写它并尽可能防止递归调用，这将导致更有效的解决方案。检查这个->递归与迭代
之后，您必须在外部进行不必要的重复调用U()和V()内部for循环。

第一步：临时变量

function U(t::Float64, agef::Float64, eduf::Float64, ?::Vector{Float64}) 
    ret::Float64 = 0.0
    if agef >= Tage
        ret = 0.0
    elseif t >= Ttime
        ret = 1.5
    else
        suma::Float64 = 0.0
        U1 = U(t+1, agef+1, eduf, ?) # nothing to do with following loop, so I take it outside
        wagef1 = wagef(t, agef, eduf, ?) # same as above comment
        for agem in 16.0:Tage-1, edum in 1.0:K
            suma += exp((1-s)*z(t, agem, edum, agef, eduf, ?) + wagef1 + ?*((1-s)*(Z(t+1, agem+1, edum, agef+1, eduf, ?, ui) - V(t+1, agem+1, edum, ?) - U1) + U1))
        end
        ret = log(
                  exp(wagef1 + ?*(U1)) + suma
                  )
    end

    return ret
end

julia> @time U(1984.0, 17.0, 3.0, t0) # a test with const Tage = 19.0
  0.869675 seconds (102.77 k allocations: 2.148 MB, 0.57% gc time)
3.785563216949393
Run Code Online (Sandbox Code Playgroud)

与原始代码相比，上述两者的编辑U() and V()使我的运行速度提高了 28 倍。

第二步：使用正确的类型

接下来是Float64数据类型的mall-usage ，IMO 只需?,?,s要是Float64类型，所有其他变量都应该声明为整数。

julia> @time U(1984, 17, 3, t0) # a test with const Tage = 19
  0.608982 seconds (125.77 k allocations: 2.530 MB, 0.89% gc time)
3.785563216949393
Run Code Online (Sandbox Code Playgroud)

现在它的速度提高了 30%。与之前的代码比较

第三步 Memoize.jl

although already some improvements have been happened for bigger problems they are not enough, but with using the right data types, in the previous step now `Memoize.jl` will be usable:  


using Memoize
const Tage = 60%Int
const initT = 1975%Int
const Ttime = 1990%Int
const K = 6%Int
const ? = 0.95
const s = 0.5

t0 = [
    0.0011
 -0.0045
 -0.0075
 -0.0013
  3.60   
  0.2587
 -0.0026  
 -0.0528
  0.0060
  1.3772
  0.2932
 -0.0030
  0.0021
  0.0044
 10.2593    
  0.7498
  1.1765
]./10000

function wagem(t::Int, agem::Int, edum::Int, ?::Vector{Float64})
    ?[5] + ?[6]*agem + ?[7]*agem^2 + ?[8]*edum + ?[9]*edum^2
end

function wagef(t::Int, agef::Int, eduf::Int, ?::Vector{Float64})
    ?[10] + ?[11]*agef + ?[12]*agef^2 + ?[13]*eduf + ?[14]*eduf^2
end

function ?(agem::Int, edum::Int, agef::Int, eduf::Int, ?::Vector{Float64})
    ?[1]*agem*agef + ?[2]*agem*eduf + ?[3]*edum*agef + ?[4]*edum*eduf
end

function z(t::Int, agem::Int, edum::Int, agef::Int, eduf::Int, ?::Vector{Float64})
    ?(agem, edum, agef, eduf, ?) + wagem(t, agem, edum, ?) + wagef(t, agef, eduf, ?)
end

function dc(?::Vector{Float64})
    ret::Float64 = 0.0
    ret = ?[15]
    return ret
end

function Z(t::Int, agem::Int, edum::Int, agef::Int, eduf::Int, ?::Vector{Float64})
    ret::Float64 = 0.0
    if agem >= Tage || agef >= Tage # terminal age conditions
        ret = 0.0
    elseif t >= Ttime # terminal time condition
        ret = 1.5
    else
        ret = log(
        exp(z(t, agem, edum, agef, eduf, ?) + ?*Z(t+1, agem+1, edum, agef+1, eduf, ?))
        + exp(-dc(?) + ?*(V(t+1, agem+1, edum, ?) + U(t+1, agef+1, eduf, ?)))
        )
    end
    return ret
end

@memoize function V(t::Int, agem::Int, edum::Int, ?::Vector{Float64})
    ret::Float64 = 0.0
    agef::Int = 0
    eduf::Int = 0
    if agem >= Tage
        ret = 0.0
    elseif t >= Ttime
        ret = 1.5
    else
        suma::Float64 = 0.0
        V1 = V(t+1, agem+1, edum, ?)
        wagem1 = wagem(t, agem, edum, ?)
        for agef in 16:Tage-1, eduf in 1:K
            suma += exp(s*z(t, agem, edum, agef, eduf, ?) + wagem1 + ?*(s*(Z(t+1, agem+1, edum, agef+1, eduf, ?) - V1 - U(t+1, agef+1, eduf, ?)) + V1))
        end
        ret = log(
                exp(wagem1 + ?*(V1)) + suma
                )
    end
    return ret
end

@memoize  function U(t::Int, agef::Int, eduf::Int, ?::Vector{Float64})
    ret::Float64 = 0.0
    agem::Int = 0
    edum::Int = 0
    if agef >= Tage
        ret = 0.0
    elseif t >= Ttime
        ret = 1.5
    else
        suma::Float64 = 0.0
        U1 = U(t+1, agef+1, eduf, ?)
        wagef1 = wagef(t, agef, eduf, ?)
        for agem in 16:Tage-1, edum in 1:K
            suma += exp((1-s)*z(t, agem, edum, agef, eduf, ?) + wagef1 + ?*((1-s)*(Z(t+1, agem+1, edum, agef+1, eduf, ?) - V(t+1, agem+1, edum, ?) - U1) + U1))
        end
        ret = log(
                  exp(wagef1 + ?*(U1)) + suma
                  )
    end  
    return ret
end
Run Code Online (Sandbox Code Playgroud)

现在为Tage = 60：

julia> @time U(1984, 17, 3, t0) # a test with const Tage = 60
  3.789108 seconds (27.08 M allocations: 494.135 MB, 15.12% gc time)
16.99266161490023
Run Code Online (Sandbox Code Playgroud)

我也测试了一个Int16版本：

julia> @time U(1984%Int16, 17%Int16, 3%Int16, t0)
  3.596871 seconds (27.05 M allocations: 413.808 MB, 11.93% gc time)
16.99266161490023
Run Code Online (Sandbox Code Playgroud)

与Int32版本相比，它的运行速度提高了 5%，内存分配减少了 20%

归档时间：

查看次数：

888 次

最近记录：

10 年 前
递归与迭代 96
更多相关链接
相关归档
如何根据喜好安排每个人？ 6
获取 Julia 中分配的视图的视图 6
关于在具有许多子目录的项目上使用make的反馈 5
在递归函数中重新排序匹配子句 5
当递归变得非常深时,理解中的递归调用有什么特别之处？ 5
如何更新正在开发的 Julia 包中的依赖项 5
在设计编译器时,处理递归是否需要特殊处理？ 2
在Julia中通过"N"由3个数组生成一个HSV图像 2
有效地使用函数作为参数 1
XSL v1.0和XPath遍历深度嵌套并递归求和所有行 0
难疑归档
什么是Python中的元类？ 5409
event.preventDefault()与return false 2891
如何在JavaScript中获取查询字符串值？ 2701
PHP:从数组中删除元素 2362
如何确定数组是否包含Java中的特定值？ 2194
什么是sleep()的JavaScript版本？ 2115
如果目录尚不存在,如何mkdir？ 1784
如何删除GitHub上的提交？ 1619
如何恢复Git中丢失的存储？ 1617
vim"用sudo写"技巧如何工作？ 1347