modified: anasen_fem/garfield_sim.py changed mesh input for garfield

This commit is contained in:
Vignesh Sitaraman 2026-04-24 14:34:50 -04:00
parent 67bb9b6fd9
commit e3bfcda35c
2 changed files with 34 additions and 7 deletions

View File

@ -4,7 +4,7 @@ import sys
# 1. FIX: Manually load the Garfield library if it's not in the ROOT namespace # 1. FIX: Manually load the Garfield library if it's not in the ROOT namespace
# Update this path to your actual installation location # Update this path to your actual installation location
garfield_lib_path = "/home/vsitaraman/garfieldpp/install/lib/libGarfield.so" garfield_lib_path = "/home/vs19g/garfieldpp/install/lib/libGarfield.so"
if os.path.exists(garfield_lib_path): if os.path.exists(garfield_lib_path):
ROOT.gSystem.Load(garfield_lib_path) ROOT.gSystem.Load(garfield_lib_path)
@ -44,12 +44,13 @@ else:
# --- 3. FIELD MAP SETUP --- # --- 3. FIELD MAP SETUP ---
fm = ROOT.Garfield.ComponentElmer() fm = ROOT.Garfield.ComponentElmer()
# Update these filenames to match your Elmer SIF output exactly
# Assuming ElmerGrid was run on 'wires2d' directory fm.Initialise("wires2d/mesh.header",
fm.Initialise("wires2d/mesh.nodes",
"wires2d/mesh.elements", "wires2d/mesh.elements",
"wires2d/mesh.boundary", "wires2d/mesh.nodes",
"wires2d/elfield_anasen.result", "mm") "wires2d/dielectrics.dat", # Dielectrics (leave as empty string)
"wires2d/elstatics.result",
"mm")
# Set the medium (Body 13 from your Gmsh script) # Set the medium (Body 13 from your Gmsh script)
fm.SetMedium(0, gas) fm.SetMedium(0, gas)
@ -73,8 +74,33 @@ x0, y0, z0, t0 = 35.0, 0.0, 0.0, 0.0
print(f"Simulating heavy ion drift from r={x0}...") print(f"Simulating heavy ion drift from r={x0}...")
drift.DriftIon(x0, y0, z0, t0) drift.DriftIon(x0, y0, z0, t0)
# Create a file to store the heavy ion track
with open("heavy_ion_track.csv", "w") as f:
f.write("x,y,z,t\n")
# After running drift.DriftIon(x0, y0, z0, t0):
n_points = drift.GetNumberOfDriftLinePoints()
for i in range(n_points):
xi, yi, zi, ti = ROOT.double(0), ROOT.double(0), ROOT.double(0), ROOT.double(0)
drift.GetDriftLinePoint(i, xi, yi, zi, ti)
f.write(f"{xi},{yi},{zi},{ti}\n")
print(f"Simulating electron avalanche from r={x0}...") print(f"Simulating electron avalanche from r={x0}...")
# AvalancheElectron(x, y, z, t, energy, dx, dy, dz) # AvalancheElectron(x, y, z, t, energy, dx, dy, dz)
aval.AvalancheElectron(x0, y0, z0, t0, 0.1, 0.0, 0.0, 0.0) aval.AvalancheElectron(x0, y0, z0, t0, 0.1, 0.0, 0.0, 0.0)
with open("avalanche_endpoints.csv", "w") as f:
f.write("x,y,z,t\n")
# After aval.AvalancheElectron(...)
n_endpoints = aval.GetNumberOfEndpoints()
for i in range(n_endpoints):
# Get start and end points of each electron in the avalanche
x1, y1, z1, t1, e1 = ROOT.double(0), ROOT.double(0), ROOT.double(0), ROOT.double(0), ROOT.double(0)
x2, y2, z2, t2, e2, status = ROOT.double(0), ROOT.double(0), ROOT.double(0), ROOT.double(0), ROOT.double(0), ROOT.int(0)
aval.GetEndpoint(i, x1, y1, z1, t1, e1, x2, y2, z2, t2, e2, status)
# We save the endpoint (x2, y2, z2) where the electron was collected or attached
f.write(f"{x2},{y2},{z2},{t2}\n")
print("Simulation complete.") print("Simulation complete.")

View File

@ -7,9 +7,10 @@ count=11
while val<178.3+0.1: while val<178.3+0.1:
print(val) print(val)
os.system("python3 wires_gmsh2d_bc.py "+str(val)) os.system("python3 wires_gmsh2d_bc.py "+str(val))
os.system("ElmerGrid 14 2 wires2d.msh") os.system("ElmerGrid 14 2 wires2d.msh -2d")
os.system("ElmerSolver wires2d.sif") os.system("ElmerSolver wires2d.sif")
os.system("./paraview_plotter.py") os.system("./paraview_plotter.py")
# os.system("python3 garfield_sim.py")
os.system("cp wires2d.msh wires2d/mesh_files/wires2d%02d_%1.4f.msh"%(count,val)) os.system("cp wires2d.msh wires2d/mesh_files/wires2d%02d_%1.4f.msh"%(count,val))
os.system("cp wires2d.sif wires2d/sif_files/wires2d_%02d_%1.4f.sif"%(count,val)) os.system("cp wires2d.sif wires2d/sif_files/wires2d_%02d_%1.4f.sif"%(count,val))
os.system("cp wires2d/elfield_anasen_t0001.vtu wires2d/vtu_files/elfield_anasen_%02d_%1.4f.vtu"%(count,val)) os.system("cp wires2d/elfield_anasen_t0001.vtu wires2d/vtu_files/elfield_anasen_%02d_%1.4f.vtu"%(count,val))