When solving system of linear equations by Tridiagonal matrix algorithm in Haskell I met following problem.
We have three vectors: a
, b
and c
, and we want to make a third vector c'
which is a combination of them:
c'[i] = c[i] / b[i], i = 0 c'[i] = c[i] / (b[i] - a[i] * c'[i-1]), 0 < i < n - 1 c'[i] = undefined, i = n - 1
Naive implementation of the formula above in Haskell is a开发者_运维问答s follows:
calcC' a b c = Data.Vector.generate n f where n = Data.Vector.length a f i = | i == 0 = c!0 / b!0 | i == n - 1 = 0 | otherwise = c!i / (b!i - a!i * f (i - 1))
It looks like this function calcC'
has complexity O(n2) due to recurrence. But all we actualy need is to pass to inner function f
one more parameter with previously generated value.
I wrote my own version of generate
with complexity O(n) and helper function mapP
:
mapP f xs = mapP' xs Nothing where mapP' [] _ = [] mapP' (x:xs) xp = xn : mapP' xs (Just xn) where xn = f x xp generateP n f = Data.Vector.fromList $ mapP f [0 .. n-1]
As one can see, mapP
acts like a standard map
, but also passes to mapping function previously generated value or Nothing
for first call.
My question: is there any pretty standard ways to do this in Haskell? Don't I reinvent the weel?
Thanks.
There are two standard function called mapAccumL
and mapAccumR
that do precisely what you want.
mapAccumL :: (acc -> x -> (acc, y)) -> acc -> [x] -> (acc, [y])
mapAccumR :: (acc -> x -> (acc, y)) -> acc -> [x] -> (acc, [y])
Basically, they behave like a combination of fold
and map
.
map f = snd . mapAccumL (\_ x -> (() , f x) ()
foldl f b = fst . mapAccumL (\b x -> (f b x, () ) b
If you use Data.Array, which is lazy, you can express the recurrence directly by referring to c' while defining c'.
Following code seems to be the simplest implementation of formula above in my case:
import qualified Data.Vector.Generic as V calcC' a b c = V.postscanl' f 0.0 $ V.zip3 a b c where f c' (a, b, c) = c / (b - a * c')
Thanks to the authors of Vector
who added helpfull postscanl'
method.
精彩评论