_tests_simulate.jl 6.71 KB
Newer Older
1
2
3

# File for benchmarking simulation and memory access of the package.

4
5
# Trajectories

6
_get_values_col(σ::AbstractTrajectory, var::String) = 
7
@view σ.values[(σ.m)._map_obs_var_idx[var],:] 
8
_get_values_row(σ::AbstractTrajectory, var::String) = 
9
@view σ.values[:,(σ.m)._map_obs_var_idx[var]] 
10

11
_get_state_col(σ::AbstractTrajectory, idx::Int) = 
12
@view σ.values[:,idx]
13
_get_state_row(σ::AbstractTrajectory, idx::Int) = 
14
@view σ.values[idx,:]
15
16
17
18
19
20
21

_get_value_col(σ::AbstractTrajectory, var::String, idx::Int) = 
σ.values[(σ.m)._map_obs_var_idx[var],idx] 
_get_value_row(σ::AbstractTrajectory, var::String, idx::Int) = 
σ.values[idx,(σ.m)._map_obs_var_idx[var]] 

# Model
22
23
24
25
26
27
28
29

function _simulate_col(m::ContinuousTimeModel)
    # trajectory fields
    full_values = zeros(m.d, 0)
    times = zeros(0)
    transitions = Vector{Union{String,Nothing}}(undef,0)
    # values at time n
    n = 0
30
    xn = @view m.x0[:]
31
32
33
34
35
    tn = m.t0 
    tr = [""]
    # at time n+1
    xnplus1 = zeros(Int, m.d)
    tnplus1 = zeros(Float64, 1)
36
37
    isabsorbing = (m.isabsorbing(m.p,xn))::Bool
    while !isabsorbing && (tn <= m.time_bound)
38
39
40
41
42
43
        m.f!(xnplus1, tnplus1, tr, xn, tn, m.p)
        full_values = hcat(full_values, xnplus1)
        push!(times, tnplus1[1])
        push!(transitions, tr[1])
        xn, tn = xnplus1, tnplus1[1]
        n += 1
44
        isabsorbing = m.isabsorbing(m.p,xn)::Bool
45
46
    end
    values = @view full_values[m._g_idx,:] 
47
    if isbounded(m)
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
        if times[end] > m.time_bound
            values[:, end] = values[:,end-1]
            times[end] = m.time_bound
            transitions[end] = nothing
        end
    end
    return Trajectory(m, values, times, transitions)
end

function _simulate_row(m::ContinuousTimeModel)
    # trajectory fields
    full_values = zeros(m.d, 0)
    times = zeros(0)
    transitions = Vector{Union{String,Nothing}}(undef,0)
    # values at time n
    n = 0
64
    xn = @view m.x0[:]
65
66
67
68
69
    tn = m.t0 
    tr = [""]
    # at time n+1
    xnplus1 = zeros(Int, m.d)
    tnplus1 = zeros(Float64, 1)
70
71
    isabsorbing = (m.isabsorbing(m.p,xn))::Bool
    while !isabsorbing && (tn <= m.time_bound)
72
73
74
75
76
77
        m.f!(xnplus1, tnplus1, tr, xn, tn, m.p)
        full_values = vcat(full_values, xnplus1)
        push!(times, tnplus1[1])
        push!(transitions, tr[1])
        xn, tn = xnplus1, tnplus1[1]
        n += 1
78
        isabsorbing = m.isabsorbing(m.p,xn)::Bool
79
80
    end
    values = @view full_values[m._g_idx,:] 
81
    if isbounded(m)
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
        if times[end] > m.time_bound
            values[:, end] = values[:,end-1]
            times[end] = m.time_bound
            transitions[end] = nothing
        end
    end
    return Trajectory(m, values, times, transitions)
end


function _simulate_col_buffer(m::ContinuousTimeModel; buffer_size::Int = 5)
    # trajectory fields
    full_values = zeros(m.d, 0)
    times = zeros(0)
    transitions = Vector{Union{String,Nothing}}(undef,0)
    # values at time n
    n = 0
99
    xn = @view m.x0[:]
100
101
102
103
104
    tn = m.t0 
    # at time n+1
    mat_x = zeros(Int, m.d, buffer_size)
    l_t = zeros(Float64, buffer_size)
    l_tr = Vector{String}(undef, buffer_size)
105
106
    isabsorbing = m.isabsorbing(m.p,xn)::Bool
    while !isabsorbing && (tn <= m.time_bound)
107
        i = 0
108
        while i < buffer_size && !isabsorbing && (tn <= m.time_bound)
