6  Basisfunktionen

Die Roothaan-Gleichungen aus Kapitel 03 benötigen eine endliche Menge von Basisfunktionen \(\{\chi_\mu\}\). Die Wahl dieser Funktionen ist keine technische Nebensache – sie bestimmt massgeblich, wie genau die Rechnung ist und wie teuer sie wird. Dieses Kapitel erklärt, wie Basisfunktionen konstruiert sind, welche Typen es gibt, und wie man für eine konkrete Rechnung einen geeigneten Basissatz wählt.


6.1 Slater-Typ-Orbitale (STO)

Die natürlichste Wahl für Basisfunktionen sind Atomorbitale. Die exakten Wasserstoff-Atomorbitale haben eine radiale Abhängigkeit \(\sim e^{-\zeta r}\) – einen exponentiellen Abfall mit dem Kernabstand \(r\). Slater-Typ-Orbitale (STO) verallgemeinern das:

\[\chi^\text{STO}_{nlm}(\mathbf{r}) = N\, r^{n-1}\, e^{-\zeta r}\, Y_l^m(\theta,\phi) \tag{4.1}\]

STOs beschreiben die Elektronendichte nahe am Kern und weit entfernt davon sehr gut. Das Problem: Die Zweielectronenintegrale \((\nu\mu|\sigma\lambda)\) aus Gl. (3.11) lassen sich mit STOs nicht analytisch berechnen – und numerisch ist das für grössere Moleküle zu teuer.


6.2 Gausstyp-Orbitale (GTO)

6.2.1 Form und Klassifikation

Gausstyp-Orbitale (GTO) ersetzen den exponentiellen Abfall durch eine Gaussfunktion:

\[\chi^\text{GTO}_{ijk}(\mathbf{r}) = N\, x^i y^j z^k\, e^{-\alpha r^2} \tag{4.2}\]

Die ganzzahligen Potenzen \(i, j, k \geq 0\) kodieren die Winkelabhängigkeit: ihre Summe \(l = i+j+k\) entspricht der Drehimpulsquantenzahl des Orbitaltyps.

\(l = i+j+k\) Typ Kombinationen Anzahl
0 s \((0,0,0)\) 1
1 p \((1,0,0),\,(0,1,0),\,(0,0,1)\) 3
2 d \((2,0,0),\,(0,2,0),\ldots,(1,1,0)\) 6

Ein s-GTO (\(i=j=k=0\)) hat keinen Vorfaktor ausser \(N\) – er ist kugelsymmetrisch. Ein p-GTO z.B. mit \((i,j,k)=(1,0,0)\) enthält den Faktor \(x\), der einen Knoten in der \(yz\)-Ebene erzwingt – genau das Knotenverhalten eines \(p_x\)-Orbitals.

6.2.2 Der entscheidende Rechenvorteil

Das Gaussian-Produkttheorem besagt: Das Produkt zweier Gaussfunktionen auf verschiedenen Zentren \(\mathbf{R}_A\) und \(\mathbf{R}_B\) ist wieder eine Gaussfunktion, zentriert auf einem Punkt zwischen den beiden:

\[e^{-\alpha|\mathbf{r}-\mathbf{R}_A|^2} \cdot e^{-\beta|\mathbf{r}-\mathbf{R}_B|^2} = K\, e^{-\gamma|\mathbf{r}-\mathbf{R}_P|^2}\]

Das macht alle Zweielectronenintegrale \((\nu\mu|\sigma\lambda)\) analytisch berechenbar – der entscheidende Vorteil gegenüber STOs, der ab initio Rechnungen für Moleküle überhaupt erst praktikabel macht.

6.2.3 GTO vs. STO: eine visuelle Gegenüberstellung

GTOs haben eine bekannte Schwäche: Sie zeigen am Kern keinen Cusp (d.h. keine endliche Steigung bei \(r=0\), obwohl das physikalisch korrekt wäre), und sie fallen für grosse \(r\) zu schnell ab. Beide Probleme lassen sich durch Linearkombinationen mehrerer GTOs mildern.

Code
r = np.linspace(0, 5, 500)

# STO-1s (ζ=1, analytisch exakt für H)
zeta = 1.0
sto = 2 * zeta**1.5 * np.exp(-zeta * r)

