Commit e15e34d606a7f455b27f44cb45ce7eb6d94c2b15

Authored by Jalabert
1 parent b6361456
Exists in master

ignored

Showing 1 changed file with 0 additions and 345 deletions   Show diff stats
test_plot.py deleted
... ... @@ -1,345 +0,0 @@
1   -# -*- coding: utf-8 -*-
2   -"""
3   -Created on Sun Jul 17 14:31:57 2022
4   -
5   -@author: frede
6   -"""
7   -
8   -import main_program as mp
9   -from astropy.table import Table
10   -import four_levels_model as four_lvl
11   -import matplotlib.pyplot as plt
12   -import numpy as np
13   -
14   -'''================================ initialization of parameters ================================'''
15   -
16   -filename = './converted/kp00_30000.txt'
17   -star_radius = 5
18   -t_gas = 500
19   -n_e = 1.6
20   -n_c = 54
21   -parsec = 6.5e-8
22   -fc_pah = 0.1
23   -ISRF = False
24   -
25   -# test3 = mp.HeatingGas(filename, star_radius, t_gas, n_e, n_c, parsec, fc_pah, ISRF)
26   -# test3.parameters()
27   -# test4 = mp.HeatingGas(filename, star_radius, t_gas, n_e, n_c, parsec, fc_pah, ISRF)
28   -# test4.parameters()
29   -
30   -''' parameter '''
31   -distance_list = np.array([parsec])
32   -for j in range(2,1002,1):
33   - d = distance_list[0]
34   - distance_list = np.append(distance_list, j * d) #in parsec
35   -
36   -test = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, n_c, parsec, fc_pah, ISRF)
37   -test.parameters()
38   -
39   -# calculus=0
40   -# for i in distance_list:
41   -# test = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, n_c, i, fc_pah, ISRF)
42   -# test.parameters()
43   -# calculus += test.total_gas_heating*1e-7/(test.g_0 * test.heating_efficiency * 0.1)
44   -
45   -# print('value = ',calculus/len(distance_list))
46   -
47   -'''================================ plots ================================'''
48   -
49   -# ''' intensity per wavelength '''
50   -# plt.figure()
51   -# plt.xlim([0,1500])
52   -# plt.plot(test.wavelength,test.wavelength_intensity*1e-3,'k')
53   -# plt.ylabel('intensity spectrum (W m$^{-2}$ sr$^{-1}$ nm$^{-1}$)')
54   -# plt.xlabel('wavelength (nm)')
55   -
56   -''' average of 3 sizes categories of pah of photoabsorption sigma
57   -and photoionization sigma_ion cross sections per C atom of pahs
58   -as functions of photon energy in eV'''
59   -
60   -''' conversion from mb/C atom to cmยฒ/C atom '''
61   -test.pah_cross_a = test.pah_cross_a*1e-18
62   -test.pah_cross_n = test.pah_cross_n*1e-18
63   -test.pah_cross_c = test.pah_cross_c*1e-18
64   -test.pah_cross_dc = test.pah_cross_dc*1e-18
65   -
66   -''' small n_c '''
67   -test_small_n_c_1 = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, 32, parsec, fc_pah, ISRF)
68   -test_small_n_c_2 = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, 36, parsec, fc_pah, ISRF)
69   -test_small_n_c_3 = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, 38, parsec, fc_pah, ISRF)
70   -
71   -test_small_n_c_1.parameters()
72   -test_small_n_c_2.parameters()
73   -test_small_n_c_3.parameters()
74   -
75   -''' medium n_c '''
76   -test_medium_n_c_1 = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, 40, parsec, fc_pah, ISRF)
77   -test_medium_n_c_2 = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, 42, parsec, fc_pah, ISRF)
78   -
79   -test_medium_n_c_1.parameters()
80   -test_medium_n_c_2.parameters()
81   -
82   -''' large n_c '''
83   -test_large_n_c_1 = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, 48, parsec, fc_pah, ISRF)
84   -test_large_n_c_2 = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, 54, parsec, fc_pah, ISRF)
85   -test_large_n_c_3 = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, 66, parsec, fc_pah, ISRF)
86   -
87   -test_large_n_c_1.parameters()
88   -test_large_n_c_2.parameters()
89   -test_large_n_c_3.parameters()
90   -
91   -''' average of first and second photoionization yield '''
92   -yield_of_first_photoionization = (test_small_n_c_1.yield_of_first_photoionization +\
93   - test_small_n_c_2.yield_of_first_photoionization +\
94   - test_small_n_c_3.yield_of_first_photoionization +\
95   - test_medium_n_c_1.yield_of_first_photoionization +\
96   - 2*test_medium_n_c_2.yield_of_first_photoionization +\
97   - test_large_n_c_1.yield_of_first_photoionization +\
98   - test_large_n_c_2.yield_of_first_photoionization +\
99   - test_large_n_c_3.yield_of_first_photoionization )/9
100   -yield_of_second_photoionization = (test_small_n_c_1.yield_of_second_photoionization +\
101   - test_small_n_c_2.yield_of_second_photoionization +\
102   - test_small_n_c_3.yield_of_second_photoionization +\
103   - test_medium_n_c_1.yield_of_second_photoionization +\
104   - 2*test_medium_n_c_2.yield_of_second_photoionization +\
105   - test_large_n_c_1.yield_of_second_photoionization +\
106   - test_large_n_c_2.yield_of_second_photoionization +\
107   - test_large_n_c_3.yield_of_second_photoionization )/9
108   -
109   -''' average of ionization cross-section for the cases of neutral molecules, and
110   -cations '''
111   -ionization_cross_a = test.pah_cross_a * test.detachment_yield
112   -ionization_cross_n = test.pah_cross_n * yield_of_first_photoionization
113   -ionization_cross_c = test.pah_cross_c * yield_of_second_photoionization
114   -#ionization_cross_dc is too small to be traced
115   -
116   -''' average of energy per cross-section for the cases of neutral molecules,
117   -cations and dications '''
118   -energy = (test.energy_negative_charged + test.energy_neutral + test.energy_charged + test.energy_double_charged)/4
119   -
120   -''' plot sigma and sigma_ion in cmยฒ/(C atom) as a function of photon energy in eV
121   -for the cases of neutral, charged and double charged molecules '''
122   -plt.figure(figsize=(3.57,9.0))
123   -
124   -plt.subplot(4,1,1)
125   -plt.plot(test.energy_negative_charged,test.pah_cross_a,'k', label='$\sigma(E,Z)$')
126   -plt.plot(energy,ionization_cross_a,'k',ls='--', label='$\sigma_{ion}(E,Z)$')
127   -plt.title('Z=-1', fontsize = 10)
128   -plt.legend(fontsize = 8)
129   -
130   -plt.subplot(4,1,2)
131   -plt.plot(test.energy_neutral,test.pah_cross_n,'k')
132   -plt.plot(energy,ionization_cross_n,'k',ls='--')
133   -plt.title('Z=0', fontsize = 10)
134   -
135   -plt.subplot(4,1,3)
136   -plt.plot(test.energy_neutral,test.pah_cross_c,'k')
137   -plt.plot(energy,ionization_cross_c,'k',ls='--')
138   -plt.title('Z=1', fontsize = 10)
139   -plt.ylabel('Cross sections (cm$^2$/C)')
140   -
141   -plt.subplot(4,1,4)
142   -plt.plot(test.energy_neutral,test.pah_cross_dc,'k')
143   -plt.title('Z=2', fontsize = 10)
144   -plt.xlabel('Photon energy (eV)')
145   -
146   -plt.subplots_adjust(hspace = 0.35)
147   -plt.savefig('cross_ioni_section_categories.pdf', bbox_inches='tight')
148   -if ISRF == False:
149   -
150   - ''' evolution of the different populations of neutral pah molecules,
151   - # cations and dications as a function of the distance to the star considered '''
152   -
153   - ''' parameter '''
154   - # population_frac3 = np.zeros([1, len(distance_list)])
155   - # population_frac4 = np.zeros([4, len(distance_list)])
156   -
157   - ''' indentation '''
158   - # j=0
159   - # gamma_list3 = []
160   - # gamma_list4 = []
161   - # for i in distance_list:
162   - # population_plot_test3 = mp.HeatingGas(filename, star_radius, t_gas, n_e, n_c, i, fc_pah, ISRF)
163   - # population_plot_test3.parameters()
164   - # population_frac3[0,j] =(population_plot_test3.frac_charged +\
165   - # population_plot_test3.frac_double_charged)/\
166   - # (population_plot_test3.frac_neutral +\
167   - # population_plot_test3.frac_charged +\
168   - # population_plot_test3.frac_double_charged)
169   - # # population_frac3[0,j] = population_plot_test3.frac_neutral
170   - # # population_frac3[1,j] = population_plot_test3.frac_charged
171   - # # population_frac3[2,j] = population_plot_test3.frac_double_charged
172   - # gamma_list3.append(population_plot_test3.gamma)
173   -
174   - # population_plot_test4 = four_lvl.HeatingGas(filename, star_radius, t_gas, n_e, n_c, i, fc_pah, ISRF)
175   - # population_plot_test4.parameters()
176   - # population_frac4[0,para] =(population_plot_test4.frac_anion +\
177   - # population_plot_test4.frac_charged +\
178   - # population_plot_test4.frac_double_charged)/\
179   - # (population_plot_test4.frac_anion +\
180   - # population_plot_test4.frac_neutral +\
181   - # population_plot_test4.frac_charged +\
182   - # population_plot_test4.frac_double_charged)
183   - # population_frac4[0,j] = population_plot_test4.frac_anion
184   - # population_frac4[1,j] = population_plot_test4.frac_neutral
185   - # population_frac4[2,j] = population_plot_test4.frac_charged
186   - # population_frac4[3,j] = population_plot_test4.frac_double_charged
187   - # gamma_list4.append(population_plot_test4.gamma)
188   - # j=j+1
189   -
190   - ''' population fractions as function of the distance '''
191   - # plt.figure(figsize=(5.90,3.93))
192   - # plt.plot(gamma_list3,population_frac3[0,:],c='k',label='3 levels model')
193   - # plt.plot(gamma_list4,population_frac4[0,:],ls='--', c='k',label='4 levels model')
194   - # plt.xlim([1e1, 1e7])
195   - # plt.plot(gamma_list3, population_frac3[0,:], color = 'green', label='Z=0, 3 levels model')
196   - # plt.plot(gamma_list3, population_frac3[1,:], color = 'blue', label='Z=1, 3 levels model')
197   - # plt.plot(gamma_list3, population_frac3[2,:], color = 'red', label='Z=2, 3 levels model')
198   - # plt.plot(gamma_list4, population_frac4[0,:], ls='--', color = 'k', label='Z=-1, 4 levels model')
199   - # plt.plot(gamma_list4, population_frac4[1,:], ls='--', color = 'green', label='Z=0, 4 levels model')
200   - # plt.plot(gamma_list4, population_frac4[2,:], ls='--', color = 'k', label='Z=1, 4 levels model')
201   - # plt.plot(gamma_list4, population_frac4[3,:], ls='--', color = 'red', label='Z=2, 4 levels model')
202   - # plt.plot(gamma_list3, population_frac3[0,:]+population_frac3[1,:]+population_frac3[2,:],'blue')
203   - # plt.plot(gamma_list4, population_frac4[0,:]+population_frac4[1,:]+population_frac4[2,:]+population_frac4[3,:],'darkblue')
204   - # plt.xlabel('$\gamma$ ($G_0\sqrt{T}/n_e$)')
205   - # plt.ylabel('population fraction R$_i$')
206   - # plt.legend()
207   -
208   - ''' photoelectric efficiencies with the pah heating model '''
209   -
210   - # gamma_list_T_one3 = []
211   - # heating_efficiency_list_T_one3 = []
212   -
213   - # gamma_list_T_two3 = []
214   - # heating_efficiency_list_T_two3 = []
215   -
216   - # gamma_list_T_one4 = []
217   - # heating_efficiency_list_T_one4 = []
218   -
219   - # gamma_list_T_two4 = []
220   - # heating_efficiency_list_T_two4 = []
221   -
222   - # ratio1_1 = []
223   - # ratio2_1 = []
224   - # gamma_mean1_1 = []
225   - # gamma_mean2_1 = []
226   - # gamma_list_1 = []
227   - # j = 0
228   - # for i in distance_list:
229   - # j = j+1
230   - # if j%100 == 0:
231   - # print('Progression ร  {}'.format(j/10),'%')
232   - # print('')
233   - # test_for_efficiencies_T_one3 = mp.HeatingGas(filename, star_radius, 100, n_e, 54, i, fc_pah, ISRF)
234   - # test_for_efficiencies_T_one3.parameters()
235   - # gamma_list_T_one3.append(test_for_efficiencies_T_one3.gamma)
236   - # heating_efficiency_list_T_one3.append(test_for_efficiencies_T_one3.heating_efficiency)
237   -
238   - # test_for_efficiencies_T_two3 = mp.HeatingGas(filename, star_radius, 1000, n_e, 54, i, fc_pah, ISRF)
239   - # test_for_efficiencies_T_two3.parameters()
240   - # gamma_list_T_two3.append(test_for_efficiencies_T_two3.gamma)
241   - # heating_efficiency_list_T_two3.append(test_for_efficiencies_T_two3.heating_efficiency)
242   -
243   - # test_for_efficiencies_T_one4 = four_lvl.HeatingGas(filename, star_radius, 100, n_e, n_c, i, fc_pah, ISRF)
244   - # test_for_efficiencies_T_one4.parameters()
245   - # gamma_list_T_one4.append(test_for_efficiencies_T_one4.gamma)
246   - # heating_efficiency_list_T_one4.append(test_for_efficiencies_T_one4.heating_efficiency)
247   -
248   - # test_for_efficiencies_T_two4 = four_lvl.HeatingGas(filename, star_radius, 1000, n_e, n_c, i, fc_pah, ISRF)
249   - # test_for_efficiencies_T_two4.parameters()
250   - # gamma_list_T_two4.append(test_for_efficiencies_T_two4.gamma)
251   - # heating_efficiency_list_T_two4.append(test_for_efficiencies_T_two4.heating_efficiency)
252   -
253   - # ratio1_1.append( (test_for_efficiencies_T_one3.heating_efficiency/test_for_efficiencies_T_one4.heating_efficiency)*100)
254   - # ratio2_1.append( (test_for_efficiencies_T_two3.heating_efficiency/test_for_efficiencies_T_two4.heating_efficiency)*100)
255   - # gamma_mean1_1.append( (test_for_efficiencies_T_one3.gamma + test_for_efficiencies_T_one4.gamma)/2 )
256   - # gamma_mean2_1.append( (test_for_efficiencies_T_two3.gamma + test_for_efficiencies_T_two4.gamma)/2 )
257   -
258   - # plt.figure()
259   - # plt.xlim([1e1, 1e6])
260   - # plt.loglog(gamma_mean1_1, ratio1_1, 'blue', label='ratio of 3 and 4 levels model $\epsilon_{PAH}$ (T = 100K) ')
261   - # plt.loglog(gamma_mean2_1, ratio2_1, 'red', ls='--' , label='ratio of 3 and 4 levels model $\epsilon_{PAH}$ (T = 1000K)')
262   - # plt.loglog(gamma_list_T_one3, heating_efficiency_list_T_one3, 'blue', label='$\epsilon_{PAH}$ (T = 100K) 3 levels model')
263   - # plt.loglog(gamma_list_T_two3, heating_efficiency_list_T_two3, 'red', ls='--' , label='$\epsilon_{PAH}$ (T = 1000K) 3 levels model')
264   - # plt.loglog(gamma_list_T_one4, heating_efficiency_list_T_one4, 'darkblue', label='$\epsilon_{PAH}$ (T = 100K) 4 levels model')
265   - # plt.loglog(gamma_list_T_two4, heating_efficiency_list_T_two4, 'darkred', ls='--' , label='$\epsilon_{PAH}$ (T = 1000K) 4 levels model')
266   - # plt.xlabel('$\gamma(G_0\sqrt{T}/n_e)$')
267   - # plt.ylabel('$\epsilon_{PAH}$')
268   - # plt.ylabel('$\epsilon_{PAH, 3lvlmod}$ / $\epsilon_{PAH, 4lvlmod}$ (%)')
269   - # plt.legend()
270   -
271   - # ''' Total photoelectric heating rates of the gas '''
272   -
273   - # gamma_list_fc_one3 = []
274   - # total_gas_heating_per_rf_list_fc_one3 = []
275   -
276   - # gamma_list_fc_two3 = []
277   - # total_gas_heating_per_rf_list_fc_two3 = []
278   -
279   - # gamma_list_fc_one4 = []
280   - # total_gas_heating_per_rf_list_fc_one4 = []
281   -
282   - # gamma_list_fc_two4 = []
283   - # total_gas_heating_per_rf_list_fc_two4 = []
284   -
285   - # ratio1_2 = []
286   - # ratio2_2 = []
287   - # gamma_mean1_2 = []
288   - # gamma_mean2_2 = []
289   - # Table_list_g_fc5 = []
290   - # Table_list_Gpg0_fc5 = []
291   -
292   - # Table_list_g_fc10 = []
293   - # Table_list_Gpg0_fc10 = []
294   -
295   - # for i in distance_list:
296   - # test_for_total_heating_rates_fc_one3 = mp.HeatingGas(filename, star_radius, t_gas, n_e, n_c, i, 0.05, ISRF)
297   - # test_for_total_heating_rates_fc_one3.parameters()
298   - # gamma_list_fc_one3.append(test_for_total_heating_rates_fc_one3.gamma)
299   - # total_gas_heating_per_rf_list_fc_one3.append( (test_for_total_heating_rates_fc_one3.total_gas_heating*1e-7) / test_for_total_heating_rates_fc_one3.g_0)
300   -
301   - # test_for_total_heating_rates_fc_two3 = mp.HeatingGas(filename, star_radius, t_gas, n_e, n_c, i, 0.1, ISRF)
302   - # test_for_total_heating_rates_fc_two3.parameters()
303   - # gamma_list_fc_two3.append(test_for_total_heating_rates_fc_two3.gamma)
304   - # total_gas_heating_per_rf_list_fc_two3.append( (test_for_total_heating_rates_fc_two3.total_gas_heating*1e-7) / test_for_total_heating_rates_fc_two3.g_0)
305   -
306   - # test_for_total_heating_rates_fc_one4 = four_lvl.HeatingGas(filename, star_radius, 500, n_e, 54, i, 0.05, ISRF)
307   - # test_for_total_heating_rates_fc_one4.parameters()
308   - # gamma_list_fc_one4.append(test_for_total_heating_rates_fc_one4.gamma)
309   - # total_gas_heating_per_rf_list_fc_one4.append( (test_for_total_heating_rates_fc_one4.total_gas_heating*1e-7) / test_for_total_heating_rates_fc_one4.g_0)
310   -
311   - # test_for_total_heating_rates_fc_two4 = four_lvl.HeatingGas(filename, star_radius, 500, n_e, 54, i, 0.1, ISRF)
312   - # test_for_total_heating_rates_fc_two4.parameters()
313   - # gamma_list_fc_two4.append(test_for_total_heating_rates_fc_two4.gamma)
314   - # total_gas_heating_per_rf_list_fc_two4.append( (test_for_total_heating_rates_fc_two4.total_gas_heating*1e-7) / test_for_total_heating_rates_fc_two4.g_0)
315   -
316   - # ratio1_2.append( ( ( (test_for_total_heating_rates_fc_one3.total_gas_heating*1e-7) /\
317   - # test_for_total_heating_rates_fc_one3.g_0) /\
318   - # ( (test_for_total_heating_rates_fc_one4.total_gas_heating*1e-7) /\
319   - # test_for_total_heating_rates_fc_one4.g_0) )*100 )
320   - # ratio2_2.append( ( ( (test_for_total_heating_rates_fc_two3.total_gas_heating*1e-7) /\
321   - # test_for_total_heating_rates_fc_two3.g_0) /\
322   - # ( (test_for_total_heating_rates_fc_two4.total_gas_heating*1e-7) /\
323   - # test_for_total_heating_rates_fc_two4.g_0) )*100 )
324   - # gamma_mean1_2.append( (test_for_total_heating_rates_fc_one3.gamma+test_for_total_heating_rates_fc_one4.gamma)/2 )
325   - # gamma_mean2_2.append( (test_for_total_heating_rates_fc_two3.gamma+test_for_total_heating_rates_fc_two4.gamma)/2 )
326   -
327   - # plt.figure()
328   - # plt.xlim([1e1, 1e7])
329   - # plt.loglog(gamma_mean1_2,ratio1_2, 'red', label='ratio of 3 and 4 levels model $\Gamma_{tot}$ (5% of C in PAHs)')
330   - # plt.loglog(gamma_mean2_2,ratio2_2, 'blue', label='ratio of 3 and 4 levels model $\Gamma_{tot}$ (10% of C in PAHs)')
331   - # plt.loglog(gamma_list_fc_one3, total_gas_heating_per_rf_list_fc_one3, c='red', label='This work (5% of C in PAHs) 3 levels model')
332   - # plt.loglog(gamma_list_fc_two3, total_gas_heating_per_rf_list_fc_two3,ls='--', c='blue', label='This work (10% of C in PAHs) 3 levels model')
333   - # plt.loglog(gamma_list_fc_one4, total_gas_heating_per_rf_list_fc_one4, c='darkred', label='This work (5% of C in PAHs) 4 levels model')
334   - # plt.loglog(gamma_list_fc_two4, total_gas_heating_per_rf_list_fc_two4,ls='--', c='darkblue', label='This work (10% of C in PAHs) 4 levels model')
335   - # plt.xlabel('$\gamma(G_0\sqrt{T}/n_e)$')
336   - # # plt.ylabel('$\Gamma_{tot}$ (W H$^{-1}$ G$_0^{-1}$)')
337   - # plt.ylabel('$\Gamma_{tot, 3lvlmod}$ / $\Gamma_{tot, 4lvlmod}$ (%)')
338   - # plt.legend()
339   -
340   - # Table_values.append(heating_efficiency_list_T_one4)
341   - # Table_list_g_fc5 = gamma_list_fc_one4[::-1]
342   - # Table_list_Gpg0_fc5 = total_gas_heating_per_rf_list_fc_one4[::-1]
343   -
344   - # Table_list_g_fc10 = gamma_list_fc_two4[::-1]
345   - # Table_list_Gpg0_fc10 = total_gas_heating_per_rf_list_fc_two4[::-1]
346 0 \ No newline at end of file