PLaSK library

Loading...
Searching...
No Matches
fields2d.py
Go to the documentation of this file.
1#!/usr/bin/env python3
2# This file is part of PLaSK (https://plask.app) by Photonics Group at TUL
3# Copyright (c) 2022 Lodz University of Technology
4#
5# This program is free software: you can redistribute it and/or modify
6# it under the terms of the GNU General Public License as published by
7# the Free Software Foundation, version 3.
8#
9# This program is distributed in the hope that it will be useful,
10# but WITHOUT ANY WARRANTY; without even the implied warranty of
11# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12# GNU General Public License for more details.
13
14 import unittest
15
16 from numpy import *
17
18 from plask import *
19 from plask import material, geometry, mesh
20 from optical import modal
21
22 config.axes = 'zxy'
23
24@material.simple()
25 class Glass(material.Material):
26 def Nr(self, wl, T=300., n=0.): return 1.3
27
28
29solver = modal.Fourier2D("solver")
30 solver.ft = 'analytic'
31wire_stack = geometry.Stack2D()
32 wire_stack.append(geometry.Rectangle(0.75, 0.125, Glass()))
33rect = geometry.Rectangle(0.75, 0.375, Glass())
34 rect.role = 'interface'
35 wire_stack.append(rect)
36
37shelf = geometry.Shelf2D()
38 shelf.append(wire_stack)
39 shelf.append(geometry.Rectangle(0.15, 0.50, 'air'))
40
41 solver.geometry = geometry.Cartesian2D(shelf, left="mirror", right="periodic")
42
43 solver.size = 12
44
45det = False
46
47start = 1.1216 # 1.1196
48
49
50
51 def plotFields(symmetry):
52 solver.lam = 1000.
53 solver.symmetry = symmetry
54 solver.invalidate()
55
56 if det:
57 figure()
58 neffs = linspace(1.0, 1.2, 1001)
59 dets = solver.get_determinant(neff=neffs)
60 plot(neffs, abs(dets))
61 window_title(f"det {symmetry}")
62
63 mn = solver.find_mode(neff=start)
64
65 msh = mesh.Rectangular2D(mesh.Regular(-2., 2., 401), mesh.Regular(-0.4, 1.0, 1501))
66 figure(figsize=(14, 6))
67
68 E = solver.outLightE(mn, msh)
69 m = max(abs(E))
70 E /= m
71
72 ax_mag = subplot2grid((3,3), (0,1), rowspan=3)
73 plot_field(Data(sum(real(E*conj(E)),-1), E.mesh))
74 plot_geometry(solver.geometry, color='w', mirror=True, periods=2, lw=0.5)
75 #colorbar(use_gridspec=True)
76
77 c = 2
78
79 subplot2grid((3,3), (0,2), rowspan=3, sharey=ax_mag)
80 msh1 = mesh.Rectangular2D([0], msh.axis1)
81 E1 = solver.outLightE(mn, msh1)
82 E1 /= m
83 mr = max(abs(E.array[:,:,c].real).ravel())
84 mi = max(abs(E.array[:,:,c].imag).ravel())
85 if mi > mr: E1 = E1.imag
86 else: E1 = E1.real
87 axhline(0.0, ls='--', color='k', lw=0.5)
88 axhline(0.5, ls='--', color='k', lw=0.5)
89 plot_profile(E1, swap_axes=True, comp=c)
90 Ec = E1.array[0,:,c].copy()
91 y = array(msh.axis1)
92 Ec[(0.0 <= y) & (y <= 0.5)] *= 1.3**2
93 plot(Ec, y)
94 ylim(-0.4, 1.0)
95 window_title(f"Field {symmetry}")
96
97 m = max(abs(E))
98 E /= m
99 levels = linspace(-1, 1, 16)
100 for c in range(3):
101 subplot2grid((3,3), (c, 0), sharey=ax_mag, sharex=ax_mag)
102 mr = max(abs(E.array[:,:,c].real).ravel())
103 mi = max(abs(E.array[:,:,c].imag).ravel())
104 if mi > mr: Ec = E.imag
105 else: Ec = E.real
106 plot_field(Ec, levels, comp=c, cmap='bwr')
107 plot_geometry(solver.geometry, color='k', mirror=True, periods=2, lw=0.5)
108 xlim(msh.axis0[0], msh.axis0[-1])
109
110 tight_layout()
111
112
113 plotFields('Htran')
114 plotFields(None)
115 show()