# Einzelne GTO (optimaler Exponent für 1s-Annäherung)
alpha_single = 0.4166
N_single = 2 * (alpha_single/np.pi)**0.75
gto_single = N_single * np.exp(-alpha_single * r**2)

# STO-3G Kontraktionskoeffizienten und Exponenten für H 1s
# (Hehre, Stewart, Pople 1969)
d_coeff = np.array([0.4446345, 0.5353281, 0.1543290])
alphas   = np.array([0.1688554, 0.6239137, 3.4252509])
gto3g = sum(d * 2*(a/np.pi)**0.75 * np.exp(-a*r**2)
            for d, a in zip(d_coeff, alphas))

fig, axes = plt.subplots(1, 2, figsize=(11, 4))

# Linkes Bild: Wellenfunktionen
ax = axes[0]
ax.plot(r, sto,        '#333333', lw=2.5, label='STO (exakt)')
ax.plot(r, gto_single, '#d6604d', lw=2,   ls='--', label='einzelne GTO')
ax.plot(r, gto3g,      '#2166ac', lw=2,   ls='-.',  label='STO-3G (3 GTOs)')
ax.set_xlabel('$r$ / Bohr', fontsize=12)
ax.set_ylabel('$\\chi(r)$', fontsize=12)
ax.set_title('Wellenfunktionen', fontsize=12)
ax.legend(fontsize=10)
ax.set_xlim(0, 4)
ax.grid(True, alpha=0.3)

# Rechtes Bild: Differenz zur STO
ax2 = axes[1]
ax2.plot(r, gto_single - sto, '#d6604d', lw=2, ls='--', label='einzelne GTO − STO')
ax2.plot(r, gto3g - sto,      '#2166ac', lw=2, ls='-.',  label='STO-3G − STO')
ax2.axhline(0, color='gray', lw=0.8)
ax2.set_xlabel('$r$ / Bohr', fontsize=12)
ax2.set_ylabel('Fehler $\\Delta\\chi(r)$', fontsize=12)
ax2.set_title('Abweichung von der STO', fontsize=12)
ax2.legend(fontsize=10)
ax2.set_xlim(0, 4)
ax2.set_ylim(-0.15, 0.15)
ax2.grid(True, alpha=0.3)

plt.suptitle('STO vs. GTO für H 1s', fontsize=13)
plt.tight_layout()
plt.show()
Abbildung 6.1: STO-1s vs. einzelne GTO vs. STO-3G (Linearkombination von 3 GTOs) für ein wasserstoffartiges 1s-Orbital. Die STO-3G-Kontraktion approximiert die STO sehr gut, obwohl jede einzelne GTO deutliche Abweichungen zeigt.

6.3 Primitive und kontrahierte Gaussfunktionen

6.3.1 Das Problem mit einzelnen GTOs

Eine einzelne GTO kann das Verhalten eines Atomorbitals nicht gut beschreiben – weder nahe am Kern noch weit entfernt. Man braucht mehrere.

Würde man alle GTOs einzeln (als primitive Gaussfunktionen) in der LCAO-Entwicklung verwenden, stiege die Zahl der Unbekannten \(c_{\mu m}\) und damit der Aufwand für die Zweielectronenintegrale (\(\sim N_b^4\)) stark an.

6.3.2 Kontraktion: feste Linearkombinationen

Die Lösung: Man fasst mehrere primitive GTOs zu einer kontrahierten Gaussfunktion zusammen, deren Koeffizienten \(d_i\) und Exponenten \(\alpha_i\) einmal (aus atomaren SCF-Rechnungen) bestimmt werden und dann fest bleiben:

\[\chi_\mu^\text{kontr}(\mathbf{r}) = \sum_{i=1}^{L} d_i\, g_i(\mathbf{r}) \tag{4.3}\]

In der LCAO-Entwicklung sind \(\chi_\mu^\text{kontr}\) die Basisfunktionen. Die Koeffizienten \(c_{\mu m}\) (die variationell optimiert werden) sind damit viel weniger zahlreich als wenn man alle Primitiven einzeln verwenden würde.

Code
r = np.linspace(0, 4, 400)

# gto3g auf dem neuen Gitter neu berechnen (d_coeff und alphas aus vorigem Block)
gto3g = sum(d * 2*(a/np.pi)**0.75 * np.exp(-a*r**2)
            for d, a in zip(d_coeff, alphas))

