Difference between revisions of "Simulations of Particle Interactions with Matter"

From New IAC Wiki
Jump to navigation Jump to search
 
(995 intermediate revisions by 4 users not shown)
Line 1: Line 1:
−
== Overview ==
+
==Class Admin==
  
−
=== Particle Detection ===
 
−
A device detects a particle only after the particle transfers energy to the device.
 
  
−
Energy intrinsic to a device depends on the material used in a device
+
[[TF_SPIM_ClassAdmin]]
  
−
Some device of material with an average atomic number (<math>Z</math>) is at some temperature (<math>T</math>).  The materials atoms are in constant thermal motion (unless T = zero degrees Klevin).
+
== Homework Problems==
 +
[[HomeWork_Simulations_of_Particle_Interactions_with_Matter]]
  
−
Statistical Thermodynamics tells us that the canonical energy distribution of the atoms is given by the Maxwell-Boltzmann statistics such that
+
=Introduction=
  
−
<math>P(E) = \frac{1}{kT} e^{-\frac{E}{kT}}</math>
+
[[TF_SPIM_Intro]]
  
−
<math>P(E)</math> represents the probability of any atom in the system having an energy <math>E</math> where
+
= Energy Loss =
  
−
<math>k= 1.38 \times 10^{-23} \frac{J}{mole \cdot K}</math>
+
[[TF_SPIM_StoppingPower]]
  
−
Note:  You may be more familiar with the Maxwell-Boltzmann distribution in the form
+
Ann. Phys. vol. 5, 325, (1930)
  
−
<math>N(\nu) = 4 \pi N \left ( \frac{m}{2\pi k T} \right ) ^{3/2} v^2  e^{-mv^2/2kT}</math>
+
=Interactions of Electrons and Photons with Matter=
  
−
where <math>N(v) \Delta v</math> would represent the molesules in the gas sample with speeds between <math>v</math> and <math>v + \Delta v</math>
+
[[TF_SPIM_e-gamma]]
  
−
==== Example 1: P(E=5 eV) ====
+
Physics Reference
  
−
;What is the probability that an atom in a 12.011 gram block of carbon would have and energy of 5 eV?
+
https://indico.cern.ch/event/679723/contributions/2792554/attachments/1559217/2454299/PhysicsReferenceManual.pdf
  
−
First lets check that the probability distribution is Normailized; ie: does <math>\int_0^{\infty} P(E) dE =1</math>?
+
Physics lists
 +
https://geant4.web.cern.ch/documentation/dev/plg_html/PhysicsListGuide/physicslistguide.html
  
 +
Livermore is the default model
  
−
<math>\int_0^{\infty} P(E) dE = \int_0^{\infty} \frac{1}{kT} e^{-\frac{E}{kT}} dE = \frac{1}{kT} \frac{1}{\frac{1}{-kT}} e^{-\frac{E}{kT}} \mid_0^{\infty} = - [e^{-\infty} - e^0]= 1</math>
+
https://www.epj-conferences.org/articles/epjconf/pdf/2019/19/epjconf_chep2018_02046.pdf
  
−
<math>P(E=5eV)</math> is calculated by integrating P(E) over some energy interval ( ie:<math> N(v) dv</math>).  I will arbitrarily choose 4.9 eV to 5.1 eV as a starting point.
 
  
 +
New PW models
  
−
<math>\int_{4.9 eV}^{5.1 eV} P(E) dE =  - [e^{-5.1 eV/kT} - e^{4.9 eV/kT}]</math>
+
https://indico.cern.ch/event/629841/contributions/2712690/attachments/1518832/2371867/PE_models_G4GV.pdf
  
−
<math>k= (1.38 \times 10^{-23} \frac{J}{mole \cdot K} )  =  (1.38 \times 10^{-23} \frac{J}{mole \cdot K} )(6.42 \times 10^{18} \frac{eV}{J})= 8.614 \times 10^{-5} \frac {eV}{mole \cdot K}</math>
+
https://geant4-userdoc.web.cern.ch/UsersGuides/PhysicsReferenceManual/html/index.html
  
−
assuming a room empterature of <math>T=300 K</math>
+
https://opengate.readthedocs.io/en/latest/introduction.html
  
−
then<math>kT = 0.0258 \frac{eV}{mole}</math>
+
= Hadronic Interactions =
  
−
and
+
[[TF_SPIM_HadronicInteractions]]
  
−
<math>\int_{4.9 eV}^{5.1 eV} P(E) dE = - [e^{-5.1/0.0258} - e^{4.9/0.0258}] = 4.48 \times 10^{-83} - 1.9 \times 10^{-86} \approx 4.48 \times 10^{-83}</math>
+
= Final Project=
  
−
or in other words the precise mathematical calculation of the probability may be approximated by just using the distribution function alone
+
A final project will be submitted that will be graded with the following metrics:
  
