Beispiele

Beispiel 1: Exponentielles Wachstum

Exponentielles Wachstum ist eines der einfachsten und zugleich wichtigsten Phänomene, die sich mit Differentialgleichungen beschreiben lassen. Es tritt immer dann auf, wenn die Änderungsrate einer Größe proportional zu ihrem aktuellen Wert ist. Typische Beispiele sind das Wachstum von Bakterienpopulationen unter idealen Bedingungen oder die Verzinsung eines Kapitals bei kontinuierlicher Zinsrechnung.

Problemstellung

Eine Bakterienkultur wächst unter optimalen Bedingungen, wobei die Wachstumsrate proportional zur aktuellen Populationsgröße ist. Zum Zeitpunkt \(t = 0\) beträgt die Populationsgröße \(N_0 = 1000\) Bakterien. Die Wachstumsrate beträgt \(r = 0.5\) pro Stunde. Wir möchten die Populationsgröße nach drei Stunden bestimmen und das Wachstumsverhalten grafisch darstellen.

Mathematische Modellierung

Die Annahme, dass die Änderungsrate proportional zur Population ist, führt auf die Differentialgleichung: \[ \frac{\mathrm{d}N}{\mathrm{d}t} = rN \tag{1}\]

Die Proportionalitätskonstante \(r\) nennt man die Wachstumsrate. Sie gibt an, um welchen Bruchteil die Population pro Zeiteinheit wächst. Ein positiver Wert von \(r\) entspricht Wachstum, während ein negativer Wert zu exponentiellem Zerfall führen würde.

Diese Differentialgleichung ist separabel, denn die rechte Seite lässt sich als Produkt \(r \cdot N\) schreiben, wobei \(r\) nur von \(t\) abhängt (als Konstante) und \(N\) nur von sich selbst. Wir trennen die Variablen: \[ \frac{1}{N} \, \mathrm{d}N = r \, \mathrm{d}t \]

Integration beider Seiten ergibt: \[ \int \frac{1}{N} \, \mathrm{d}N = \int r \, \mathrm{d}t \] \[ \ln|N| = rt + C_1 \]

Durch Anwendung der Exponentialfunktion auf beide Seiten erhalten wir: \[ N(t) = Ce^{rt} \]

wobei \(C = e^{C_1}\) eine positive Konstante ist. Die Anfangsbedingung \(N(0) = N_0\) liefert: \[ N_0 = Ce^0 = C \]

Damit lautet die spezielle Lösung: \[ N(t) = N_0 e^{rt} \tag{2}\]

Für die gegebenen Werte \(N_0 = 1000\), \(r = 0.5\) und \(t = 3\) Stunden ergibt sich: \[ N(3) = 1000 \cdot e^{0.5 \cdot 3} = 1000 \cdot e^{1.5} \approx 4482 \]

Nach drei Stunden ist die Population auf etwa 4482 Bakterien angewachsen, also mehr als das Vierfache der Ausgangspopulation.

Visualisierung

Die folgende Grafik zeigt den exponentiellen Verlauf der Populationsgröße über einen Zeitraum von zehn Stunden. Die charakteristische Eigenschaft des exponentiellen Wachstums ist, dass die Steigung der Kurve ständig zunimmt. Je größer die Population wird, desto schneller wächst sie weiter.

import numpy as np
import matplotlib.pyplot as plt

# Parameter
N0 = 1000
r = 0.5
t = np.linspace(0, 10, 100)

# Lösung
N = N0 * np.exp(r * t)

# Plot
plt.figure(figsize=(7, 4))
plt.plot(t, N, 'b-', linewidth=2, label=f'$N(t) = {N0}e^{{{r}t}}$')
plt.axhline(y=N0, color='r', linestyle='--', alpha=0.5, label=f'Anfangswert $N_0 = {N0}$')
plt.xlabel('Zeit $t$ (Stunden)', fontsize=12)
plt.ylabel('Population $N(t)$', fontsize=12)
plt.title('Exponentielles Bakterienwachstum', fontsize=14)
plt.grid(True, alpha=0.3)
plt.legend(fontsize=11)
plt.show()

Abbildung 1: Exponentielles Wachstum einer Bakterienpopulation

Numerische Lösung mit solve_ivp

Die analytische Lösung ist hier bequem, weil sich die Gleichung separieren und das Integral geschlossen auswerten lässt. Für die meisten Differentialgleichungen aus der Praxis gelingt das nicht, und man greift auf numerische Verfahren zurück. Das Anfangswertproblem aus Gleichung 1 eignet sich gut, um den Ablauf einer solchen Rechnung mit scipy.integrate.solve_ivp (vgl. die Einführung in Kapitel 1) zu zeigen, denn die exakte Lösung Gleichung 2 steht zur Kontrolle bereit.

Der Löser benötigt drei Angaben: eine Funktion, die zu gegebenem Zeitpunkt \(t\) und Zustand \(N\) die rechte Seite \(rN\) der Differentialgleichung zurückgibt, das Zeitintervall, über das integriert werden soll, und den Anfangswert \(N_0\). Die Funktion erhält den Zustand als Array y, auch wenn er wie hier nur aus einer einzigen Größe besteht, und gibt die Ableitung ebenfalls als Array zurück. Die Schrittweite wählt solve_ivp selbst so, dass eine vorgegebene Fehlertoleranz eingehalten wird; mit t_eval legen wir lediglich fest, an welchen Zeitpunkten die Lösung ausgegeben wird.

import numpy as np
import matplotlib.pyplot as plt
from scipy.integrate import solve_ivp

