Skip to content

Commit fbd58de

Browse files
ChrisRackauckasChrisRackauckas
committed
hermite nirvanna
1 parent 87e13a5 commit fbd58de

File tree

1 file changed

+177
-1
lines changed

1 file changed

+177
-1
lines changed

‎src/dense/generic_dense.jl

Lines changed: 177 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -403,6 +403,69 @@ Herimte Interpolation, chosen if no other dispatch for ode_interpolant
403403
out
404404
end
405405

406+
"""
407+
Herimte Interpolation, chosen if no other dispatch for ode_interpolant
408+
"""
409+
@inline function hermite_interpolant(Θ,dt,y₀,y₁,k,cache,idxs,T::Type{Val{1}}) # Default interpolant is Hermite
410+
if typeof(y₀) <: AbstractArray
411+
if typeof(idxs) <: Tuple
412+
out = similar(y₀,idxs)
413+
iter_idxs = enumerate(idxs)
414+
else
415+
out = similar(y₀,indices(idxs))
416+
iter_idxs = enumerate(idxs)
417+
end
418+
@inbounds for (j,i) in iter_idxs
419+
out[j] = k[1][i] + Θ*(-4*dt*k[1][i] - 2*dt*k[2][i] - 6*y₀[i] + Θ*(3*dt*k[1][i] + 3*dt*k[2][i] + 6*y₀[i] - 6*y₁[i]) + 6*y₁[i])/dt
420+
end
421+
else
422+
out = k[1] + Θ*(-4*dt*k[1] - 2*dt*k[2] - 6*y₀ + Θ*(3*dt*k[1] + 3*dt*k[2] + 6*y₀ - 6*y₁) + 6*y₁)/dt
423+
end
424+
out
425+
end
426+
427+
"""
428+
Herimte Interpolation, chosen if no other dispatch for ode_interpolant
429+
"""
430+
@inline function hermite_interpolant(Θ,dt,y₀,y₁,k,cache,idxs,T::Type{Val{2}}) # Default interpolant is Hermite
431+
if typeof(y₀) <: AbstractArray
432+
if typeof(idxs) <: Tuple
433+
out = similar(y₀,idxs)
434+
iter_idxs = enumerate(idxs)
435+
else
436+
out = similar(y₀,indices(idxs))
437+
iter_idxs = enumerate(idxs)
438+
end
439+
@inbounds for (j,i) in iter_idxs
440+
out[j] = (-4*dt*k[1][i] - 2*dt*k[2][i] - 6*y₀[i] + Θ*(6*dt*k[1][i] + 6*dt*k[2][i] + 12*y₀[i] - 12*y₁[i]) + 6*y₁[i])/(dt*dt)
441+
end
442+
else
443+
out = (-4*dt*k[1] - 2*dt*k[2] - 6*y₀ + Θ*(6*dt*k[1] + 6*dt*k[2] + 12*y₀ - 12*y₁) + 6*y₁)/(dt*dt)
444+
end
445+
out
446+
end
447+
448+
"""
449+
Herimte Interpolation, chosen if no other dispatch for ode_interpolant
450+
"""
451+
@inline function hermite_interpolant(Θ,dt,y₀,y₁,k,cache,idxs,T::Type{Val{3}}) # Default interpolant is Hermite
452+
if typeof(y₀) <: AbstractArray
453+
if typeof(idxs) <: Tuple
454+
out = similar(y₀,idxs)
455+
iter_idxs = enumerate(idxs)
456+
else
457+
out = similar(y₀,indices(idxs))
458+
iter_idxs = enumerate(idxs)
459+
end
460+
@inbounds for (j,i) in iter_idxs
461+
out[j] = (6*dt*k[1][i] + 6*dt*k[2][i] + 12*y₀[i] - 12*y₁[i])/(dt*dt*dt)
462+
end
463+
else
464+
out = (6*dt*k[1] + 6*dt*k[2] + 12*y₀ - 12*y₁)/(dt*dt*dt)
465+
end
466+
out
467+
end
468+
406469
"""
407470
Hairer Norsett Wanner Solving Ordinary Differential Euations I - Nonstiff Problems Page 190
408471
@@ -418,6 +481,45 @@ Herimte Interpolation, chosen if no other dispatch for ode_interpolant
418481
end
419482
end
420483

484+
"""
485+
Herimte Interpolation, chosen if no other dispatch for ode_interpolant
486+
"""
487+
@inline function hermite_interpolant!(out,Θ,dt,y₀,y₁,k,cache,idxs,T::Type{Val{1}}) # Default interpolant is Hermite
488+
if out == nothing
489+
return k[1][idxs] + Θ*(-4*dt*k[1][idxs] - 2*dt*k[2][idxs] - 6*y₀[idxs] + Θ*(3*dt*k[1][idxs] + 3*dt*k[2][idxs] + 6*y₀[idxs] - 6*y₁[idxs]) + 6*y₁[idxs])/dt
490+
else
491+
@inbounds for (j,i) in enumerate(idxs)
492+
out[j] = k[1][i] + Θ*(-4*dt*k[1][i] - 2*dt*k[2][i] - 6*y₀[i] + Θ*(3*dt*k[1][i] + 3*dt*k[2][i] + 6*y₀[i] - 6*y₁[i]) + 6*y₁[i])/dt
493+
end
494+
end
495+
end
496+
497+
"""
498+
Herimte Interpolation, chosen if no other dispatch for ode_interpolant
499+
"""
500+
@inline function hermite_interpolant!(out,Θ,dt,y₀,y₁,k,cache,idxs,T::Type{Val{2}}) # Default interpolant is Hermite
501+
if out == nothing
502+
return (-4*dt*k[1][idxs] - 2*dt*k[2][idxs] - 6*y₀[idxs] + Θ*(6*dt*k[1][idxs] + 6*dt*k[2][idxs] + 12*y₀[idxs] - 12*y₁[idxs]) + 6*y₁[idxs])/(dt*dt)
503+
else
504+
@inbounds for (j,i) in enumerate(idxs)
505+
out[j] = (-4*dt*k[1][i] - 2*dt*k[2][i] - 6*y₀[i] + Θ*(6*dt*k[1][i] + 6*dt*k[2][i] + 12*y₀[i] - 12*y₁[i]) + 6*y₁[i])/(dt*dt)
506+
end
507+
end
508+
end
509+
510+
"""
511+
Herimte Interpolation, chosen if no other dispatch for ode_interpolant
512+
"""
513+
@inline function hermite_interpolant!(out,Θ,dt,y₀,y₁,k,cache,idxs,T::Type{Val{3}}) # Default interpolant is Hermite
514+
if out == nothing
515+
return (6*dt*k[1][idxs] + 6*dt*k[2][idxs] + 12*y₀[idxs] - 12*y₁[idxs])/(dt*dt*dt)
516+
else
517+
@inbounds for (j,i) in enumerate(idxs)
518+
out[j] = (6*dt*k[1][i] + 6*dt*k[2][i] + 12*y₀[i] - 12*y₁[i])/(dt*dt*dt)
519+
end
520+
end
521+
end
522+
421523
"""
422524
Hairer Norsett Wanner Solving Ordinary Differential Euations I - Nonstiff Problems Page 190
423525
@@ -431,10 +533,42 @@ Herimte Interpolation, chosen if no other dispatch for ode_interpolant
431533
end
432534
end
433535