fig, ax = plt.subplots(figsize=(7, 4))
colors_prim = ['#f4a582', '#d6604d', '#b2182b']

# Einzelne normierte primitive GTOs (ungewichtet)
for i, (d, a, col) in enumerate(zip(d_coeff, alphas, colors_prim)):
    prim = 2*(a/np.pi)**0.75 * np.exp(-a*r**2)
    ax.plot(r, prim, color=col, lw=1.5, ls='--',
            label=f'Primitiv $g_{i+1}$ ($\\alpha={a:.3f}$)', alpha=0.8)
    # Gewichteter Beitrag (klein)
    ax.fill_between(r, 0, d*prim, alpha=0.08, color=col)

# Kontrahierte Funktion
ax.plot(r, gto3g, '#2166ac', lw=3, label='Kontraktion $\\chi^\\mathrm{STO-3G}$')

ax.set_xlabel('$r$ / Bohr', fontsize=12)
ax.set_ylabel('Amplitude', fontsize=12)
ax.set_title('Primitive GTOs → kontrahierte Basisfunktion (STO-3G, H 1s)', fontsize=12)
ax.legend(fontsize=9, loc='upper right')
ax.set_xlim(0, 3.5)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
Abbildung 6.2: Kontraktion: Drei primitive GTOs (dünn, gestrichelt) werden mit festen Koeffizienten zur kontrahierten STO-3G Basisfunktion (dick, blau) kombiniert.

6.4 Die Basissatzhierarchie

6.4.1 Minimaler Basissatz

Im minimalen Basissatz wird jedes Atomorbital der Valenztheorie durch genau eine (kontrahierte) Basisfunktion beschrieben:

  • H, He: 1 Funktion (1s)
  • C, N, O, F: 5 Funktionen (1s, 2s, 2p\(_x\), 2p\(_y\), 2p\(_z\))

Für H₂O ergibt das \(2 \times 1 + 5 = 7\) Basisfunktionen – das STO-3G-Set ist der Standardvertreter. Es ist schnell, aber für quantitative Vorhersagen zu unflexibel.

6.4.2 Split-Valence: der praktische Kompromiss

Die physikalische Motivation: Core-Elektronen (z.B. O 1s) sind chemisch inaktiv und kaum polarisierbar – ein Satz Basisfunktionen genügt. Valenzelectronen bilden Bindungen, reagieren auf Nachbaratome und sind polarisierbar – sie brauchen mehr Flexibilität.

Ein Split-Valence (SV) Basissatz gibt jedem Valenzorbital zwei statt einer Basisfunktion: eine enge (kernnahe) und eine diffuse (ausgedehnte).

Die Notation \(m\)-\(np\)G bedeutet:

  • \(m\): Core-AOs = 1 kontrahierte GTO aus \(m\) Primitiven
  • \(n\): äussere Valenz-GTO aus \(n\) Primitiven
  • \(p\): innere Valenz-GTO = \(p\) Primitive (oft \(p=1\): eine einzelne GTO)

Das bekannteste Beispiel ist 6-31G: Core aus 6 Primitiven, Valenz aufgeteilt in eine 3-Primitiven-Kontraktion (aussen) und eine einzelne Primitive (innen).

Code
# Echte 6-31G Exponenten und Koeffizienten für O
# Quelle: Hehre, Ditchfield, Pople (1972)
# Valenz-2s: enger Anteil (3 Primitive) und diffuser Anteil (1 Primitiv)
alpha_tight = np.array([5.0331513, 1.1695961, 0.3803890])
d_tight_2s  = np.array([-0.1117496, -0.1480263,  1.1308439])  # 2s contracted

alpha_diff  = np.array([0.1212778])
d_diff_2s   = np.array([1.0])  # unnormiert, zur Illustration

r = np.linspace(0, 6, 500)

def gto_contracted(r, alphas, d_coeffs):
    result = np.zeros_like(r)
    for a, d in zip(alphas, d_coeffs):
        N = (2*a/np.pi)**0.75
        result += d * N * np.exp(-a * r**2)
    return result

chi_tight = gto_contracted(r, alpha_tight, d_tight_2s)
chi_diff  = gto_contracted(r, alpha_diff,  d_diff_2s)

# Normieren zur Darstellung
chi_tight /= np.max(np.abs(chi_tight))
chi_diff  /= np.max(np.abs(chi_diff))

fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(r, chi_tight, '#d6604d', lw=2.5, label='enge Funktion (3 Primitive, $\\alpha$ gross)')
ax.plot(r, chi_diff,  '#2166ac', lw=2.5, label='diffuse Funktion (1 Primitiv, $\\alpha$ klein)')
ax.axhline(0, color='gray', lw=0.8)
ax.set_xlabel('$r$ / Bohr', fontsize=12)
ax.set_ylabel('Amplitude (normiert)', fontsize=12)
ax.set_title('Split-Valence: O 2s im 6-31G-Basissatz', fontsize=12)
ax.legend(fontsize=10)
ax.set_xlim(0, 5)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
Abbildung 6.3: Split-Valence-Idee für O 2s im 6-31G-Basissatz: zwei Basisfunktionen (eng + diffus) statt einer. Die Überlagerung beider beschreibt die Polarisierbarkeit des Orbitals in der Molekülumgebung viel besser.

6.4.3 s- vs. p-Funktionen: der geometrische Unterschied

Im 6-31G-Set für Sauerstoff verwenden 2s- und 2p-Funktionen dieselben Exponenten \(\alpha\) – das sogenannte SP-Kontraktionsschema spart Integrale. Was sich unterscheidet, ist der präexponentielle Faktor:

\[\chi_s = N\, e^{-\alpha r^2} \quad \text{vs.} \quad \chi_{p_x} = N\, x\, e^{-\alpha r^2}\]

Der Faktor \(x\) in \(\chi_{p_x}\) erzwingt einen Knoten in der \(yz\)-Ebene (\(x=0\)) – das p-Orbital hat also einen geometrischen Knoten, der s-Orbital keinen.

Code
alpha = 0.3803890
lim, n = 3.5, 200
x1d = np.linspace(-lim, lim, n)
X, Z = np.meshgrid(x1d, x1d)
r2 = X**2 + Z**2

chi_s  =          np.exp(-alpha * r2)
chi_px = X      * np.exp(-alpha * r2)

fig, axes = plt.subplots(1, 2, figsize=(10, 4))
vmax_s  = np.max(np.abs(chi_s))
vmax_px = np.max(np.abs(chi_px))

for ax, chi, vmax, title in zip(
    axes,
    [chi_s, chi_px],
    [vmax_s, vmax_px],
    ['s-GTO: $N e^{-\\alpha r^2}$', '$p_x$-GTO: $N\\,x\\,e^{-\\alpha r^2}$']
):
    im = ax.contourf(X, Z, chi, levels=20,
                     cmap='RdBu_r', vmin=-vmax, vmax=vmax)
    ax.contour(X, Z, chi, levels=[0], colors='k', linewidths=1.5)
    plt.colorbar(im, ax=ax)
    ax.set_xlabel('$x$ / Bohr', fontsize=11)
    ax.set_ylabel('$z$ / Bohr', fontsize=11)
    ax.set_title(title, fontsize=12)
    ax.set_aspect('equal')

plt.suptitle(f'Gleicher Exponent $\\alpha = {alpha}$, Schnitt durch $y=0$',
             fontsize=12)
plt.tight_layout()
plt.show()
Abbildung 6.4: s-GTO vs. p-GTO mit demselben Exponenten α=0.38. Links: Schnitt durch die xz-Ebene (y=0). Das s-GTO ist rotationssymmetrisch; das p_x-GTO hat einen Knoten bei x=0 und wechselt das Vorzeichen.

6.4.4 Polarisationsfunktionen

Eine weitere Verbesserung: Polarisationsfunktionen fügen Basisfunktionen mit höherer Drehimpulsquantenzahl hinzu als für das Atom nötig wäre. Sie erlauben es, Atomorbitale in der Molekülumgebung zu verzerren (polarisieren):

  • d-Funktionen auf C, N, O: beschreiben, wie p-Orbitale sich in Bindungen asymmetrisch verformen
  • p-Funktionen auf H: beschreiben, wie das 1s-Orbital von H in einer Bindung seine Kugelsymmetrie verliert

Die Notationen * und ** in Basissatznamen kodieren das:

Basissatz Bedeutung
6-31G nur Split-Valence
6-31G* + d auf C,N,O,F (nicht auf H)
6-31G** + d auf C,N,O,F und p auf H

6.5 Überblick: Basissätze im Vergleich für H₂O