# Parameter und Anfangswert
r = 0.5      # Wachstumsrate in 1/h
N0 = 1000    # Population zum Zeitpunkt t = 0

# Rechte Seite der Differentialgleichung dN/dt = r N.
# solve_ivp übergibt die Zeit t und den Zustand als Array y = [N].
def growth(t, y):
    N = y[0]
    return [r * N]

# Zeitintervall und Zeitpunkte, an denen die Lösung ausgegeben wird
t_span = (0, 10)
t_eval = np.linspace(0, 10, 101)

# Numerische Integration des Anfangswertproblems
sol = solve_ivp(growth, t_span, [N0], t_eval=t_eval, dense_output=True)

# Analytische Lösung zum Vergleich
N_exact = N0 * np.exp(r * t_eval)

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(7, 4))

ax1.plot(t_eval, N_exact, 'b-', linewidth=2, label='Analytisch $N_0 e^{rt}$')
ax1.plot(sol.t[::4], sol.y[0][::4], 'o', color='red', markersize=5, alpha=0.8,
         label='Numerisch (solve_ivp)')
ax1.set_xlabel('Zeit $t$ (Stunden)', fontsize=12)
ax1.set_ylabel('Population $N(t)$', fontsize=12)
ax1.set_title('(a) Lösung', fontsize=13)
ax1.grid(True, alpha=0.3)
ax1.legend(fontsize=10)

rel_error = np.abs(sol.y[0] - N_exact) / N_exact
ax2.semilogy(t_eval, rel_error, 'o-', color='purple', markersize=3, linewidth=1,
             label='Relativer Fehler')
ax2.axhline(1e-3, color='gray', linestyle='--', label='Standardtoleranz rtol')
ax2.set_xlabel('Zeit $t$ (Stunden)', fontsize=12)
ax2.set_ylabel('Relativer Fehler $|N_\\mathrm{num} - N| / N$', fontsize=12)
ax2.set_title('(b) Fehler', fontsize=13)
ax2.set_ylim(1e-7, 3e-3)
ax2.grid(True, alpha=0.3)
ax2.legend(fontsize=10)

plt.tight_layout()
plt.show()

Abbildung 2: Numerische Lösung des Bakterienwachstums mit solve_ivp. (a) Die numerische Lösung (Punkte) liegt auf der analytischen Lösung (Linie). (b) Der relative Fehler der numerischen Lösung bleibt unterhalb der Standardtoleranz \(\mathrm{rtol} = 10^{-3}\) des Lösers.

Das Ergebnis sol enthält im Feld sol.t die Zeitpunkte und im Feld sol.y die zugehörigen Werte der Lösung; sol.y[0] ist die erste (und hier einzige) Komponente des Zustands. Weil wir dense_output=True gesetzt haben, lässt sich die Lösung über sol.sol zusätzlich an beliebigen Zwischenzeitpunkten auswerten. Damit beantworten wir die Frage nach der Population nach drei Stunden:

N3_num = sol.sol(3.0)[0]
N3_exact = N0 * np.exp(r * 3.0)
print(f"N(3) numerisch:  {N3_num:.1f}")
print(f"N(3) analytisch: {N3_exact:.1f}")

N(3) numerisch:  4481.2
N(3) analytisch: 4481.7

Die numerische Lösung weicht um weniger als ein Bakterium von der analytischen Lösung ab. Diese Abweichung ist der Diskretisierungsfehler des Verfahrens. Er ist kleiner als die relative Standardtoleranz \(\mathrm{rtol} = 10^{-3}\), mit der solve_ivp seine Schrittweite steuert, und lässt sich durch die Argumente rtol und atol weiter verkleinern. Abbildung 2 (b) zeigt, dass der relative Fehler über das gesamte Intervall in dieser Größenordnung bleibt, obwohl die Lösung selbst um mehr als zwei Größenordnungen wächst. Die Fehlerkontrolle adaptiver Verfahren ist Gegenstand von Kapitel 13.

Beispiel 2: Radioaktiver Zerfall

Der radioaktive Zerfall ist das natürliche Gegenstück zum exponentiellen Wachstum. Bei radioaktiven Substanzen nimmt die Anzahl der noch nicht zerfallenen Atomkerne mit der Zeit ab. Die Zerfallsrate ist dabei proportional zur Anzahl der vorhandenen Kerne, was zu einem exponentiellen Abklingen führt.

Problemstellung

Kohlenstoff-14 ist ein radioaktives Isotop, das in der Archäologie zur Altersbestimmung verwendet wird. Es zerfällt mit einer Halbwertszeit von \(T_{1/2} = 5730\) Jahren. Wir möchten bestimmen, wie viel von einer ursprünglich vorhandenen Menge nach 10000 Jahren noch übrig ist.

Mathematische Modellierung

Die Annahme, dass die Zerfallsrate proportional zur vorhandenen Menge ist, führt auf die Differentialgleichung: \[ \frac{\mathrm{d}N}{\mathrm{d}t} = -\lambda N \tag{3}\]

Das Minuszeichen zeigt an, dass die Menge abnimmt. Die Zerfallskonstante \(\lambda > 0\) charakterisiert die Geschwindigkeit des Zerfalls. Je größer \(\lambda\), desto schneller zerfällt die Substanz.

Diese Differentialgleichung ist strukturell identisch mit der Wachstumsgleichung, nur dass das Vorzeichen negativ ist. Die Lösung erhält man analog durch Separation der Variablen: \[ N(t) = N_0 e^{-\lambda t} \tag{4}\]

