PLaSK library

Loading...
Searching...
No Matches
xance.hpp
Go to the documentation of this file.
1/*
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#ifndef PLASK__SOLVER_OPTICAL_MODAL_XANCE_H
15#define PLASK__SOLVER_OPTICAL_MODAL_XANCE_H
16
17#include "matrices.hpp"
18#include "transfer.hpp"
19#include "solver.hpp"
20
21
22 namespace plask { namespace optical { namespace modal {
23
27 struct PLASK_SOLVER_API XanceTransfer: public Transfer {
28
30 struct FieldsDiagonalized {
31 cvector E0, H0, Ed, Hd;
32 };
33
34 protected:
35
36 cmatrix Y;
37
38 bool needAllY;
39
40 std::vector<FieldsDiagonalized> fields;
41
42 std::vector<cmatrix> memY;
43
45 double get_d(size_t n, double& z, PropagationDirection& part) {
46 double d = (n == 0 || std::size_t(n) == solver->vbounds->size())?
47 solver->vpml.dist :
48 solver->vbounds->at(n) - solver->vbounds->at(n-1);
49 if (std::ptrdiff_t(n) >= solver->interface) {
50 z = d - z;
51 if (part == PROPAGATION_DOWNWARDS) part = PROPAGATION_UPWARDS;
52 else if (part == PROPAGATION_UPWARDS) part = PROPAGATION_DOWNWARDS;
53 }
54 else if (n == 0) z += d;
55 return d;
56 }
57
59 double get_d(size_t n, double& z1, double& z2) {
60 double d = (n == 0 || std::size_t(n) == solver->vbounds->size())?
61 solver->vpml.dist :
62 solver->vbounds->at(n) - solver->vbounds->at(n-1);
63 if (std::ptrdiff_t(n) >= solver->interface) {
64 const double zl = z1;
65 z1 = d - z2; z2 = d - zl;
66 }
67 else if (n == 0) {
68 z1 += d;
69 z2 += d;
70 }
71 return d;
72 }
73
74 public:
75
76 cvector getTransmissionVector(const cvector& incident, IncidentDirection side) override;
77
78 XanceTransfer(ModalBase* solver, Expansion& expansion);
79
80 protected:
81
82 cvector getFieldVectorE(double z, std::size_t n, PropagationDirection part = PROPAGATION_TOTAL) override;
83
84 cvector getFieldVectorH(double z, std::size_t n, PropagationDirection part = PROPAGATION_TOTAL) override;
85
90 void storeY(size_t n);
91
96 const cmatrix& getY(std::size_t n) {
97 if (memY.size() == solver->stack.size() && needAllY)
98 return memY[n];
99 else
100 throw CriticalException("{0}: Y matrices are not stored", solver->getId());
101 }
102
104 inline void get_y1(const cdiagonal& gamma, double d, cdiagonal& y1) const {
105 const std::size_t N = gamma.size();
106 assert(y1.size() == N);
107 for (std::size_t i = 0; i < N; i++) {
108 dcomplex t = tanh(I*gamma[i]*d);
109 if (isinf(real(t)) || isinf(imag(t))) y1[i] = 0.;
110 else if (abs(t) < SMALL)
111 throw ComputationError(solver->getId(),
112 "Matrix y1 has some infinite value (try changing wavelength or layer thickness a bit)");
113 else y1[i] = 1. / t;
114 }
115 }
116
118 inline void get_y2(const cdiagonal& gamma, double d, cdiagonal& y2) const {
119 const std::size_t N = gamma.size();
120 assert(y2.size() == N);
121 for (std::size_t i = 0; i < N; i++) {
122 dcomplex s = sinh(I*gamma[i]*d);
123 if (isinf(real(s)) || isinf(imag(s))) y2[i] = 0.;
124 else if (abs(s) < SMALL)
125 throw ComputationError(solver->getId(),
126 "Matrix y2 has some infinite value (try changing wavelength or layer thickness a bit)");
127 else y2[i] = - 1. / s;
128 }
129 }
130
131 double integrateField(WhichField field, size_t n, double z1, double z2) override;
132};
133
134
135}}} // namespace plask::optical::modal
136
137#endif // PLASK__SOLVER_OPTICAL_MODAL_XANCE_H