Skip to content

Commit

Permalink
Merge branch 'FCI' of github.com:gustavojra/Fermi.jl into FCI
Browse files Browse the repository at this point in the history
  • Loading branch information
@gustavojra
gustavojra committed Jul 22, 2022
2 parents b228b8b + 90be2b0 commit fadf046
Show file tree
Hide file tree
Showing 5 changed files with 15 additions and 22 deletions.
11 changes: 1 addition & 10 deletions src/Core/Molecule.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -23,8 +23,6 @@ Object storing information about a molecule (group of atoms).
multiplicity Multiplicity ``(2Mₛ + 1)``
Nα Number of Alpha electrons
Nβ Number of Beta electrons
Vnuc Nuclear repulsion energy
# Examples:
Expand All @@ -46,7 +44,6 @@ H 0.919788188200 2.458018557000 0.629793883200
Charge: 0 Multiplicity: 1
Nuclear repulsion: 8.8880641737
julia> Molecule(charge=2, multiplicity=3)
Molecule:
Expand All @@ -57,7 +54,6 @@ H 0.919788188200 2.458018557000 0.629793883200
Charge: 2 Multiplicity: 3
Nuclear repulsion: 8.8880641737
```
"""
struct Molecule
Expand All @@ -66,7 +62,6 @@ struct Molecule
multiplicity::Int
Nα::Int
Nβ::Int
Vnuc::Float64
end

function Molecule(;
Expand All @@ -81,9 +76,6 @@ function Molecule(;
end

function Molecule(atoms::Vector{T}, charge::Int, multiplicity::Int) where T <: Atom

# Compute Nuclear repulsion
Vnuc = Molecules.nuclear_repulsion(atoms)

# Compute number of electrons
nelec = -charge
Expand All @@ -107,7 +99,7 @@ function Molecule(atoms::Vector{T}, charge::Int, multiplicity::Int) where T <: A

Nβ = (nelec - αexcess)/2
Nα = nelec - Nβ
out = Molecule(atoms, charge, multiplicity, Nα, Nβ, Vnuc)
out = Molecule(atoms, charge, multiplicity, Nα, Nβ)
return out
end

Expand All @@ -123,7 +115,6 @@ function string_repr(M::Molecule)
out = out*format("\n")
out = out*format("\nCharge: {} ", M.charge)
out = out*format("Multiplicity: {} \n", M.multiplicity)
out = out*format("Nuclear repulsion: {:15.10f}", M.Vnuc)
return out
end

Expand Down
6 changes: 4 additions & 2 deletions src/Methods/HartreeFock/RHF/Direct.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -54,15 +54,17 @@ function RHF(ints::IntegralHelper{Float64, <:AbstractERI, AtomicOrbitals}, C::Ar
diis_start = Options.get("diis_start")
end

# Grab ndocc,nvir
# Grab ndocc,nvir,Vnuc
ndocc = try
Int((molecule.Nα + molecule.Nβ)/2)
catch InexactError
throw(FermiException("Invalid number of electrons $(molecule.Nα + molecule.Nβ) for RHF method."))
end
nvir = size(C,2) - ndocc
nao = size(C,1)
Vnuc = Molecules.nuclear_repulsion(molecule.atoms)

output("Nuclear repulsion: {:15.10f}", Vnuc)
output(" Number of AOs: {:5.0d}", nao)
output(" Number of Doubly Occupied Orbitals: {:5.0d}", ndocc)
output(" Number of Virtual Spatial Orbitals: {:5.0d}", nvir)
Expand Down Expand Up @@ -117,7 +119,7 @@ function RHF(ints::IntegralHelper{Float64, <:AbstractERI, AtomicOrbitals}, C::Ar
Eelec = RHFEnergy(D, T + V, F)

# Compute Energy
Enew = Eelec + molecule.Vnuc
Enew = Eelec + Vnuc

# Store vectors for DIIS
if do_diis
Expand Down
6 changes: 4 additions & 2 deletions src/Methods/HartreeFock/RHF/RHFa.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -93,15 +93,17 @@ function RHF(ints::IntegralHelper{Float64, <:AbstractERI, AtomicOrbitals}, C::Ab
diis_start = Options.get("diis_start")
end

# Grab ndocc,nvir
# Grab ndocc,nvir,Vnuc
ndocc = try
Int((molecule.Nα + molecule.Nβ)/2)
catch InexactError
throw(FermiException("Invalid number of electrons $(molecule.Nα + molecule.Nβ) for RHF method."))
end
nvir = size(C,2) - ndocc
nao = size(C,1)
Vnuc = Molecules.nuclear_repulsion(molecule.atoms)

output("Nuclear repulsion: {:15.10f}", Vnuc)
output(" Number of AOs: {:5.0d}", nao)
output(" Number of Doubly Occupied Orbitals: {:5.0d}", ndocc)
output(" Number of Virtual Spatial Orbitals: {:5.0d}", nvir)
Expand Down Expand Up @@ -149,7 +151,7 @@ function RHF(ints::IntegralHelper{Float64, <:AbstractERI, AtomicOrbitals}, C::Ab
Eelec = RHFEnergy(D, T + V, F)

# Compute Energy
Enew = Eelec + molecule.Vnuc
Enew = Eelec + Vnuc

# Store vectors for DIIS
if do_diis
Expand Down
5 changes: 3 additions & 2 deletions src/Methods/HartreeFock/UHF/UHFa.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -50,6 +50,7 @@ function UHF(ints::IntegralHelper{Float64, <:AbstractERI, AtomicOrbitals}, Cα::

Nα = molecule.Nα
Nβ = molecule.Nβ
Vnuc = Molecules.nuclear_repulsion(molecule.atoms)
S = ints["S"]
T = ints["T"]
V = ints["V"]
Expand Down Expand Up @@ -80,7 +81,7 @@ function UHF(ints::IntegralHelper{Float64, <:AbstractERI, AtomicOrbitals}, Cα::


build_fock!(Fα, Fβ, Jα, Jβ, Kα, Kβ, Dα, Dβ, ints)
output(" Guess Energy {:20.14f}", UHFEnergy(H, Dα, Dβ, Fα, Fβ, molecule.Vnuc))
output(" Guess Energy {:20.14f}", UHFEnergy(H, Dα, Dβ, Fα, Fβ, Vnuc))

output("\n Iter. {:>15} {:>10} {:>10} {:>8} {:>8} {:>8}", "E[UHF]", "ΔE", "Dᵣₘₛ", "t", "DIIS", "damp")
output(repeat("-",80))
Expand Down Expand Up @@ -112,7 +113,7 @@ function UHF(ints::IntegralHelper{Float64, <:AbstractERI, AtomicOrbitals}, Cα::
build_fock!(Fα, Fβ, Jα, Jβ, Kα, Kβ, Dα, Dβ, ints)

# Calculate energy
E = UHFEnergy(H, Dα, Dβ, Fα, Fβ, molecule.Vnuc)
E = UHFEnergy(H, Dα, Dβ, Fα, Fβ, Vnuc)

# Store vectors for DIIS
if do_diis
Expand Down
9 changes: 3 additions & 6 deletions test/test_molecule.jl
View file
Original file line number Diff line number Diff line change
Expand Up @@ -12,8 +12,7 @@
length(mol.atoms) == 5 &&
mol.charge == 0 &&
mol.multiplicity == 1 &&
mol.Nα == mol.Nβ == 5 &&
mol.Vnuc ≈ 13.6865180253
mol.Nα == mol.Nβ == 5
end

molstringB = """
Expand All @@ -27,8 +26,7 @@
mol.charge == 0 &&
mol.multiplicity == 3 &&
mol.Nα == 5 &&
mol.Nβ == 3 &&
mol.Vnuc ≈ 6.0160249492
mol.Nβ == 3
end

# Test wrong charge, mult and unit
Expand All @@ -45,7 +43,6 @@
occursin(r"C\s+?[+-]{0,1}\d+?\.\d+?\s+?[+-]{0,1}\d+?\.\d+?\s+[+-]{0,1}\d+?\.\d+?",x) &&
occursin(r"H\s+?[+-]{0,1}\d+?\.\d+?\s+?[+-]{0,1}\d+?\.\d+?\s+[+-]{0,1}\d+?\.\d+?",x) &&
occursin(r"Charge:\s+?0", x) &&
occursin(r"Multiplicity:\s+?", x) &&
occursin(r"Nuclear\s+?repulsion:\s+?\d+?\.\d+?", x)
occursin(r"Multiplicity:\s+?", x)
end
end

0 comments on commit fadf046

Please sign in to comment.