Welcome to RadiationDamageCalculations’s documentation!¶
Introduction¶
Radiation Damage Simulation and Calculations with multi-scale tools in the Calculations and Molecular Dynamics.

Tools to complentent my research works

Basic Radiation Damage model¶
Overview¶
def some_function():
interesting = False
print 'This line is highlighted.'
print 'This one is not...'
print '...but this one is.'
Columb scattering with target atom¶
A universal one-parameter differential scattering cross section equation in reduced notation is expressed by Lindhard et al. as:
where t is a dimensionless collision parameter defined by:
T is the transferred energy to the target,and $T_{max}$ is the maximum transferred energy as:
MCCM algorithm¶
The possible damage induced by radiation in solid materials has been a topic of great interest for the scientific community.The effects that may be induced by this way depend on the nature of the study material. To obtain a good effectiveness inside the material volume,inclusive at high depth,it is well-known that the gamma radiation is the best candiate. One of researchs of our group is to research the influence of materials properties under defects.Different defects have different impacts on the optical、electrical and magnetic material properties,However it is difficult to study the single defect’s effect which can give a more clear point of material properties interpretation.Because when the heavy-ions and neutron irradiated the materials,it is a complicated process with various defects and difficult to distinguish different defects effect on material properties.But the good news is that gamma-ray irradiation can separate the defects on the basis of material atom displacement threshold energy(\(E_d\)) when our group reasearch the gamma-ray influence different dioxides.
Basic Models and Methods¶

MCCM algorithm¶
The number of atom displacements in the voxel volume \(N_{dpa}\) from the MCCM algorithm is an extended sum over the calculated discrete electron flux values at a given depth \(z\) ,generalized for any k-atom in the materials:
where \(n_k\) is the relative fraction of the k-atom in its crystalline sublattice and \(N^e_{dpa,k}\) is obtained for each one of these atoms.The number of atoms displaced per electron was calculated following the expression proposed by Oen and Holmes:
where
is the cutoff kinetic energy of electrons in order to displace an atom from its crystalline site,:math:N_a is the number of atoms in a unit of volume in the sample and (-dE/dx) is the electron stopping power,calculated following the expression given by Bethe and Ashkin:
where \(Z^2_m\) is the mass square atomic number of the sample material and \(I\) is the mean excitation potential of the atom.
The incident gamma radiation produces secondary electrons inside material. These electrons can remove an atom from its lattice position through elastic scattering processes.The removed atom is known as a primary knock-in atom(PKA) and the corresponding cross section, \(\sigma_{PKA}\),is obtained starting from the McKinley-Feshbach approximation.If any of these recoil atoms has a kinetic energy above the displacement threshold energy \(T_d\),the secondary atoms can be knocked on by the PKA and removed from the lattice.The number of secondary displaced atoms can be caclulated by introducing the damage function \(v(T)\). Then,the total number of displaced atoms per target atom can be obtained by writing the displacement per atom cross section as follows:
with
For the damage function we use the Kinchin-Pease model
where
Notation | Meanings |
---|---|
\(Z_a\) | atomic number of target atom |
\(r_0\) | electron classical radius |
\(\alpha\) | \(Z_a/137\) |
\(\beta\) | ration of electron velocity to light velocity |
\(\gamma\) | \(1/(1-\beta^2)\) |
\(T_m\) | \(2E(E+2mc^2)/Mc^2\) (maximum kinetic energy of recoil atoms) |
\(E\) | electron kinetic energy |
\(mc^2\) | electron rest energy |
\(M\) | atomic mass |
In order to evaluate the damage function,the average atom recoil kinetic energy \(T_{av}\) was used obtained by the expression
MCNP¶
MCNP is a general purpose Monte Carlo N-Particle code that can be used for neutron,photon,electron,or coupled neutron/photon/electron transport,including the capbility to calculate eigenvalues for critical systems.For beginners,its tutorials are unfriendly and obscure, and for this purpose,I rewrite a simple and intelligible tutorial for beginners which can reduce time and vigour for the research works.Okay,let’s begin our tour.

Surface Cards¶
Surface can be defined by equations, points ,or macrobodies.In this simple tutorial,I only introduce the equations section which is most popurlar use in the simulation works.If your work involes other section,you can sql official documentation. Okay,You can define your surfaces with the following form and relative explainations are in the table below.
def surface_card():
/*
Form: j k a list
*/
return 0;

