module Diagrams.CubicSpline.Internal
(
solveCubicSplineDerivatives
, solveCubicSplineDerivativesClosed
, solveCubicSplineCoefficients
) where
import Diagrams.Solve.Tridiagonal
import Data.List
solveCubicSplineDerivatives :: Fractional a => [a] -> [a]
solveCubicSplineDerivatives (x:xs) = solveTriDiagonal as bs as ds
where
as = replicate (l - 1) 1
bs = 2 : replicate (l - 2) 4 ++ [2]
l = length ds
ds = zipWith f (xs ++ [last xs]) (x:x:xs)
f a b = 3*(a - b)
solveCubicSplineDerivatives _ = error "argument to solveCubicSplineDerivatives must be nonempty"
solveCubicSplineDerivativesClosed :: Fractional a => [a] -> [a]
solveCubicSplineDerivativesClosed xs = solveCyclicTriDiagonal as bs as ds 1 1
where
as = replicate (l - 1) 1
bs = replicate l 4
l = length xs
xs' = cycle xs
ds = take l $ zipWith f (drop 1 xs') (drop (l - 1) xs')
f a b = 3*(a - b)
solveCubicSplineCoefficients :: Fractional a => Bool -> [a] -> [[a]]
solveCubicSplineCoefficients closed xs =
[ [x,d,3*(x1-x)-2*d-d1,2*(x-x1)+d+d1]
| (x,x1,d,d1) <- zip4 xs' (tail xs') ds' (tail ds')
]
where
ds | closed = solveCubicSplineDerivativesClosed xs
| otherwise = solveCubicSplineDerivatives xs
close as | closed = as ++ [head as]
| otherwise = as
xs' = close xs
ds' = close ds