1+ from libRustBCA import *
2+ import numpy as np
3+ import matplotlib .pyplot as plt
4+ import sys
5+ import os
6+ #This should allow the script to find materials and formulas from anywhere
7+ sys .path .append (os .path .dirname (__file__ )+ '/../scripts' )
8+ sys .path .append ('scripts' )
9+ from materials import *
10+ from formulas import *
11+
12+ R_T_Be = np .array (
13+ [
14+ [9.67742 , 0.016754618 ],
15+ [20.32258 , 0.09511873 ],
16+ [30.0 , 0.13153034 ],
17+ [49.677418 , 0.102242745 ],
18+ [74.83871 , 0.08403694 ],
19+ [100.32258 , 0.031002639 ],
20+ ]
21+ )
22+
23+ R_D_Be = np .array (
24+ [
25+ [10.32258 , 0.11649077 ],
26+ [20.32258 , 0.19168866 ],
27+ [30.32258 , 0.19722955 ],
28+ [49.677418 , 0.17269129 ],
29+ [75.16129 , 0.11174142 ],
30+ [100.0 , 0.044459105 ],
31+ ]
32+ )
33+
34+ R_H_Be = np .array (
35+ [
36+ [9.67742 , 0.19485489 ],
37+ [19.67742 , 0.32150397 ],
38+ [30.0 , 0.33575198 ],
39+ [50.322582 , 0.25105542 ],
40+ [75.16129 , 0.1378628 ],
41+ [100.645164 , 0.044459105 ],
42+ ]
43+ )
44+
45+ #Plotting the MD data points.
46+ energies = R_H_Be [:, 0 ]
47+ R = R_H_Be [:, 1 ]
48+ plt .plot (energies , R , marker = 'o' , linestyle = '' , label = 'H on Be MD+NNP' )
49+
50+ #Plotting the MD data points.
51+ energies = R_D_Be [:, 0 ]
52+ R = R_D_Be [:, 1 ]
53+ plt .plot (energies , R , marker = 'x' , linestyle = '' , label = 'D on Be MD+NNP' )
54+
55+ #Plotting the MD data points.
56+ energies = R_T_Be [:, 0 ]
57+ R = R_T_Be [:, 1 ]
58+ plt .plot (energies , R , marker = '^' , linestyle = '' , label = 'T on Be MD+NNP' )
59+
60+
61+ for ion in [hydrogen , deuterium , tritium ]:
62+ energies = np .logspace (1 , 2 , 25 )
63+ r_rustbca = [reflection_coefficient (ion , beryllium , energy , 0.0 , 10000 )[0 ] for energy in energies ]
64+ r_thomas = [thomas_reflection (ion , beryllium , energy ) for energy in energies ]
65+ line = plt .plot (energies , r_rustbca , label = f'{ ion ["symbol" ]} on Be RustBCA' )[0 ]
66+ plt .plot (energies , r_thomas , label = f'{ ion ["symbol" ]} on Be Thomas' , linestyle = '--' , color = line .get_color ())
67+
68+ plt .xlabel ('E [eV]' )
69+ plt .ylabel ('R_N' )
70+ plt .title ('Reflection Coefficients' )
71+ plt .legend ()
72+ plt .savefig ('reflection_md_nnp_rustbca.png' )
73+
74+ #Data digitized from https://doi.org/10.1088/1741-4326/ac592a
75+ #MD+NNP Sputtering Yields
76+ y_H_Be = np .array ([
77+ [20.004726 , 0.0029361406 ],
78+ [29.880848 , 0.013341782 ],
79+ [49.791046 , 0.024575423 ],
80+ [75.00529 , 0.016207783 ],
81+ [99.68844 , 0.017530797 ],
82+ ])
83+
84+ y_D_Be = np .array ([
85+ [19.9699 , 0.005002439 ],
86+ [29.601622 , 0.021064475 ],
87+ [49.76866 , 0.03461476 ],
88+ [74.71549 , 0.030082114 ],
89+ [99.954605 , 0.013559438 ],
90+ ])
91+
92+ y_T_Be = np .array ([
93+ [20.013433 , 0.0025699555 ],
94+ [29.85846 , 0.01879205 ],
95+ [49.756844 , 0.04147379 ],
96+ [74.707405 , 0.034042966 ],
97+ [99.68098 , 0.01965117 ],
98+ ])
99+
100+
101+ plt .figure ()
102+
103+ #Plotting the MD data points.
104+ energies = y_H_Be [:, 0 ]
105+ y = y_H_Be [:, 1 ]
106+ plt .semilogy (energies , y , marker = 'o' , linestyle = '' , label = 'H on Be MD+NNP' )
107+
108+ #Plotting the MD data points.
109+ energies = y_D_Be [:, 0 ]
110+ y = y_D_Be [:, 1 ]
111+ plt .semilogy (energies , y , marker = 'x' , linestyle = '' , label = 'D on Be MD+NNP' )
112+
113+ #Plotting the MD data points.
114+ energies = y_T_Be [:, 0 ]
115+ y = y_T_Be [:, 1 ]
116+ plt .semilogy (energies , y , marker = '^' , linestyle = '' , label = 'T on Be MD+NNP' )
117+
118+ for ion in [hydrogen , deuterium , tritium ]:
119+ energies = np .logspace (1 , 2 , 25 )
120+ y_rustbca = [sputtering_yield (ion , beryllium , energy , 0.0 , 100000 ) for energy in energies ]
121+ plt .semilogy (energies , y_rustbca , label = f'{ ion ["symbol" ]} on Be RustBCA' )
122+
123+ yamamura_yield = [yamamura (hydrogen , beryllium , energy ) for energy in energies ]
124+ plt .semilogy (energies , yamamura_yield , label = 'Yamamura H on Be' , linestyle = '--' )
125+
126+ plt .title ('Sputtering Yields of Hydrogenic Species on Be' )
127+ plt .xlabel ('E [eV]' )
128+ plt .ylabel ('Y [at/ion]' )
129+ plt .gca ().set_ylim ([1e-4 , 1e-1 ])
130+ plt .legend (loc = 'lower right' )
131+
132+ plt .savefig ('sputtering_md_nnp_rustbca.png' )
0 commit comments