Code
import pandas as pd

bases = ['sto-3g', '6-31g', '6-31g*', '6-31g**', 'cc-pvdz', 'cc-pvtz']
rows  = []

for basis in bases:
    try:
        m = gto.Mole(verbose=0)
        m.atom = '''
            O  0.000000  0.000000  0.117176
            H  0.000000  0.756950 -0.468703
            H  0.000000 -0.756950 -0.468703
        '''
        m.basis = basis
        m.build()
        mf = scf.RHF(m)
        mf.verbose = 0
        E = mf.kernel()
        rows.append({'Basissatz': basis, 'N_b': m.nao_nr(),
                     'E_HF / Hartree': round(E, 6)})
    except Exception:
        rows.append({'Basissatz': basis, 'N_b': '—', 'E_HF / Hartree': '—'})

df = pd.DataFrame(rows)
df['ΔE vs. STO-3G / kcal/mol'] = [
    round((r - rows[0]['E_HF / Hartree']) * 627.509, 2)
    if isinstance(r, float) else '—'
    for r in df['E_HF / Hartree']
]
print(df.to_string(index=False))
Tabelle 6.1: Anzahl Basisfunktionen und HF-Energie für H₂O mit verschiedenen Basissätzen.
Basissatz  N_b  E_HF / Hartree  ΔE vs. STO-3G / kcal/mol
   sto-3g    7      -74.962928                      0.00
    6-31g   13      -75.983998                   -640.73
   6-31g*   18      -76.009132                   -656.50
  6-31g**   24      -76.022648                   -664.98
  cc-pvdz   24      -76.026799                   -667.59
  cc-pvtz   58      -76.057169                   -686.65
Code
numeric_rows = [r for r in rows if isinstance(r['E_HF / Hartree'], float)]
labels = [r['Basissatz'] for r in numeric_rows]
energies = [r['E_HF / Hartree'] for r in numeric_rows]
nb_vals  = [r['N_b'] for r in numeric_rows]

fig, ax = plt.subplots(figsize=(8, 4))
ax.plot(range(len(labels)), energies, 'o-', color='#2166ac',
        lw=2.5, ms=9, markerfacecolor='white', markeredgewidth=2.5)

for i, (e, nb) in enumerate(zip(energies, nb_vals)):
    ax.annotate(f'$N_b={nb}$', xy=(i, e), xytext=(0, 12),
                textcoords='offset points', ha='center',
                fontsize=9, color='#555')

ax.set_xticks(range(len(labels)))
ax.set_xticklabels(labels, rotation=15, fontsize=11)
ax.set_ylabel('$E_\\mathrm{HF}$ / Hartree', fontsize=12)
ax.set_title('HF-Energie von H₂O als Funktion des Basissatzes', fontsize=12)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()
Abbildung 6.5: Konvergenz der HF-Energie von H₂O mit wachsendem Basissatz. Die Energie sinkt monoton gegen das HF-Limit.
WichtigBasissatz-Trunkierungsfehler

Mit einem endlichen Basissatz erhält man immer eine Energie oberhalb des HF-Limits (Variationsprinzip). Die Differenz heisst Basissatz-Trunkierungsfehler. Er ist verschieden vom Korrelationsfehler (Kapitel 03): Selbst mit einem unendlich grossen Basissatz würde HF die Korrelationsenergie nicht erfassen.


6.6 Zusammenfassung

STO  →  physikalisch korrekt, aber Integrale nicht analytisch
GTO  →  Integrale analytisch (Gauss-Produkttheorem), aber einzeln ungenau
       ↓
Kontraktion: feste Linearkomb. von Primitiven → eine Basisfunktion
       ↓
Minimaler BS  →  1 Funktion/AO          (schnell, ungenau)
Split-Valence →  2 Funktionen/Valenz-AO (guter Kompromiss)
+ Polarisation → höheres l              (wichtig für Geometrien, Energien)
TippVorschau

Kapitel 05 verlässt die Wellenfunktions-Welt und formuliert das Problem neu: Statt die 30-dimensionale Wellenfunktion von H₂O zu suchen, genügt es – so das Hohenberg-Kohn-Theorem – die dreidimensionale Elektronendichte \(\rho(\mathbf{r})\) zu kennen. Das ist der Ausgangspunkt der Dichtefunktionaltheorie.