PLaSK library

Loading...
Searching...
No Matches
electr3d.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__MODULE_ELECTRICAL_ELECTR3D_H
15#define PLASK__MODULE_ELECTRICAL_ELECTR3D_H
16
17#include "common.hpp"
18
19 namespace plask { namespace electrical { namespace shockley {
20
24 struct PLASK_SOLVER_API ElectricalFem3DSolver : public FemSolverWithMaskedMesh<Geometry3D, RectangularMesh<3>> {
25 protected:
26
28 struct Active {
29 struct Region {
30 size_t bottom, top, left, right, back, front, lon, tra;
31 bool warn;
32 Region() {}
33 Region(size_t b, size_t t, size_t x, size_t y)
34 : bottom(b),
35 top(t),
36 left(std::numeric_limits<size_t>::max()),
37 right(0),
38 back(std::numeric_limits<size_t>::max()),
39 front(0),
40 lon(x),
41 tra(y),
42 warn(true) {}
43 };
44 size_t bottom, top, left, right, back, front;
45 size_t ld;
46 ptrdiff_t offset;
47 double height;
48 Active()
49 : bottom(0),
50 top(0),
51 left(0),
52 right(0),
53 back(0),
54 front(0),
55 ld(0),
56 offset(0),
57 height(0.) {}
58 Active(size_t tot, size_t v0, size_t v1, size_t t0, size_t t1, size_t l0, size_t l1, double h)
59 : bottom(v0),
60 top(v1),
61 left(t0),
62 right(t1),
63 back(l0),
64 front(l1),
65 ld(l1 - l0),
66 offset(tot - (l1 - l0) * t0 - l0),
67 height(h) {}
68 Active(size_t tot, Region r, double h)
69 : bottom(r.bottom),
70 top(r.top),
71 left(r.left),
72 right(r.right),
73 back(r.back),
74 front(r.front),
75 ld(r.front - r.back),
76 offset(tot - (r.front - r.back) * r.left - r.back),
77 height(h) {}
78 };
79
80 double pcond;
81 double ncond;
82
83 int loopno;
84 double toterr;
86
87 DataVector<Tensor2<double>> junction_conductivity;
88 Tensor2<double> default_junction_conductivity;
89
90 DataVector<Tensor2<double>> conds;
91
92 DataVector<double> potential;
93 DataVector<Vec<3, double>> current;
94 DataVector<double> heat;
95
96 std::vector<Active> active;
97
104 void setMatrix(FemMatrix& A,
105 DataVector<double>& B,
106 const BoundaryConditionsWithMesh<RectangularMesh<3>::Boundary, double>& bvoltage,
107 const LazyData<double>& temperature);
108
115 virtual Tensor2<double> activeCond(size_t n, double U, double jy, double T) = 0;
116
120 LazyData<double> loadConductivity();
121
123 void saveConductivity();
124
126 void saveHeatDensity();
127
129 void onInitialize() override;
130
132 void onInvalidate() override;
133
134 // void onMeshChange(const typename RectangularMesh<3>::Event& evt) override {
135 // SolverWithMesh<Geometry3D, RectangularMesh<3>>::onMeshChange(evt);
136 // setupActiveRegions();
137 // }
138
139 // void onGeometryChange(const Geometry::Event& evt) override {
140 // SolverWithMesh<Geometry3D, RectangularMesh<3>>::onGeometryChange(evt);
141 // setupActiveRegions();
142 // }
143
145 void setupActiveRegions();
146
148 size_t isActive(const Vec<3>& point) const {
149 size_t no(0);
150 auto roles = this->geometry->getRolesAt(point);
151 for (auto role : roles) {
152 size_t l = 0;
153 if (role.substr(0, 6) == "active")
154 l = 6;
155 else if (role.substr(0, 8) == "junction")
156 l = 8;
157 else
158 continue;
159 if (no != 0) throw BadInput(this->getId(), "multiple 'active'/'junction' roles specified");
160 if (role.size() == l)
161 no = 1;
162 else {
163 try {
164 no = boost::lexical_cast<size_t>(role.substr(l)) + 1;
165 } catch (boost::bad_lexical_cast&) {
166 throw BadInput(this->getId(), "bad junction number in role '{0}'", role);
167 }
168 }
169 }
170 return no;
171 }
172
174 size_t isActive(const RectangularMesh<3>::Element& element) const { return isActive(element.getMidpoint()); }
175
177 size_t isActive(const RectangularMaskedMesh<3>::Element& element) const { return isActive(element.getMidpoint()); }
178
179 public:
180 Convergence convergence;
181
182 double maxerr;
183 Vec<3, double> maxcur;
184
185 // Boundary conditions
186 BoundaryConditions<RectangularMesh<3>::Boundary, double> voltage_boundary;
187
188 typename ProviderFor<Voltage, Geometry3D>::Delegate outVoltage;
189
190 typename ProviderFor<CurrentDensity, Geometry3D>::Delegate outCurrentDensity;
191
192 typename ProviderFor<Heat, Geometry3D>::Delegate outHeat;
193
194 typename ProviderFor<Conductivity, Geometry3D>::Delegate outConductivity;
195
196 ReceiverFor<Temperature, Geometry3D> inTemperature;
197
198 ElectricalFem3DSolver(const std::string& name = "");
199
200 void loadConfiguration(XMLReader& source, Manager& manager) override;
201
202 void parseConfiguration(XMLReader& source, Manager& manager);
203
204 ~ElectricalFem3DSolver();
205
211 double compute(unsigned loops = 1);
212
219 double integrateCurrent(size_t vindex, bool onlyactive = false);
225 double getTotalCurrent(size_t nact = 0);
226
231 double getTotalEnergy();
232
237 double getCapacitance();
238
243 double getTotalHeat();
244
246 double getErr() const { return toterr; }
247
249 double getCondPcontact() const { return pcond; }
251 void setCondPcontact(double cond) {
252 pcond = cond;
253 this->invalidate();
254 }
255
257 double getCondNcontact() const { return ncond; }
259 void setCondNcontact(double cond) {
260 ncond = cond;
261 this->invalidate();
262 }
263
265 DataVector<const Tensor2<double>> getCondJunc() const { return junction_conductivity; }
267 void setCondJunc(double cond) {
268 junction_conductivity.reset(max(junction_conductivity.size(), size_t(1)), cond);
269 default_junction_conductivity = Tensor2<double>(0., cond);
270 }
272 void setCondJunc(Tensor2<double> cond) {
273 junction_conductivity.reset(max(junction_conductivity.size(), size_t(1)), cond);
274 default_junction_conductivity = cond;
275 }
277 void setCondJunc(const DataVector<Tensor2<double>>& cond) {
278 size_t condsize = 0;
279 for (const auto& act : active) condsize += (act.right - act.left) * act.ld;
280 condsize = max(condsize, size_t(1));
281 if (!this->mesh || cond.size() != condsize)
282 throw BadInput(this->getId(), "provided junction conductivity vector has wrong size");
283 junction_conductivity = cond.claim();
284 }
285
286 protected:
287 const LazyData<double> getVoltage(shared_ptr<const MeshD<3>> dest_mesh, InterpolationMethod method) const;
288
289 const LazyData<Vec<3>> getCurrentDensity(shared_ptr<const MeshD<3>> dest_mesh, InterpolationMethod method);
290
291 const LazyData<double> getHeatDensity(shared_ptr<const MeshD<3>> dest_mesh, InterpolationMethod method);
292
293 const LazyData<Tensor2<double>> getConductivity(shared_ptr<const MeshD<3>> dest_mesh, InterpolationMethod method);
294};
295
296}}} // namespace plask::electrical::shockley
297
298#endif