Die Halbwertszeit \(T_{1/2}\) ist definiert als die Zeit, nach der die Hälfte der ursprünglichen Menge zerfallen ist. Mathematisch bedeutet dies \(N(T_{1/2}) = N_0/2\). Einsetzen in die Lösung ergibt: \[ \frac{N_0}{2} = N_0 e^{-\lambda T_{1/2}} \]

Division durch \(N_0\) und Logarithmieren führt auf: \[ \ln\left(\frac{1}{2}\right) = -\lambda T_{1/2} \]

Daraus folgt der Zusammenhang zwischen Halbwertszeit und Zerfallskonstante: \[ \lambda = \frac{\ln 2}{T_{1/2}} \tag{5}\]

Für Kohlenstoff-14 mit \(T_{1/2} = 5730\) Jahren ergibt sich: \[ \lambda = \frac{0.693}{5730} \approx 1.21 \times 10^{-4} \text{ Jahr}^{-1} \]

Nach 10000 Jahren beträgt die verbleibende Menge: \[ N(10000) = N_0 e^{-\lambda \cdot 10000} = N_0 e^{-1.21} \approx 0.298 N_0 \]

Etwa 30 Prozent der ursprünglichen Menge sind noch vorhanden, während 70 Prozent bereits zerfallen sind.

Visualisierung

Die Grafik zeigt den charakteristischen exponentiellen Abfall der radioaktiven Substanz. Die vertikalen Linien markieren die Halbwertszeiten. Nach jeder Halbwertszeit ist genau die Hälfte der zu Beginn dieser Periode vorhandenen Menge zerfallen. Dieses Halbierungsmuster wiederholt sich unbegrenzt.

import numpy as np
import matplotlib.pyplot as plt

# Parameter
T_half = 5730  # Halbwertszeit in Jahren
lambda_decay = np.log(2) / T_half
t = np.linspace(0, 20000, 1000)

# Lösung (normalisiert auf N0 = 1)
N = np.exp(-lambda_decay * t)

# Plot
plt.figure(figsize=(7, 4))
plt.plot(t, N, 'b-', linewidth=2, label='$N(t)/N_0 = e^{-\\lambda t}$')
plt.axhline(y=0.5, color='r', linestyle='--', alpha=0.5, label=f'Halbwertszeit $T_{{1/2}} = {T_half}$ Jahre')
plt.axvline(x=T_half, color='r', linestyle='--', alpha=0.5)
plt.axvline(x=2*T_half, color='r', linestyle='--', alpha=0.3)
plt.axvline(x=3*T_half, color='r', linestyle='--', alpha=0.3)
plt.xlabel('Zeit $t$ (Jahre)', fontsize=12)
plt.ylabel('Verbleibendes Material $N(t)/N_0$', fontsize=12)
plt.title('Radioaktiver Zerfall: Kohlenstoff-14', fontsize=14)
plt.grid(True, alpha=0.3)
plt.legend(fontsize=11)
plt.show()

Abbildung 3: Radioaktiver Zerfall von Kohlenstoff-14

Beispiel 3: Freier Fall mit Luftwiderstand

Ein fallender Körper in der Erdatmosphäre unterliegt nicht nur der Schwerkraft, sondern auch dem bremsenden Einfluss des Luftwiderstands. Für moderate Geschwindigkeiten kann der Luftwiderstand als proportional zur Geschwindigkeit angenommen werden. Dies führt zu einer linearen Differentialgleichung erster Ordnung.

Problemstellung

Ein Objekt der Masse \(m\) fällt aus der Ruhe unter dem Einfluss der Schwerkraft. Der Luftwiderstand erzeugt eine der Bewegung entgegengesetzte Kraft, die proportional zur Geschwindigkeit ist. Wir möchten die Geschwindigkeit als Funktion der Zeit bestimmen und insbesondere die Endgeschwindigkeit untersuchen, die der Körper asymptotisch erreicht.

Mathematische Modellierung

Die auf den Körper wirkenden Kräfte sind die Gewichtskraft \(F_\mathrm{g} = mg\) in Fallrichtung und die Luftwiderstandskraft \(F_\mathrm{d} = -kv\) entgegen der Fallrichtung. Dabei bezeichnet \(k\) den Luftwiderstandskoeffizienten und \(v\) die Geschwindigkeit. Das Minuszeichen zeigt an, dass die Kraft der Bewegung entgegenwirkt. Das Kräftediagramm in Abbildung 4 veranschaulicht diese Situation.

import matplotlib.pyplot as plt
from matplotlib.patches import Circle

fig, ax = plt.subplots(figsize=(7, 4))

# Kugel
ball = Circle((0, 0), 0.3, facecolor='lightsteelblue', edgecolor='black',
              linewidth=1.5, zorder=3)
ax.add_patch(ball)
ax.text(-0.55, 0.35, '$m$', fontsize=13)

# Gewichtskraft mg nach unten
ax.annotate('', xy=(0, -1.5), xytext=(0, -0.3),
            arrowprops=dict(arrowstyle='-|>', color='red', linewidth=2,
                            mutation_scale=20))
ax.text(0.12, -1.05, '$mg$', fontsize=14, color='red')

# Luftwiderstand kv nach oben
ax.annotate('', xy=(0, 1.3), xytext=(0, 0.3),
            arrowprops=dict(arrowstyle='-|>', color='blue', linewidth=2,
                            mutation_scale=20))
ax.text(0.12, 0.85, '$kv$', fontsize=14, color='blue')

# Bewegungsrichtung
ax.annotate('', xy=(1.1, -0.7), xytext=(1.1, 0.1),
            arrowprops=dict(arrowstyle='-|>', color='gray', linewidth=1.5,
                            mutation_scale=15))
