Skip to content
Snippets Groups Projects
Commit ca91349d authored by moritz.buchhorn's avatar moritz.buchhorn
Browse files

First draft on week 4

parent c1b85bd4
Branches
No related tags found
No related merge requests found
%% Cell type:markdown id: tags:
# B.CTC: DFT
In der Vorlesung haben sie den Hartree-Fock (HF)-Ansatz und die Dichtefunktionaltheorie (DFT) kennengelernt. Die Implementation als Kohn-Sham-DFT ist dem HF-Ansatz aus technischer Perspektive sehr ähnlich, weshalb jede (uns bekannte) Quantenchemie-Software beide Methoden zur Verfügung stellt.
In dieser Demonstration wollen wir uns die Ergebnisse einer Geometrieoptimierung für verschiedene Methoden genauer anschauen. Da DFT nicht systematisch verbesserbar ist, wird zu Beginn einer theoretischen Studie immer ein sogenanntes _Benchmarking_ durchgeführt, ähnlich einer Kalibration in der Analytik. Mit verschiedenen Funktionalen wird die geforderte Größe (hier Bindungslängen) einer bekannten Verbindung berechnet. Mit dem Funktional, das die experimentelle Größe am besten wiedergibt, wird dann die entsprechende Größe für eine andere Verbindung berechnet. Dieses Verfahren funktioniert umso besser, je ähnlicher sich die Referenzverbindungen und die untersuchten Verbindungen sind.
In diesem Notebook verwenden wir Formaldehyd, das hier beste Funktional könnte also zum Beispiel sehr gut für Acetaldehyd verwendet werden.
%% Cell type:markdown id: tags:
%% Cell type:code id: tags:
``` python
import subprocess as sub
import support
```
%% Cell type:code id: tags:
``` python
for basis_set in ("sto-3g", "def2-svp", "def2-tzvp"):
with open(f"hcho-hf-{basis_set}-opt.inp", "w") as f:
f.write(f"""\
! RHF {basis_set} Opt
* xyzfile 0 1 hcho-guess.xyz
""")
sub.run(f"orca hcho-hf-{basis_set}-opt.inp > hcho-hf-{basis_set}-opt.out",
shell=True)
```
%% Cell type:markdown id: tags:
Für HF-Berechnungen ist ersichtlich, dass die Bindungslänge für größere Basissätze als zu kurz berechnet wird. Da ein größerer Basissatz _immer_ ein besseres Ergebnis liefert, muss der Fehler in der Vorhersage durch den HF-Ansatz verursacht sein.
Es mag antiintuitiv erscheinen, dass die Bindungslänge mit einem angeblich schlechteren (kleineren) Basissatz scheinbar besser vorhergesagt wird. Das ist ein klassisches Beispiel für ein "richtig aus dem falschen Grund", ein in der Quantenchemie recht häufig auftretendes Ergebnis. Durch Fehlerausgleich werden scheinbar korrekte Größen bestimmt, und es liegt an ihnen als Wissenschaftler, die Plausibilität dieser Größen zu beurteilen.
%% Cell type:code id: tags:
``` python
for functional in ("bp", "pbe", "tpss", "tpssh"):
with open(f"hcho-{functional}-opt.inp", "w") as f:
f.write(f"""\
! RKS {functional} def2-TZVP Opt
* xyzfile 0 1 hcho-guess.xyz
""")
sub.run(f"orca hcho-{functional}-opt.inp > hcho-{functional}-opt.out",
shell=True)
```
%% Cell type:code id: tags:
``` python
exp = 1.21
for basis_set in ("sto-3g", "def2-svp", "def2-tzvp"):
coords = support.XYZ(f"hcho-hf-{basis_set}-opt.xyz")
print(coords.get_bond_length(0, 1))
for functional in ("bp", "pbe", "tpss", "tpssh"):
coords = support.XYZ(f"hcho-{functional}-opt.xyz")
print(coords.get_bond_length(0, 1))
```
%% Cell type:markdown id: tags:
Im Gegensatz zur Basissatzgröße sind Dichtefunktionale nicht systematisch verbesserbar, das heißt es gibt keinen Parameter, der mathematisch beweisbar das Ergebnis verbessert. Es ist daher ein plausibles Ergebnis, dass das "anspruchsvollere" Funktional TPSSh die schlechteste Bindungslänge liefert. In diesem Fall dürfen wir uns einfach freuen, dass die Funktionale mit geringen Rechenzeiten so gute Ergebnisse liefern und diese als bevorzugte Funktionale für unsere Rechnung verbuchen.
import numpy as np
class XYZ:
def __init__(self, file: str):
"""
Reads in files with the structure 'element x y z' and returns
a numpy array with the stored coordinates.
@param file: String. Path to the file to be read in.
@return: 2D numpy array. Contains the coordinates from the file at the
corresponding indices.
"""
with open(file) as f:
lines = f.readlines()
if len(lines[0].split()) == 1:
lines = lines[2:]
self.xyz = np.array([line.split()[1:] for line in lines]).astype(float)
self.elements = [line.split()[0] for line in lines]
def get_bond_length(self, i: int, j: int) -> float:
"""
Returns the bond length between atoms i and j.
@param i: Integer. Index of the first atom.
@param j: Integer. Index of the second atom.
@return: Float. The bond length between the two atoms.
"""
return np.linalg.norm(self.xyz[i] - self.xyz[j])
\ No newline at end of file
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment