module InitialGuessModule
export computeDensityGuessSAD, InitialGuess, ZeroGuess, SADGuess
using ..AtomModule

abstract type InitialGuess end
type ZeroGuess_ <: InitialGuess end
ZeroGuess = ZeroGuess_()
type SADGuess_ <: InitialGuess end
SADGuess = SADGuess_()

function computeDensityMatrixGuess(mat::Matrix,dim::Int64)
  if (size(mat) != (dim,dim))
    error("Guess matrix has wrong dimensions: $(size(mat)) (expected $((dim,dim)))")
  end
  return mat
end

function computeDensityMatrixGuess(::ZeroGuess_,dim::Int64)
  return zeros(dim,dim)
end

type DensityGuessSAD
  α
  β
end

# These are constructed by running a unrestricted SCF calculation on the single atom with the corresponding method and basis set
densityGuessSADHartreeFock = DensityGuessSAD(
#alpha spin density
Dict("STO-3G" => Dict(
	Element("H") => [ 1.000000 ],
	Element("O") => [ 1.065739   -0.264689    0.000000    0.000000    0.000000;
                         -0.264689    1.065739    0.000000    0.000000    0.000000;
                          0.000000    0.000000    1.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    1.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000 ]
	),
     "def2-SVP" => Dict(
        Element("H") => [ 0.345487    0.294542    0.000000    0.000000    0.000000;
                          0.294542    0.251108    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000 ],
	Element("O") => [ 1.083840    0.188236    0.143345    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.188236    0.470190    0.284077    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.143345    0.284077    0.172502    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.500580    0.000000    0.000000    0.308265   -0.000000   -0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.473992    0.000000    0.000000    0.315137    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.500580    0.000000   -0.000000    0.308265    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.308265    0.000000    0.000000    0.189835    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000   -0.000000    0.315137   -0.000000    0.000000    0.209522   -0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000   -0.000000    0.000000    0.308265    0.000000   -0.000000    0.189835    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000 ]
	)
),
#beta spin density
Dict("STO-3G" => Dict(
	Element("H") => [ 0.000000 ],
	Element("O") => [ 1.065739   -0.264689    0.000000    0.000000    0.000000;
                         -0.264689    1.065739    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000 ]
	),
     "def2-SVP" => Dict(
        Element("H") => [ 0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000 ],
	Element("O") => [ 1.075394    0.156592    0.157593    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.156592    0.367558    0.298568    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.157593    0.298568    0.243442    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000   -0.000000   -0.000000    0.000000   -0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000   -0.000000    0.425423    0.000000   -0.000000    0.324962   -0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000   -0.000000    0.000000    0.000000   -0.000000    0.000000   -0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000   -0.000000   -0.000000    0.000000   -0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000   -0.000000    0.324962    0.000000   -0.000000    0.248224   -0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000   -0.000000   -0.000000    0.000000   -0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000;
                          0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000    0.000000 ]
	)
)
)

function computeDensityGuessSAD(methodName, basisSetName, geometry)
  if methodName=="HF" || methodName=="HartreeFock"
    αDensities = densityGuessSADHartreeFock.α[basisSetName]
    βDensities = densityGuessSADHartreeFock.β[basisSetName]
    αResult = cat([1,2],[αDensities[atom.element] for atom in geometry.atoms]...)
    βResult = cat([1,2],[βDensities[atom.element] for atom in geometry.atoms]...)
    return (αResult, βResult)
  else
    error("No SAD guess implemented for method "*methodName)
  end
end
end #module