ax.text(1.25, -0.4, '$v$', fontsize=12, color='gray')

ax.set_xlim(-2.5, 2.5)
ax.set_ylim(-2.0, 1.8)
ax.set_aspect('equal')
ax.axis('off')
plt.show()

Abbildung 4: Kräftediagramm des fallenden Körpers: Die Gewichtskraft \(mg\) wirkt nach unten, der Luftwiderstand \(kv\) entgegen der Bewegungsrichtung nach oben

Nach dem zweiten Newtonschen Gesetz gilt: \[ m\frac{\mathrm{d}v}{\mathrm{d}t} = mg - kv \tag{6}\]

Division durch \(m\) führt auf die Standardform: \[ \frac{\mathrm{d}v}{\mathrm{d}t} = g - \frac{k}{m}v \tag{7}\]

Dies ist eine inhomogene lineare Differentialgleichung erster Ordnung. Die Umformung in die Standardform \(v' + p(t)v = q(t)\) ergibt: \[ v' + \frac{k}{m}v = g \]

Hier ist \(p(t) = k/m\) konstant und \(q(t) = g\) ebenfalls konstant.

Lösung mit integrierendem Faktor

Der integrierende Faktor ist: \[ \mu(t) = e^{\int \frac{k}{m} \, \mathrm{d}t} = e^{\frac{k}{m}t} \]

Multiplikation der Differentialgleichung mit \(\mu(t)\) führt auf: \[ e^{\frac{k}{m}t} \frac{\mathrm{d}v}{\mathrm{d}t} + e^{\frac{k}{m}t} \frac{k}{m}v = ge^{\frac{k}{m}t} \]

Die linke Seite ist die Ableitung von \(e^{\frac{k}{m}t} v\): \[ \frac{\mathrm{d}}{\mathrm{d}t}\left[e^{\frac{k}{m}t} v\right] = ge^{\frac{k}{m}t} \]

Integration liefert: \[ e^{\frac{k}{m}t} v = g \cdot \frac{m}{k} e^{\frac{k}{m}t} + C \]

Auflösen nach \(v\) ergibt: \[ v(t) = \frac{mg}{k} + Ce^{-\frac{k}{m}t} \tag{8}\]

Die Anfangsbedingung \(v(0) = 0\) (Körper fällt aus der Ruhe) bestimmt die Konstante: \[ 0 = \frac{mg}{k} + C \quad \Rightarrow \quad C = -\frac{mg}{k} \]

Die spezielle Lösung lautet: \[ v(t) = \frac{mg}{k}\left(1 - e^{-\frac{k}{m}t}\right) \tag{9}\]

Physikalische Interpretation

Für große Zeiten strebt der Exponentialterm \(e^{-\frac{k}{m}t}\) gegen null. Die Geschwindigkeit nähert sich asymptotisch dem Grenzwert: \[ v_{\infty} = \lim_{t \to \infty} v(t) = \frac{mg}{k} \tag{10}\]

Dieser Grenzwert heißt Endgeschwindigkeit oder Terminalgeschwindigkeit. Bei dieser Geschwindigkeit kompensiert der Luftwiderstand genau die Gewichtskraft, sodass die Nettokraft verschwindet und keine weitere Beschleunigung mehr stattfindet.

Die Endgeschwindigkeit ist umso größer, je schwerer der Körper und je kleiner der Luftwiderstandskoeffizient ist. Ein schwerer, stromlinienförmiger Körper erreicht eine höhere Endgeschwindigkeit als ein leichter Körper mit großer Querschnittsfläche.

Visualisierung

Die Grafik zeigt die Geschwindigkeit als Funktion der Zeit für verschiedene Luftwiderstandskoeffizienten. Je größer der Luftwiderstand, desto kleiner ist die Endgeschwindigkeit und desto schneller wird sie erreicht. Je kleiner der Luftwiderstandskoeffizient, desto länger folgt die Geschwindigkeit anfangs dem linearen Anstieg \(v \approx gt\) des freien Falls ohne Luftwiderstand.

import numpy as np
import matplotlib.pyplot as plt

# Parameter
g = 9.81  # m/s^2
m = 1.0   # kg
k_values = [0.1, 0.5, 1.0, 2.0]  # verschiedene Luftwiderstände

t = np.linspace(0, 10, 1000)

plt.figure(figsize=(7, 4))
for k in k_values:
    v_terminal = m * g / k
    v = v_terminal * (1 - np.exp(-k * t / m))
    plt.plot(t, v, linewidth=2, label=f'$k = {k}$ kg/s (${v_terminal:.1f}$ m/s)')

plt.xlabel('Zeit $t$ (s)', fontsize=12)
plt.ylabel('Geschwindigkeit $v(t)$ (m/s)', fontsize=12)
plt.title('Freier Fall mit Luftwiderstand', fontsize=14)
plt.grid(True, alpha=0.3)
plt.legend(fontsize=10, title='Luftwiderstand (Endgeschwindigkeit)')
plt.show()

Abbildung 5: Fallende Objekte mit verschiedenen Luftwiderstandskoeffizienten

Beispiel 4: RC-Schaltkreis

Der RC-Schaltkreis ist ein grundlegendes Beispiel aus der Elektrotechnik, bei dem die Spannung über einem Kondensator einer Differentialgleichung erster Ordnung genügt. Das Verhalten dieses einfachen Schaltkreises ist prototypisch für viele elektrische und elektronische Systeme.

Problemstellung

