PLaSK library

Loading...
Searching...
No Matches
3d.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__VECTORCART3D_H
15#define PLASK__VECTORCART3D_H
16
21#include <iostream>
22
23#include "../math.hpp"
24#include <plask/exceptions.hpp>
25
26#include "common.hpp"
27
28#include "../utils/metaprog.hpp" // for is_callable
29#include "../utils/warnings.hpp"
30
31 namespace plask {
32
36 template <typename T>
37 struct Vec<3,T> {
38
39 static const int DIMS = 3;
40
42 T c0, c1, c2;
43
44 T& lon() { return c0; }
45 constexpr const T& lon() const { return c0; }
46
47 T& tran() { return c1; }
48 constexpr const T& tran() const { return c1; }
49
50 T& vert() { return c2; }
51 constexpr const T& vert() const { return c2; }
52
53 // radial coordinates
54 T& rad_p() { return c0; }
55 constexpr const T& rad_p() const { return c0; }
56
57 T& rad_r() { return c1; }
58 constexpr const T& rad_r() const { return c1; }
59
60 T& rad_z() { return c2; }
61 constexpr const T& rad_z() const { return c2; }
62
63 // for surface-emitting lasers (z-axis up)
64 T& se_x() { return c0; }
65 constexpr const T& se_x() const { return c0; }
66
67 T& se_y() { return c1; }
68 constexpr const T& se_y() const { return c1; }
69
70 T& se_z() { return c2; }
71 constexpr const T& se_z() const { return c2; }
72
73 // for surface-emitting lasers (z-axis up)
74 T& zup_x() { return c0; }
75 constexpr const T& z_up_x() const { return c0; }
76
77 T& zup_y() { return c1; }
78 constexpr const T& z_up_y() const { return c1; }
79
80 T& zup_z() { return c2; }
81 constexpr const T& z_up_z() const { return c2; }
82
83 // for edge emitting lasers (y-axis up), we keep the coordinates right-handed
84 T& ee_z() { return c0; }
85 constexpr const T& ee_z() const { return c0; }
86
87 T& ee_x() { return c1; }
88 constexpr const T& ee_x() const { return c1; }
89
90 T& ee_y() { return c2; }
91 constexpr const T& ee_y() const { return c2; }
92
93 // for edge emitting lasers (y-axis up), we keep the coordinates right-handed
94 T& yup_z() { return c0; }
95 constexpr const T& y_up_z() const { return c0; }
96
97 T& yup_x() { return c1; }
98 constexpr const T& y_up_x() const { return c1; }
99
100 T& yup_y() { return c2; }
101 constexpr const T& y_up_y() const { return c2; }
102
106 typedef T* iterator;
107
111 typedef const T* const_iterator;
112
114 Vec() {}
115
120 template <typename OtherT>
121 constexpr Vec(const Vec<3,OtherT>& p): c0(p.c0), c1(p.c1), c2(p.c2) {}
122
127 constexpr Vec(const T& c0__lon, const T& c1__tran, const T& c2__up): c0(c0__lon), c1(c1__tran), c2(c2__up) {}
128
133 template <typename T0, typename T1, typename T2>
134 constexpr Vec(const std::tuple<T0,T1,T2>& comp): c0(std::get<0>(comp)), c1(std::get<1>(comp)), c2(std::get<2>(comp)) {}
135
141 template <typename InputIteratorType>
142 static inline Vec<3,T> fromIterator(InputIteratorType inputIt) {
143 Vec<3,T> result;
144 result.c0 = T(*inputIt);
145 result.c1 = T(*++inputIt);
146 result.c2 = T(*++inputIt);
147 return result;
148 }
149
154 iterator begin() { return &c0; }
155
160 const_iterator begin() const { return &c0; }
161
166 iterator end() { return &c0 + 3; }
167
172 const_iterator end() const { return &c0 + 3; }
173
179 template <typename OtherT>
180 constexpr bool operator==(const Vec<3,OtherT>& p) const { return p.c0 == c0 && p.c1 == c1 && p.c2 == c2; }
181
188 template <typename OtherT, typename SuprType>
189 constexpr bool equals(const Vec<3, OtherT>& p, const SuprType& abs_supremum) const {
190 return is_zero(p.c0 - c0, abs_supremum) && is_zero(p.c1 - c1, abs_supremum) && is_zero(p.c2 - c2, abs_supremum); }
191
197 template <typename OtherT>
198 constexpr bool equals(const Vec<3, OtherT>& p) const {
199 return is_zero(p.c0 - c0) && is_zero(p.c1 - c1) && is_zero(p.c2 - c2);
200 }
201
207 template <typename OtherT>
208 constexpr bool operator!=(const Vec<3,OtherT>& p) const { return p.c0 != c0 || p.c1 != c1 || p.c2 != c2; }
209
216 inline T& operator[](size_t i) {
217 assert(i < 3);
218 return *(&c0 + i);
219 }
220
227 inline const T& operator[](size_t i) const {
228 assert(i < 3);
229 return *(&c0 + i);
230 }
231
237 template <typename OtherT>
238 constexpr auto operator+(const Vec<3,OtherT>& other) const -> Vec<3,decltype(c0 + other.c0)> {
239 return Vec<3,decltype(this->c0 + other.c0)>(c0 + other.c0, c1 + other.c1, c2 + other.c2);
240 }
241
247 Vec<3,T>& operator+=(const Vec<3,T>& other) {
248 c0 += other.c0;
249 c1 += other.c1;
250 c2 += other.c2;
251 return *this;
252 }
253
259 template <typename OtherT>
260 constexpr auto operator-(const Vec<3,OtherT>& other) const -> Vec<3,decltype(c0 - other.c0)> {
261 return Vec<3, decltype(this->c0 - other.c0)>(c0 - other. c0, c1 - other. c1, c2 - other.c2);
262 }
263
269 Vec<3,T>& operator-=(const Vec<3,T>& other) {
270 c0 -= other.c0;
271 c1 -= other.c1;
272 c2 -= other.c2;
273 return *this;
274 }
275
281 template <typename OtherT>
282 constexpr auto operator*(const OtherT scale) const -> Vec<3,decltype(c0*scale)> {
283 PLASK_NO_CONVERSION_WARNING_BEGIN
284 return Vec<3,decltype(c0*scale)>(c0 * scale, c1 * scale, c2 * scale);
285 PLASK_NO_WARNING_END
286 }
287
293 Vec<3,T>& operator*=(const T scalar) {
294 c0 *= scalar;
295 c1 *= scalar;
296 c2 *= scalar;
297 return *this;
298 }
299
305 constexpr Vec<3,T> operator/(const T scalar) const { return Vec<3,T>(c0 / scalar, c1 / scalar, c2 / scalar); }
306
312 Vec<3,T>& operator/=(const T scalar) {
313 c0 /= scalar;
314 c1 /= scalar;
315 c2 /= scalar;
316 return *this;
317 }
318
323 constexpr Vec<3,T> operator-() const {
324 return Vec<3,T>(-c0, -c1, -c2);
325 }
326
331 Vec<3,T> sqr() const {
332 return Vec<3,T>(c0*c0, c1*c1, c2*c2);
333 }
334
339 Vec<3,T>& sqr_inplace() {
340 c0 *= c0; c1 *= c1; c2 *= c2;
341 return *this;
342 }
343
348 Vec<3,T> sqrt() const {
349 return Vec<3,T>(std::sqrt(c0), std::sqrt(c1), std::sqrt(c2));
350 }
351
356 Vec<3,T>& sqrt_inplace() {
357 c0 = std::sqrt(c0); c1 = std::sqrt(c1); c2 = std::sqrt(c2);
358 return *this;
359 }
360
365 template <typename OtherT>
366 Vec<3,T> pow(OtherT a) const {
367 return Vec<3,T>(std::pow(c0, a), std::pow(c1, a), std::pow(c2, a));
368 }
369
375 inline void flip(size_t i) {
376 assert(i < 3);
377 operator[](i) = -operator[](i);
378 }
379
386 inline Vec<3,T> flipped(size_t i) {
387 Vec<3,T> res = *this;
388 res.flip(i);
389 return res;
390 }
391
398 friend inline std::ostream& operator<<(std::ostream& out, const Vec<3,T>& to_print) {
399 return out << '[' << to_print.c0 << ", " << to_print.c1 << ", " << to_print.c2 << ']';
400 }
401
409 template<class OT>
410 bool operator< (Vec<3, OT> const& v) const {
411 if (dbl_compare_lt(this->c0, v.c0)) return true;
412 if (dbl_compare_gt(this->c0, v.c0)) return false;
413 if (dbl_compare_lt(this->c1, v.c1)) return true;
414 if (dbl_compare_gt(this->c1, v.c1)) return false;
415 return dbl_compare_lt(this->c2, v.c2);
416 }
417
418
419};
420
426 template <typename T>
427 inline constexpr Vec<3,T> conj(const Vec<3,T>& v) { return Vec<3,T>(conj(v.c0), conj(v.c1), conj(v.c2)); }
428
435 template <typename T1, typename T2>
436 inline auto dot(const Vec<3,T1>& v1, const Vec<3,T2>& v2) -> decltype(v1.c0*v2.c0) {
437 return ::plask::fma(v1.c0, v2.c0, ::plask::fma(v1.c1, v2.c1, v1.c2 * v2.c2)); //MSVC needs ::plask::
438}
439
446//template <> //MSVC2015 doesn't support this specialization, and using overloding shouldn't have negative consequences
447 inline auto dot(const Vec<3,dcomplex>& v1, const Vec<3,double>& v2) -> decltype(v1.c0*v2.c0) {
448 return ::plask::fma(conj(v1.c0), v2.c0, ::plask::fma(conj(v1.c1), v2.c1, conj(v1.c2) * v2.c2)); //MSVC needs ::plask::
449}
450
457//template <> //MSVC2015 doesn't support this specialization, and using overloding shouldn't have negative consequences
458 inline auto dot(const Vec<3,dcomplex>& v1, const Vec<3,dcomplex>& v2) -> decltype(v1.c0*v2.c0) {
459 return ::plask::fma(conj(v1.c0), v2.c0, ::plask::fma(conj(v1.c1), v2.c1, conj(v1.c2) * v2.c2)); //MSVC needs ::plask::
460}
461
467 template <typename T>
468 inline constexpr Vec<3,T> vec(const T c0__lon, const T c1__tran, const T c2__up) {
469 return Vec<3,T>(c0__lon, c1__tran, c2__up);
470}
471
473 template <typename T>
474 struct NaNImpl<Vec<3,T>> {
475 static constexpr Vec<3,T> get() { return Vec<3,T>(NaN<T>(), NaN<T>(), NaN<T>()); }
476};
477
479 template <typename T>
480 struct ZeroImpl<Vec<3,T>> {
481 static constexpr Vec<3,T> get() { return Vec<3,T>(0., 0., 0.); }
482};
483
486 template <typename T>
487 inline bool is_zero(const Vec<3,T>& v) {
488 return is_zero(v.c0) && is_zero(v.c1) && v.c2;
489}
490
491 PLASK_API_EXTERN_TEMPLATE_SPECIALIZATION_STRUCT(Vec<3, double>)
492 PLASK_API_EXTERN_TEMPLATE_SPECIALIZATION_STRUCT(Vec<3, std::complex<double> >)
493
494} //namespace plask
495
496 namespace std {
497 template <typename T>
498 plask::Vec<3,T> sqrt(plask::Vec<3,T> vec) {
499 return vec.sqrt();
500 }
501
502 template <typename T, typename OtherT>
503 plask::Vec<3,T> pow(plask::Vec<3,T> vec, OtherT a) {
504 return vec.pow(a);
505 }
506}
507
508#if FMT_VERSION >= 90000
509 template <typename T> struct fmt::formatter<plask::Vec<3, T>> : ostream_formatter {};
510#endif
511
512
513#endif // PLASK__VECTORCART3D_H