Merge branch 'devel_vignesh' of https://fsunuc.physics.fsu.edu/git/rtang/ANASEN_analysis into devel_vignesh

modified:   anasen_fem/paraview_plotter.py
	modified:   anasen_fem/wires2d.sif
	modified:   anasen_fem/wires_gmsh2d_bc.py

.gitignore

@@ -33,3 +33,4 @@ QQQStage2.C
 anasen_fem/scalars.dat.names
 myenv/
 eloss_calculations/
+anasen_fem/He96_CO2_4_260Torr.gas

anasen_fem/garfield_sim.py

@@ -4,7 +4,7 @@ import sys
 
 # 1. FIX: Manually load the Garfield library if it's not in the ROOT namespace
 # 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):
     ROOT.gSystem.Load(garfield_lib_path)
@@ -44,12 +44,13 @@ else:
 # --- 3. FIELD MAP SETUP ---
 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.nodes", 
+
+fm.Initialise("wires2d/mesh.header", 
               "wires2d/mesh.elements", 
-              "wires2d/mesh.boundary", 
-              "wires2d/elfield_anasen.result", "mm")
+              "wires2d/mesh.nodes", 
+              "wires2d/dielectrics.dat", # Dielectrics (leave as empty string)
+              "wires2d/elstatics.result", 
+              "mm")
 
 # Set the medium (Body 13 from your Gmsh script)
 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}...")
 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}...")
 # AvalancheElectron(x, y, z, t, energy, dx, dy, dz)
 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.")

anasen_fem/paraview_plotter.py

@@ -1,4 +1,4 @@
-#!/home/vsitaraman/ParaView-6.1.0-RC1-MPI-Linux-Python3.12-x86_64/bin/pvbatch
+#!/home/vs19g/ParaView-6.1.0-MPI-Linux-Python3.12-x86_64/bin/pvbatch
 import numpy as np
 import sys
 from paraview.simple import *
@@ -35,15 +35,14 @@ contour_display.SetScalarBarVisibility(renderView, True)
 view = GetActiveView()
 
 # 2. Define your desired coordinate ranges (x_min, x_max, y_min, y_max, z_min, z_max)
-# Example: Look at a box from -10 to 10 in all dimensions
-x_min, x_max = -50.0, 50.0
-y_min, y_max = -50.0, 50.0
-z_min, z_max = -50.0, 50.0
+x_min, x_max = -0.05, 0.05
+y_min, y_max = -0.05, 0.05
+z_min, z_max = -0.05, 0.05
 
 # 3. Calculate Center, Position, and Parallel Scale
 center = [(x_min + x_max) / 2.0, (y_min + y_max) / 2.0, (z_min + z_max) / 2.0]
 # Position the camera far away along Z to look at the center
-position = [center[0], center[1], z_min - 30.0] 
+position = [center[0], center[1], 1.0] 
 # Parallel scale defines how much of the scene is visible. 
 # It is usually half the height of the viewed area.
 view.CameraParallelScale = max((x_max - x_min), (y_max - y_min))/1.6
@@ -69,15 +68,13 @@ contour_display.RenderLinesAsTubes = 0    # Makes lines look smoother at high re
 # 1. Get the active view
 view = GetActiveView()
 
-# 4. Apply settings
-# 1. Set the Focal Point to the middle of the quadrant
-zoom_center = [-25, 25, 0.0] 
+# 1. Set the Focal Point to the middle of the quadrant in metres
+zoom_center = [-0.025, 0.025, 0.0] 
 
 # 2. Tighten the Parallel Scale 
-view.CameraParallelScale = 15 
+view.CameraParallelScale = 0.015 
 
-# 3. Position the Camera
-# Keep it 0.5m away looking "down" at the new center
+# 3. Position the Camera (0.5m away is fine)
 view.CameraPosition = [zoom_center[0], zoom_center[1], 0.5]
 view.CameraFocalPoint = zoom_center
 view.CameraViewUp = [0.0, 1.0, 0.0]