Ein RC-Schaltkreis besteht aus einem Widerstand \(R\) und einem Kondensator \(C\), die in Reihe mit einer konstanten Spannungsquelle \(V_0\) geschaltet sind. Zum Zeitpunkt \(t = 0\) ist der Kondensator ungeladen. Wir möchten die Spannung über dem Kondensator als Funktion der Zeit bestimmen. Den Aufbau des Schaltkreises skizziert Abbildung 6.

import numpy as np
import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(7, 4))

# Untere Leitung
ax.plot([0, 6], [0, 0], 'k', linewidth=1.5)

# Linke Leitung mit Spannungsquelle (Batteriesymbol)
ax.plot([0, 0], [0, 1.2], 'k', linewidth=1.5)
ax.plot([0, 0], [1.8, 3], 'k', linewidth=1.5)
ax.plot([-0.45, 0.45], [1.8, 1.8], 'k', linewidth=2)   # langer Strich (+)
ax.plot([-0.22, 0.22], [1.2, 1.2], 'k', linewidth=3.5)  # kurzer Strich (−)
ax.text(-1.2, 1.4, '$V_0$', fontsize=14)
ax.text(0.6, 1.85, '$+$', fontsize=11)
ax.text(0.6, 0.95, '$-$', fontsize=11)

# Obere Leitung mit Widerstand (Zickzack-Symbol)
ax.plot([0, 2], [3, 3], 'k', linewidth=1.5)
x_zz = np.linspace(2, 4, 9)
y_zz = 3 + 0.25 * np.array([0, 1, -1, 1, -1, 1, -1, 1, 0])
ax.plot(x_zz, y_zz, 'k', linewidth=1.5)
ax.plot([4, 6], [3, 3], 'k', linewidth=1.5)
ax.text(3, 3.55, '$R$', fontsize=14, ha='center')

# Rechte Leitung mit Kondensator (zwei parallele Platten)
ax.plot([6, 6], [3, 1.65], 'k', linewidth=1.5)
ax.plot([6, 6], [1.35, 0], 'k', linewidth=1.5)
ax.plot([5.5, 6.5], [1.65, 1.65], 'k', linewidth=2.5)
ax.plot([5.5, 6.5], [1.35, 1.35], 'k', linewidth=2.5)
ax.text(6.8, 1.4, '$C$', fontsize=14)

# Strompfeil I auf der oberen Leitung
ax.annotate('', xy=(1.7, 3.3), xytext=(0.6, 3.3),
            arrowprops=dict(arrowstyle='-|>', color='red', linewidth=1.5,
                            mutation_scale=15))
ax.text(1.1, 3.5, '$I$', fontsize=13, color='red')

ax.set_xlim(-2, 8)
ax.set_ylim(-0.8, 4.3)
ax.set_aspect('equal')
ax.axis('off')
plt.show()

Abbildung 6: Schaltplan des RC-Schaltkreises: Spannungsquelle \(V_0\), Widerstand \(R\) und Kondensator \(C\) in Serie; der Strom \(I\) lädt den Kondensator auf

Mathematische Modellierung

Nach dem Kirchhoffschen Spannungsgesetz entspricht die Quellenspannung der Summe der Spannungen über Widerstand und Kondensator: \[ V_0 = V_\mathrm{R} + V_\mathrm{C} \tag{11}\]

Die Spannung über dem Widerstand ist nach dem Ohmschen Gesetz \(V_\mathrm{R} = RI\), wobei \(I\) der Strom ist. Die Spannung über dem Kondensator ist \(V_\mathrm{C} = Q/C\), wobei \(Q\) die Ladung auf dem Kondensator bezeichnet. Der Strom ist die zeitliche Änderungsrate der Ladung: \(I = \mathrm{d}Q/\mathrm{d}t\).

Mit \(Q = CV_\mathrm{C}\) folgt \(I = C \, \mathrm{d}V_\mathrm{C}/\mathrm{d}t\). Einsetzen in das Kirchhoffsche Gesetz ergibt: \[ V_0 = RC\frac{\mathrm{d}V_\mathrm{C}}{\mathrm{d}t} + V_\mathrm{C} \tag{12}\]

Umformung in Standardform liefert: \[ \frac{\mathrm{d}V_\mathrm{C}}{\mathrm{d}t} + \frac{1}{RC}V_\mathrm{C} = \frac{V_0}{RC} \tag{13}\]

Dies ist eine inhomogene lineare Differentialgleichung erster Ordnung mit konstanten Koeffizienten.

Lösung durch Variation der Konstanten

Die Beispiele 3 und 5 lösen ihre Gleichungen mit dem integrierenden Faktor. Hier führen wir stattdessen die Variation der Konstanten vor. Der erste Schritt ist die Lösung der homogenen Gleichung \(\dot{V}_\mathrm{C} + V_\mathrm{C}/(RC) = 0\). Trennung der Variablen liefert: \[ V_{\mathrm{C},\mathrm{h}}(t) = K e^{-\frac{t}{RC}} \tag{14}\]

Im zweiten Schritt ersetzen wir die Konstante \(K\) durch eine Funktion \(K(t)\) und setzen den Ansatz \(V_\mathrm{C}(t) = K(t)\,e^{-t/(RC)}\) in Gleichung 13 ein. Die Produktregel ergibt: \[ \dot{K} e^{-\frac{t}{RC}} - \frac{K}{RC} e^{-\frac{t}{RC}} + \frac{K}{RC} e^{-\frac{t}{RC}} = \frac{V_0}{RC} \]