434-
######################## Linear Interpolants
536+
"""
537+
Herimte Interpolation, chosen if no other dispatch for ode_interpolant
538+
"""
539+
@inline function hermite_interpolant!(all_out::ArrayPartition,Θ,dt,all_y₀,all_y₁,all_k,cache,all_idxs,T::Type{Val{1}}) # Default interpolant is Hermite
540+
for (out,y₀,y₁,idxs,k1,k2) in zip(all_out.x,all_y₀.x,all_y₁.x,all_idxs,all_k[1].x,all_k[2].x)
541+
@inbounds for (j,i) in enumerate(idxs...)
542+
out[j] = k[1][i] + Θ*(-4*dt*k[1][i] - 2*dt*k[2][i] - 6*y₀[i] + Θ*(3*dt*k[1][i] + 3*dt*k[2][i] + 6*y₀[i] - 6*y₁[i]) + 6*y₁[i])/dt
543+
end
544+
end
545+
end
435546

547+
"""
548+
Herimte Interpolation, chosen if no other dispatch for ode_interpolant
549+
"""
550+
@inline function hermite_interpolant!(all_out::ArrayPartition,Θ,dt,all_y₀,all_y₁,all_k,cache,all_idxs,T::Type{Val{2}}) # Default interpolant is Hermite
551+
for (out,y₀,y₁,idxs,k1,k2) in zip(all_out.x,all_y₀.x,all_y₁.x,all_idxs,all_k[1].x,all_k[2].x)
552+
@inbounds for (j,i) in enumerate(idxs...)
553+
out[j] = (-4*dt*k[1][i] - 2*dt*k[2][i] - 6*y₀[i] + Θ*(6*dt*k[1][i] + 6*dt*k[2][i] + 12*y₀[i] - 12*y₁[i]) + 6*y₁[i])/(dt*dt)
554+
end
555+
end
556+
end
557+
558+
"""
559+
Herimte Interpolation, chosen if no other dispatch for ode_interpolant
560+
"""
561+
@inline function hermite_interpolant!(all_out::ArrayPartition,Θ,dt,all_y₀,all_y₁,all_k,cache,all_idxs,T::Type{Val{3}}) # Default interpolant is Hermite
562+
for (out,y₀,y₁,idxs,k1,k2) in zip(all_out.x,all_y₀.x,all_y₁.x,all_idxs,all_k[1].x,all_k[2].x)
563+
@inbounds for (j,i) in enumerate(idxs...)
564+
out[j] = (6*dt*k[1][i] + 6*dt*k[2][i] + 12*y₀[i] - 12*y₁[i])/(dt*dt*dt)
565+
end
566+
end
567+
end
436568

