@@ -427,34 +427,32 @@ end
427
427
let
428
428
@variables X (t) Y (t)
429
429
@parameters k1 k2
430
- vrj1 = VariableRateJump (k1 * X, [X ~ X - 1 ])
431
- vrj2 = VariableRateJump (k2 , [Y ~ Y + 1 ])
430
+ vrj1 = VariableRateJump (k1 * X, [X ~ X - 1 ]; save_positions = ( false , false ) )
431
+ vrj2 = VariableRateJump (k1 , [Y ~ Y + 1 ]; save_positions = ( false , false ) )
432
432
eqs = [D (X) ~ k2, D (Y) ~ - k2/ 100 * Y]
433
433
@named jsys = JumpSystem ([vrj1, vrj2, eqs[1 ], eqs[2 ]], t, [X, Y], [k1, k2])
434
434
jsys = complete (jsys)
435
435
X0 = 0.0 ; Y0 = 0.0
436
436
u0 = [X => X0, Y => Y0]
437
437
k1val = 1.0 ; k2val = 20.0
438
438
p = [k1 => k1val, k2 => k2val]
439
- tspan = (0.0 , 20 .0 )
439
+ tspan = (0.0 , 10 .0 )
440
440
oprob = ODEProblem (jsys, u0, tspan, p)
441
- jprob = JumpProblem (jsys, oprob; rng)
441
+ jprob = JumpProblem (jsys, oprob; rng, save_positions = ( false , false ) )
442
442
443
- times = range (0.0 , 20.0 , length = 100 )
443
+ times = range (0.0 , tspan[ 2 ] , length = 100 )
444
444
Nsims = 2000
445
- X = zeros (length (times))
446
- Y = similar (X )
445
+ Xv = zeros (length (times))
446
+ Yv = copy (Xv )
447
447
for n in 1 : Nsims
448
- sol = solve (jprob, Tsit5 ())
449
- X .+ = Array ( sol (times; idxs = X) )
450
- Y .+ = Array ( sol (times; idxs = Y) )
448
+ sol = solve (jprob, Tsit5 (); saveat = times )
449
+ Xv .+ = sol[ 1 ,:] # sol (times; idxs = X)
450
+ Yv .+ = sol[ 2 ,:] # sol (times; idxs = Y)
451
451
end
452
- X ./= Nsims; Y ./= Nsims;
452
+ Xv ./= Nsims; Yv ./= Nsims;
453
453
454
454
Xact (t) = X0 * exp (- k1val * t) + (k2val / k1val) * (1 - exp (- k1val * t))
455
- Yact (t) = Y0 * exp (- k2val/ 100 * t) + (k1 / (k2val/ 100 )) * (1 - exp (- k2val/ 100 * t))
456
-
457
- @test all (abs .(X .- Xact .(times)) .<= 0.05 .* X)
458
- @test all (abs .(Y .- Yact .(times)) .<= 0.05 .* Y)
459
-
455
+ Yact (t) = Y0 * exp (- k2val/ 100 * t) + (k1val / (k2val/ 100 )) * (1 - exp (- k2val/ 100 * t))
456
+ @test all (abs .(Xv .- Xact .(times)) .<= 0.05 .* Xv)
457
+ @test all (abs .(Yv .- Yact .(times)) .<= 0.05 .* Yv)
460
458
end
0 commit comments