Die beiden Terme mit \(K\) heben sich auf, wie es der Ansatz vorsieht. Übrig bleibt eine Gleichung, die \(\dot{K}\) direkt angibt und sich durch Integration lösen lässt: \[ \dot{K} = \frac{V_0}{RC} e^{\frac{t}{RC}} \quad \Rightarrow \quad K(t) = V_0 e^{\frac{t}{RC}} + K_0 \]

Im dritten Schritt setzen wir \(K(t)\) in den Ansatz ein: \[ V_\mathrm{C}(t) = V_0 + K_0 e^{-\frac{t}{RC}} \tag{15}\]

Die Lösung hat die im Konzepte-Teil beschriebene Struktur: \(V_0\) ist eine partikuläre Lösung, nämlich die konstante Kondensatorspannung im stationären Zustand, und \(K_0 e^{-t/(RC)}\) ist die homogene Lösung Gleichung 14, die den Übergang in diesen Zustand beschreibt. Mit dem integrierenden Faktor \(\mu(t) = e^{t/(RC)}\) erhält man dieselbe Lösung in derselben Anzahl von Schritten.

Die Anfangsbedingung \(V_\mathrm{C}(0) = 0\) (ungeladener Kondensator) bestimmt die Konstante: \[ 0 = V_0 + K_0 \quad \Rightarrow \quad K_0 = -V_0 \]

Die spezielle Lösung lautet: \[ V_\mathrm{C}(t) = V_0\left(1 - e^{-\frac{t}{RC}}\right) \tag{16}\]

Physikalische Interpretation

Die Kondensatorspannung steigt exponentiell von null auf den Endwert \(V_0\). Die Geschwindigkeit dieses Anstiegs wird durch die Zeitkonstante \(\tau = RC\) charakterisiert. Nach einer Zeitkonstante hat die Spannung etwa 63 Prozent ihres Endwerts erreicht: \[ V_\mathrm{C}(\tau) = V_0(1 - e^{-1}) \approx 0.632 V_0 \]

Nach etwa fünf Zeitkonstanten ist der Kondensator praktisch vollständig aufgeladen. Ein großer Widerstand oder eine große Kapazität führen zu einer langen Zeitkonstante und damit zu einem langsamen Ladevorgang. Diese RC-Zeitkonstante spielt in der Elektronik eine fundamentale Rolle für Zeitverzögerungen, Filter und viele andere Anwendungen.

Visualisierung

Die Grafik zeigt den charakteristischen exponentiellen Verlauf der Kondensatorspannung. Die Tangente am Ursprung würde die Endspannung \(V_0\) genau nach einer Zeitkonstante \(\tau\) erreichen. Tatsächlich erreicht die Kurve die Endspannung aber nur asymptotisch.

import numpy as np
import matplotlib.pyplot as plt

# Parameter
V0 = 10  # Volt
RC = 2   # Sekunden (Zeitkonstante)

t = np.linspace(0, 10, 1000)

# Lösung
V_C = V0 * (1 - np.exp(-t / RC))

plt.figure(figsize=(7, 4))
plt.plot(t, V_C, 'b-', linewidth=2, label='$V_\mathrm{C}(t) = V_0(1 - e^{-t/(RC)})$')
plt.axhline(y=V0, color='r', linestyle='--', alpha=0.5, label=f'$V_0 = {V0}$ V')
plt.axvline(x=RC, color='g', linestyle='--', alpha=0.5, label=f'$\\tau = RC = {RC}$ s')
plt.axhline(y=V0*0.632, color='g', linestyle=':', alpha=0.5)
plt.xlabel('Zeit $t$ (s)', fontsize=12)
plt.ylabel('Spannung $V_\mathrm{C}(t)$ (V)', fontsize=12)
plt.title('Kondensator-Aufladung im RC-Schaltkreis', fontsize=14)
plt.grid(True, alpha=0.3)
plt.legend(fontsize=11)
plt.show()

Abbildung 7: Aufladen eines Kondensators im RC-Schaltkreis

Beispiel 5: Mischungsproblem mit veränderlichem Volumen

Mischungsprobleme treten in vielen praktischen Situationen auf, etwa in der chemischen Verfahrenstechnik, der Umwelttechnik oder der Pharmazie. Ein typisches Problem ist ein Behälter mit einer Lösung, in den kontinuierlich frische Lösung einfließt, während gleichzeitig gemischte Lösung abfließt. In den Beispielen 3 und 4 war der Koeffizient \(p(t)\) der linearen Differentialgleichung eine Konstante, und der integrierende Faktor war deshalb stets eine Exponentialfunktion. Fließt mehr Lösung zu als ab, ändert sich das Volumen im Tank mit der Zeit, und \(p(t)\) wird zeitabhängig. Dieses Beispiel zeigt, dass die Methode des integrierenden Faktors auch diesen Fall ohne zusätzlichen Aufwand erfasst.

Problemstellung

Ein Tank mit einem Fassungsvermögen von 200 Litern enthält zu Beginn 100 Liter Wasser, in dem 20 kg Salz gelöst sind. Eine Salzlösung mit einer Konzentration von 0.5 kg pro Liter fließt mit einer Rate von 4 Litern pro Minute in den Tank. Die Lösung im Tank wird kontinuierlich gut durchmischt. Gleichzeitig fließt die Mischung mit einer Rate von nur 2 Litern pro Minute aus dem Tank ab. Wir möchten die Salzmenge im Tank nach 30 Minuten bestimmen und wissen, wie viel Salz der Tank in dem Moment enthält, in dem er überläuft. Abbildung 8 zeigt den Aufbau dieses Mischungsproblems schematisch.

import matplotlib.pyplot as plt
from matplotlib.patches import Rectangle

fig, ax = plt.subplots(figsize=(7, 4))