Cell Cards¶
If you have read the surface cards section,it is cheerful to use these surfaces to define the cell or geometry for your model. The geometry of MCNP treats an arbitrary three-dimensional configuration of user-fined materials in geometric cells bounded by first- and second-degree surfaces. The calls are defined by the intersections, unions, and complements of the regions bounded by the surfaces.Each surface divides all space into two regions,one with positive sense with respect to the surface and the other with negative sense. Define \(S=f(x,y,z)=0\) as the equation of a surface in the problem.For any set of points \((x,y,z)\),if \(S=0\) the points are on the surface. If \(S\) is negative,the points are said to have a negative sense with respect to that surface and,conversely,a positive sense if \(S\) is positive.
Data Cards¶
Some Examples¶
Fmn examples:
F4:N 10 FM4 0.04789 999 102 M999 92238.13 1
The F4 neutron tally is the track length estimate of the average fluence in cell 10.Material 999 is 238U with an atomic fraction of 100% where C=0.04787 normalization factor(such as atom/barn.cm),which you can regard it as atom density and calculate with the following formula:
Conclusions¶
Neutron Damage¶
My another research work is the neutron radiation damage in materials and mainly focus on the displacement cross section and displacement per atom(dpa) calculations.Quantification of displacement damage is an essential parameter to assess the strength and lifetime of functional materials in the ITER and upcoming fusion devices such as European Demo fusion reactor.Fusion neutrons produce different species of primary knocked on atoms(PKA) based on the different reaction channels such as (n,n’),(n,2n),(n, \(\alpha\) ), (n,p) and (n, \(\gamma\) )etc.Damage cascading of PKA and other knocked on atoms produce vacancies and interstitials,both collectively known as Frenkel pair.Displacement of atoms has adverse effects on the liftime of reactor materials and quantized with displacement per atom(dpa).Prediction of dpa values is required to accurately predict the liftime of functional materialsto be placed in fusion reactor. And my research route shows the following picture:

