from hfscf import *
def SCF (r = 1.4632, Z=[1,1], b1 = GTO["H"], b2 = GTO["H"], b = 3, vbs=True):
R = [0, r]
if vbs: print("1) Генерируем матрицу перекрываний S")
s_scf = S(R, b1, b2, b)
if vbs: print("\n 2) Генерируем гамильтониан H")
h_scf = H(R, Z, b1, b2, b)
# Диагонализируем матрицу S и находим матрицу X
if vbs: print("\n 3) Диагонализируем матрицу S и находим матрицу X")
X = diagon(m=s_scf)
Xa = X.getH()
# Рассчитываем матрицу расчетной плотности P
if vbs: print("\n 4) Создаем матрицу расчетной плотности P")
p_scf = np.matrix([[0,0],[0,0]], dtype=np.float64) # Referencia (7) p. 148
# Начинаем итерационный процесс
if vbs: print("\n 5) Начинаем с функции SCF")
for iteracion in range(50):
# Конструируем матрицу Фока
# F = H + G
if vbs: print("\n 6) Генерируем матрицу Фока: рассчитываем интегралы для 2х электронов")
g_scf = G(r, p_scf, b1, b2, b)
f_scf = h_scf + g_scf
# Генерируем матрицу F'
# F' = X_adj * F * X
if vbs: print("n\ 7) Смена базиса F")
f_tra = Xa * f_scf * X
# Диагонализируем F' и генерируем C'
if vbs: print("n\ 8) Диагонализируем F' и генерируем C'")
c_tra = diagon2(m=f_tra)
# Конструируем матрицу C
# C = X * C'
if vbs: print("n\ 9) Конструируем матрицу коэффициентов C")
c_scf = X * c_tra
# Конструируем матрицу Р на основе матрицы С
if vbs: print("n\ 10) Пересчитываем матрицу расчетной плотности P")
p_temp = P(C=c_scf)
# рассчитайтываем электронную энергию молекулы и полную энергию молекулы
energy_electron = ener_elec(p_temp, h_scf, f_scf)
energy_sum = ener_tot(r, Z, 0)
# строим графики 1D и 2D
orbital(c_tra, r, b1, b2, b)
orbital2D(c_scf, X, f_scf, r, b1, b2, b)
print("\n Номер цикла " + str(iteracion + 1) + "\n")
# Проверяем сходимость
if np.linalg.norm(p_temp - p_scf) < 1E-4:
print("\n\n-->Расчитанное поле сходится")
return {"S":s_scf,"H":h_scf,"X": X,"F":f_scf,"C":c_scf,"P":p_temp, "EE":energy_electron, "ES":energy_sum}
else:
p_scf = p_temp
print("\n\n-->Расчитанное поле НЕ сходится!\n")
return {"S":s_scf,"H":h_scf,"X": X,"F":f_scf,"C":c_scf,"P":p_temp, "EE":energy_electron, "ES":energy_sum}
SCF()
1) Генерируем матрицу перекрываний S 2) Генерируем гамильтониан H 3) Диагонализируем матрицу S и находим матрицу X 4) Создаем матрицу расчетной плотности P 5) Начинаем с функции SCF 6) Генерируем матрицу Фока: рассчитываем интегралы для 2х электронов n\ 7) Смена базиса F n\ 8) Диагонализируем F' и генерируем C' n\ 9) Конструируем матрицу коэффициентов C n\ 10) Пересчитываем матрицу расчетной плотности P
Номер цикла 1 6) Генерируем матрицу Фока: рассчитываем интегралы для 2х электронов n\ 7) Смена базиса F n\ 8) Диагонализируем F' и генерируем C' n\ 9) Конструируем матрицу коэффициентов C n\ 10) Пересчитываем матрицу расчетной плотности P
Номер цикла 2 6) Генерируем матрицу Фока: рассчитываем интегралы для 2х электронов n\ 7) Смена базиса F n\ 8) Диагонализируем F' и генерируем C' n\ 9) Конструируем матрицу коэффициентов C n\ 10) Пересчитываем матрицу расчетной плотности P
Номер цикла 3 6) Генерируем матрицу Фока: рассчитываем интегралы для 2х электронов n\ 7) Смена базиса F n\ 8) Диагонализируем F' и генерируем C' n\ 9) Конструируем матрицу коэффициентов C n\ 10) Пересчитываем матрицу расчетной плотности P
Номер цикла 4 6) Генерируем матрицу Фока: рассчитываем интегралы для 2х электронов n\ 7) Смена базиса F n\ 8) Диагонализируем F' и генерируем C' n\ 9) Конструируем матрицу коэффициентов C n\ 10) Пересчитываем матрицу расчетной плотности P
Номер цикла 5 6) Генерируем матрицу Фока: рассчитываем интегралы для 2х электронов n\ 7) Смена базиса F n\ 8) Диагонализируем F' и генерируем C' n\ 9) Конструируем матрицу коэффициентов C n\ 10) Пересчитываем матрицу расчетной плотности P
Номер цикла 6 6) Генерируем матрицу Фока: рассчитываем интегралы для 2х электронов n\ 7) Смена базиса F n\ 8) Диагонализируем F' и генерируем C' n\ 9) Конструируем матрицу коэффициентов C n\ 10) Пересчитываем матрицу расчетной плотности P
Номер цикла 7 -->Расчитанное поле сходится
{'S': matrix([[0.99999999, 0.63749012], [0.63749012, 0.99999999]]), 'H': matrix([[-1.09920375, -0.91954841], [-0.91954841, -1.09920375]]), 'X': matrix([[ 0.55258063, 1.17442445], [ 0.55258063, -1.17442445]]), 'F': matrix([[-0.34684107, -0.57771139], [-0.57771139, -0.34684028]]), 'C': matrix([[ 0.55258114, -1.17442421], [ 0.55258013, 1.17442468]]), 'P': matrix([[0.61069182, 0.61069071], [0.61069071, 0.6106896 ]]), 'EE': -1.7974485548087507, 'ES': 0.683433570256971}