import numpy as np import pandas as pd import matplotlib as mpl import matplotlib.pyplot as plt # Konfiguriere matplotlib um schönere Plots zu erhalten plt.rcParams.update({ "interactive": True, "figure.dpi": 200, "figure.figsize": (10, 6), "font.family": "serif", "mathtext.fontset": "cm", "axes.grid": True, "axes.labelsize": 14, "grid.linewidth": .5, "grid.alpha": .5, "lines.linewidth": 1, "lines.dashed_pattern": (6, 4), "axes.prop_cycle": mpl.cycler(color=['k', 'r', 'b', 'g', 'y']) })

h, delta = np.loadtxt("abweichung.dat", unpack=True) plt.scatter(h, delta, s=5, label="Messwerte") plt.xlabel("$h$") plt.ylabel("$\\Delta_g(h)$") plt.xscale("log") plt.yscale("log") plt.legend();

def fit_linear_region(i1, i2): lin_h = h[i1:i2] lin_delta = delta[i1:i2] m, b = np.polyfit(np.log(lin_h), np.log(lin_delta), 1) print(f"m = {m}, b = {b}") fit_delta = np.exp(b) * lin_h ** m plt.scatter(lin_h, lin_delta, s=20, label="Messwerte") plt.plot(lin_h, fit_delta, c="r", label="Fit der linearen Region") plt.xlabel("$h$") plt.ylabel("$\\Delta_g(h)$") plt.xscale("log") plt.yscale("log") plt.legend() return m, b m_V, b_V = fit_linear_region(110, 150)

m_M, b_M = fit_linear_region(0, 90)

delta_M = np.sqrt(6) * np.exp(b_M) print(f"delta_M = {delta_M}")

def f(x, omega): return np.cos(2 * np.pi * omega * x) def f_deriv(x, omega): return -2 * np.pi * omega * np.sin(2 * np.pi * omega * x)

def g2(x, h, omega): return (f(x + h, omega) - f(x, omega)) / h def g3(x, h, omega): return (f(x + h, omega) - f(x - h, omega)) / (2 * h)

def find_h_max(delta_g: np.ndarray, h: np.ndarray, start_h: float = -np.inf): """Finds the first local maximum of delta_g for h >= start_h.""" # Iteriere über alle relevanten h- und Delta-Werte for i in range(len(h) - 1): # Überspringe Werte vor dem angegebenen h-Startwert. # Dies wird benötigt, um die Maximums-Suche erst hinter dem Maschinenfehler-Bereich zu starten. if h[i] < start_h: continue if delta_g[i+1] < delta_g[i]: # Wenn der nächste Wert kleiner ist, als der aktuelle, hat man das Maximum erreicht. # Der zugehörige h-Wert wird zurückgegeben. return h[i] raise "Could not find local maximum."

def analyze_delta_g( x: float, omega: float, h: np.ndarray, start_h: float = -np.inf, plot: bool = True, ): # Berechne die Approximationsfehler delta_g2 = np.abs(g2(x, h, omega) - f_deriv(x, omega)) delta_g3 = np.abs(g3(x, h, omega) - f_deriv(x, omega)) # Bestimme h_max für beide Approximationen h_max2 = find_h_max(delta_g2, h, start_h) h_max3 = find_h_max(delta_g3, h, start_h) if plot: # Trage Approximationsfehler gegen h auf plt.scatter(h, delta_g2, s=5, label="2-Punkt-Formel", c="r", alpha=.7) plt.scatter(h, delta_g3, s=5, label="3-Punkt-Formel", c="g", alpha=.7) plt.xlabel("$h$") plt.ylabel("$\\Delta_g(h)$") plt.xscale("log") plt.yscale("log") # Trage h_max im Plot ein plt.axvline(h_max2, c="r", label="$h_\\mathrm{max, 2P} = " + f"{h_max2:.4f}$") plt.axvline(h_max3, c="g", label="$h_\\mathrm{max, 3P} = " + f"{h_max3:.4f}$") plt.legend() return h_max2, h_max3 analyze_delta_g(x=.1, omega=1, h=np.logspace(-14, 3, 400), start_h=1e-3)

def plot_for_many_omegas(omega_arr, x_arr): plt.figure(figsize=(10, 20)) for i, (x, omega) in enumerate(zip(x_arr, omega_arr)): # Erstelle Subplot für jeden omega-Wert plt.subplot(len(omega_arr), 1, i + 1) plt.title(f"$\\omega = {omega:.2f}$") analyze_delta_g(x, omega, h=np.logspace(-14, 3, 400), start_h=1e-5) plt.tight_layout() # omega_arr = np.linspace(1, 20, num=6) plot_for_many_omegas(omega_arr, x_arr=.1 / omega_arr)

def find_all_h_max(omega_arr, x_arr): h_max2_list = [] h_max3_list = [] for x, omega in zip(x_arr, omega_arr): h_max2, h_max3 = analyze_delta_g(x, omega, h=np.logspace(-3, 1, 2000), plot=False) h_max2_list.append(h_max2) h_max3_list.append(h_max3) return np.array(h_max2_list), np.array(h_max3_list)

omega_arr = np.linspace(1, 20, num=40) h_max2_arr, h_max3_arr = find_all_h_max(omega_arr, x_arr=.1 / omega_arr) plt.scatter(omega_arr, h_max2_arr, s=20, label="2-Punkt-Formel") plt.scatter(omega_arr, h_max3_arr, s=20, label="3-Punkt-Formel") plt.xlabel("$\\omega$") plt.ylabel("$h_\\mathrm{max}$") plt.legend();

def plot_h_max_against_omega(omega_arr, x_arr): h_max2_arr, h_max3_arr = find_all_h_max(omega_arr, x_arr) m2, b2 = np.polyfit(np.log(omega_arr), np.log(h_max2_arr), 1) m3, b3 = np.polyfit(np.log(omega_arr), np.log(h_max3_arr), 1) print(f"m2 = {m2}, m3 = {m3}") plt.scatter(omega_arr, h_max2_arr, s=20, label="2-Punkt-Formel", alpha=.5) plt.scatter(omega_arr, h_max3_arr, s=20, label="3-Punkt-Formel", alpha=.5) plt.xlabel("$\\omega$") plt.ylabel("$h_\\mathrm{max}$") plt.xscale("log") plt.yscale("log") plt.legend() omega_arr = np.logspace(0, 1, base=20, num=40) plot_h_max_against_omega(omega_arr, x_arr=.1/omega_arr)

omega_arr = np.linspace(1, 20, num=6) plot_for_many_omegas(omega_arr, x_arr=.249 / omega_arr)

omega_arr = np.logspace(0, 1, base=20, num=40) plot_h_max_against_omega(omega_arr, x_arr=.249/omega_arr)