article-magnetic-field-circular/A_z_surf_grad_A_z_vec_top.py at main · shredEngineer/article-magnetic-field-circular · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
import numpy as np
import pyvista as pv

from matplotlib import cm
from matplotlib.colors import ListedColormap

# --- Eigene dunklere Version von 'cool' ---
base_cmap = cm.get_cmap("cool", 256)
dark_cmap_colors = base_cmap(np.linspace(0, 1, 256))
dark_cmap_colors[:, :3] *= 0.6  # <-- RGB abdunkeln (0.0 - 1.0)
dark_cmap = ListedColormap(dark_cmap_colors)

# --- Parameter ---
z_offset = 0.05  # ← Höhe der Pfeile über der Fläche (0.0 = direkt auf Fläche)

# --- Konstanten ---
I = 1
epsilon_0 = 8.854e-12
c = 3e8
mu_0 = 4 * np.pi * 1e-7

# --- Az-Funktion ---
def Az(x, y):
	r = np.sqrt(x**2 + y**2)
	r[r == 0] = 1e-10
	return I / (2 * np.pi * epsilon_0 * c**2) * np.log(1 / r)

# === Fläche (flach!) ===
x_fine = np.linspace(-1, 1, 100)
y_fine = np.linspace(-1, 1, 100)
Xf, Yf = np.meshgrid(x_fine, y_fine)
Zf = np.zeros_like(Xf)  # <-- Fläche bleibt flach
Az_vals = Az(Xf, Yf)

grid = pv.StructuredGrid()
grid.points = np.c_[Xf.ravel(), Yf.ravel(), Zf.ravel()]
grid.dimensions = Xf.shape[0], Xf.shape[1], 1
grid["Az"] = Az_vals.ravel()

# --- Isolinien ---
contours = grid.contour(isosurfaces=25, scalars="Az")

# === Vektoren (in leicht erhobener Ebene, steuerbar mit z_offset) ===
x_coarse = np.linspace(-1, 1, 8)
y_coarse = np.linspace(-1, 1, 8)
Xc, Yc = np.meshgrid(x_coarse, y_coarse)
Zc = np.full_like(Xc, z_offset)  # <-- leicht über der Fläche

R = np.sqrt(Xc**2 + Yc**2)
R[R == 0] = 1e-10
grad_A_mag = -mu_0 * I / (2 * np.pi * R)
grad_Ax = grad_A_mag * Xc / R
grad_Ay = grad_A_mag * Yc / R
grad_Az = np.zeros_like(grad_Ax)

points = np.c_[Xc.ravel(), Yc.ravel(), Zc.ravel()]
vectors = np.c_[grad_Ax.ravel(), grad_Ay.ravel(), grad_Az.ravel()]
magnitudes = np.linalg.norm(vectors, axis=1)
directions = vectors / magnitudes[:, np.newaxis]

pdata = pv.PolyData(points)
pdata["gradA"] = directions
pdata["magnitude"] = np.log10(magnitudes)

factor = 0.25
pdata.points -= 0.5 * directions * factor
glyphs = pdata.glyph(orient="gradA", scale=False, factor=factor)

# --- Draht als Zylinder ---
cylinder = pv.Cylinder(center=(0, 0, 0), direction=(0, 0, 1), radius=0.02, height=2.0)

# === Plotten ===
plotter = pv.Plotter(window_size=(1000, 600))
plotter.set_background("white")
plotter.renderer.SetBackgroundAlpha(0)

plotter.add_mesh(grid, scalars="Az", cmap="cool", show_scalar_bar=False, opacity=1.0)
plotter.add_mesh(contours, color="black", line_width=1)
plotter.add_mesh(glyphs, scalars="magnitude", cmap=dark_cmap, show_scalar_bar=False)
plotter.add_mesh(cylinder, color="black")

plotter.view_vector([0, 0, 1], viewup=[0, 1, 0])
plotter.camera.zoom(1.6)
plotter.show_axes()

# === Anzeige + Screenshot ===
plotter.show(auto_close=False, interactive=False)
plotter.screenshot("A_z_surf_grad_A_z_vec_top.png", transparent_background=True)
plotter.close()