437569

570+
######################## Linear Interpolants
571+
438572
@inline function linear_interpolant(Θ,dt,y₀,y₁,idxs,T::Type{Val{0}})
439573
if typeof(y₀) <: AbstractArray
440574
if typeof(idxs) <: Tuple
@@ -454,6 +588,24 @@ end
454588
out
455589
end
456590

591+
@inline function linear_interpolant(Θ,dt,y₀,y₁,idxs,T::Type{Val{1}})
592+
if typeof(y₀) <: AbstractArray
593+
if typeof(idxs) <: Tuple
594+
out = similar(y₀,idxs)
595+
iter_idxs = enumerate(eachindex(y₀))
596+
else
597+
out = similar(y₀,indices(idxs))
598+
iter_idxs = enumerate(idxs)
599+
end
600+
@inbounds for (j,i) in iter_idxs
601+
out[j] = (y₁[i] - y₀[i])/dt
602+
end
603+
else
604+
out = (1-Θ)*y₀ + Θ*y₁
605+
end
606+
out
607+
end
608+
457609
"""
458610
Linear Interpolation
459611
"""
@@ -468,6 +620,19 @@ Linear Interpolation
468620
end
469621
end
470622

623+
"""
624+
Linear Interpolation
625+
"""
626+
@inline function linear_interpolant!(out,Θ,dt,y₀,y₁,idxs,T::Type{Val{1}})
627+
if out == nothing
628+
return (y₁[idxs] - y₀[idxs])/dt
629+
else
630+
@inbounds for (j,i) in enumerate(idxs)
631+
out[j] = (y₁[i] - y₀[i])/dt
632+
end
633+
end
634+
end
635+
471636
"""
472637
Linear Interpolation
473638
"""
@@ -479,3 +644,14 @@ Linear Interpolation
479644
end
480645
end
481646
end
647+
648+
"""
649+
Linear Interpolation
650+
"""
651+
@inline function linear_interpolant!(all_out::ArrayPartition,Θ,dt,all_y₀,all_y₁,cache,all_idxs,T::Type{Val{1}})
652+
for (out,y₀,y₁,idxs) in zip(all_out.x,all_y₀.x,all_y₁.x,all_idxs)
653+
@inbounds for (j,i) in enumerate(idxs...)
654+
out[j] = (y₁[i] - y₀[i])/dt
655+
end
656+
end
657+
end

0 commit comments

Comments
 (0)