{-|
Module      : TBit.Topological.Spillage
Description : calculate spin-orbit and other spillages
Copyright   : (c) Matthew Daniels, 2014
License     : New BSD
Maintainer  : danielsmw@gmail.com
Stability   : experimental

Provide routines for determining the spillage, and in particular
Vanderbilt and Liu's spin-orbit spillage, given a Hamiltonian whose
spin-orbit coupling term can be turned on and off.
-}
module TBit.Topological.Spillage (spillage) where

import TBit.Types
import TBit.Parameterization (crunch, getScalar)
import TBit.Hamiltonian.Eigenstates (eigenbras,eigenkets)
import Debug.Trace

import Numeric.LinearAlgebra.HMatrix ((<>), fromBlocks, size, (!), Matrix)

import Control.Monad (liftM)
import Control.Monad.State

import Data.Map (adjust)
import Data.Complex

                        
spillage :: String -> Filling -> Parameterized Hamiltonian -> Wavevector -> Parameterized Double
spillage ps n m k = do h' <- (modify (setScalar ps 0) >> m)
                       h  <- m
                       s  <- spillMatrix h h' k
                       return $ (fromIntegral n) 
                              - sum [ (realPart $ abs $ s!j!k)^2 
                                      | j <- [0 .. pred n]
                                      , k <- [0 .. pred n]]
    where setScalar str val ps = ps { scalarParams = adjust (\_ -> val :+ 0.0) str 
                                                   $ scalarParams ps }
                                                   
spillMatrix :: Hamiltonian 
            -> Hamiltonian 
            -> Wavevector
            -> Parameterized (Matrix (Complex Double))
spillMatrix h h' k = do es  <- eigenkets h  k
                        es' <- eigenbras h' k
                        let dim = fst . size . h $ k
                        return $ fromBlocks
                               $ [[ (es'!!m) <> (es!!n) | m <- [0 .. pred dim]]
                                                        | n <- [0 .. pred dim]]