−
<math>P(E=5eV) = e^{-5/0.0258} \approx 10^{-85}</math>
+
1.) The document must be less than 15 pages.
  
−
This approximation breaks down as <math>E \rightarrow 0.0258 eV</math>
+
2.) The document must contain references in a bibliography (5 points) .
  
−
Since we have 12.011 grams of carbon and 1 mole of carbon = 12.011 g = <math>6 \times 10^{23} </math>carbon atoms
+
3.) A comparison must be made between GEANT4's prediction and either the prediction of someone else or an experimental result(30 points).
  
−
We do not expect to see a 5 eV carbon atom in a sample size of <math>6 \times 10^{23} </math> carbon atoms when the probability of observing such an atom is <math>\approx 10^{-85}</math>
+
4.) The graphs must be of publication quality with font sizes similar or larger than the 12 point font (10 points).
  
−
The energy we expect to see would be calculated by
+
5.) The document must be grammatically correct (5 points).
  
−
<math><E> = \int_{0}^{\infty} E \cdot P(E) dE</math>
+
6.) The document format must contain the following sections: An abstract of 5 sentences (5 points) , an Introduction(10 points), a Theory section (20 points) , if applicable a section describing the experiment that was simulated, a section delineating the comparisons that were made, and a conclusion( 15 points).
  
−
If you used this block of carbon as a detector you would easily notice an event in which a carbon atom absorbed 5 eV of energy as compared to the energy of a typical atom in the carbon block.
+
=Resources=
  