109
110
111
112
            i += 1
            m.f!(mat_x, l_t, l_tr, i, xn, tn, m.p)
            xn = @view mat_x[:,i]
            tn = l_t[i]
113
            isabsorbing = m.isabsorbing(m.p,xn)::Bool
114
115
116
117
118
        end
        full_values = hcat(full_values, @view mat_x[:,1:i])
        append!(times, @view l_t[1:i])
        append!(transitions,  @view l_tr[1:i])
        n += i
119
        isabsorbing = m.isabsorbing(m.p,xn)::Bool
120
121
    end
    values = @view full_values[m._g_idx,:] 
122
    if isbounded(m)
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
        if times[end] > m.time_bound
            values[:, end] = values[:,end-1]
            times[end] = m.time_bound
            transitions[end] = nothing
        end
    end
    return Trajectory(m, values, times, transitions)
end

function _simulate_row_buffer(m::ContinuousTimeModel; buffer_size::Int = 5)
    # trajectory fields
    full_values = zeros(0, m.d)
    times = zeros(0)
    transitions = Vector{Union{String,Nothing}}(undef,0)
    # values at time n
    n = 0
139
    xn = @view m.x0[:]
140
141
142
143
144
    tn = m.t0 
    # at time n+1
    mat_x = zeros(Int, buffer_size, m.d)
    l_t = zeros(Float64, buffer_size)
    l_tr = Vector{String}(undef, buffer_size)
145
146
    isabsorbing = m.isabsorbing(m.p,xn)::Bool
    while !isabsorbing && (tn <= m.time_bound)
147
        i = 0
148
        while i < buffer_size && !isabsorbing && (tn <= m.time_bound)
149
150
            i += 1
            m.f!(mat_x, l_t, l_tr, i, xn, tn, m.p)
151
            xn = @view mat_x[i,:]
152
            tn = l_t[i]
153
            isabsorbing = m.isabsorbing(m.p,xn)::Bool
154
155
156
157
158
        end
        full_values = vcat(full_values, @view mat_x[1:i,:])
        append!(times, @view l_t[1:i])
        append!(transitions,  @view l_tr[1:i])
        n += i
159
        isabsorbing = m.isabsorbing(m.p,xn)::Bool
160
161
    end
    values = @view full_values[:,m._g_idx] 
162
    if isbounded(m)
163
164
165
166
167
168
169
170
171
        if times[end] > m.time_bound
            values[end,:] = values[end-1,:]
            times[end] = m.time_bound
            transitions[end] = nothing
        end
    end
    return Trajectory(m, values, times, transitions)
end

172
173
174
175
176
177
178
179
function _simulate_without_view(m::ContinuousTimeModel)
    # trajectory fields
    full_values = Matrix{Int}(undef, 1, m.d)
    full_values[1,:] = m.x0
    times = Float64[m.t0]
    transitions = Union{String,Nothing}[nothing]
    # values at time n
    n = 0
180
    xn = @view m.x0[:]
181
182
183
184
185
    tn = m.t0 
    # at time n+1
    mat_x = zeros(Int, m.buffer_size, m.d)
    l_t = zeros(Float64, m.buffer_size)
    l_tr = Vector{String}(undef, m.buffer_size)
186
187
    isabsorbing = m.isabsorbing(m.p,xn)::Bool
    while !isabsorbing && (tn <= m.time_bound)
188
        i = 0
189
        while i < m.buffer_size && !isabsorbing && (tn <= m.time_bound)
190
            i += 1
191
            m.f!(mat_x, l_t, l_tr, i, @view(xn[:]), tn, m.p)
192
193
            xn = mat_x[i,:]
            tn = l_t[i]
194
            isabsorbing = m.isabsorbing(m.p,@view(xn[:]))::Bool
195
196
197
198
199
        end
        full_values = vcat(full_values, mat_x[1:i,:])
        append!(times, l_t[1:i])
        append!(transitions,  l_tr[1:i])
        n += i
200
        isabsorbing = m.isabsorbing(m.p,@view(xn[:]))::Bool
201
    end
202
    if isbounded(m)
203
204
205
206
207
208
209
210
211
212
213
214
215
216
        if times[end] > m.time_bound
            full_values[end,:] = full_values[end-1,:]
            times[end] = m.time_bound
            transitions[end] = nothing
        else
            full_values = vcat(full_values, reshape(full_values[end,:], 1, m.d))
            push!(times, m.time_bound)
            push!(transitions, nothing)
        end
    end
    values = full_values[:,m._g_idx] 
    return Trajectory(m, values, times, transitions)
end