Cross sections for ionization-assisted electron knock-on damage of hBN
Author
Toma Susi
Title
Cross sections for ionization-assisted electron knock-on damage of hBN
Description
Supplemental material for: T. An Bui et al., "Creation of single vacancies in hBN with electron irradiation." https://doi.org/10.1002/smll.202301926
Category
Academic Articles & Supplements
Keywords
electron irradiation, displacement cross section, hexagonal boron nitride, scanning transmission electron microscopy
URL
http://www.notebookarchive.org/2023-05-4lw178w/
DOI
https://notebookarchive.org/2023-05-4lw178w
Date Added
2023-05-10
Date Last Modified
2023-05-10
File Size
1.2 megabytes
Supplements
Rights
CC BY 4.0
Download
Open in Wolfram Cloud
Cross sections for ionization-assisted electron knock-on damage of hBN
Cross sections for ionization-assisted electron knock-on damage of hBN
Toma Susi
(*========================================================================Supplementalmaterialfor:T.AnBuietal.,"Creation of single vacancies in hBN with electron irradiation",publishedinSmall(2023).
https://doi.org/10.1002/smll.202301926
.========================================================================*)(*Settingthedirectoryofthenotebookfileastheworkingdirectory.*)setDir:=Quiet@Check[SetDirectory@DirectoryName@$InputFileName,SetDirectory@NotebookDirectory[]]In[]:=
(*CONSTANTSDEFINITION*)c=Entity["PhysicalConstant","SpeedOfLight"]["Value"][[1]];(*speedoflight*)e=Entity["PhysicalConstant","ElementaryCharge"]["Value"][[1]];(*elementarycharge*)m=Entity["PhysicalConstant","ElectronMass"]["Value"][[1]];(*massofelectron*)m1=m*c^2./e;(*electronmassineV*)ℏ=Entity["PhysicalConstant","ReducedPlanckConstant"]["Value"][[1]];(*reducedPlanckconstant*)eps=Entity["PhysicalConstant","ElectricConstant"]["Value"][[1]];(*vacuumpermittivity*)a0=Entity["PhysicalConstant","BohrRadius"]["Value"][[1]];(*Bohrradius.*)ER=Entity["PhysicalConstant","RydbergEnergy"]["Value"][[1]]/e;(*Rydbergenergy.*)NA=Entity["PhysicalConstant","AvogadroConstant"]["Value"][[1]];(*Avogadro'sconstant.*)kB=Entity["PhysicalConstant","BoltzmannConstant"]["Value"][[1]];(*Boltzmann'sconstant.*)
In[]:=
(*Re-definingconstantsdirectlytospeedupcomputationofnumericalintegrals.*)c=299792458.;e=;m=9.10938370150000079003439766`9.211300968617545*^-31;m1=510998.9499961641`;ℏ=;eps=8.85418781762039`*^-12;a0=5.291772109029999813118170167`9.519481150304392*^-11;ER=13.60569312299423437751198959631745574385`11.42021640338319;NA=602214076000000000000000.;kB=;ZB=5;(*Boronatomicnumber*)ZC=6;(*Carbonatomicnumber*)ZN=7;(*Nitrogenatomicnumber*)ZS=16;(*Sulfuratomicnumber*)
801088317
5000000000000000000000000000
132521403
400000000000000000000000000000000000000000π
1380649
100000000000000000000000000000
In[]:=
(*FunctionsforcalculatingionizationenergiesandshellpopulationsfortheBethecrosssection.*)β[U_?NumericQ]:=Sqrt[1.-1./((U*1000./m1)+1.)^2.];(*RelativisticbetaforkineticenergyinkeV.*)Ec[Z_,i_]:=QuantityMagnitude[UnitConvert[ElementData[Z]["IonizationEnergies"][[i]],"eV/mol"]/NA];(*Nthionizationenergy*)ns[Z_,s_]:=ElementData[Z]["ElectronShellConfiguration"][[s]];(*Numberofelectronsonshells.*)bethe[U_,Z_,ne_,b_,cc_,Ec_]:=Piecewise*Log(cc*(m)*)(2*(Ec*e))-Log[1-β[U]^2]-β[U]^2,(cc*m)*(2*(Ec*e))>=1,0;(*Betheionizationcrosssection.*)
28
10
2π*ne*b
4
e
(4*Pi*eps)^2(m)(Ec*e)
2
c
2
β[U]
2
c
2
β[U]
2
c
2
β[U]
In[]:=
(*========================================================================ATOMICFITTINGFORBORON========================================================================*)(*AtomicionizationcrosssectionsforB:"Total counting"valuesfromTableIIofdoi:10.1103/PhysRevA.64.052707.*)tableB={{9,0.094`},{10,0.865`},{12,1.503`},{15,2.193`},{20,3.058`},{25,3.558`},{30,3.842`},{40,4.066`},{50,4.071`},{60,3.984`},{70,3.859`},{80,3.721`},{90,3.581`},{100,3.444`},{125,3.13`},{150,2.862`},{175,2.635`},{200,2.441`},{225,2.276`},{250,2.132`},{275,2.007`},{300,1.896`},{400,1.559`},{500,1.329`},{600,1.162`},{700,1.034`},{800,0.934`},{900,0.852`},{1000,0.784`}};tableB[[All,1]]*=1/1000;(*ConverttokeV.*)tableB[[All,2]]*=1/10^20*10^28;(*Converttobarn.*)
In[]:=
(*Pre-assigningionizationenergiesandshellpopulationsforspeed.*)EcB1=Ec[ZB,1];EcB2=Ec[ZB,2];EcB3=Ec[ZB,3];neB=ns[5,2];
In[]:=
(*FittingBatomicdatawithtwoionizedstates.*)dB=5;(*Numberofpointstodropfromthestart;ensuresbetterfitofthecrucialhigh-energytail.*)nlmB2=NonlinearModelFit[Drop[tableB,dB],{bethe[U,ZB,neB,B,CC,EcB1]+bethe[U,ZB,neB,B,CC,EcB2],{{0<B<1.0},{0<CC<1.0}}},{{B,0.5},{CC,0.5}},U,Method->"NMinimize",MaxIterations1000];BoptB2=nlmB2["BestFitParameters"][[1,2]];CoptB2=nlmB2["BestFitParameters"][[2,2]];Show[{Plot[{bethe[U,ZB,neB,BoptB2,CoptB2,EcB1],bethe[U,ZB,2,BoptB2,CoptB2,EcB2],bethe[U,ZB,neB,BoptB2,CoptB2,EcB1]+bethe[U,ZB,2,BoptB2,CoptB2,EcB2]},{U,0,1},PlotRangeFull,PlotStyle{Directive[Dashed,RGBColor[0.75,0.0,0.75]],Directive[Dashed,RGBColor[0.9,0.0,0.9]],Purple},FrameTrue,FrameStyle->Thickness[0.003],FrameLabel{"Electron energy (keV)","Cross section (barn)"(*,"Boron (Bethe ="<>ToString[neB]<>", b and c constants for this shell)"*)}],ListPlot[tableB,PlotStyleBlack]},PlotRange{{0,1},{0,5*10^8}},ImageSizeLarge,EpilogInset[Column[PointLegend[{Black},{"Kim & Stone (2001)"}],LineLegend[{Purple},{"Total (1st + 2nd ionized state)"}],LineLegend[{Directive[Dashed,Purple]},"="<>ToString[EcB1]<>" eV, ="<>ToString[BoptB2]<>", ="<>ToString[CoptB2]],LineLegend[{Directive[Dashed,Lighter[Purple]]},"="<>ToString[EcB2]<>" eV"]],Scaled[{0.65,0.7}]],BaseStyle14]
n
2
ion
E
1
b
2
c
2
ion
E
2
Out[]=
Export["Bethe_boron_valence_2ionized_shellconst.pdf",%]
Out[]=
Bethe_boron_valence_2ionized_shellconst.pdf
In[]:=
(*========================================================================NITROGENATOMICFITTING========================================================================*)(*AtomicionizationcrosssectionsforN:"70-30 mix"valuesfromTableIVofdoi:10.1103/PhysRevA.66.012708.*)tableN={{14,0.0822`},{16,0.177`},{18,0.298`},{20,0.419`},{25,0.675`},{30,0.89`},{35,1.068`},{40,1.205`},{45,1.309`},{50,1.389`},{60,1.493`},{70,1.55`},{80,1.576`},{90,1.583`},{100,1.577`},{120,1.544`},{150,1.47`},{180,1.389`},{200,1.336`},{250,1.215`},{300,1.111`},{350,1.022`},{400,0.947`},{450,0.882`},{500,0.826`},{600,0.734`},{700,0.661`},{800,0.602`},{900,0.553`},{1000,0.512`},{2000,0.3`},{3000,0.216`},{4000,0.17`},{5000,0.141`}};tableN[[All,1]]*=1/1000;(*ConverttokeV.*)tableN[[All,2]]*=1/10^20*10^28;(*Converttobarn.*)
In[]:=
(*Pre-calculatingionizationenergiesandshellpopulationsforspeed.*)EcN1=Ec[ZN,1];EcN2=Ec[ZN,2];EcN3=Ec[ZN,3];neN=ns[7,2];
In[]:=
(*FittingNatomicdatawiththreeionizedstates.*)dN=0;(*Numberofpointstodropfromthestart.*)nlmN3=NonlinearModelFit[Drop[tableN,dN],{bethe[U,ZN,neN,B,CC,EcN1]+bethe[U,ZN,neN,B,CC,EcN2]+bethe[U,ZN,neN,B,CC,EcN3],{{B>0},{0<CC<1.0}}},{{B,0.5},{CC,0.5}},U,Method->"NMinimize",MaxIterations1000];BoptN3=nlmN3["BestFitParameters"][[1,2]];CoptN3=nlmN3["BestFitParameters"][[2,2]];Show{Plot[{bethe[U,ZN,neN,BoptN3,CoptN3,EcN1],bethe[U,ZN,neN,BoptN3,CoptN3,EcN2],bethe[U,ZN,neN,BoptN3,CoptN3,EcN3],bethe[U,ZN,neN,BoptN3,CoptN3,EcN1]+bethe[U,ZN,neN,BoptN3,CoptN3,EcN2]+bethe[U,ZN,neN,BoptN3,CoptN3,EcN3]},{U,0,5},PlotRangeFull,PlotStyle{Directive[Dashed,RGBColor[0,0.5,1]],Directive[Dashed,RGBColor[0,0.7,1]],Directive[Dashed,RGBColor[0,0.9,1]],Darker[Blue]},FrameTrue,BaseStyle->13,FrameStyle->Thickness[0.003],FrameLabel{"Electron energy (keV)","Cross section (barn)"(*,"Nitrogen (Bethe ="<>ToString[ns[7,2]]<>", b and c constant for this shell)"*)}],ListPlot[tableN,PlotStyleBlack]},PlotRange{{0,5},{0,2*10^8}},ImageSizeLarge,EpilogInsetColumnPointLegend[{Black},{"Kim & Desclaux (2002)"}],LineLegend[{Darker[Blue]},{"Total (1st + 2nd + 3rd ionized state)"}],LineLegend[{Directive[Dashed,RGBColor[0,0.5,1]]},"="<>ToString[EcN1]<>" eV, ="<>ToString[BoptN3]<>", ="<>ToString[CoptN3]],LineLegend[{Directive[Dashed,RGBColor[0,0.7,1]]},"="<>ToString[EcN2]<>" eV"],LineLegend{Directive[Dashed,RGBColor[0,0.9,1]]},"="<>ToString[EcN3]<>" eV",Scaled[{0.6,0.65}],BaseStyle14
n
2
ion
E
1
b
2
c
2
ion
E
2
ion
E
3
Out[]=
Export["Bethe_nitrogen_valence_3ionized_shellconst.pdf",%]
Out[]=
Bethe_nitrogen_valence_3ionized_shellconst.pdf
In[]:=
(*NvsBatomiccrosssectionsextrapolatedtoourexperimentalenergies.*)Plot[{bethe[U,ZB,neB,BoptB2,CoptB2,EcB1]+bethe[U,ZB,neB,BoptB2,CoptB2,EcB2],bethe[U,ZB,neB,BoptB2,CoptB2,EcB1],bethe[U,ZB,neB,BoptB2,CoptB2,EcB2],bethe[U,ZN,neN,BoptN3,CoptN3,EcN1]+bethe[U,ZN,neN,BoptN3,CoptN3,EcN2]+bethe[U,ZN,neN,BoptN3,CoptN3,EcN3],bethe[U,ZN,neN,BoptN3,CoptN3,EcN1],bethe[U,ZN,neN,BoptN3,CoptN3,EcN2],bethe[U,ZN,neN,BoptN3,CoptN3,EcN3]},{U,35,105},PlotStyle{Purple,Directive[Dashed,RGBColor[0.75,0.0,0.75]],Directive[Dashed,RGBColor[0.9,0.0,0.9]],Darker[Blue],Directive[Dashed,RGBColor[0,0.5,1]],Directive[Dashed,RGBColor[0,0.7,1]],Directive[Dashed,RGBColor[0,0.9,1]]},FrameTrue,FrameLabel{"Electron energy (keV)","Ionization cross section (barn)"},ImageSizeLarge,BaseStyle13,FrameStyle->Thickness[0.003],PlotRange->{{40,100},{0,5*10^6}},PlotLegends->Placed[{"B (total)","B (1st ionized)","B (2nd ionized)","N (total)","N (1st ionized)","N (2nd ionized)","N (3rd ionized)"},{0.67,0.74}]]
Out[]=
In[]:=
Export["Bethe_NvsB_valence_all.pdf",%]
Out[]=
Bethe_NvsB_valence_all.pdf
In[]:=
(*========================================================================CROSSSECTIONSFORHBNFROMHEREONWARDS========================================================================*)(*AdditionalformulasforelasticMcKinley-FeshbachcrosssectionwithvibrationalenhancementofenergytransferbasedonDOI:10.1038/ncomms13040.*)rho=3.683245906802317*10^19;(*ArealdensityofhBN*)MB=1.79504374`*^-26;(*MassofBinamu.*)MN=2.325918378`*^-26;(*MassofNinamu.*)M[Z_]:=Switch[Z,5,MB,7,MN];(*FunctiontoreturnmassdependingonatomicnumberZ.*)velo[v_,VRMS_]:=1./Sqrt[2.*Pi*VRMS]*Exp[-v^2./(2.*VRMS)];X[Z_,v_]:=0.5*M[Z]*v^2.;(*ConvertedtokeVinsteadofeVforconvenience.*)r[Z_,v_,U_]:=(1./c)*Sqrt[(e*U*1000.)^2.+2.*e*U*1000.*m*c^2.]+M[Z]*v;t[Z_,v_,U_]:=Sqrt[(e*U*1000.+X[Z,v])*(e*U*1000.+2.*m*c^2.+X[Z,v])];TV[Z_,v_,U_]:=r[Z,v,U]*(r[Z,v,U]+2.*t[Z,v,U]/c)/(2.*M[Z])+t[Z,v,U]^2./(2.*M[Z]*c^2.);σiso[Z_,v_?NumericQ,U_?NumericQ,Thres_]:=4.Pi-1.-Log+Piβ[U]2.-Log-2.;csiso[Z_,U_,Thres_,Vfit_]:=NIntegrate[Boole[TV[Z,v,U]/e≥Thres]*σiso[Z,v,U,Thres]*velo[v,Vfit]*10.^28.,{v,-vlim,vlim},Method->{"GlobalAdaptive"},MinRecursion3,PrecisionGoal2];
2
Z
2
e
4.πeps*2.m*
2
c
(1.-)
2.
β[U]
4.
β[U]
TV[Z,v,U]/e
Thres
2
β[U]
TV[Z,v,U]/e
Thres
Z
2
e
cℏ
TV[Z,v,U]/e
Thres
TV[Z,v,U]/e
Thres
In[]:=
(*Calculatingthemean-squaredout-of-planevelocityfromthephonondispersion.*)filepathN=ToString[setDir]<>"/hBN_zDOS_N.csv";filepathB=ToString[setDir]<>"/hBN_zDOS_B.csv";zdosN=Import[filepathN,HeaderLines->1];(*Z-componentofthephononDOSforNincludedassupplementaryfile;updatepathtomatchyourfilesystem.*)zdosB=Import[filepathB,HeaderLines->1];(*Z-componentofthephononDOSforNincludedassupplementaryfile;updatepathtomatchyourfilesystem.*)gN=Interpolation[zdosN,InterpolationOrder1];gB=Interpolation[zdosB,InterpolationOrder1];VN[T_]:=2**NIntegrate0.5+*x*gN[x*(kB*T)/ℏ],{x,ℏ*zdosN[[1,1]]/(kB*T),ℏ*zdosN[[-1,1]]/(kB*T)};VB[T_]:=2**NIntegrate0.5+*x*gB[x*(kB*T)/ℏ],{x,ℏ*zdosB[[1,1]]/(kB*T),ℏ*zdosB[[-1,1]]/(kB*T)};vrmB=Sqrt[Quiet[VB[300]]];(*Root-mean-squarevelocityforout-of-planevibrationofBinhBN.*)vrmN=Sqrt[Quiet[VN[300]]];(*Root-mean-squarevelocityforout-of-planevibrationofNinhBN.*)vlim=vrmB*6.;(*Integrationlimitforvelocitiestocovertheirvariation.*)
(kB*T)^2
2MN*ℏ
1
Exp[x]-1
(kB*T)^2
2MB*ℏ
1
Exp[x]-1
In[]:=
(*ExperimentalcrosssectionsforBejectionfromhBNwithuncertaintylimits.*)Bwitherr={{50,Around[0.03068534`,{0.003776,0.004003}]},{55,Around[0.04324756`,{0.010491`,0.011109}]},{60,Around[0.02340882`,{0.005737`,0.006041`}]},{65,Around[0.02647839`,{0.010512`,0.010870`}]},{70,Around[0.04503451`,{0.016962`,0.018028`}]},{75,Around[0.02483343`,{0.004831,0.005043}]},{80,Around[0.09476548`,{0.018229,0.024104}]},{85,Around[0.22598338`,{0.049477`,0.068482}]},{90,Around[0.32592229`,{0.10527982`,0.1497560}]}};(*ExperimentalcrosssectionsforNejectionfromhBNwithuncertaintylimits.*)Nwitherr={{50,Around[0.01716163`,{0.000459,0.000539}]},{55,Around[0.01593338`,{0.000404`,0.000505}]},{60,Around[0.01835762`,{0.000402`,0.000503`}]},{65,Around[0.02554704`,{0.001601`,0.001774`}]},{70,Around[0.03212862`,{0.001215`,0.001511`}]},{75,Around[0.02120466`,{0.001343,0.001486}]},{80,Around[0.08379026`,{0.009860,0.011957}]},{85,Around[0.16173439`,{0.030753`,0.040135}]},{90,Around[0.35773825`,{0.153565`,0.208970}]}};
In[]:=
(*Omit75keVpointforfitting;restoreallvaluesforplotting.*)(*Separatevariablesforvaluesanduncertaintylimitsforeaseofuseinfittingbelow.*)Bwitherrfit=Drop[Bwitherr,{6}];Bcsfit=Bwitherrfit[[All,2]][[All,1]];Basymmerrfit=Bwitherrfit[[All,2]][[All,2]][[All,1]];Nwitherrfit=Drop[Nwitherr,{6}];Ncsfit=Nwitherrfit[[All,2]][[All,1]];Nasymmerrfit=Nwitherrfit[[All,2]][[All,2]][[All,1]];
In[]:=
(*Optimalfitvalues;thesewerefoundviaextensivesearch,usedbelowasabasisforexamplesofhowthetotalerrorisoptimizedwithrespecttoacoupleofparametersatatime.*)TgsB=21.5;Tex1B=16.45;Tex2B=6.0;eta1B=2.95;eta2B=0.09;
In[]:=
(*Calcutinguncertainty-weightedmean-squaredsumoferrorbetweenexperimentandtheory.Duetocomplicatednumericalnatureofthefittingfunctionandthecomputationalcostofcalculatingvalues,needtobeiterativelyoptimizedbyvaryingthedifferentfittingparametersoneafteranother.WARNING:thesewilltakeseveralminutestoevaluateevenon10parallelkernels.Firstoptimizethesecondexcitedstate.*)Clear[Tex2B,eta2B](*Clearsthevaluesoftheparametersthataretobeoptimized.*)errorsB2=Table[{Tex2B,eta2B,Total[(Flatten[ParallelTable[{csiso[ZB,U,TgsB,vrmB^2.]+10^(-28)*rho*(eta1B*csiso[ZB,U,Tex1B,vrmB^2.]*bethe[U,ZB,neB,BoptB2,CoptB2,EcB1]+eta2B*csiso[ZB,U,Tex2B,vrmB^2.]*bethe[U,ZB,neB,BoptB2,CoptB2,EcB2])},{U,{50,55,60,65,70,80,85,90}}]]-Bcsfit)^2./Basymmerrfit]},{Tex2B,5.9,6.1,0.05},{eta2B,0.08,0.10,0.01}];posB2=Position[errorsB2,Min[errorsB2[[All,All,3]]]];(*Findsthepositionoftheerror-minimizingparametersinthelist.*)errorsB2[[posB2[[1]][[1]],posB2[[1]][[2]]]](*Printstheoptimalparameters.*)(*Setthevaluestotheoptimalones.*)Tex2B=errorsB2[[posB2[[1]][[1]],posB2[[1]][[2]]]][[1]];eta2B=errorsB2[[posB2[[1]][[1]],posB2[[1]][[2]]]][[2]];
Out[]=
{6.05,0.09,0.09133}
In[]:=
(*Thenoptimizetheground-statethresholdenergyandthefirstexcitedstate.*)Clear[TgsB,Tex1B,eta1B]errorsB=Table[{TgsB,Tex1B,eta1B,Total[(Flatten[ParallelTable[{csiso[ZB,U,TgsB,vrmB^2.]+10^(-28)*rho*(eta1B*csiso[ZB,U,Tex1B,vrmB^2.]*bethe[U,ZB,neB,BoptB2,CoptB2,EcB1]+eta2B*csiso[ZB,U,Tex2B,vrmB^2.]*bethe[U,ZB,neB,BoptB2,CoptB2,EcB2])},{U,{50,55,60,65,70,80,85,90}}]]-Bcsfit)^2./Basymmerrfit]},{TgsB,21.4,21.6,0.05},{Tex1B,16.4,16.5,0.05},{eta1B,2.93,2.95,0.01}];posB=Position[errorsB,Min[errorsB[[All,All,3]]]];errorsB[[posB[[1]][[1]],posB[[1]][[2]],posB[[1]][[3]]]]TgsB=errorsB[[posB[[1]][[1]],posB[[1]][[2]],posB[[1]][[3]]]][[1]];Tex1B=errorsB[[posB[[1]][[1]],posB[[1]][[2]],posB[[1]][[3]]]][[2]];eta1B=errorsB[[posB[[1]][[1]],posB[[1]][[2]],posB[[1]][[3]]]][[3]];
Out[]=
{21.5,16.45,2.95,0.09133}
In[]:=
(*PlottingtheoreticalcrosssectionvaluesforBwithoptimalfitparameters;pre-calculatingvaluesfirst.*)Ene=Table[U,{U,40,100,1}];csB=ParallelTable[{csiso[ZB,U,TgsB,vrmB^2],csiso[ZB,U,Tex1B,vrmB^2]*bethe[U,ZB,neB,BoptB2,CoptB2,EcB1]*10^(-28)*rho,csiso[ZB,U,Tex2B,vrmB^2]*bethe[U,ZB,neB,BoptB2,CoptB2,EcB2]*10^(-28)*rho},{U,40,100,1}];gs0B=Partition[Riffle[Ene,csB[[All,1]]],2];ex1B=Partition[Riffle[Ene,eta1B*csB[[All,2]]],2];ex2B=Partition[Riffle[Ene,eta2B*csB[[All,3]]],2];totB=Partition[Riffle[Ene,csB[[All,1]]+eta1B*csB[[All,2]]+eta2B*csB[[All,3]]],2];Show[{ListLinePlot[{gs0B,ex1B,ex2B,totB},PlotRange->{{40,100},{0,0.6}},InterpolationOrder->2,BaseStyle->13,ImageSize->Large,FrameStyle->Thickness[0.003],Joined->True,Filling->None,Frame->True,(*PlotLabel->"Boron",*)FrameLabel->{"Primary beam energy (keV)","Cross section (barn)"},PlotLegends->Placed[{"Ground-state (="<>ToString[TgsB]<>" eV)","1st ionized (="<>ToString[Tex1B]<>" eV, ="<>ToString[eta1B]<>")","2nd ionized (="<>ToString[NumberForm[Tex2B,{1,1}]]<>" eV, ="<>ToString[eta2B]<>")","Total"},{0.35,0.6}],PlotStyle->{Directive[DotDashed,Purple],Directive[Dashed,RGBColor[0.75,0.0,0.75]],Directive[Dashed,RGBColor[0.9,0.0,0.9]],Purple}],ListPlot[Bwitherr,IntervalMarkers"Fences",PlotStyle->{Thickness[0.002],Purple},PlotLegends->Placed[{"Experiment"},{0.236,0.81}]]}]
gs
T
ex
T
1
ξ
1
ex
T
2
ξ
2
Out[]=
In[]:=
Export["boron_cs_twoexcited_revised.pdf",%]
Out[]=
boron_cs_twoexcited_revised.pdf
In[]:=
(*========================================================================NITROGENFROMHEREONWARDS========================================================================*)(*Optimalfitvalues;thesewerefoundviaextensivesearch,usedbelowasabasisforexamplesofhowthetotalerrorisoptimizedwithrespecttoacoupleofparametersatatime.*)TgsN=19.0;Tex1N=14.6;Tex2N=10.55;Tex3N=5.2;eta1N=17.5;eta2N=0.63;eta3N=0.07;
In[]:=
(*Calcutinguncertainty-weightedmean-squaredsumoferrorbetweenexperimentandtheory.Duetocomplicatednumericalnatureofthefittingfunctionandthecomputationalcostofcalculatingvalues,needtobeiterativelyoptimizedbyvaryingthedifferentfittingparametersoneafteranother.WARNING:willtakeseveralminutestoevaluateevenon10parallelkernels.Firstoptimizethethirdexcitedstate.*)Clear[Tex3N,eta3N]errorsN3=Table[{Tex3N,eta3N,Total[(Flatten[ParallelTable[{csiso[ZN,U,TgsN,vrmN^2.]+10^(-28)*rho*(eta1N*csiso[ZN,U,Tex1N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN1]+eta2N*csiso[ZN,U,Tex2N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN2]+eta3N*csiso[ZN,U,Tex3N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN3])},{U,{50,55,60,65,70,80,85,90}}]]-Ncsfit)^2./Nasymmerrfit]},{Tex3N,5.1,5.3,0.05},{eta3N,0.06,0.08,0.01}];posN3=Position[errorsN3,Min[errorsN3[[All,All,3]]]];errorsN3[[posN3[[1]][[1]],posN3[[1]][[2]]]]Tex3N=errorsN3[[posN3[[1]][[1]],posN3[[1]][[2]]]][[1]];eta3N=errorsN3[[posN3[[1]][[1]],posN3[[1]][[2]]]][[2]];
Out[]=
{5.2,0.07,0.0160077}
In[]:=
(*Thenoptimizethesecondexcitedstate.*)Clear[Tex2N,eta2N]errorsN2=Table[{Tex2N,eta2N,Total[(Flatten[ParallelTable[{csiso[ZN,U,TgsN,vrmN^2.]+10^(-28)*rho*(eta1N*csiso[ZN,U,Tex1N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN1]+eta2N*csiso[ZN,U,Tex2N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN2]+eta3N*csiso[ZN,U,Tex3N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN3])},{U,{50,55,60,65,70,80,85,90}}]]-Ncsfit)^2./Nasymmerrfit]},{Tex2N,10.5,10.6,0.05},{eta2N,0.62,0.64,0.01}];posN2=Position[errorsN2,Min[errorsN2[[All,All,3]]]];errorsN2[[posN2[[1]][[1]],posN2[[1]][[2]]]]Tex2N=errorsN2[[posN2[[1]][[1]],posN2[[1]][[2]]]][[1]];eta2N=errorsN2[[posN2[[1]][[1]],posN2[[1]][[2]]]][[2]];
Out[]=
{10.55,0.64,0.0159147}
In[]:=
(*Andfinallythefirstexcitedstate.Notethatanyground-statethresholdenergyhigherthanaverylowvalueof19eVwillnotaffectthefittingatall,andthusisleftout.*)Clear[Tex1N,eta1N]errorsN=Table[{Tex1N,eta1N,Total[(Flatten[ParallelTable[{csiso[ZN,U,TgsN,vrmN^2.]+10^(-28)*rho*(eta1N*csiso[ZN,U,Tex1N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN1]+eta2N*csiso[ZN,U,Tex2N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN2]+eta3N*csiso[ZN,U,Tex3N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN3])},{U,{50,50,60,65,70,80,85,90}}]]-Ncsfit)^2./Nasymmerrfit]},{Tex1N,14.5,14.7,0.05},{eta1N,17.4,17.6,0.1}];posN=Position[errorsN,Min[errorsN[[All,All,3]]]];errorsN[[posN[[1]][[1]],posN[[1]][[2]]]]Tex1N=errorsN[[posN[[1]][[1]],posN[[1]][[2]]]][[1]];eta1N=errorsN[[posN[[1]][[1]],posN[[1]][[2]]]][[2]];
Out[]=
{14.6,17.4,0.0153857}
In[]:=
(*PlottingtheoreticalcrosssectionvaluesforNwithoptimalfitparameters;pre-calculatingvaluesfirst.*)Ene=Table[U,{U,40,100,1}];crossN=ParallelTable[{csiso[ZN,U,TgsN,vrmN^2],csiso[ZN,U,Tex1N,vrmN^2]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN1]*10^(-28)*rho,csiso[ZN,U,Tex2N,vrmN^2]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN2]*10^(-28)*rho,csiso[ZN,U,Tex3N,vrmN^2]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN3]*10^(-28)*rho},{U,40,100,1}];gs0N=Partition[Riffle[Ene,crossN[[All,1]]],2];ex1N=Partition[Riffle[Ene,eta1N*crossN[[All,2]]],2];ex2N=Partition[Riffle[Ene,eta2N*crossN[[All,3]]],2];ex3N=Partition[Riffle[Ene,eta3N*crossN[[All,4]]],2];totN=Partition[Riffle[Ene,crossN[[All,1]]+eta1N*crossN[[All,2]]+eta2N*crossN[[All,3]]+eta3N*crossN[[All,4]]],2];Show[ListLinePlot[{gs0N,ex1N,ex2N,ex3N,totN},PlotRange->{{40,100},{0,0.6}},BaseStyle->13,InterpolationOrder->2,ImageSize->Large,FrameStyle->Thickness[0.003],Joined->True,Filling->None,Frame->True,(*PlotLabel->"Nitrogen",*)FrameLabel->{"Primary beam energy (keV)","Cross section (barn)"},PlotLegends->Placed["Ground-state (>"<>ToString[NumberForm[TgsN,{2,1}]]<>" eV)","1st ionized (="<>ToString[Tex1N]<>" eV, ="<>ToString[eta1N]<>")","2nd ionized (="<>ToString[Tex2N]<>" eV, ="<>ToString[eta2N]<>")","3rd ionized ="<>ToString[Tex3N]<>" eV, ="<>ToString[eta3N]<>")","Total",{0.35,0.55}],PlotStyle->{Directive[DotDashed,Darker[Blue]],Directive[Dashed,RGBColor[0,0.5,1]],Directive[Dashed,RGBColor[0,0.7,1]],Directive[Dashed,RGBColor[0,0.9,1]],Darker[Blue]}],ListPlot[Nwitherr,IntervalMarkers"Fences",PlotStyle->{Thickness[0.002],Darker[Blue]},PlotLegends->Placed[{"Experiment"},{0.23,0.8}]]]
gs
T
ex
T
1
ξ
1
ex
T
1
ξ
2
ex
T
3
ξ
3
Out[]=
In[]:=
Export["nitrogen_cs_threeexcited_revised.pdf",%]
Out[]=
nitrogen_cs_threeexcited_revised.pdf
In[]:=
(*========================================================================RELATIVEUNCERTAINTIESANDSTATISTICSFORBANDNDATA.========================================================================*)(*Calculatingrelativeuncertaintiesfortheexperimentaldata.*)scaledNwitherr=Partition[Riffle[Nwitherr[[All,1]],Nwitherr[[All,2]]/Nwitherr[[All,2]][[All,1]]],2];scaledBwitherr=Partition[Riffle[Bwitherr[[All,1]],Bwitherr[[All,2]]/Bwitherr[[All,2]][[All,1]]+1],2];normerrors=ListPlot[{scaledBwitherr,scaledNwitherr},IntervalMarkers"Fences",PlotStyle->{Directive[Thickness[0.002],Purple],Directive[Thickness[0.002],Darker[Blue]]},PlotRange->{{47,93},{0,3}},BaseStyle->13,ImageSize->600,FrameStyle->Thickness[0.003],Filling->None,Frame->True,(*PlotLabel->"Nitrogen",*)FrameLabel->{None,"Normalized cross section (a.u.)"},FrameTicks->{Automatic,None},PlotLegends->Placed[{"B (normalized & offset)","N (normalized)"},{0.2,0.88}],Epilog->{Dashed,Line[{{40,1.0},{100,1.0}}],{Dashed,Line[{{40,2.0},{100,2.0}}]}}];(*Numberofvalidseriesforeachelementateachvoltage.*)statsN={{50,111},{55,117},{60,167},{65,35},{70,75},{75,53},{80,52},{85,31},{90,32}};statsB={{50,44},{55,31},{60,41},{65,9},{70,16},{75,24},{80,68},{85,54},{90,45}};stats=ListPlot[{statsB,statsN},PlotStyle->{Directive[Thickness[0.003],Purple],Directive[Thickness[0.003],Darker[Blue]]},PlotRange->{{47,93},{0,180}},BaseStyle->13,ImageSize->632,FrameStyle->Thickness[0.003],Frame->True,Joined->True,Mesh->All,(*PlotLabel->"Nitrogen",*)FrameLabel->{"Primary beam energy (keV)","Number of valid data series"},PlotLegends->Placed[{"B","N"},{0.08,0.87}]];Show[GraphicsColumn[{normerrors,stats}],ImageSize->Large]
Out[]=
In[]:=
Export["nitrogenboron_cs_errors_stats.pdf",%]
Out[]=
nitrogenboron_cs_errors_stats.pdf
In[]:=
(*========================================================================EXTRAPOLATIONOFTOTALCROSSSECTIONSFORhBNDAMAGETOHIGHERENERGIES.========================================================================*)(*Estimatingtherelativechangeintotaluncertainty-weightederrorduetochangeofground-statethresholdenergywithrespecttooptimalfitvaluesforB.*)Clear[TgsB]errorsrelB=Table[{TgsB,Total[(Flatten[ParallelTable[{csiso[ZB,U,TgsB,vrmB^2.]+10^(-28)*rho*(eta1B*csiso[ZB,U,Tex1B,vrmB^2.]*bethe[U,ZB,neB,BoptB2,CoptB2,EcB1]+eta2B*csiso[ZB,U,Tex2B,vrmB^2.]*bethe[U,ZB,neB,BoptB2,CoptB2,EcB2])},{U,{50,55,60,65,70,80,85,90}}]]-Bcsfit)^2./Basymmerrfit]},{TgsB,21.2,21.8,0.05}];(*Printsouttherelativeerrors;valuesforplottingbelowarechosentobethosethatarewithin25%oftheminimum.*)errorsrelB[[All,2]]/Min[errorsrelB[[All,2]]]
Out[]=
{1.75933,1.50219,1.30914,1.17046,1.07757,1.02295,1.,1.00301,1.027,1.06769,1.12139,1.18497,1.25574}
In[]:=
(*Estimatingtherelativechangeintotaluncertainty-weightederrorduetochangeofground-statethresholdenergywithrespecttooptimalfitvaluesforN.*)Clear[TgsN]errorsrelN=Table[{TgsN,Total[(Flatten[ParallelTable[{csiso[ZN,U,TgsN,vrmN^2.]+10^(-28)*rho*(eta1N*csiso[ZN,U,Tex1N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN1]+eta2N*csiso[ZN,U,Tex2N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN2]+eta3N*csiso[ZN,U,Tex3N,vrmN^2.]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN3])},{U,{50,55,60,65,70,80,85,90}}]]-Ncsfit)^2./Nasymmerrfit]},{TgsN,17.3,21.0,0.1}];(*Printsouttherelativeerrors:notethatforNthevaluestendto1sincehighervaluesofthethresholdenergydonotaffectthefit.Valuesforplottingbelowarechosentobethosethatarewithin25%oftheminimumfromthebottomand22.25eVfromDFTfromabove.*)errorsrelN[[All,2]]/errorsrelN[[-1,2]]
Out[]=
{2.19007,1.65663,1.34083,1.15963,1.0599,1.00827,0.984151,0.975129,0.973897,0.976292,0.980038,0.983963,0.987526,0.990524,0.992934,0.994808,0.996233,0.997296,0.998078,0.998646,0.999053,0.999343,0.999548,0.999691,0.99979,0.999859,0.999906,0.999937,0.999959,0.999973,0.999983,0.999989,0.999993,0.999996,0.999998,0.999999,1.,1.}
In[]:=
(*PlottingextrapolatedtotalcrosssectionsforBandNwiththeabove-estimateduncertainties.*)TgsN=19.0;TgsNlow=17.6;(*LowestthresholdenergyvalueforNwithin5%oftotalweightederrorsumcomparedtooptimalfit.*)TgsNhigh=22.25;(*UpperlimitforNthresholdbasedonDFT/MD.*)TgsB=21.5;TgsBhigh=21.8;(*HighestthresholdenergyvalueforBwithin5%oftotalweightederrorsumcomparedtooptimalfit.*)TgsBlow=21.35;(*LowestthresholdenergyvalueforBwithin5%oftotalweightederrorsumcomparedtooptimalfit.*)Ene2=Table[U,{U,60,305,5}];csextB=ParallelTable[{csiso[ZB,U,TgsBlow,vrmB^2],csiso[ZB,U,TgsB,vrmB^2],csiso[ZB,U,TgsBhigh,vrmB^2],csiso[ZB,U,Tex1B,vrmB^2]*bethe[U,ZB,neB,BoptB2,CoptB2,EcB1]*10^(-28)*rho,csiso[ZB,U,Tex2B,vrmB^2]*bethe[U,ZB,neB,BoptB2,CoptB2,EcB2]*10^(-28)*rho},{U,60,305,5}];totextBlow=Partition[Riffle[Ene2,csextB[[All,1]]+eta1B*csextB[[All,4]]+eta2B*csextB[[All,5]]],2];totextB=Partition[Riffle[Ene2,csextB[[All,2]]+eta1B*csextB[[All,4]]+eta2B*csextB[[All,5]]],2];totextBhigh=Partition[Riffle[Ene2,csextB[[All,3]]+eta1B*csextB[[All,4]]+eta2B*csextB[[All,5]]],2];csextN=ParallelTable[{csiso[ZN,U,TgsNlow,vrmN^2],csiso[ZN,U,TgsN,vrmN^2],csiso[ZN,U,TgsNhigh,vrmN^2],csiso[ZN,U,Tex1N,vrmN^2]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN1]*10^(-28)*rho,csiso[ZN,U,Tex2N,vrmN^2]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN2]*10^(-28)*rho,csiso[ZN,U,Tex3N,vrmN^2]*bethe[U,ZN,neN,BoptN3,CoptN3,EcN3]*10^(-28)*rho},{U,60,305,5}];totextNlow=Partition[Riffle[Ene2,csextN[[All,1]]+eta1N*csextN[[All,4]]+eta2N*csextN[[All,5]]+eta3N*csextN[[All,6]]],2];totextN=Partition[Riffle[Ene2,csextN[[All,2]]+eta1N*csextN[[All,4]]+eta2N*csextN[[All,5]]+eta3N*csextN[[All,6]]],2];totextNhigh=Partition[Riffle[Ene2,csextN[[All,3]]+eta1N*csextN[[All,4]]+eta2N*csextN[[All,5]]+eta3N*csextN[[All,6]]],2];ListLinePlot{totextNlow,totextN,totextNhigh,totextBlow,totextB,totextBhigh},InterpolationOrder->2,Filling{{3{{1},Directive[Blue,Opacity[0.1]]}},6{{4},Directive[Purple,Opacity[0.1]]}},PlotStyle{None,Darker[Blue],None,None,Purple,None},FrameTrue,FrameLabel{"Electron energy (keV)","Extrapolated total displacement cross section (barn)"},ImageSizeLarge,BaseStyle13,FrameStyle->Thickness[0.003],PlotRange->{{45,305},{0,25}},PlotLegends->PlacedNone,"N = eV",None,None,"B = eV",None,{0.2,0.85}
gs
T
+3.25
19.0
-1.40
gs
T
+0.30
21.5
-0.15
Out[]=
In[]:=
Export["nitrogen_boron_higher_energy_revised.pdf",%]
Out[]=
nitrogen_boron_higher_energy_revised.pdf
In[]:=
(*========================================================================CDFTTHRESHOLDENERGIESPLOT.========================================================================*)hole={0.0,0.25,0.5,0.75,1.0};TNcDFT={22.25,21.1,17.9,12.9,9.15};TBcDFT={20.15,18.7,14.3,9.5,5.05};NcDFT=Partition[Riffle[hole,TNcDFT],2];BcDFT=Partition[Riffle[hole,TBcDFT],2];ListLinePlot[{NcDFT,BcDFT},PlotStyle->{Directive[Dashed,Darker[Blue]],Directive[Dashed,Purple]},PlotLegends->Placed[{Style["Nitrogen",18],Style["Boron",18]},{0.8,0.85}],BaseStyle->18,ImageSize->Large,AspectRatio->0.8,PlotMarkers->Automatic,FrameStyle->Thickness[0.003],Frame->True,FrameTicks->{{All,None},{{0,0.25,0.5,0.75,1.0},None}},FrameLabel->{"Constrained charge (e)","Displacement threshold energy (eV)"}]
Out[]=
In[]:=
Export["nitrogen_boron_cdft_thresholds.pdf",%]
Out[]=
nitrogen_boron_cdft_thresholds.pdf
Cite this as: Toma Susi, "Cross sections for ionization-assisted electron knock-on damage of hBN" from the Notebook Archive (2023), https://notebookarchive.org/2023-05-4lw178w
Download
