field-SiAgSi-flow.py 3.5 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127
  1. #!/usr/bin/env python
  2. # -*- coding: UTF-8 -*-
  3. #
  4. # Copyright (C) 2009-2015 Ovidio Peña Rodríguez <ovidio@bytesfall.com>
  5. # Copyright (C) 2013-2015 Konstantin Ladutenko <kostyfisik@gmail.com>
  6. #
  7. # This file is part of python-scattnlay
  8. #
  9. # This program is free software: you can redistribute it and/or modify
  10. # it under the terms of the GNU General Public License as published by
  11. # the Free Software Foundation, either version 3 of the License, or
  12. # (at your option) any later version.
  13. #
  14. # This program is distributed in the hope that it will be useful,
  15. # but WITHOUT ANY WARRANTY; without even the implied warranty of
  16. # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
  17. # GNU General Public License for more details.
  18. #
  19. # The only additional remark is that we expect that all publications
  20. # describing work using this software, or all commercial products
  21. # using it, cite the following reference:
  22. # [1] O. Pena and U. Pal, "Scattering of electromagnetic radiation by
  23. # a multilayered sphere," Computer Physics Communications,
  24. # vol. 180, Nov. 2009, pp. 2348-2354.
  25. #
  26. # You should have received a copy of the GNU General Public License
  27. # along with this program. If not, see <http://www.gnu.org/licenses/>.
  28. # This test case calculates the electric field in the
  29. # E-k plane, for an spherical Si-Ag-Si nanoparticle. Core radius is 17.74 nm,
  30. # inner layer 23.31nm, outer layer 22.95nm. Working wavelength is 800nm, we use
  31. # silicon epsilon=13.64+i0.047, silver epsilon= -28.05+i1.525
  32. import scattnlay
  33. from scattnlay import fieldnlay
  34. from scattnlay import scattnlay
  35. from fieldplot import fieldplot
  36. import numpy as np
  37. import cmath
  38. epsilon_Si = 13.64 + 0.047j
  39. epsilon_Ag = -28.05 + 1.525j
  40. # epsilon_Si = 2.0 + 0.047j
  41. # epsilon_Ag = -2.0 + 1.525j
  42. # air = 1
  43. # epsilon_Si = air*2
  44. # epsilon_Ag = air*2
  45. index_Si = np.sqrt(epsilon_Si)
  46. index_Ag = np.sqrt(epsilon_Ag)
  47. print(index_Si)
  48. print(index_Ag)
  49. # # Values for 800 nm, taken from http://refractiveindex.info/
  50. # index_Si = 3.69410 + 0.0065435j
  51. # index_Ag = 0.18599 + 4.9886j
  52. WL=800 #nm
  53. core_width = 17.74 #nm Si
  54. inner_width = 23.31 #nm Ag
  55. outer_width = 22.95 #nm Si
  56. core_r = core_width
  57. inner_r = core_r+inner_width
  58. outer_r = inner_r+outer_width
  59. # n1 = 1.53413
  60. # n2 = 0.565838 + 7.23262j
  61. nm = 1.0
  62. x = np.ones((3), dtype = np.float64)
  63. x[0] = 2.0*np.pi*core_r/WL
  64. x[1] = 2.0*np.pi*inner_r/WL
  65. x[2] = 2.0*np.pi*outer_r/WL
  66. m = np.ones((3), dtype = np.complex128)
  67. m[0] = index_Si/nm
  68. m[1] = index_Ag/nm
  69. m[2] = index_Si/nm
  70. print "x =", x
  71. print "m =", m
  72. npts = 501
  73. factor=2.2
  74. flow_total = 21
  75. crossplane='XZ'
  76. #crossplane='YZ'
  77. #crossplane='XY'
  78. # Options to plot: Eabs, Habs, Pabs, angleEx, angleHy
  79. field_to_plot='Eabs'
  80. #field_to_plot='angleEx'
  81. comment='SiAgSi-absorber-flow'
  82. WL_units='nm'
  83. import matplotlib.pyplot as plt
  84. fig, axs = plt.subplots(1,1)#, sharey=True, sharex=True)
  85. fig.tight_layout()
  86. fieldplot(fig, axs, x,m, WL, comment, WL_units, crossplane, field_to_plot, npts, factor, flow_total,
  87. subplot_label=' ',is_flow_extend=False, outline_width=1.5)
  88. #fieldplot(x,m, WL, comment, WL_units, crossplane, field_to_plot, npts, factor, flow_total, is_flow_extend=False)
  89. # for ax in axs:
  90. # ax.locator_params(axis='x',nbins=5)
  91. # ax.locator_params(axis='y',nbins=5)
  92. fig.subplots_adjust(hspace=0.3, wspace=-0.1)
  93. plt.savefig(comment+"-R"+str(int(round(x[-1]*WL/2.0/np.pi)))+"-"+crossplane+"-"
  94. +field_to_plot+".pdf",pad_inches=0.02, bbox_inches='tight')
  95. plt.draw()
  96. # plt.show()
  97. plt.clf()
  98. plt.close()