--------------------------------------------------------------------------------
module AC where
import Circuit
import Array
import List
import Matrix.Complex.LU
import Complex
---- Y node admittance matrix
---- v vector of unknown voltages
---- i equivalent current source vector
--- Build matrix Y of Yv = i
mk_Y :: [([Int], [Double])] -> [Double] -> Array (Int, Int) (Complex Double)
mk_Y (e:es) param = array ((1, 1), (n, n))
([((i, j), z i j) | i <- [1..n], j <- [1..n]])
where z i j | i <= m && j <= m
= yy_1 i j (admList (e:es))
| i > m || j > m
= yy_2 i j (srcList (e:es))
| otherwise = 0
yy_2 _ _ [] = 0
yy_2 i j (e:es) | typeDev (fst e) == tVAC
= vsrc i j (e:es) + yy_2 i j es
| otherwise = yy_2 i j es
yy_1 _ _ [] = 0
yy_1 i j (e:es) | typeDev (fst e) == tRES
= adm i j (e:es) + yy_1 i j es
| typeDev (fst e) == tCAP
= cap i j (e:es) + yy_1 i j es
| typeDev (fst e) == tIND
= ind i j (e:es) + yy_1 i j es
| otherwise = yy_1 i j es
adm m n (e:es) | m == nodez 1 (fst e) && n == nodez 2 (fst e)
= -1 * (((snd e)!!0) :+ 0)
| n == nodez 1 (fst e) && m == nodez 2 (fst e)
= -1 * (((snd e)!!0) :+ 0)
| n == m && (n == nodez 1 (fst e) ||
m == nodez 2 (fst e))
= (((snd e)!!0) :+ 0)
| otherwise = 0
ind m n (e:es) | m == nodez 1 (fst e) && n == nodez 2 (fst e)
= -1.0 / (0 :+ (xL ((snd e)!!0)))
| n == nodez 1 (fst e) && m == nodez 2 (fst e)
= -1.0 / (0 :+ (xL ((snd e)!!0)))
| n == m && (n == nodez 1 (fst e) ||
m == nodez 2 (fst e))
= 1.0 / (0 :+ (xL ((snd e)!!0)))
| otherwise = 0
cap m n (e:es) | m == nodez 1 (fst e) && n == nodez 2 (fst e)
= -1.0 / (0 :+ (xC ((snd e)!!0)))
| n == nodez 1 (fst e) && m == nodez 2 (fst e)
= -1.0 / (0 :+ (xC ((snd e)!!0)))
| n == m && (n == nodez 1 (fst e) ||
m == nodez 2 (fst e))
= 1.0 / (0 :+ (xC ((snd e)!!0)))
| otherwise = 0
vsrc i j (e:es) | length es == (u - i) &&
j == nodez 1 (fst e) = 1
| length es == (u - i) &&
j == nodez 2 (fst e) = 1
| length es == (u - j) &&
i == nodez 1 (fst e) = 1
| length es == (u - j) &&
i == nodez 2 (fst e) = -1
| otherwise = 0
m = ncol (e:es)
n = nrow (e:es)
u = m + length (srcList (e:es))
freq = param!!pFREQ
xC value = -1.0 / (value * 2 * pi * freq)
xL value = value * 2 * pi * freq
--------------------------------------------------------------------------------
--- Build vector i of Yv = i
mk_i :: [([Int], [Double])] -> [Double] -> Array Int (Complex Double)
mk_i (e:es) param = array (1, n) ([(j, z j) | j <- [1..n]])
where z j | j <= m = ii_1 j (srcList (e:es))
| otherwise = ii_2 j (srcList (e:es))
ii_1 _ [] = 0
ii_1 j (e:es) | typeDev (fst e) == tIAC
= isrc j (e:es) + ii_1 j es
| otherwise = ii_1 j es
ii_2 _ [] = 0
ii_2 j (e:es) | typeDev (fst e) == tVAC
= vsrc j (e:es) + ii_2 j es
| otherwise = ii_2 j es
vsrc j (v:vs) | length vs == (u - j) = ((snd v)!!0):+0
| otherwise = 0
isrc j (i:is) | j == nodez 1 (fst i) = (((snd i)!!0) :+ 0)
| j == nodez 2 (fst i) = -1 * (((snd i)!!0) :+ 0)
| otherwise = 0
m = ncol (e:es)
n = nrow (e:es)
u = m + length (srcList (e:es))
--- Utilities
typeDev [] = 0
typeDev parm1 = parm1!!0
nodez _ [] = 0
nodez n parm1 = parm1!!n
nodeList [] = []
nodeList (n:ns) = nub (union (tail $ fst n) (nodeList ns))
ncol netlist = length (nodeList netlist) - 1
nrow netlist = ncol netlist + length (srcList netlist)
- length (isrcList netlist)
srcList [] = []
srcList netlist = vsrcList netlist ++ isrcList netlist
admList [] = []
admList netlist = resList netlist ++ indList netlist ++ capList netlist
resList :: [([Int], [Double])] -> [([Int], [Double])]
resList [] = []
resList (n:ns) | typeDev (fst n) == tRES = n : resList ns
| otherwise = resList ns
capList :: [([Int], [Double])] -> [([Int], [Double])]
capList [] = []
capList (n:ns) | typeDev (fst n) == tCAP = n : capList ns
| otherwise = capList ns
indList :: [([Int], [Double])] -> [([Int], [Double])]
indList [] = []
indList (n:ns) | typeDev (fst n) == tIND = n : indList ns
| otherwise = indList ns
vsrcList :: [([Int], [Double])] -> [([Int], [Double])]
vsrcList [] = []
vsrcList (n:ns) | typeDev (fst n) == tVAC = n : vsrcList ns
| otherwise = vsrcList ns
isrcList :: [([Int], [Double])] -> [([Int], [Double])]
isrcList [] = []
isrcList (n:ns) | typeDev (fst n) == tISRC = n : isrcList ns
| otherwise = isrcList ns
-------------------------------------------------------------------------------
--- Solve v of Yv = i
pFREQ :: Int
pFREQ = 0
pPHASE :: Int
pPHASE = 1
acCalc :: [([Int], [Double])] -> [Double] -> Array Int (Complex Double)
acCalc netlist param = psolve (mk_Y netlist param) (mk_i netlist param)
|