In order to calculate the dpa,you should calculate the two physical quantity neutron fluxes and displacement cross section where their formulas is:
where \((\frac{d\sigma}{dE})_i\) is the energy spectra of PKA or recoil atom from \(i^{th}\) reaction channel, \(\nu(T)_i\) is the number of produced Frenkel pair from the dynamics of PKA atom of T energy, \(E_n\) is the energy of incident neutron, \(E_d\) is the displacement threshold energy of target lattice, \(T_{max}\) is the maximum energy available for displacement damage and \(\Phi_n\) is the neutron flux on the material.
Basic Models and Methods¶
SPECTRA and SPECOMP_Distribution¶
SPECTER is a neutron damage calculations for materials irradiations programm.
SPECTER input is as follows(all free format):
Title (80 characters)
ITYP,ISIG,IGP,IPKA,ACNM,TIME(Free Format)
ITYP = 0 If no flux uncertainties,= 2 Otherwise
ISIG = 1 To print collapsed cross sections, =0 To stop,
= 2 For 1-group Thermal Neutron(2200 m/s)
IGP = 0 If group differential flux(n/cm2-MeV, = OUTPUT)
= 1 If true differential flux
Talys1.9¶
NJOY Nuclear Data Processing¶
The NJOY nuclear data processing system is a comprehensive computer code package for producing pointwise and multigroup nuclear cross sections and related quantities from evaluated nuclear data in the ENDF format.The U.S. Evaluated Nuclear Data Files(ENDF) have progressed through a number of versions,notably ENDF/B-III,ENDF/B-IV,ENDF/B-V,ENDF/B-VI,ENDF/B-VII.The ENDF format has also evolved through many versions.Variations of the format called “ENDF-6” were used for ENDF/B-VI and ENDF/B-VII,and will be used in ENDF/B-VIII.The ENDF format is also used in other nuclear data libraries such as the JEFF libraries in Europe and the JENDL libraries in Iapan,or in specialized libraries distributed through the Nuclear Data Section of International Atomic Energy Agency(IAEA0.These libraries represent the underlying nuclear data from a physics viewpoint,but practical calculations usually require special libraries for particle transport codes or reactor core physics codes.This is the mission of NJOY-to take the basic data from the nuclear data library and convert it into the forms needed for applications.The NJOY code consists of a set of main modules,each performing a well-defined processing task.Each of these main modules is essentially a separate computer program.They are linked to one another by input and output files.And in this simple tutorial,I will introduce some main modules which is important to radiation damage calculations.Okay,let’s start our tour.
Before our learns,you should read the examples test 1 seriously and I will regard it as a main case to our learning. And what I write here is to let you are familar with NJOY as soon as possible,in contrast to this,the theory or complete tutorials need you to read.However I believe that what this tutorial will make you love it and do your research better.

MODER¶
MODER converts ENDF “tapes” back and forth between formatted(that is,ASCII) and blocked binary modes.
#---input specifications(free) format)--------------------
**card 1** unit numbers
nin input unit
nout output unit
a positive unit is coded (mode 3)
a negative unit is blocked binary(njoy mode)
#note:abs(min) ge 1 and le 19 is a flag to select various materials fro one or more
#input tapes,with or without mode conversion,the kind of data to be processed is
#keyed to min as follows:
nin = 1,for endf or pendf input and output
2,for gendf input and output
3, for errorr-format input and output
cards 2 and 3 for abs(min) ge 1 and le 19 only
**card** 2
tpid tapeid for nout,66 characters allowed
(delimited with ',ended with /)
**card** 3
nin input unit
terminate moder by setting nin=0
matd material on this tape to add to nout
RECONR¶
RECONR reconstructs pointwise(energy-dependent) cross sections from ENDF resonance parameters and interpolation schemes.
#---input specifications (free format)-------------
**card** 1
nendf unit for endf tape
npend unit for pendf tape
**card** 2
tlabel 66 character label for new pendf tape
delimited with quotes,ended with /.
**card** 3
mat material to be reconstructed
ncards number of cards of descriptive data for new mf1(default=0)
ngrid number of user energy grid points to be added.(default=0)
**card** 4
err fractional reconstruction tolerance used when resonance-integral
error criterion(see errint) is not satisfied.
tempr reconstruction temperature(deg kelvin)(default=0)
errmax fractional reconstruction tolerance used when resonance-integral
error criterion is satisfied(errmax.ge.err,default=10*err)
errint maximum resonance-integral error(in barns) per grid point (default=err/20000)
**card** 5
cards ncards of descriptive comments for mt451 each card
delimited with quotes,ended with /.
**card** 6
enode users energy grid points
cards 3,4,5,6 must be input for each material desired
mat=0/ terminates execution of reconr
BROADR¶
BROADR Doppler-broadens and thins pointwise cross sections.
#---input specification(free format)--------
**card** 1
nendf input endf tape(for thermal nubar only)
nin input pendf tape
nout output pendf tape
**card** 2
mat1 material to broadened and thinned
ntemp2 number of final temperature(default=1)
istart restart(0 no,1 yes,default 0)
istrap bootstrap(0 no,1 yes,default 0)
temp1 starting temperature from nin(default=OK)
**card** 3
errthn fractional tolerance for thinning
thnmax max. energy for broadening and thinning(default=1 MeV)
errmax fractional tolerance used when integral criterion is satisfied
(same usage as in reconr)
(errmax.ge.errthn, default=10*errthn)
errint parameter to control integral thinning(usage as in reconr)
(default=errthn/20000) set very small to turn off integral thinning
#(A good choice for the convergence parameters errthn,errmax,and errint is
# the same set of values used in reconr)
**card** 4
temp2 final temperatures (deg Kelvin)
**card** 5
mat1 next MAT number to be processed with these parameters.Terminate with mat1=0.
HEATR¶
HEATER generative pointwise heat production cross section (neutron KERMA factors) and radiation damage production cross sections.
#---input specifications(free format)------------------
**card** 1
nendf unit for endf tape
nin unit for input pendf tape
nout unit for output pendf tape
nplot unit for graphical check output
**card** 2
matd material to be processed
npk number of partial kermax desired(default=0)
nqa number of user q values(default=0)
ntemp number of temperatures to process(default=0,meaning all on pendf)
local 0/1=gamma rays transported/deposited locally(default=0)
iprint print(0 min, 1 max, 2 deck)(default=0)
ed displacement energy for damage(default from built-in table)
**card** 3 for npk gt 0 only
mtk mt numbers for partial kermas desired total (mt301) will be provided
automatically.partical kerma for reaction mt is mt+300 and may not be
properly defined unless a gamma file for mt is on endf tape.
special values allowed--
303 non-elastic (all but mt2)
304 inelastic (mt51 thru 91)
318 fission(mt18 or mt19,20,21,38)
401 disappearance(mt102 thru 120)
442 total ev-barns
443 total kinematic kerma(high limit)
damage energy production values--
444 total
445 elastic(mt2)
446 inelastic (mt51 thru 91)
447 disappearance (mt102 thru 120)
cards and 5 for nqa gt 0 only
**card** 4
mta mt number for users q values
**card** 5
qa user specified q values (ev)
(if qa.ge.99.e6,read in variable qbar for this reaction)
**card** 5a variable qbar (for reactions with qa flag only)
qbar tab1 record giving qbar versus e (1000 words max)
GROUPR¶
GROUPR generates self-shielded multigroup cross sections,group-to-group scattering matrices,photon production matrices,and charged-particle multi-group cross sections from pointwise input.
#---input specifications(free format)-----------
**card** 1
nendf unit for endf tape
npend unit for pendf tape
ngout1 unit for input gout
ACER¶
ACER prepares libraries in ACE format for the Los Alamos continuous-energy Monte Carlo MCNP and MCNPX codes.The ACER module is supported by subsidiary modules for the different classes of the ACE format.
PLOTR and VIEWR¶
PLOTR makes plots of cross sections and perspective plots of distributions for both pointwise and multigroup data by generating input for the VIEWR module and VIEWR converts plotting files produced by the other modules into high-quality color Postscript plots.
Conclusions¶
Molecular Dynamics¶
Radiation Damage Simulations with Molecular Dynamics
Before introducing the simulation,many people may think why we need simulations,what can we do through simulations
Basic Concepts and Methods¶
Molecular dynamics simulations compute the motions of individual molecules in models of solids,liquids and gases and it is a technique that allows a numerical “thought experiment” to be carried out using a model that,to a limited extent,approximates a real physical or chemical system.Such a “virtual laboratory” approach has the advantage that many such “experiments” can be easily set up and carried out in succession by simply varying the control parameters.Morever,extreme conditions,such high temperature and obvious downside is that the results are only as good as the numerical model.In addition,the results can be artificially biased if the molecular dynamics calculation is unable to sample an adequate number of microstates over the time it is allowed to run.
Lammps Tutortials¶
LAMMPS is a classical molecular dynamics(MD) code that models ensembles of particles in a liquid,solid,or gaseous state.It can model atomic,polymeric,biological,solid-state(metals,ceramics,oxides),granular,coarse-grained,or macroscopic systems using a variety of interatomic potentials(force fields) and boundary conditions.It can model 2d or 3d systems with only a few particles up to millions or billions.
Install LAMMPS¶
In this section,I will show you the most direct method to compile and buid the parallel lammps in you Linux(Ubuntu) system. Before install the lammps,you should install the MPICH and FFTW in your PC machine.
Prerequisites¶
Note
- MPICH should be installed firstly
Before you install the mpich,you should execute the following command in order to meet some of the dependencies:
sudo apt install gfortran
sudo apt install gcc
sudo apt install g++
Then you can go to it’s source website and download the package,like mpich-3.3.1.tar.gz,when you unzip the package,read the README file,you can can implement the following commands:
tar -zxvf mpich-3.3.1.tar.gz
cd mpich-3.2
./configure --prefix=/home/day/computer/mpich3/( = you should remmember this path)
make & make install
Okay you have install the mpich in your system,next you should configure enviroment
cd ~
vim .bashrc
In the .bashrc file,add the following statements
export PATH=/home/day/computer/mpich3/bin:$PATH
export MANPATH=//home/day/computer/mpich3/man:$MANPATH
Exit the .bashrc file and execute the following command.Of course you can test your installation after execute the command with the help of mpirun -np k ./examples/cpi statement.
source .bashrc
- FFTW
Second Let’s start the FFTW installation,you can download the source code from the official website and execute the following commands
tar -zxvf fftw
cd fftw-3.3.8
./configure --prefix=/home/day/material/lammps/fftw3/ --enable-shared=yes ( = you should remmember this path)
make -j n
make install
Okay It’s done!
Obtaining and building the Lammps¶
Next download the lammps source code and execute the following commands
tar -zxvf lammps
cd lammps-22Aug18/src
cd STUBS/
make
cd ../MAKE
vim Makefile.mpi(modify details see below)
cd ..
make mpi
When you finish the commands,you can see a executable file lmp_mpi in the src/ folder,you can make a directory lammps sibling directory and copy the lmp_mpi in this folder.Remmember to add the its path to .bashrc in order to you can call lammps command in any folders.
Note
If you want install the optional packages that extend LAMMPS functionality which enable a specific set of features.You can use the following command to check relative packages:
cd lammps-22Aug18/src
make package-status
And use the following commands to install and uninstall the packages.
make yes-package
make no-package
DL_POLY Tutorials¶
DL_POLY is a package of subroutines,programs and data files,designed to facilitate molecular dynamics simulations of macromolecules, polymers,ionic systems and solutions on a distributed memory parallel computer.
Install DL_POLY¶
In this section,I will show you the simplest way to install the DL_POLY in you Linux(Ubuntu) system.What you should install gfortran and MPICH in your computer before install the DL_POLY which are introduced in the Lammps Install section. Consequently you only need to download the source code from the official website and execute the following commands
unzip dl_poly_4.09.zip
cd dl_poly_4.09
cd source
ln -s ../build/Makefile-MPI Makefile
vim Makefile (**modify the Makefile**)
make hpc
If you want a JAVA-GUI to visualize your simulation system and examine basic properties,you should execute
cd ../java
./build
Okay you have install the DL_POLY in your PC and you can see a executable file DLPOLY.Z in the execute folder.Hope you have installed successfully!
Help Turitorials¶
This section will provide some useful website and resources in the future!