@@ -94,9 +91,8 @@ glyph = Glyph(Input=contour_filter, GlyphType='Arrow') #
 # Orientation Array: Use the 'electric field' vector from Elmer
 glyph.OrientationArray = ['POINTS', 'electric field']
 glyph.ScaleArray = ['POINTS', 'No scale array']
-glyph.ScaleFactor = 1
+glyph.ScaleFactor = 0.001  
 
-# Sampling: Every nth point (Stride 16)
 glyph.GlyphMode = 'Every Nth Point'
 glyph.Stride = 24
 

anasen_fem/run.py

@@ -7,9 +7,10 @@ count=11
 while val<178.3+0.1:
     print(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("./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.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))

anasen_fem/wires2d.sif

@@ -10,6 +10,7 @@ Simulation
   Steady State Max Iterations = 1
   Output File = "elstatics.result"
   Post File = "elstatics.ep"
+  Coordinate Scaling = 0.001  ! Converts mm from Gmsh to meters for Elmer
 End
 
 Constants

anasen_fem/wires_gmsh2d_bc.py

@@ -188,17 +188,17 @@ i2wire_surfs = get_surfs(ic2_wires) if include_ic_wires else []
 all_active_wire_surfs = needle_surfs + gwire_surfs + awire_surfs + cwire_surfs + i1wire_surfs + i2wire_surfs
 gmsh.model.mesh.embed(1, all_active_wire_surfs, 2, anasen_barrel)
     
-# f1 = gmsh.model.mesh.field.add("Distance")
-# gmsh.model.mesh.field.setNumbers(f1, "CurvesList", all_active_wire_surfs)
+f1 = gmsh.model.mesh.field.add("Distance")
+gmsh.model.mesh.field.setNumbers(f1, "CurvesList", all_active_wire_surfs)
 
-# f2 = gmsh.model.mesh.field.add("Threshold")
-# gmsh.model.mesh.field.setNumber(f2, "InField", f1)
-# gmsh.model.mesh.field.setNumber(f2, "SizeMin", 0.1)   # Fine mesh near wires
-# gmsh.model.mesh.field.setNumber(f2, "SizeMax", 10.0)  # Large mesh in empty space
-# gmsh.model.mesh.field.setNumber(f2, "DistMin", 1.0)   # Apply SizeMin within 1mm
-# gmsh.model.mesh.field.setNumber(f2, "DistMax", 20.0)  # Transition to SizeMax by 20mm
+f2 = gmsh.model.mesh.field.add("Threshold")
+gmsh.model.mesh.field.setNumber(f2, "InField", f1)
+gmsh.model.mesh.field.setNumber(f2, "SizeMin", 0.05)   # Fine mesh near wires
+gmsh.model.mesh.field.setNumber(f2, "SizeMax", 5.0)  # Large mesh in empty space
+gmsh.model.mesh.field.setNumber(f2, "DistMin", 0.5)   # Apply SizeMin within 1mm
+gmsh.model.mesh.field.setNumber(f2, "DistMax", 15.0)  # Transition to SizeMax by 20mm
 
-# gmsh.model.mesh.field.setAsBackgroundMesh(f2)
+gmsh.model.mesh.field.setAsBackgroundMesh(f2)
 
 
 # --- Physical Groups ---
@@ -222,8 +222,9 @@ gmsh.model.addPhysicalGroup(2, [anasen_barrel], tag=13, name="gas")
 gmsh.option.setNumber("Mesh.Algorithm", 6)
 
 gmsh.model.mesh.generate(dim=2)
-gmsh.model.mesh.refine()
-gmsh.model.mesh.refine()
+# gmsh.model.mesh.refine()
+# gmsh.model.mesh.refine()
 gmsh.write("wires2d.msh")
+gmsh.model.mesh.setOrder(1)
 #gmsh.fltk.run()
 gmsh.finalize()