# Tank (offener Behälter)
ax.plot([2, 2, 6, 6], [4.4, 1, 1, 4.4], color='black', linewidth=2)

# Flüssigkeit im Tank
water = Rectangle((2.05, 1.05), 3.9, 2.3, facecolor='lightblue',
                  edgecolor='none', zorder=0)
ax.add_patch(water)
ax.plot([2.05, 5.95], [3.35, 3.35], color='steelblue', linewidth=1.5)
ax.annotate('', xy=(5.6, 3.9), xytext=(5.6, 3.4),
            arrowprops=dict(arrowstyle='-|>', color='steelblue', linewidth=1.5,
                            mutation_scale=14))
ax.text(4, 2.2, '$V(t) = 100\\,\\mathrm{L} + 2\\,\\frac{\\mathrm{L}}{\\mathrm{min}} \\cdot t$\ngut durchmischt,\n$S(t)$ kg Salz',
        fontsize=11, ha='center', va='center')

# Fassungsvermögen
ax.plot([2.05, 5.95], [4.2, 4.2], color='gray', linestyle='--', linewidth=1)
ax.text(6.1, 4.2, '200 L', fontsize=10, color='gray', va='center')

# Zufluss-Pfeil
ax.annotate('', xy=(2.7, 4.1), xytext=(1.0, 5.5),
            arrowprops=dict(arrowstyle='-|>', color='blue', linewidth=2,
                            mutation_scale=20))
ax.text(0.2, 5.7, 'Zufluss: $4$ L/min\nKonzentration $0{,}5$ kg/L',
        fontsize=11, color='blue')

# Abfluss-Pfeil
ax.annotate('', xy=(8.0, 0.3), xytext=(6.1, 1.2),
            arrowprops=dict(arrowstyle='-|>', color='blue', linewidth=2,
                            mutation_scale=20))
ax.text(6.4, 1.6, 'Abfluss: $2$ L/min\nKonzentration $S(t)/V(t)$',
        fontsize=11, color='blue')

ax.set_xlim(-0.5, 10.5)
ax.set_ylim(-0.3, 6.7)
ax.set_aspect('equal')
ax.axis('off')
plt.show()

Abbildung 8: Schema des Mischungsproblems: Zufluss mit fester Konzentration, gut durchmischter Tank mit wachsendem Volumen und Abfluss mit der aktuellen Konzentration

Mathematische Modellierung

Weil pro Minute 4 Liter zufließen, aber nur 2 Liter abfließen, nimmt das Volumen im Tank um 2 Liter pro Minute zu. Mit \(t\) in Minuten und \(V\) in Litern gilt: \[ V(t) = 100 + 2t \tag{17}\]

Der Tank ist voll, sobald \(V(t) = 200\) L erreicht ist, also nach 50 Minuten. Das folgende Modell gilt bis zu diesem Zeitpunkt.

Sei \(S(t)\) die Salzmenge im Tank zur Zeit \(t\) in Kilogramm. Die Änderungsrate von \(S\) setzt sich aus zwei Beiträgen zusammen: \[ \frac{\mathrm{d}S}{\mathrm{d}t} = \text{(Salzeinstrom)} - \text{(Salzausstrom)} \tag{18}\]

Der Salzeinstrom ist wie zuvor konstant: \[ \text{Einstrom} = 0.5 \, \frac{\text{kg}}{\text{L}} \times 4 \, \frac{\text{L}}{\text{min}} = 2 \, \frac{\text{kg}}{\text{min}} \]

Der Salzausstrom hängt von der aktuellen Konzentration \(S(t)/V(t)\) im Tank ab. Da das Volumen wächst, ist diese Konzentration nicht mehr proportional zu \(S\) mit festem Faktor, sondern der Faktor ändert sich mit der Zeit: \[ \text{Ausstrom} = \frac{S(t)}{100 + 2t} \, \frac{\text{kg}}{\text{L}} \times 2 \, \frac{\text{L}}{\text{min}} = \frac{S(t)}{50 + t} \, \frac{\text{kg}}{\text{min}} \]

Die Differentialgleichung lautet: \[ \frac{\mathrm{d}S}{\mathrm{d}t} = 2 - \frac{S}{50 + t} \tag{19}\]

mit der Anfangsbedingung \(S(0) = 20\) kg. In der Standardform \(S' + p(t) S = q(t)\) ist \[ S' + \frac{1}{50 + t} S = 2 . \] Der Koeffizient \(p(t) = 1/(50 + t)\) hängt von der Zeit ab, die Inhomogenität \(q(t) = 2\) ist konstant. Die Gleichung ist nicht separabel, denn die rechte Seite von Gleichung 19 lässt sich nicht als Produkt einer Funktion von \(t\) und einer Funktion von \(S\) schreiben.

Lösung mit integrierendem Faktor

Der integrierende Faktor \(\mu = e^{\int p \, \mathrm{d}t}\) aus dem Konzepte-Teil erfordert jetzt eine echte Integration: \[ \mu(t) = e^{\int \frac{1}{50 + t} \, \mathrm{d}t} = e^{\ln(50 + t)} = 50 + t \]

Das Betragszeichen aus \(\int \mathrm{d}t/(50+t) = \ln|50 + t|\) entfällt, weil \(50 + t > 0\) für alle \(t \geq 0\) gilt. Der integrierende Faktor ist keine Exponentialfunktion, sondern eine lineare Funktion der Zeit; bis auf den Faktor \(1/2\) ist er das Volumen \(V(t)\).

