fmm_balance.m 2.4 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162
  1. %{
  2. Copyright © 2020 Alexey A. Shcherbakov. All rights reserved.
  3. This file is part of GratingFMM.
  4. GratingFMM is free software: you can redistribute it and/or modify
  5. it under the terms of the GNU General Public License as published by
  6. the Free Software Foundation, either version 2 of the License, or
  7. (at your option) any later version.
  8. GratingFMM is distributed in the hope that it will be useful,
  9. but WITHOUT ANY WARRANTY; without even the implied warranty of
  10. MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
  11. GNU General Public License for more details.
  12. You should have received a copy of the GNU General Public License
  13. along with GratingFMM. If not, see <https://www.gnu.org/licenses/>.
  14. %}
  15. %% description:
  16. % calculate the power balance (check the energy conservation law) in case
  17. % of the collinear diffraction by 1D gratings
  18. % the returned value should be close to zero for pure dielectric structures
  19. % with no losses
  20. %% input:
  21. % no: number of Fourier harmonics
  22. % V_inc: incident field amplitude matrix of size (no, 2)
  23. % V_dif: diffracted field amplitude matrix of size (no, 2)
  24. % first index of V_inc, V_dif indicates diffraction harmonics
  25. % (0-th order index is ind_0 = ceil(no/2))
  26. % second index of V_inc, V_dif, V_eff indicates whether the diffraction orders
  27. % are in the substrate (V(:,1)) or in the superstrate (V(:,2))
  28. % kx0: incident plane wave wavevector x-projection (Bloch wavevector)
  29. % kg: wavelength-to-period ratio (grating vector)
  30. % eps1, eps2: substrate and superstrate permittivities
  31. % pol: polarization (either "TE" or "TM")
  32. %% output:
  33. % if the incident field has propagating harmonics the function returns the
  34. % normalized difference between the incident and diffractied field total
  35. % power, otherwise (if the incident field is purely evanescent) it returns
  36. % the total power carried by propagating diffraction orders
  37. %% implementation
  38. function [b] = fmm_balance(no, V_inc, V_dif, kx0, kg, eps1, eps2, pol)
  39. [kz1, kz2] = fmm_kxz(no, kx0, 0, kg, eps1, eps2);
  40. kz1 = transpose(kz1);
  41. kz2 = transpose(kz2);
  42. if (strcmp(pol,'TM'))
  43. kz1 = kz1/eps1;
  44. kz2 = kz2/eps2;
  45. end
  46. P_inc = sum( abs(V_inc(:,1).^2).*real(kz1) + abs(V_inc(:,2).^2).*real(kz2) );
  47. P_dif = sum( abs(V_dif(:,1).^2).*real(kz1) + abs(V_dif(:,2).^2).*real(kz2) );
  48. if (abs(P_inc) > 1e-15)
  49. b = abs(P_dif/P_inc-1);
  50. else
  51. b = 0.5*P_dif;
  52. end
  53. end
  54. %
  55. % END
  56. %