−
----
+
[http://geant4.web.cern.ch/geant4/  GEANT4 Home Page]
  
−
;Silicon detectors and Ionization chambers are two commonly used devices for detecting radiation.
+
[http://root.cern.ch ROOT Home page]
  
−
approximately 1 eV of energy is all that you need to create an electron-ion pair in Silicon
+
[http://conferences.fnal.gov/g4tutorial/g4cd/Documentation/WorkshopExercises/  Fermi Lab Example]
  
−
<math>P(E=1 eV) = e^{-1/0.0258} \approx 10^{-17}</math>
 
  
−
approximately 10 eV of energy is needed to ionize an atom in a gas chamber
+
[http://physics.nist.gov/PhysRefData/Star/Text/ESTAR.html NIST Range Tables]
  
−
<math>P(E=10 eV) = e^{-10/0.0258} \approx 10^{-169}</math>
+
[http://ie.lbl.gov/xray/  X-ray specturm]
  
 +
[[Installing_GEANT4.9.3_Fsim]]
  
 +
== Saving/restoring Random number seed==
  
−
The low probability of having an atom with 10 eV of energy means that an ionization chamber would have a better Signal to Noise ratio (SNR) for detecting 10 eV radiation than a silicon detector
+
You save the current state of the random number generator with the command
  
−
But if you cool the silicon detector to 200 degrees Kelvin (200 K) then
+
/random/setSavingFlag 1
  
−
<math>P(E=1 eV) = e^{-1/0.0172} \approx 10^{-26}  10^{-17}<< </math>
+
/run/beamOn 100
  
−
So cooling your detector will slow the atoms down making it more noticable when one of the atoms absorbs energy.
+
/random/saveThisRun
  
−
also, if the radiation flux is large, more electron-hole pairs are created and you get a more noticeable signal.
+
A file is created called
  
−
Unfortunately, with some detectore, like silicon, you can cause radiation damage that diminishes it's quantum efficiency for absorbing energy.
+
currentEvent.rndm
  
−
=== The Monte Carlo method ===
+
/control/shell mv currentEvent.rndm currentEvent10.rndm
−
; Stochastic
 
−
: from the greek word "stachos"
 
−
: a means of, relating to, or characterized by conjecture and randomness.
 
  
  
−
A stochastic process is one whose behavior is non-deterministic in that the next state of the process is partially determined.
+
You can restore the random number generator and begin generating random number from the last save time
  
−
Physics has many such non-deterministic systems:
+
/random/resetEngineFrom currentEvent.rndm
  
−
*Quantum Mechanics
 
−
*Thermodynamics
 
  
 +
== Creating Template==
  
−
Basically the monte-carlo method uses a random number generator (RNG) to generate a distribution (gaussian, uniform, Poission,...) which is used to solve a stochastic process based on an astochastic description.
+
[[TForest_G4_Template]]
  
 +
==Building GEANT4.11==
  
−
==== Example 2 Calculation of <math>\pi</math>====
+
===4.11.2===
 +
[[TF_GEANT4.11]]
  
−
;Astochastic description:
+
==Building GEANT4.10==
−
: <math>\pi</math> may be measured as the ratio of the area of a circle of radius <math>r</math> divided by the area of a square of length <math>2r</math>
 
  
−
[[Image:PI_from_AreaRatio.jpg]]<math>\frac{A_{circle}}{A_{square}} = \frac{\pi r^2}{4r^2} = \frac{\pi}{4}</math>
 
  
−
You can measure the value of <math>\pi</math> if you physically measure the above ratios.
+
===4.10.02===
  
−
; Stochastic description:
+
[[TF_GEANT4.10.2]]
−
: Construct a dart board representing the above geometry, throw several darts at it, and look at a ratio of the number of darts in the circle to the total number of darts thrown (assuming you always hit the dart board).
+
===4.10.01===
  
−
; Monte-Carlo Method
+
[[TF_GEANT4.10.1]]
−
:Here is an outline of a program to calulate <math>\pi</math> using the Monte-Carlo method with the above Stochastic description
 
−
[[Image:MC_PI_fromAreaRatio.jpg]]
 
−
begin loop
 
−
  x=rnd
 
−
  y=rnd
 
−
  dist=sqrt(x*x+y*y)
 
−
  if dist <= 1.0 then numbCircHits+=1.0
 
−
  numbSquareHist += 1.0
 
−
end loop
 
−
  print PI = 4*numbCircHits/numbSquareHits
 
  
−
=== A Unix Primer ===
+
==Building GEANT4.9.6==
−
To get our feet wet using the UNIX operating system, we will try to solve example 2 above using a RNG under UNIX
 
  
−
==== List of important Commands====
+
[[TF_GEANT4.9.6]]
  
−
# ls
+
==Building GEANT4.9.5==
−
# pwd
 
−
# cd
 
−
# df
 
−
# ssh
 
−
# scp
 
−
# mkdir
 
−
# printenv
 
−
# emacs, vi, vim
 
−
# make, gcc
 
−
# man
 
−
# less
 
−
# rm
 
  
−
----
+
[[TF_GEANT4.9.5]]
−
Most of the commands executed within a shell under UNIX have command line arguments (switches) which tell the command to print information about using the command to the screen. The common forms of these switches are "-h", "--h", or "--help"
 
  
−
ls --help
+
An old version of Installation notes for versions prior to 9.5
−
ssh -h
 
  
−
'' the switch deponds on your flavor of UNIX''
+
[http://brems.iac.isu.edu/~tforest/NucSim/Day3/ Old Install Notes]
  
−
if using the switch doesn;t help you can try the "man" (sort for manual) pages (if they were installed). 
 
−
Try
 
−
man -k pwd
 
  
−
the above command will search the manual for the key word "pwd"
+
Visualization Libraries:
  
−
==== Example 3: using UNIX ====
+
[http://www.opengl.org/ OpenGL]
  
−
Step
+
[http://geant4.kek.jp/~tanaka/DAWN/About_DAWN.html DAWN]
−
# login to inca.<br> [http://physics.isu.edu/~tforest/Classes/NucSim/Day1/XwindowsOnWindows.html click here for a description of logging in if using windows]
 
−
# mkdir src
 
−
# cd src
 
−
# cp -R ~tforest/NucSim/Day1 ./
 
−
# ls
 
−
# cd Day1
 
−
# make
 
−
#./rndtest
 
  
−
[http://physics.isu.edu/~tforest/Classes/NucSim/Day1/RNG/Marsaglia/noviceExample/ Here is a web link to the source files you can copy in case the above doesn't work]
 
  
−
=== A Root Primer ===
+
[http://doc.coin3d.org/Coin/ Coin3D]
−
==== Example 1: Create Ntuple and Draw Histogram====
 
  
−
=== Cross Sections ===
+
==Compiling G4 with ROOT==
−
==== Definitions ====
 
−
<math>\sigma(\theta)</math> = scattering cross-section <math>\equiv \frac{\frac{\# particles\; scattered}{solid \; angle}} {\frac{ \# incident \; particles}{Area}}</math>
 
  
−
; Solid Angle
+
These instruction describe how you can create a tree within ExN02SteppingVerbose to store tracking info in an array (max number of steps in a track is set to 100 for the desired particle)
−
:[[Image:SolidAngleDefinition.jpg]]
 
−
: <math>\Omega</math>= surface area of a sphere covered by the detector
 
−
: ie;the detectors area projected onto the surface of a sphere
 
−
:A= surface area of detector
 
−
:r=distance from interaction point to detector
 
−
:<math>\Omega = \frac{A}{r^2} </math>sterradians
 
−
: <math>A_{sphere} = 4 \pi r^2</math> if your detector was a hollow ball
 
−
:<math>\Omega_{max} = \frac{4 \pi r^2}{r^2} = 4\pi</math>sterradians
 
  
−
;Units
+
[[G4CompileWRootforTracks]]
−
:Cross-sections have the units of Area
 
−
:1 barn = <math>10^{-28} m^2</math>
 
−
; [units of <math>\sigma(\theta)</math>] =<math>\frac{\frac{[particles]}{[sterradian]}} {\frac{ [ particles]}{[m^2]}} = m^2</math>
 
  
 +
==Using SLURM==
  
−
[[Image:FixedTargetScatteringCrossSection.jpg]]
+
http://slurm.schedmd.com/quickstart.html
−
; Fixed target scattering
 
−
: <math>N_{in}</math>= # of particles in = <math>I \cdot A_{in}</math>
 
−
:: <math>A_{in}</math> is the area of the ring of incident particles
 
−
:<math>dN_{in} = I \cdot dA = I (2\pi b) db</math>= # particles in a ring of radius <math>b</math> and thickness <math>db</math>
 
  
−
You can measure <math>\sigma(\theta)</math> if you measure the # of particles detected <math>d N</math> in a known detector solid angle <math>d \Omega</math> from a know incident particle Flux (<math>I</math>)  as
 
  
−
<math>\sigma(\theta) = \frac{\frac{d N}{ d \Omega}}{I}</math>
+
https://rc.fas.harvard.edu/resources/documentation/convenient-slurm-commands/
  
−
Alternatively if you have a theory which tells you <math>\sigma(\theta)</math> which you want to test experimentally with a beam of flux <math>I</math> then you would measure counts (particles)
+
===simple batch script for one process job===
  
−
<math>dN = I \sigma(\theta) d \Omega = I \sigma(\theta)  \frac{d A}{r^2} = I \sigma(\theta) \frac{r^2 \sin(\theta) d \theta d \phi}{r^2}</math>
+
create the file submit.sbatch below
  
−
;Units
+
<pre>
−
: <math>[d N] = [\frac {particles}{m^2}][m^2] [sterradian] </math> = # of particles
+
#!/bin/sh
−
: or for a count rate divide both sides by time and you get beam current on the RHS
+
#SBATCH --time=1
−
: integrate and you have the total number of counts
+
cd src/PI
 +
./PI_MC 100000000000000
 +
</pre>
  
−
;Classical Scattering
+
the execute
−
: In classical scattering you get the same number of particle out that you put in (no capture, conversion,..)
 
−
: <math>d N_{in} = dN</math>
 
−
:<math>d N_{in} = I dA = I (2\pi b) db</math>
 
−
: <math>d N = I \sigma(\theta) d \Omega =  I \sigma(\theta) \sin(\theta) d \theta d \phi = I \sigma(\theta) \sin(\theta) d \theta (2 \pi )</math>
 
−
:<math>  I (2\pi b) db =  I \sigma(\theta) \sin(\theta) d \theta (2 \pi )</math>
 
−
:<math>    b  db =  \sigma(\theta) \sin(\theta) d \theta </math>
 
−
:<math>\sigma(\theta) =  \frac{b}{\sin(\theta)}\frac{db}{d \theta}</math>
 
−
:<math>\frac{db}{d \theta}</math> tells you how the impact parameter <math>b</math> changes with scattering angle <math>\theta</math>
 
  
−
==== Example : Elastic Scattering ====
+
:sbatch submit.sbatch
−
This example is an example of classical scattering.
 
  
−
Our goal is to find <math>\sigma(\theta)</math> for an elastic collision of 2 impenetrable spheres of diameter <math>a</math>.  To solve this elastic scattering problem we will describe the collision in the Center of Mass (C.M.) frame.  As we shall see, in the C.M. fram the 2-body collision becomes a 1-body problem.
+
check if its running with
  
−
[[Image:SPIM_ElasCollis_Lab_CM_Frame.jpg]]
+
:squeue
  
−
; Variable definitions
+
to kill a batch job
−
:<math>b</math>= impact parameter ; distance of closest approach
 
−
:<math>m_1</math>= mass of incoming ball
 
−
:<math>m_2</math>= mass of target ball
 
−
:<math>u_1</math>= iniital velocity of  incoming ball in Lab Frame
 
−
:<math>v_1</math>= final velocity of  <math>m_1</math> in Lab Frame
 
−
:<math>\psi_1</math>= scattering angle of <math>m_1</math> in lab frame after collision
 
−
:<math>u_1^{\prime}</math>= iniital velocity of  <math>m_1</math> in C.M.  Frame
 
−
:<math>v_1^{\prime}</math>= final velocity of  <math>m_1</math> in C.M. Frame
 
−
:<math>u_2^{\prime}</math>= iniital velocity of  <math>m_2</math> in C.M.  Frame
 
−
:<math>v_2^{\prime}</math>= final velocity of  <math>m_2</math> in C.M. Frame
 
−
:<math>\theta</math>= scattering angle of <math>m_1</math> in C.M. frame after collision
 
  
 +
:scancel JOBID
  
−
We can reduct the 2-body problem to a 1-body problem using the following coordinates
+
===On minerve===
  
−
[[Image:SPIM_2Body-1BodyCoordSystem.jpg]]
+
Sample script to submit 10 batch jobs.
  
−
; vector definitions
+
the filename is minervesubmit and you run like
−
:<math>\vec{r_1}</math> = a position vector pointing to the location of <math>m_1</math>
+
source minervesubmit
−
:<math>\vec{r_2}</math> = a position vector pointing to the location of <math>m_2</math>
 
−
:<math>\vec{R}</math> = a position vector pointing to the center of mass of the two ball system
 
−
:<math>\vec{r} \equiv \vec{r_1} - \vec{r_2}</math> = the magnitude of this vector is the distance between the two masses
 
  
−
In the C.M. reference frame the above vectors have the following relationships
+
<pre>
 +
cd /home/foretony/src/GEANT4/geant4.9.5/Simulations/N02wROOT/batch
 +
qsub submit10mil
 +
qsub submit20mil
 +
qsub submit30mil
 +
qsub submit40mil
 +
qsub submit50mil
 +
qsub submit60mil
 +
qsub submit70mil
 +
qsub submit80mil
 +
qsub submit90mil
 +
qsub submit100mil
 +
</pre>
  
−
# <math>\vec{R} = 0 = \frac{m1 \vec{r_1} + m2 \vec{r_2}}{m_1 + m_2} \Rightarrow m_2 \vec{r_1} = -m_2 \vec{r_2}</math>
+
The file submit10mil looks like this
−
# <math>\vec{r_1} - \vec{r_2} = \vec{r}</math>
+
<pre>
 +
#!/bin/sh
 +
#PBS -l nodes=1
 +
#PBS -A FIAC
 +
#PBS -M foretony@isu.edu
 +
#PBS -m abe
 +
#
 +
source /home/foretony/src/GEANT4/geant4.9.5/geant4.9.6-install/bin/geant4.sh
 +
cd /home/foretony/src/GEANT4/geant4.9.5/Simulations/N02wROOT/batch/10mil
 +
../../exampleN02 run1.mac > /dev/null
 +
</pre>
  
−
solving the above equations for <math>\vec{r_1}</math> and <math>\vec{r_2}</math> and defining the reduced mass <math>\mu</math> as
 
  
−
:<math>\mu = \frac{m_1 \cdot m_2}{m_1 + m_2} \equiv</math> reduced mass
+
use
  
−
leads to
+
qstat
  
−
: <math>\vec{r_1} = \frac{\mu}{m_1} \vec{r}</math>
+
to check that the process is still running
−
: <math>\vec{r_2} = \frac{\mu}{m_2} \vec{r}</math>
 
  
−
We can use the above relationships to construct a Hamilton in terms of <math>\vec{r}</math> instead of <math>\vec{r_1}</math> and <math> \vec{r_2}</math> thereby reducing the problem from a 2-body problem to a 1-body problem.
+
use
  
−
; Construct the Hamiltonian
+
qdel jobID#
  
−
To construct the Hamiltonian for this problem we will start with the Lagrangian.
+
if you want to kill the batch job, the jobID number shows up when you do stat.
  
 +
for example
 +
<pre>
 +
[foretony@minerve HW10]$ qstat
 +
Job id                    Name            User            Time Use S Queue
 +
------------------------- ---------------- --------------- -------- - -----
 +
27033.minerve            submit          foretony        00:41:55 R default       
 +
[foretony@minerve HW10]$ qdel 27033
 +
[foretony@minerve HW10]$ qstat
 +
</pre>
  
−
<math>\mathcal{L} = T - U</math>
+
==Definitions of Materials==
  
−
where
+
[[File:MCNP_Compendium_of_Material_Composition.pdf]]
  
−
<math>T \equiv</math> kinetic energy of the system
+
==Minerve2 GEANT 4.10.1 Xterm error==
  
−
<math>U \equiv</math> Potential energy of the system which describes the interaction
 
  
 +
On OS X El Capitan V 10.11.4 using XQuartz
  
−
<math>\mathcal{L} = \frac{1}{2} |\dot{\vec{r_1}}|^2 + \frac{1}{2} |\dot{\vec{r_2}}|^2 - U</math>
+
~/src/GEANT4/geant4.10.1/Simulations/B2/B2a/exsmpleB2a
−
:= <math>\frac{1}{2} m_1 \left (\frac{m_2}{m1+m_2} \right )^2  |\dot{\vec{r}}|^2 + \frac{1}{2} m_2 \left (\frac{m_1}{m1+m_2} \right )^2  |\dot{\vec{r}}|^2 -U(\vec{r})</math>
 
  
−
after substituting derivative of the expressions for <math>\vec{r_1}</math> and <math>\vec{r_2}</math>
+
<pre>
  
−
: = <math>\frac{1}{2} \mu |\dot{\vec{r}}|^2 -U(\vec{r})</math> The 2-body problem is now described by a 1-body Lagrangian
+
# Use this open statement to create an OpenGL view:
 +
/vis/open OGL 600x600-0+0
 +
/vis/sceneHandler/create OGL
 +
/vis/viewer/create ! ! 600x600-0+0
 +
libGL error: failed to load driver: swrast
 +
X Error of failed request:  BadValue (integer parameter out of range for operation)
 +
  Major opcode of failed request:  150 (GLX)
 +
  Minor opcode of failed request:  3 (X_GLXCreateContext)
 +
  Value in failed request:  0x0
 +
  Serial number of failed request:  25
 +
  Current serial number in output stream:  26
 +
</pre>
  
−
Lagranges equations of motion are given by
 
−
: <math>\frac{\partial \mathcal{L}}{\partial q} = \frac{d}{dt} \frac{\partial \mathcal{L}}{\dot{q}}</math>
 
−
where <math>q</math> represents on of the coordinate (cannonical variables).
 
  
−
To get the classical scattering cross section we are interested in finding an expression for the dependence of the impact parameter on the scattering angle,<math>\frac{d b}{d \theta}</math>.  In this case <math>\theta</math> is the C.M. scattering angle so the cross section is calcuated in the C.M. reference Frame.
 
  
−
Now lets redraw the collision in terms of the reduced mass in the Lab frame.
 
  
−
[[Image:SPIM_ElasCollis_ReducedMass_CM_Frame_1.jpg]]
 
  
−
The C.M. Frame rides along the center of mass <math>\mu</math>.  If <math>b > a</math> then there is no collision (<math>\theta=0</math>),  otherwise a collision happens when r=a (the distance between the balls is equal to their diameter).  A head on collision is defined as <math>b=0</math> (<math>\theta=\pi</math>).
+
[[TF_SPIM_OLD]]
−
 
 
−
;Observation
 
−
: as <math>\theta</math> gets smaller, <math>b</math> gets bigger
 
−
: <math>\frac{d b}{d \theta} < 0</math>
 
−
 
 
−
Using plane polar coordinates (<math>R, \phi</math>) we can describe the problem in the lab frame as:
 
−
 
 
−
<math>\vec{v} = \dot{R} \hat{e}_R + R \dot{\phi} \hat{e}_{\phi}</math>
 
−
 
 
−
<math>T = \frac{1}{2} \mu ( \dot{R}^2 + R^2 \dot{\phi}^2)</math>
 
−
 
 
−
 
 
−
<math>U(R) = \left \{  {0 \; R > a \atop \infty \; R \le a} \right .</math>
 
−
 
 
−
<math>\mathcal{L} = T -U = \frac{1}{2} \mu ( \dot{R}^2 + R^2 \dot{\phi}^2) - U(R)</math>
 
−
 
 
−
Lagranges Equation of Motion:
 
−
 
 
−
<math>\frac{\partial \mathcal {L}}{\partial \phi} = \frac{d}{d t} \frac{\partial \mathcal{L}}{\partial \dot{\phi}}</math>
 
−
<math>0 = \frac{d}{d t} [ \mu R^2 \dot{\phi}] \Rightarrow</math>  there is a constant of motion ( Constant angular momentum)
 
−
 
 
−
<math>\ell \equiv \mu R^2 \dot{\phi} = \vec{R} \times \vec{p} = \vec{R} \times \mu \vec{v} = R^2 \mu \dot{\phi}</math>
 
−
 
 
−
substitute <math>\ell</math> into <math>\mathcal{L}</math>
 
−
 
 
−
<math>\mathcal{L} = \frac{1}{2} ( \mu  \dot{R}^2 + \frac{\ell}{\mu R^2} ) - U(R)</math>
 
−
 
 
−
The two equations above are in terms of <math>R</math> and <math>\phi</math> whereas our goal is to find an expression for <math>\frac{ d b}{ d \theta}</math>.  Since <math>R</math> is related to <math>b</math> and <math>\phi</math> is related to<math> \theta</math> (<math>\theta = \pi - 2\phi</math>; see figure above) we should try and find expressions for <math>d \phi</math> in terms of <math>R(b)</math>
 
−
 
 
−
;Trick
 
−
: <math>\dot{\phi} = \frac{d \phi}{d t} = \frac{d \phi}{d R} \frac{d R}{d t}</math>
 
−
:<math>\Rightarrow \ell = \mu R^2 \frac{d \phi}{d R} \dot{R}</math>
 
−
:or
 
−
: <math>d \phi = \frac{\ell}{\mu R^2 \dot{R}} dR</math>
 
−
 
 
−
We now need an expression for <math>\dot{R}</math> in order to integrate the above equation to determine the functional dependence of <math>\phi</math> and hence<math> \theta</math>.
 
−
 
 
−
Since Energy is conserved (Elastic Scattering), we may define the Hamiltonian as
 
−
 
 
−
<math>H = T + U = \frac{1}{2} (mu \dot{R}^2 + \frac{\ell}{\mu R^2}) + U(R) = constant \equiv E</math>
 
−
 
 
−
solving for <math>\dot{R}</math>
 
−
 
 
−
<math>\dot{R} = \pm \sqrt{\frac{2(E-U(R))}{\mu} - \frac{\ell^2}{\mu^2 R^2}}</math>
 
−
 
 
−
substituting the above into the equation for <math>d \phi</math> and integrating:
 
−
 
 
−
<math>\phi = \int d \phi = \int_{r_{min}}^{r_{max}} \frac{\ell}{\mu R^2 \dot{R}} dR</math>
 
−
 
 
−
<math>r_{min} = a  \; \; \;  r_{max}= \infty  \; \; \;  U(R) = 0 : a \le R \le \infty</math>
 
−
 
 
−
<math>\phi = \int_a^{\infty} \frac{\ell} {R^2 \sqrt{2 \mu E - \frac{\ell^2}{R^2}} }dR</math>
 
−
 
 
−
For <math>a \le R \le \infty</math> : <math>E = \frac{1}{2} \mu v^2_{cm} \Rightarrow v_{cm} = \sqrt{\frac{2E}{\mu}}</math>
 
−
 
 
−
<math>\vec{\ell} = \vec{R} \times \vec{p} \Rightarrow |\vec{\ell}| = |\vec{R}| |\vec{p}| \sin(\phi) = R \mu v_{cm} \sin(\phi) = R \mu \left ( \sqrt{\frac{2E}{\mu}} \right) \sin(\phi) = \sqrt{2 \mu E} R\sin(\phi) =\sqrt{2 \mu E} b</math>
 
−
 
 
−
substituting this expression for \ell into the last expression for \phi above :
 
−
 
 
−
<math>\phi =\int_a^{\infty}  \frac{b dR}{R\sqrt{(R^2-b^2)}}</math>
 
−
 
 
−
;Integral Table
 
−
: <math>\int  \frac{dx}{x\sqrt{(\alpha x^2+\beta x+\gamma)}} = \frac{-1}{\sqrt{-\gamma}} \sin^{-1} \left (\frac{\beta x+2\gamma}{|x|\sqrt{\beta^2-4\alpha \gamma}} \right )</math>
 
−
 
 
−
let <math>x=R \;\; \alpha=1 \;\; \beta=0 \;\; \gamma=-b^2</math>
 
−
 
 
−
then
 
−
 
 
−
<math>\phi = \left . b \frac{1}{\sqrt{-(-b^2)}} \sin^{-1} \left (\frac{-2b^2}{R\sqrt{0-4(1)(-b^2) } }\right ) \right |_a^{\infty} = \sin^{-1} (0)- \sin^{-1}(-\frac{b}{a})</math>
 
−
 
 
−
or
 
−
<math>\sin(\phi) = \frac{b}{a} = \sin \left ( \frac{\pi}{2} - \frac{\theta}{2} \right ) = \cos \left ( \frac{\theta}{2} \right )</math>
 
−
 
 
−
<math>\Rightarrow b = a \cos \left( \frac{\theta}{2} \right)</math>
 
−
 
 
−
==== Lab Frame Cross Sections ====
 
−
 
 
−
== Stopping Power ==
 
−
=== Bethe Equation ===
 
−
====Classical Energy Loss ====
 
−
====Bethe-Bloch Equation ====
 
−
=== Energy Straggling ===
 
−
==== Thick Absorber ====
 
−
====Thin Absorbers====
 
−
=== Range Straggling===
 
−
=== Electron Capture and Loss ===
 
−
=== Multiple Scattering ===
 
−
 
 
−
== Interactions of Electrons and Photons with Matter==
 
−
=== Bremsstrahlung===
 
−
=== Photo-electric effect===
 
−
=== Compton Scattering ===
 
−
=== Pair Production ===
 
−
 
 
−
== Hadronic Interactions ==
 
−
=== Neutron Interactions ===
 
−
==== Elastic scattering====
 
−
 
 
−
==== Inelasstic Scattering====
 

Latest revision as of 00:15, 27 January 2025

Class Admin

TF_SPIM_ClassAdmin

Homework Problems

HomeWork_Simulations_of_Particle_Interactions_with_Matter

Introduction

TF_SPIM_Intro

Energy Loss

TF_SPIM_StoppingPower

Ann. Phys. vol. 5, 325, (1930)

Interactions of Electrons and Photons with Matter

TF_SPIM_e-gamma

Physics Reference

https://indico.cern.ch/event/679723/contributions/2792554/attachments/1559217/2454299/PhysicsReferenceManual.pdf

Physics lists https://geant4.web.cern.ch/documentation/dev/plg_html/PhysicsListGuide/physicslistguide.html

Livermore is the default model

https://www.epj-conferences.org/articles/epjconf/pdf/2019/19/epjconf_chep2018_02046.pdf


New PW models

https://indico.cern.ch/event/629841/contributions/2712690/attachments/1518832/2371867/PE_models_G4GV.pdf

https://geant4-userdoc.web.cern.ch/UsersGuides/PhysicsReferenceManual/html/index.html

https://opengate.readthedocs.io/en/latest/introduction.html

Hadronic Interactions

TF_SPIM_HadronicInteractions

Final Project

A final project will be submitted that will be graded with the following metrics:

1.) The document must be less than 15 pages.

2.) The document must contain references in a bibliography (5 points) .

3.) A comparison must be made between GEANT4's prediction and either the prediction of someone else or an experimental result(30 points).

4.) The graphs must be of publication quality with font sizes similar or larger than the 12 point font (10 points).

5.) The document must be grammatically correct (5 points).

6.) The document format must contain the following sections: An abstract of 5 sentences (5 points) , an Introduction(10 points), a Theory section (20 points) , if applicable a section describing the experiment that was simulated, a section delineating the comparisons that were made, and a conclusion( 15 points).

Resources

GEANT4 Home Page

ROOT Home page

Fermi Lab Example


NIST Range Tables

X-ray specturm

Installing_GEANT4.9.3_Fsim

Saving/restoring Random number seed

You save the current state of the random number generator with the command

/random/setSavingFlag 1

/run/beamOn 100

/random/saveThisRun

A file is created called

currentEvent.rndm

/control/shell mv currentEvent.rndm currentEvent10.rndm


You can restore the random number generator and begin generating random number from the last save time

/random/resetEngineFrom currentEvent.rndm


Creating Template

TForest_G4_Template

Building GEANT4.11

4.11.2

TF_GEANT4.11

Building GEANT4.10

4.10.02

TF_GEANT4.10.2

4.10.01

TF_GEANT4.10.1

Building GEANT4.9.6

TF_GEANT4.9.6

Building GEANT4.9.5

TF_GEANT4.9.5

An old version of Installation notes for versions prior to 9.5

Old Install Notes


Visualization Libraries:

OpenGL

DAWN


Coin3D

Compiling G4 with ROOT

These instruction describe how you can create a tree within ExN02SteppingVerbose to store tracking info in an array (max number of steps in a track is set to 100 for the desired particle)

G4CompileWRootforTracks

Using SLURM

http://slurm.schedmd.com/quickstart.html


https://rc.fas.harvard.edu/resources/documentation/convenient-slurm-commands/

simple batch script for one process job

create the file submit.sbatch below 

#!/bin/sh
#SBATCH --time=1
cd src/PI
./PI_MC 100000000000000

the execute

sbatch submit.sbatch

check if its running with

squeue

to kill a batch job

scancel JOBID

On minerve

Sample script to submit 10 batch jobs.

the filename is minervesubmit and you run like 

source minervesubmit

cd /home/foretony/src/GEANT4/geant4.9.5/Simulations/N02wROOT/batch
qsub submit10mil
qsub submit20mil
qsub submit30mil
qsub submit40mil
qsub submit50mil
qsub submit60mil
qsub submit70mil
qsub submit80mil
qsub submit90mil
qsub submit100mil

The file submit10mil looks like this 

#!/bin/sh
#PBS -l nodes=1
#PBS -A FIAC
#PBS -M foretony@isu.edu
#PBS -m abe
#
source /home/foretony/src/GEANT4/geant4.9.5/geant4.9.6-install/bin/geant4.sh
cd /home/foretony/src/GEANT4/geant4.9.5/Simulations/N02wROOT/batch/10mil
../../exampleN02 run1.mac > /dev/null


use 

qstat

to check that the process is still running

use 

qdel jobID#

if you want to kill the batch job, the jobID number shows up when you do stat.

for example 

[foretony@minerve HW10]$ qstat
Job id                    Name             User            Time Use S Queue
------------------------- ---------------- --------------- -------- - -----
27033.minerve             submit           foretony        00:41:55 R default        
[foretony@minerve HW10]$ qdel 27033
[foretony@minerve HW10]$ qstat

Definitions of Materials

File:MCNP Compendium of Material Composition.pdf

Minerve2 GEANT 4.10.1 Xterm error

On OS X El Capitan V 10.11.4 using XQuartz

~/src/GEANT4/geant4.10.1/Simulations/B2/B2a/exsmpleB2a 


# Use this open statement to create an OpenGL view:
/vis/open OGL 600x600-0+0
/vis/sceneHandler/create OGL
/vis/viewer/create ! ! 600x600-0+0
libGL error: failed to load driver: swrast
X Error of failed request:  BadValue (integer parameter out of range for operation)
  Major opcode of failed request:  150 (GLX)
  Minor opcode of failed request:  3 (X_GLXCreateContext)
  Value in failed request:  0x0
  Serial number of failed request:  25
  Current serial number in output stream:  26



TF_SPIM_OLD

Retrieved from "https://wiki.iac.isu.edu/index.php?title=Simulations_of_Particle_Interactions_with_Matter&oldid=130029"