Multiplikation der Differentialgleichung mit \(\mu(t)\) führt auf: \[ (50 + t) \frac{\mathrm{d}S}{\mathrm{d}t} + S = 2(50 + t) \]

Die linke Seite ist nach der Produktregel die Ableitung von \((50 + t) S\): \[ \frac{\mathrm{d}}{\mathrm{d}t}\left[(50 + t) S\right] = 2(50 + t) \]

Integration beider Seiten ergibt: \[ (50 + t) S = (50 + t)^2 + C \]

Division durch \(50 + t\) liefert die allgemeine Lösung: \[ S(t) = 50 + t + \frac{C}{50 + t} \tag{20}\]

Die Anfangsbedingung \(S(0) = 20\) bestimmt die Konstante: \[ 20 = 50 + \frac{C}{50} \quad \Rightarrow \quad C = -1500 \]

Die spezielle Lösung lautet: \[ S(t) = 50 + t - \frac{1500}{50 + t} \tag{21}\]

Nach 30 Minuten ergibt sich: \[ S(30) = 80 - \frac{1500}{80} = 80 - 18.75 \approx 61 \text{ kg} \]

Beim Überlaufen nach 50 Minuten enthält der Tank: \[ S(50) = 100 - \frac{1500}{100} = 85 \text{ kg} \]

Die Konzentration im Tank folgt aus \(c(t) = S(t)/V(t)\) mit \(V(t) = 2(50 + t)\): \[ c(t) = \frac{1}{2} - \frac{750}{(50 + t)^2} \tag{22}\]

Für \(t = 0\) liefert diese Formel die Anfangskonzentration \(0.2\) kg/L, beim Überlaufen beträgt die Konzentration \(0.425\) kg/L.

Physikalische Interpretation

Anders als in den Beispielen 3 und 4 strebt die Salzmenge keinem Gleichgewichtswert zu. Für große Zeiten wächst sie linear wie \(50 + t\), denn der Tank füllt sich mit 2 Litern pro Minute, und jeder zusätzliche Liter trägt im Grenzfall 0.5 kg Salz bei. Die Annäherung an einen Grenzwert zeigt sich erst in der Konzentration Gleichung 22: sie strebt gegen die Zulaufkonzentration von 0.5 kg/L, aber die Abweichung \(750/(50 + t)^2\) klingt algebraisch ab, nicht exponentiell. Bei gleichen Zu- und Abflussraten von 4 Litern pro Minute und konstantem Volumen von 100 Litern wäre die Abweichung dagegen nach \(S(t) = 50 - 30 e^{-t/25}\) exponentiell mit der charakteristischen Zeit \(\tau = 25\) min abgeklungen. Nach 50 Minuten wäre sie auf \(e^{-2} \approx 14\,\%\) ihres Anfangswerts gesunken, im vorliegenden Fall sind es noch \(25\,\%\). Der geringere Abfluss spült das salzarme Ausgangswasser langsamer aus dem Tank.

Visualisierung

Abbildung 9 zeigt links die Salzmenge, die sich der Asymptote \(50 + t\) von unten nähert, und rechts die Konzentration im Vergleich zum Fall gleicher Raten. Beide Kurven enden beim Überlaufen des Tanks nach 50 Minuten.

import numpy as np
import matplotlib.pyplot as plt

# Parameter
S0 = 20         # Anfangssalzmenge (kg)
c_in = 0.5      # Zulaufkonzentration (kg/L)
t_full = 50     # Zeitpunkt des Überlaufens (min)
t = np.linspace(0, t_full, 500)

# Lösung bei veränderlichem Volumen
S = 50 + t - 1500 / (50 + t)
V = 100 + 2 * t
c = S / V

# Vergleich: gleiche Raten, konstantes Volumen
c_const = c_in - 0.3 * np.exp(-t / 25)

fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(7, 4))

ax1.plot(t, S, 'b-', linewidth=2, label='$S(t) = 50 + t - 1500/(50+t)$')
ax1.plot(t, 50 + t, 'r--', alpha=0.6, label='Asymptote $50 + t$')
ax1.axvline(t_full, color='gray', linestyle=':', label='Tank voll')
ax1.set_xlabel('Zeit $t$ (min)', fontsize=12)
ax1.set_ylabel('Salzmenge $S(t)$ (kg)', fontsize=12)
ax1.set_title('(a) Salzmenge', fontsize=13)
ax1.set_xlim(0, t_full + 5)
ax1.grid(True, alpha=0.3)
ax1.legend(fontsize=9, loc='lower right')

ax2.plot(t, c, 'b-', linewidth=2, label='$V(t) = 100 + 2t$')
ax2.plot(t, c_const, 'g-.', linewidth=2, label='$V = 100$ L (gleiche Raten)')
ax2.axhline(c_in, color='r', linestyle='--', alpha=0.6, label='Zulaufkonzentration')
ax2.axvline(t_full, color='gray', linestyle=':', label='Tank voll')
ax2.set_xlabel('Zeit $t$ (min)', fontsize=12)
ax2.set_ylabel('Konzentration $c(t)$ (kg/L)', fontsize=12)
ax2.set_title('(b) Konzentration', fontsize=13)
ax2.set_xlim(0, t_full + 5)
ax2.set_ylim(0.15, 0.55)
ax2.grid(True, alpha=0.3)
ax2.legend(fontsize=9, loc='lower right')

plt.tight_layout()
plt.show()

Abbildung 9: Mischungsproblem mit veränderlichem Volumen. (a) Salzmenge im Tank und ihre Asymptote. (b) Konzentration im Tank, verglichen mit dem Fall gleicher Zu- und Abflussraten bei konstantem Volumen von 100 L.