GOAT (Geometrical optics application tool) 0.1
Loading...
Searching...
No Matches
supergrid.h
Go to the documentation of this file.
1
2#pragma once
3
4#include <stdlib.h>
5#include <fstream>
6#include "vector.h"
7#include "objectshape.h"
8#include <vector>
9
10namespace GOAT
11{
12 namespace raytracing
13 {
14 constexpr int SUPERGRID_NOERRORS = -1;
22 template<class T> class SuperGrid
23 {
24 public:
25 SuperGrid();
32 SuperGrid(double r0, int nx, int ny, int nz, const int typ = IN_HOST);
33 SuperGrid(double r0, int nx, int ny, int nz, ObjectShape** Obj, int numObjs, const bool isAbsolute = false, const int typ = 0);
34 ~SuperGrid() { clear(); };
35
40 void copy(const SuperGrid& S);
47 void add(const SuperGrid& S);
54 void sub(const SuperGrid& S);
55
56 bool addPlane(maths::Vector<double> P, maths::Vector<double> n, double d1, double d2, maths::Vector<double> cellSize);
66 bool addObject(ObjectShape** Obj, int numObj, const bool isAbsolute = false); // numObjs Einschluesse hinzufuegen
67 // isAbsolute=true => Orts-/Groessenangaben in absoluten Koordinaten
76 bool addObject(ObjectShape* E, const bool isAbsolute = false);
77 maths::Vector<int> gitterpunkt(maths::Vector<double> P);
78 bool inObject(maths::Vector<double> P, int i);
79 bool inObject(maths::Vector<int> Pi, int i);
80 bool inObject(int ix, int iy, int iz, int i);
81 T & operator () (int ix, int iy, int iz);
87 T& operator () (int i, int ix, int iy, int iz, bool isObjectCoordinate = true);
88 T& operator () (maths::Vector<int> Pi);
89 T& operator () (int i, maths::Vector<int> Pi);
90 T& operator () (maths::Vector<double> P);
91 T& operator () (int i, maths::Vector<double> P);
92
93 void fill(const T& x);
94 SuperGrid& operator = (const SuperGrid& S);
95 void clear();
96 void allockugel();
97
103 maths::Vector<int> kugelindex(maths::Vector<int> Pi);
104 T kugelwert(maths::Vector<int> Pi);
105 T kugelwert(int ix, int iy, int iz);
106 friend std::ostream& operator << (std::ostream& os, const SuperGrid& S);
107
108 /*
109 Hier kommen die eigentlichen Daten:
110 Ein : Array mit den Einschluessen (Orte und Ausdehnungen der Einschluesse werden absolut angegeben)
111 numObjs : Anzahl Einschluesse
112 type: Verteilung der aktiven Molek�le
113 0: nur im Einschluss
114 1: nur im Host
115 2: �berall
116 ywerte,
117 zwerte : Arrays zur Bestimmung der Dimensionen von K
118 r0 : Radius des Host-Partikels=halbe Breite des Uebergitters
119 G : Array mit den Untergittern, wird ben�tigt falls type 0 ist; ObjectShape : G[Einschlussindex][ix][iy][iz]
120 K : Kugelgitter, wird ben�tigt falls type 1 oder 2 ist
121 Pul : Array mit den Basispunkten der Untergitter
122 n : Array mit den Ausdehnungen der " n[Einschlussindex][nx][ny][nz] (nx,ny,nz sind als Vektor zusammengefasst)
123 nges : Ausdehnung des Uebergitters als Vektor (nx,ny,nz) zusammengefasst
124 d : Gitterkonstanten als Vektor (dx,dy,dz) zusammengefasst
125 */
126
127
128 std::vector<ObjectShape*> Obj;
129 int numObjs;
130 int type;
131 int* ywerte;
132 int** zwerte;
133 // *std::vector<std:vector<std::vector<T> > > grid;
134 std::vector<std::vector<std::vector<std::vector<T> > > > G;
135 std::vector<std::vector<std::vector<T> > > K;
136 T dummy;
137 std::vector<GOAT::maths::Vector<int> >Pul;
138 std::vector<GOAT::maths::Vector<int> >n;
139 GOAT::maths::Vector<int> nges;
140 maths::Vector<double> d;
141 double r0;
142 bool isequal;
143 bool iscleared;
144 maths::Vector<std::complex<double> > pc;
145 int Error;
146 maths::Matrix<double> H, R;
147 };
148
149/* -------------------- IMPLEMENTATION ---------------------- */
150
151
152
153
154 template <class T>
155 SuperGrid<T>::SuperGrid()
156 {
157 H = maths::unity();
158 R = maths::unity();
159 Error = SUPERGRID_NOERRORS;
160 numObjs = 0;
161 isequal = false;
162
163 type = IN_HOST;
164 ywerte = 0;
165 zwerte = 0;
166 }
167
168
169 template <class T> SuperGrid<T>::SuperGrid(double r0, int nx, int ny, int nz, const int typ)
170 {
171 Error = SUPERGRID_NOERRORS;
172 double b = 2.0 * r0;
173 // G = 0;
174 isequal = false;
175 numObjs = 0;
176 // Obj = 0;
177 ywerte = 0;
178 zwerte = 0;
179 // K = 0;
180 this->r0 = r0;
181
182
183
184 nges = maths::Vector<int>(nx, ny, nz);
185 d = maths::Vector<double>(b / nx, b / ny, b / nz);
186 iscleared = true;
187 type = typ;
188 if (type & IN_HOST)
189 allockugel();
190 }
191
192 template <class T> SuperGrid<T>::SuperGrid(double r0, int nx, int ny, int nz, ObjectShape** Obj, int numObjs, const bool isAbsolute, const int typ)
193 {
194 Error = SUPERGRID_NOERRORS;
195 double b = 2.0 * r0;
196 isequal = isAbsolute;
197 numObjs = 0;
198 G = 0;
199 this->r0 = r0;
200 nges = maths::Vector<int>(nx, ny, nz);
201 type = typ;
202 d = maths::Vector<double>(b / nx, b / ny, b / nz);
203 iscleared = true;
204 if (type & IN_HOST)
205 {
206 allockugel();
207 }
208 addObject(Obj, numObjs, isAbsolute);
209 }
210
211
212 template <class T> bool SuperGrid<T>::addObject(ObjectShape** Obj, int numObjs, const bool isAbsolute)
213 {
214 bool ok = true;
215 for (int i = 0; (i < numObjs) && (ok); i++) ok = addObject(Obj[i], isAbsolute);
216 return ok;
217 }
218
219 template <class T> maths::Vector<int> SuperGrid<T>::gitterpunkt(maths::Vector<double> P)
220 {
221 int ix, iy, iz, jx, jy, jz, jxh, jyh, anzx;
222 maths::Vector<int> pi, ph;
223 maths::Vector<double> h;
224 h = H * P + maths::Vector<double>(r0, r0, r0);
225 ph = maths::Vector<int>(floor(h[0] / d[0]), floor(h[1] / d[1]), floor(h[2] / d[2]));
226
227 if (type & IN_HOST)
228 {
229 pi = ph;
230 }
231 else
232 {
233 pi = ph;
234 }
235 return pi;
236 }
237
238 template<class T> bool SuperGrid<T>::addPlane(maths::Vector<double> P, maths::Vector<double> norm, double d1, double d2, maths::Vector<double> cellSize)
239 {
240 GOAT::maths::Vector<double> e1, e2;
241 if (abs(norm % GOAT::maths::ex) > 1E-5)
242 e1 = GOAT::maths::ex - (GOAT::maths::ex * norm) * norm;
243 else
244 e1 = GOAT::maths::ey - (GOAT::maths::ey * norm) * norm;
245 e1 = e1 / abs(e1);
246 e2 = norm % e1;
247 e2 = e2 / abs(e2);
248
249 maths::trafo(e1, e2, norm, H, R);
250
251
252 if (numObjs < 1) // Es ist der erste Einschluss der hinzugef�gt wird
253 {
254 G = std::vector < std::vector < std::vector <std::vector> > >(1);
255 Pul = std::vector < GOAT::maths::Vector<int> >(1);
256 n = std::vector < GOAT::maths::Vector<int> >(1);
257 Obj = std::vector <ObjectShape*>(1);
258 }
259 else
260 {
261 G.push_back(std::vector<std::vector<std::vector<T> > >(1));
262 Pul.push_back(GOAT::maths::Vector<int> (1));
263 n.push_back(GOAT::maths::Vector<int> (1));
264 Obj.push_back(ObjectShape*);
265 }
266 n[numObjs] = GOAT::maths::Vector<int>(d1/cellSize[0], d2/cellSize[1], 1);
267
268 GOAT::maths::Vector<double> h;
269 GOAT::maths::Vector<int> hn;
270 h = ceil(ediv(E->por, d)) - floor(ediv(E->pul, d));
271 hn = maths::Vector<int>((int)h[0], (int)h[1], (int)h[2]);
272 Pul[numObjs] = maths::Vector<int>((int)h[0], (int)h[1], (int)h[2]);
273
274 G[numObjs] = std::vector<std::vector<std::vector<T> > >(n[numObjs][0]);
275 for (int ix = 0; ix < n[numObjs][0] + 1; ix++)
276 {
277 G[numObjs][ix] = std::vector<std::vector<T> >(n[numObjs][1]);
278 for (int iy = 0; iy < n[numObjs][1] + 1; iy++)
279 {
280 G[numObjs][ix][iy] = std::vector<T>(n[numObjs][2]);
281 /*for (int iz = 0; iz < n[numObjs][2] + 1; iz++)
282 G[numObjs][ix][iy][iz] = maths::Vector<std::complex<double> >(0.0, 0.0, 0.0);*/
283 G[numObjs][ix][iy][0] = maths::Vector<std::complex<double> >(0.0, 0.0, 0.0);
284 }
285 }
286 numObjs++;
287
288 }
289
290 template<class T> bool SuperGrid<T>::addObject(ObjectShape* E, const bool isAbsolute)
291 {
292 int hmem, hmem2;
293 char str[500];
294 SysMemInfo smi;
295 MemInfo mi;
296 long int allocMem;
297 maths::Vector<int> pul, por, N, hn;
298 maths::Vector<double> h, O;
299 O = maths::Vector<double>(-r0, -r0, -r0);
300
301 smi = sysmem();
302 mi = memstat();
303 h = ceil(ediv(E->por, d)) - floor(ediv(E->pul, d));
304
305 hn = maths::Vector<int>((int)h[0], (int)h[1], (int)h[2]); // Gr��e des 3D-Gitters in die drei Koordinatenrichtungen
306
307 if (numObjs < 1) // Es ist der erste Einschluss der hinzugef�gt wird
308 {
309 G = std::vector < std::vector < std:vector <std::vector> > >(1);
310 Pul = std::vector < GOAT::maths::Vector<int> >(1);
311 n = std::vector < GOAT::maths::Vector<int> >(1);
312 Obj = std::vector <ObjectShape*>(1);
313 }
314 else
315 {
316 G.push_back(std::vector<std::vector<std::vector<T> > >(1));
317 Pul.push_back(GOAT::maths::Vector<int> (1));
318 n.push_back(GOAT::maths::Vector<int> (1));
319 Obj.push_back(ObjectShape *);
320 }
321
322 n[numObjs] = hn;
323
324 Obj[numObjs] = E;
325 h = floor(ediv(Obj[numObjs]->pul + maths::Vector<double>(r0, r0, r0), d));
326 Pul[numObjs] = maths::Vector<int>((int)h[0], (int)h[1], (int)h[2]); // ???????????
327 Ellipsoid* H = (Ellipsoid*)E;
328 if (E->isActive()) // Ist der Einschluss �berhaupt inelastisch aktiv ?
329 {
330 G[numObjs]=std::vector<std::vector<std::vector<T> > > (n[numObjs][0]);
331 for (int ix = 0; ix < n[numObjs][0] + 1; ix++)
332 {
333 G[numObjs][ix]=std::vector<std::vector<T> > (n[numObjs][1]);
334 for (int iy = 0; iy < n[numObjs][1] + 1; iy++)
335 {
336 G[numObjs][ix][iy] = std::vector<T>(n[numObjs][2]);
337 for (int iz = 0; iz < n[numObjs][2] + 1; iz++)
338 G[numObjs][ix][iy][iz] = maths::Vector<std::complex<double> >(0.0, 0.0, 0.0);
339 }
340 }
341 }
342 // else G[numObjs] = NULL;
343
344 int i = 0;
345 numObjs++;
346 iscleared = false;
347 return true;
348 }
349
350 template <class T> bool SuperGrid<T>::inObject(maths::Vector<double> P, int i) // da muss noch was gemacht werden
351 {
352 return Obj[i]->isInside(P);
353 }
354
355 template <class T> bool SuperGrid<T>::inObject(maths::Vector<int> Pi, int i)
356 {
357 return (Obj[i]->isInside(emult(Pi, d)));
358 }
359
360 template <class T> bool SuperGrid<T>::inObject(int ix, int iy, int iz, int i)
361 {
362 return (Obj[i]->isInside(emult(maths::Vector<double>(ix, iy, iz), d)));
363 }
364
365 template <class T> T& SuperGrid<T>::operator () (int ix, int iy, int iz)
366 {
367 maths::Vector<int> Pi = maths::Vector<int>(ix, iy, iz);
368 int i = 0;
369 bool found = false;
370 T dummy;
371 if (type == IN_OBJECT)
372 {
373 do
374 {
375 found = inObject(ix, iy, iz, i);
376 i++;
377 } while ((i < numObjs) && (!found));
378 i--;
379 if (!found) return dummy;
380 else
381 {
382 Pi = Pi - Pul[i];
383 Pi = kugelindex(Pi);
384 if (Error == SUPERGRID_NOERRORS)
385 return G[i][Pi[0]][Pi[1]][Pi[2]];
386 else return dummy;
387 }
388 }
389 else
390 {
391 Pi = kugelindex(Pi);
392 if (Error == SUPERGRID_NOERRORS)
393 return K[Pi[0]][Pi[1]][Pi[2]];
394 else return dummy;
395 }
396 }
397
398 template<class T> T& SuperGrid<T>::operator () (maths::Vector<int> Pi)
399 {
400 T dummy;
401 if (type == IN_OBJECT)
402 {
403 int i = 0;
404 bool found = false;
405 do
406 {
407 found = inObject(Pi, i);
408 i++;
409 } while ((i < numObjs) && (!found));
410 i--;
411
412 if (!found) return dummy;
413 else
414 {
415 Pi = Pi - Pul[i];
416 Pi = kugelindex(Pi);
417 if (Error == SUPERGRID_NOERRORS)
418 return G[i][Pi[0]][Pi[1]][Pi[2]];
419 else
420 return dummy;
421 }
422 }
423 else
424 {
425 Pi = kugelindex(Pi);
426 if (Error == SUPERGRID_NOERRORS)
427 return K[Pi[0]][Pi[1]][Pi[2]];
428 else
429 return dummy;
430 }
431 }
432
433 template <class T> T& SuperGrid<T>::operator () (int i, int ix, int iy, int iz, bool isEinKoord)
434 {
435 T dummy;
436 maths::Vector<int> Pi = maths::Vector<int>(ix, iy, iz);
437 if (type == IN_OBJECT)
438 {
439 if (G[i] == NULL) return dummy;
440 if (!isEinKoord) Pi = Pi - Pul[i];
441 return G[i][Pi[0]][Pi[1]][Pi[2]];
442 }
443 else
444 {
445 Pi = kugelindex(Pi);
446 if (Error == SUPERGRID_NOERRORS)
447 return K[Pi[0]][Pi[1]][Pi[2]];
448 else
449 return dummy;
450 }
451 }
452
453 template <class T> T& SuperGrid<T>::operator () (int i, maths::Vector<int> Pi)
454 {
455 dummy = maths::Vector<std::complex<double> >(0.0, 0.0, 0.0);
456
457 if (type == IN_OBJECT)
458 {
459 Pi = Pi - Pul[i];
460 if (G[i] == NULL) return dummy;
461 if (Pi[0] < 0) return dummy; //maths::Vector<std::complex<double> > (0,0,0);
462 if (Pi[1] < 0) return dummy; //maths::Vector<std::complex<double> > (0,0,0);
463 if (Pi[2] < 0) return dummy; // maths::Vector<std::complex<double> > (0,0,0);
464 if (Pi[0] > n[i][0])
465 {
466 Error = SUPERGITTER;
467 return dummy; //maths::Vector<std::complex<double> > (0,0,0);
468 }
469
470 if (Pi[1] > n[i][1])
471 {
472 Error = SUPERGITTER;
473 return dummy; //maths::Vector<std::complex<double> > (0,0,0);
474 }
475
476 if (Pi[2] > n[i][2])
477 {
478 Error = SUPERGITTER;
479 }
480
481 Error = SUPERGRID_NOERRORS;
482 return G[i][Pi[0]][Pi[1]][Pi[2]];
483 }
484 else
485 {
486 Pi = kugelindex(Pi);
487 if (Error == SUPERGRID_NOERRORS)
488 return K[Pi[0]][Pi[1]][Pi[2]];
489 else
490 return dummy;
491 }
492 }
493
494 template<class T> T& SuperGrid<T>::operator () (maths::Vector<double> P)
495 {
496 int i;
497 maths::Vector<int> Pi;
498 maths::Vector<double> h;
499 h = H * P - maths::Vector<double>(-r0, -r0, -r0);
500 T dummy;
501 for (i = 0; i < 3; i++)
502 Pi[i] = h[i] / d[i];
503
504 if (type == IN_OBJECT)
505 {
506 bool found = false;
507 i = 0;
508 do
509 {
510 found = inObject(Pi, i);
511 i++;
512 } while ((i < numObjs) && (!found));
513 i--;
514
515 if (!found) return dummy;
516 else
517 {
518 if (G[i] == NULL) return dummy;
519 Pi = Pi - Pul[i];
520 return G[i][Pi[0]][Pi[1]][Pi[2]];
521 }
522 }
523 else
524 {
525 Pi = kugelindex(Pi);
526 if (Error == SUPERGRID_NOERRORS)
527 return K[Pi[0]][Pi[1]][Pi[2]];
528 else return dummy;
529 }
530 }
531
532 template <class T> T& SuperGrid<T>::operator () (int i, maths::Vector<double> P)
533 {
534 maths::Vector<int> Pi;
535 maths::Vector<double> h;
536 h = H * P - maths::Vector<double>(-r0, -r0, -r0);
537 T dummy;
538
539 for (int j = 0; j < 3; j++)
540 Pi[j] = h[j] / d[j];
541
542 if (type == IN_OBJECT)
543 {
544 if (G[i] == NULL) return dummy;
545 Pi = Pi - Pul[i];
546 return G[i][Pi[0]][Pi[1]][Pi[2]];
547 }
548 else
549 {
550 Pi = kugelindex(Pi);
551 if (Error == SUPERGRID_NOERRORS)
552 return K[Pi[0]][Pi[1]][Pi[2]];
553 else return dummy;
554 }
555 }
556
557 template <class T> void SuperGrid<T>::clear()
558 {
559 /* int anzx, anzx2;
560 anzx = nges[0];
561 anzx2 = nges[0] / 2;
562
563 if (type == IN_OBJECT)
564 {
565 if (numObjs > 0)
566 {
567 for (int i = numObjs - 1; i >= 0; i--)
568 {
569 if (G[i] != NULL)
570 {
571 for (int ix = n[i][0]; ix >= 0; ix--)
572 {
573 for (int iy = n[i][1]; iy >= 0; iy--)
574 if (G[i][ix][iy] != 0)
575 free(G[i][ix][iy]);
576 if (G[i][ix] != 0)
577 {
578 free(G[i][ix]);
579 } // if
580 } // for ix
581 free(G[i]);
582 } // if (G[i]!=NULL)
583 } // for i
584 free(n);
585 free(Pul);
586 free(Obj);
587 numObjs = 0;
588 iscleared = true;
589 G = 0;
590 }
591 }
592 else
593 {
594 if (K != 0)
595 {
596 for (int k = 0; k < anzx2; k++)
597 {
598 for (int l = 0; l < ywerte[k]; l++)
599 {
600 delete[] K[anzx2 - 1 - k][ywerte[k] + l];
601 delete[] K[anzx2 + k][ywerte[k] + l];
602 delete[] K[anzx2 - 1 - k][ywerte[k] - 1 - l];
603 delete[] K[anzx2 + k][ywerte[k] - 1 - l];
604 }
605 delete[] K[anzx2 - 1 - k];
606 delete[] K[anzx2 + k];
607 delete[] zwerte[k];
608 }
609
610 delete[] ywerte;
611 delete[] zwerte;
612
613 if (numObjs > 0)
614 free(Obj);
615 }
616 K = 0;
617 zwerte = 0;
618 ywerte = 0;
619 numObjs = 0;
620 iscleared = true;
621 }*/
622 }
623
624 template <class T> void SuperGrid<T>::copy(const SuperGrid<T>& S) // MUSS DRINGEND GE�NDERT WERDEN !
625 {
626 for (int i = 0; i < numObjs; i++)
627 if (S.G[i].size()>0)
628 {
629 for (int ix = 0; ix < n[i][0]; ix++)
630 for (int iy = 0; iy < n[i][1]; iy++)
631 for (int iz = 0; iz < n[i][2]; iz++)
632 G[i][ix][iy][iz] = S.G[i][ix][iy][iz];
633 }
634 }
635
636 template<class T> SuperGrid<T>& SuperGrid<T>::operator = (const SuperGrid<T>& S)
637 {
638 int i = 0;
639 int anzx, anzx2;
640
641
642 if (this == &S) return *this;
643 clear();
644 r0 = S.r0;
645 d = S.d;
646 nges = S.nges;
647 type = S.type;
648
649 anzx = nges[0]; anzx2 = nges[0] / 2;
650
651 if (type == IN_OBJECT)
652 {
653 isequal = true; addObject(S.Obj, S.numObjs, true);
654 isequal = true;
655 if (S.numObjs != 0)
656 do
657 {
658 for (int ix = 0; ix < n[i][0]; ix++)
659 for (int iy = 0; iy < n[i][1]; iy++)
660 for (int iz = 0; iz < n[i][2]; iz++)
661 G[i][ix][iy][iz] = S.G[i][ix][iy][iz];
662 i++;
663 } while (i < S.numObjs);
664 isequal = false;
665 }
666 else
667 {
668 allockugel();
669 isequal = true; addObject(S.Obj, S.numObjs, true);
670 isequal = true;
671 for (int k = 0; k < anzx2; k++)
672 {
673 ywerte[k] = S.ywerte[k];
674 for (int l = 0; l < S.ywerte[k]; l++)
675 {
676 zwerte[k][l] = S.zwerte[k][l];
677 for (int m = 0; m < S.zwerte[k][l]; m++)
678 {
679 K[anzx2 - 1 - k][ywerte[k] - 1 - l][m] = S.K[anzx2 - 1 - k][ywerte[k] - 1 - l][m];
680 K[anzx2 + k][ywerte[k] - 1 - l][m] = S.K[anzx2 + k][ywerte[k] - 1 - l][m];
681 K[anzx2 - 1 - k][ywerte[k] + l][m] = S.K[anzx2 - 1 - k][ywerte[k] + l][m];
682 K[anzx2 + k][ywerte[k] + l][m] = S.K[anzx2 + k][ywerte[k] + l][m];
683 K[anzx2 - 1 - k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m] = S.K[anzx2 - 1 - k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m];
684 K[anzx2 + k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m] = S.K[anzx2 + k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m];
685 K[anzx2 - 1 - k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m] = S.K[anzx2 - 1 - k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m];
686 K[anzx2 + k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m] = S.K[anzx2 + k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m];
687 }
688 }
689 }
690 }
691 return *this;
692 }
693
694
695
696 template <class T> void SuperGrid<T>::add(const SuperGrid<T>& S)
697 {
698 int anzx2 = nges[0] / 2;
699 isequal = false;
700 H = S.H;
701 R = S.R;
702 if (type == IN_OBJECT)
703 {
704 for (int i = 0; i < S.numObjs; i++)
705 for (int ix = 0; ix < S.n[i][0]; ix++)
706 for (int iy = 0; iy < S.n[i][1]; iy++)
707 for (int iz = 0; iz < S.n[i][2]; iz++)
708 G[i][ix][iy][iz] += S.G[i][ix][iy][iz];
709 }
710 else
711 {
712 for (int k = 0; k < anzx2; k++)
713 {
714 for (int l = 0; l < ywerte[k]; l++)
715 {
716
717 for (int m = 0; m < zwerte[k][l]; m++)
718 {
719 K[anzx2 - 1 - k][ywerte[k] - 1 - l][m] += S.K[anzx2 - 1 - k][ywerte[k] - 1 - l][m];
720 K[anzx2 + k][ywerte[k] - 1 - l][m] += S.K[anzx2 + k][ywerte[k] - 1 - l][m];
721 K[anzx2 - 1 - k][ywerte[k] + l][m] += S.K[anzx2 - 1 - k][ywerte[k] + l][m];
722 K[anzx2 + k][ywerte[k] + l][m] += S.K[anzx2 + k][ywerte[k] + l][m];
723 K[anzx2 - 1 - k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m] += S.K[anzx2 - 1 - k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m];
724 K[anzx2 + k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m] += S.K[anzx2 + k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m];
725 K[anzx2 - 1 - k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m] += S.K[anzx2 - 1 - k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m];
726 K[anzx2 + k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m] += S.K[anzx2 + k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m];
727 }
728 }
729 }
730 }
731 }
732
733 template <class T> void SuperGrid<T>::sub(const SuperGrid<T>& S)
734 {
735 int anzx2 = nges[0] / 2;
736
737 isequal = false;
738 if (type == IN_OBJECT)
739 {
740 for (int i = 0; i < S.numObjs; i++)
741 for (int ix = 0; ix < S.n[i][0]; ix++)
742 for (int iy = 0; iy < S.n[i][1]; iy++)
743 for (int iz = 0; iz < S.n[i][2]; iz++)
744 G[i][ix][iy][iz] -= S.G[i][ix][iy][iz];
745 }
746 else
747 {
748 for (int k = 0; k < anzx2; k++)
749 {
750 for (int l = 0; l < ywerte[k]; l++)
751 {
752
753 for (int m = 0; m < zwerte[k][l]; m++)
754 {
755 K[anzx2 - 1 - k][ywerte[k] - 1 - l][m] -= S.K[anzx2 - 1 - k][ywerte[k] - 1 - l][m];
756 K[anzx2 + k][ywerte[k] - 1 - l][m] -= S.K[anzx2 + k][ywerte[k] - 1 - l][m];
757 K[anzx2 - 1 - k][ywerte[k] + l][m] -= S.K[anzx2 - 1 - k][ywerte[k] + l][m];
758 K[anzx2 + k][ywerte[k] + l][m] -= S.K[anzx2 + k][ywerte[k] + l][m];
759 K[anzx2 - 1 - k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m] -= S.K[anzx2 - 1 - k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m];
760 K[anzx2 + k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m] -= S.K[anzx2 + k][ywerte[k] - 1 - l][2 * zwerte[k][l] - 1 - m];
761 K[anzx2 - 1 - k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m] -= S.K[anzx2 - 1 - k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m];
762 K[anzx2 + k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m] -= S.K[anzx2 + k][ywerte[k] + l][2 * zwerte[k][l] - 1 - m];
763 }
764 }
765 }
766 }
767 }
768
769 template<class T> std::ostream& operator << (std::ostream& os, const SuperGrid<T> & S)
770 {
771 /*os << "r0=" << S.r0 << std::endl;
772 os << "Ausdehnung: nx=" << S.nges[0] << " ny=" << S.nges[1] << " nz=" << S.nges[2] << std::endl;
773 os << S.numObjs << " Einschluesse" << std::endl;
774 if (S.numObjs > 0)
775 for (int i = 0; i < S.numObjs; i++)
776 {
777 os << "====================== Einschluss Nr. " << i << " ======================" << std::endl;
778 for (int ix = 0; ix < S.n[i][0]; ix++)
779 for (int iy = 0; iy < S.n[i][1]; iy++)
780 for (int iz = 0; iz < S.n[i][2]; iz++)
781 os << abs(S.G[i][ix][iy][iz]) << std::endl;
782 }*/
783 return os;
784 }
785
786 template <class T> void SuperGrid<T>::fill(const T& x)
787 {
788 int anzx, anzx2;
789 anzx = nges[0]; anzx2 = nges[0] / 2;
790
791 if (type == IN_OBJECT)
792 {
793 for (int i = 0; i < numObjs; i++)
794 for (int ix = 0; ix < n[i][0]; ix++)
795 for (int iy = 0; iy < n[i][1]; iy++)
796 for (int iz = 0; iz < n[i][2]; iz++)
797 G[i][ix][iy][iz] = x;
798 }
799 else
800 {
801 for (int k = 0; k < anzx2; k++)
802 {
803 for (int l = 0; l < ywerte[k]; l++)
804 {
805 for (int m = 0; m < zwerte[k][l]; m++)
806 {
807 K[anzx2 - 1 - k][ywerte[k] - 1 - l][zwerte[k][l] - 1 - m] = x;
808 K[anzx2 + k][ywerte[k] - 1 - l][zwerte[k][l] - 1 - m] = x;
809 K[anzx2 - 1 - k][ywerte[k] + l][zwerte[k][l] - 1 - m] = x;
810 K[anzx2 + k][ywerte[k] + l][zwerte[k][l] - 1 - m] = x;
811 K[anzx2 - 1 - k][ywerte[k] - 1 - l][zwerte[k][l] + m] = x;
812 K[anzx2 + k][ywerte[k] - 1 - l][zwerte[k][l] + m] = x;
813 K[anzx2 - 1 - k][ywerte[k] + l][zwerte[k][l] + m] = x;
814 K[anzx2 + k][ywerte[k] + l][zwerte[k][l] + m] = x;
815 }
816 }
817 }
818 }
819 }
820
821
822 template <class T> void SuperGrid<T>::allockugel()
823 {
824 int anzx, anzx2;
825 double dx, dy, dz, hilf;
826 anzx = nges[0];
827 anzx2 = nges[0] / 2;
828 // cout << "allockugel" << std::endl;
829 ywerte = new int[anzx2];
830 zwerte = new int* [anzx2];
831
832 /*dx=d[0];
833 dy=d[1];
834 dz=d[2];*/
835 dx = 2.0 / nges[0];
836 dy = 2.0 / nges[1];
837 dz = 2.0 / nges[2];
838
839 // cout << "dx=" << dx << std::endl;
840 for (int i = 0; i < anzx2; i++)
841 {
842 ywerte[i] = int(ceil(sqrt(1.0 - (i * dx) * (i * dx)) / dy));
843 zwerte[i] = new int[anzx2];
844 for (int j = 0; j < anzx2; j++)
845 {
846 hilf = 1.0 - (i * dx) * (i * dx) - (j * dy) * (j * dy);
847 // cout << "hilf:" << hilf << std::endl;
848 if (hilf >= 0)
849 // zwerte[i][j]=int(ceil(sqrt(1.0-(i*dx)*(i*dx)-(j*dy)*(j*dy))*nges[2]/2.0));
850 zwerte[i][j] = ceil(sqrt(hilf) / dz);
851 else
852 zwerte[i][j] = 1;
853
854 // cout << "zwerte[" << i << "][" << j << "]:" << zwerte[i][j] << std::endl;
855 }
856 }
857 }
858
859
860
861 /* template<class T> void SuperGrid<T>::makeReal()
862 {
863 for (int i = 0; i < numObjs; i++)
864 for (int ix = 0; ix < n[i][0]; ix++)
865 for (int iy = 0; iy < n[i][1]; iy++)
866 for (int iz = 0; iz < n[i][2]; iz++)
867 for (int j = 0; j < 3; j++)
868 G[i][ix][iy][iz][j] = abs(G[i][ix][iy][iz][j]);
869 }
870 */
871 }
872}
GOAT::maths::Matrix
This class represents a threedimensional (numeric) Matrix as a template.
Definition matrix.h:23
GOAT::maths::Vector
Template class for threedimensional vectors.
Definition vector.h:57
GOAT::raytracing::Ellipsoid
This class represents an ellipsoid This class represents an ellipsoid, defined by its half axis a,...
Definition ellipsoid.h:21
GOAT::raytracing::ObjectShape
Abstract base class for all volume objects This abstract class provides a template for all volume obj...
Definition objectshape.h:60
GOAT::raytracing::ObjectShape::por
maths::Vector< double > por
corners of the circumferent cuboid (lower left corner and upper right corner)
Definition objectshape.h:125
GOAT::raytracing::ObjectShape::pul
maths::Vector< double > pul
Definition objectshape.h:125
GOAT::raytracing::ObjectShape::isActive
bool isActive()
returns true if the object should be considered for inelastic calculation
Definition objectshape.h:112
GOAT::raytracing::SuperGrid::numObjs
int numObjs
Number of objects.
Definition supergrid.h:129
GOAT::raytracing::SuperGrid::operator=
SuperGrid & operator=(const SuperGrid &S)
Assignment operator.
Definition supergrid.h:636
GOAT::raytracing::SuperGrid::operator()
T & operator()(int ix, int iy, int iz)
gives the content of the cell[ix][iy][iz]
Definition supergrid.h:365
GOAT::raytracing::SuperGrid::G
std::vector< std::vector< std::vector< std::vector< T > > > > G
Definition supergrid.h:134
GOAT::raytracing::SuperGrid::fill
void fill(const T &x)
Fill the whole SuperArray with value x.
Definition supergrid.h:786
GOAT::raytracing::SuperGrid::type
int type
Type of grid (0: Field is stored only in the objects, 1: field is stored in the surroundings only,...
Definition supergrid.h:130
GOAT::raytracing::SuperGrid::R
maths::Matrix< double > R
Definition supergrid.h:146
GOAT::raytracing::SuperGrid::add
void add(const SuperGrid &S)
Adds another SuperGrid object. In this function, all array elements from S are added to the existing ...
Definition supergrid.h:696
GOAT::raytracing::SuperGrid::pc
maths::Vector< std::complex< double > > pc
Definition supergrid.h:144
GOAT::raytracing::SuperGrid::sub
void sub(const SuperGrid &S)
Subtract another SuperArray object. In this function, all array elements from S are subtracted from t...
Definition supergrid.h:733
GOAT::raytracing::SuperGrid::H
maths::Matrix< double > H
Definition supergrid.h:146
GOAT::raytracing::SuperGrid::addObject
bool addObject(ObjectShape **Obj, int numObj, const bool isAbsolute=false)
Add objects to SuperGrid Adds a list of objects to the SuperGrid. An object is added to the SuperGrid...
Definition supergrid.h:212
GOAT::raytracing::SuperGrid::d
maths::Vector< double > d
Definition supergrid.h:140
GOAT::raytracing::SuperGrid::operator<<
friend std::ostream & operator<<(std::ostream &os, const SuperGrid &S)
Definition supergrid.h:769
GOAT::raytracing::SuperGrid::kugelwert
T kugelwert(maths::Vector< int > Pi)
GOAT::raytracing::SuperGrid::K
std::vector< std::vector< std::vector< T > > > K
Definition supergrid.h:135
GOAT::raytracing::SuperGrid::Pul
std::vector< GOAT::maths::Vector< int > > Pul
Definition supergrid.h:137
GOAT::raytracing::SuperGrid::kugelwert
T kugelwert(int ix, int iy, int iz)
GOAT::raytracing::SuperGrid::clear
void clear()
Clear the SuperGrid and release the allocated memory.
Definition supergrid.h:557
GOAT::raytracing::SuperGrid::nges
GOAT::maths::Vector< int > nges
Definition supergrid.h:139
GOAT::raytracing::SuperGrid::gitterpunkt
maths::Vector< int > gitterpunkt(maths::Vector< double > P)
Definition supergrid.h:219
GOAT::raytracing::SuperGrid::n
std::vector< GOAT::maths::Vector< int > > n
Definition supergrid.h:138
GOAT::raytracing::SuperGrid::copy
void copy(const SuperGrid &S)
Copies SuperArray object Here, S will be copied into the existing SuperArray and all elements will be...
Definition supergrid.h:624
GOAT::raytracing::SuperGrid::Obj
std::vector< ObjectShape * > Obj
List of the objects.
Definition supergrid.h:128
GOAT::raytracing::SuperGrid::addPlane
bool addPlane(maths::Vector< double > P, maths::Vector< double > n, double d1, double d2, maths::Vector< double > cellSize)
Definition supergrid.h:238
GOAT::raytracing::SuperGrid::kugelindex
maths::Vector< int > kugelindex(maths::Vector< int > Pi)
Checks if the point Pi (indicated by the indices) is inside the calculation sphere The function retur...
GOAT::raytracing::SuperGrid::inObject
bool inObject(maths::Vector< double > P, int i)
checks if P is inside the i-th object (p in real coordinates)
Definition supergrid.h:350
SUPERGITTER
#define SUPERGITTER
Definition error.h:27
GOAT::maths::trafo
void trafo(const Vector< double > &e0, const Vector< double > &e1, const Vector< double > &e2, Matrix< double > &H, Matrix< double > &R)
Calculates the transformation matrices from the laboratory coordinate system into a local system and ...
GOAT::maths::unity
Matrix< double > unity()
unity matrix (double precision)
GOAT::raytracing
Raytracer used for ultrashort pulse calculation with raytracing only.
Definition asphericLens.h:6
GOAT::raytracing::sysmem
SysMemInfo sysmem()
GOAT::raytracing::memstat
MemInfo memstat()
GOAT::raytracing::SUPERGRID_NOERRORS
constexpr int SUPERGRID_NOERRORS
Definition supergrid.h:14
GOAT::raytracing::operator<<
std::ostream & operator<<(std::ostream &os, Box B)
output operator for the Box class
GOAT::raytracing::clear
void clear(GlobalParms parms, maths::Vector< std::complex< double > > **Gitter)
GOAT
This class is used for the iray class. This class is intended for internal use only....
Definition fresnel.h:7
IN_OBJECT
#define IN_OBJECT
Definition resutil.h:39
IN_HOST
#define IN_HOST
Definition resutil.h:40
vector.h
This file contains the Vector template class and some useful functions around this class.