@@ -495,17 +495,14 @@ end
495495 # problem with sparse matrices
496496 p_prototype = spdiagm(- 1 => ones(eltype(u0), N - 1 ),
497497 N - 1 => ones(eltype(u0), 1 ))
498- d_prototype = zero(u0)
499498 linear_advection_fd_upwind_PDS_sparse = PDSProblem(linear_advection_fd_upwind_P!,
500499 linear_advection_fd_upwind_D!,
501500 u0, tspan;
502- p_prototype = p_prototype,
503- d_prototype = d_prototype)
501+ p_prototype = p_prototype)
504502 linear_advection_fd_upwind_PDS_sparse_2 = PDSProblem{true }(linear_advection_fd_upwind_P!,
505503 linear_advection_fd_upwind_D!,
506504 u0, tspan;
507- p_prototype = p_prototype,
508- d_prototype = d_prototype)
505+ p_prototype = p_prototype)
509506 linear_advection_fd_upwind_ConsPDS_sparse = ConservativePDSProblem(linear_advection_fd_upwind_P!,
510507 u0, tspan;
511508 p_prototype = p_prototype)
@@ -1201,6 +1198,75 @@ end
12011198 end
12021199 end
12031200
1201+ # Here we check that the type of p_prototype actually
1202+ # defines the types of the Ps inside the algorithm caches.
1203+ # We test sparse, tridiagonal, and dense matrices.
1204+ @testset " Prototype type check" begin
1205+ # prod and dest functions
1206+ prod_inner! = (P, u, p, t) -> begin
1207+ fill!(P, zero(eltype(P)))
1208+ for i in 1 : (length(u) - 1 )
1209+ P[i, i + 1 ] = i * u[i]
1210+ end
1211+ return nothing
1212+ end
1213+ prod_sparse! = (P, u, p, t) -> begin
1214+ @test P isa SparseMatrixCSC
1215+ prod_inner!(P, u, p, t)
1216+ return nothing
1217+ end
1218+ prod_tridiagonal! = (P, u, p, t) -> begin
1219+ @test P isa Tridiagonal
1220+ prod_inner!(P, u, p, t)
1221+ return nothing
1222+ end
1223+ prod_dense! = (P, u, p, t) -> begin
1224+ @test P isa Matrix
1225+ prod_inner!(P, u, p, t)
1226+ return nothing
1227+ end
1228+ dest! = (D, u, p, t) -> begin
1229+ fill!(D, zero(eltype(D)))
1230+ end
1231+ # prototypes
1232+ P_tridiagonal = Tridiagonal([0.1 , 0.2 , 0.3 ],
1233+ [0.0 , 0.0 , 0.0 , 0.0 ],
1234+ [0.4 , 0.5 , 0.6 ])
1235+ P_dense = Matrix(P_tridiagonal)
1236+ P_sparse = sparse(P_tridiagonal)
1237+ # problem definition
1238+ u0 = [1.0 , 1.5 , 2.0 , 2.5 ]
1239+ tspan = (0.0 , 1.0 )
1240+ dt = 0.5
1241+ # # conservative PDS
1242+ prob_default = ConservativePDSProblem(prod_dense!, u0, tspan)
1243+ prob_tridiagonal = ConservativePDSProblem(prod_tridiagonal!, u0, tspan;
1244+ p_prototype = P_tridiagonal)
1245+ prob_dense = ConservativePDSProblem(prod_dense!, u0, tspan;
1246+ p_prototype = P_dense)
1247+ prob_sparse = ConservativePDSProblem(prod_sparse!, u0, tspan;
1248+ p_prototype = P_sparse)
1249+ # # nonconservative PDS
1250+ prob_default2 = PDSProblem(prod_dense!, dest!, u0, tspan)
1251+ prob_tridiagonal2 = PDSProblem(prod_tridiagonal!, dest!, u0, tspan;
1252+ p_prototype = P_tridiagonal)
1253+ prob_dense2 = PDSProblem(prod_dense!, dest!, u0, tspan;
1254+ p_prototype = P_dense)
1255+ prob_sparse2 = PDSProblem(prod_sparse!, dest!, u0, tspan;
1256+ p_prototype = P_sparse)
1257+ # solve and test
1258+ for alg in (MPE(), MPRK22(0.5 ), MPRK22(1.0 ), MPRK43I(1.0 , 0.5 ),
1259+ MPRK43I(0.5 , 0.75 ),
1260+ MPRK43II(2.0 / 3.0 ), MPRK43II(0.5 ), SSPMPRK22(0.5 , 1.0 ),
1261+ SSPMPRK43())
1262+ for prob in (prob_default, prob_tridiagonal, prob_dense, prob_sparse,
1263+ prob_default2,
1264+ prob_tridiagonal2, prob_dense2, prob_sparse2)
1265+ solve(prob, alg; dt, adaptive = false )
1266+ end
1267+ end
1268+ end
1269+
12041270 # Here we check the convergence order of pth-order schemes for which
12051271 # also an interpolation of order p is available
12061272 @testset " Convergence tests (conservative)" begin
0 commit comments