| otherwise = zipWith (-) a $ vsum . zipWith cmul (map (vdot a) d') $ d'

cmul c       = map (c*)

module Rosenbrock where
import Data.List
 
rosenbrock :: ([Double] -> Double) -> Double -> [Double] -> [Double]
rosenbrock f e x = xiter (map (+(2*e)) x) x d where
    d = take (length x) . iterate ((0:).init) $ 1:replicate (length x - 1) 0
    xiter x_ x d | (<e) . vnorm $ zipWith (-) x_ x = x
                 | otherwise = xiter x x' $ gramSchmidt e d l where
        (x',l) = yiter d x [] where
            yiter []      y l = (y, reverse l)
            yiter (di:ds) y l = yiter ds (pl li) (li:l) where
                pl x = zipWith (+) y $ cmul x di
                li   = dichotomy (f.pl) e (-100000) 100000
 
dichotomy f e a b = go a b where
    delta  = e / 10
    go a b | abs (b - a) < e = (a + b)/2
           | otherwise       = go a' b' where
               middle  = (a + b) / 2
               lambda  = middle - delta
               mu      = middle + delta
               (a',b') = if f lambda < f mu then (a,mu) else (lambda,b)
 
gramSchmidt e d l = go a [] where
    go []     d' = d'
    go (a:as) d' = go as $ d' ++ [normalize di] where
        di | null d'   = a
           | otherwise = zipWith (-) a $ vsum . zipWith cmul (map (vdot a) d') $ d'
 
    a = zipWith3 (\li di ai -> if abs li < e then di else ai) l d a'
    a' = map (vsum . map (uncurry cmul)) . tails . zip l $ d
 
---- стандартные векторные операции
vnorm        = sqrt . sum . map (^2)
normalize v  = map (/vnorm v) v
vsum         = foldr1 (zipWith (+))
vdot v1 v2   = sum $ zipWith (*) v1 v2
cmul c       = map (c*)
 
main = do
    print $ rosenbrock (\[x,y] -> (x-2)^4 + (x-2*y)^2) (1.5e-2) [0,3]
    print $ rosenbrock (\[x,y,z] -> (x-2)^4 + (x-2*y)^2 + (z-x^3-y-77)^2) (1.0e-4) [0,3,0]