Graph Framework
Loading...
Searching...
No Matches
dispersion.hpp
Go to the documentation of this file.
1//------------------------------------------------------------------------------
6//------------------------------------------------------------------------------
7//------------------------------------------------------------------------------
148//------------------------------------------------------------------------------
149
150#ifndef dispersion_h
151#define dispersion_h
152
153#include "newton.hpp"
154#include "equilibrium.hpp"
155
157namespace dispersion {
158//******************************************************************************
159// Z Function interface.
160//******************************************************************************
161//------------------------------------------------------------------------------
166//------------------------------------------------------------------------------
167 template<jit::float_scalar T, bool SAFE_MATH=false>
168 class z_function {
169 public:
170//------------------------------------------------------------------------------
175//------------------------------------------------------------------------------
176 virtual graph::shared_leaf<T, SAFE_MATH>
177 Z(graph::shared_leaf<T, SAFE_MATH> zeta) = 0;
178
180 typedef T base;
182 static constexpr bool safe_math = SAFE_MATH;
183 };
184
185//------------------------------------------------------------------------------
190//------------------------------------------------------------------------------
191 template<jit::complex_scalar T, bool SAFE_MATH=false>
192 class z_power_series final : public z_function<T, SAFE_MATH> {
193 public:
194//------------------------------------------------------------------------------
199//------------------------------------------------------------------------------
200 virtual graph::shared_leaf<T, SAFE_MATH>
201 Z(graph::shared_leaf<T, SAFE_MATH> zeta) {
202 auto zeta2 = zeta*zeta;
203 auto zeta4 = zeta2*zeta2;
204 auto zeta6 = zeta4*zeta2;
205 return graph::i<T>*std::sqrt(M_PI)/graph::exp(zeta2) -
206 2.0*(1.0 - 2.0/3.0*zeta2
207 + 4.0/15.0*zeta4
208 - 8.0/105.0*zeta6)*zeta;
209 }
210 };
211
212//------------------------------------------------------------------------------
217//------------------------------------------------------------------------------
218 template<jit::complex_scalar T, bool SAFE_MATH=false>
219 class z_erfi final : public z_function<T, SAFE_MATH> {
220 public:
221//------------------------------------------------------------------------------
226//------------------------------------------------------------------------------
227 virtual graph::shared_leaf<T, SAFE_MATH>
228 Z(graph::shared_leaf<T, SAFE_MATH> zeta) {
229 return -std::sqrt(M_PI)*graph::exp(-zeta*zeta)*(graph::erfi(zeta) -
230 graph::i<T>);
231 }
232 };
233
235 template<class Z>
236 concept z_func = std::is_base_of<z_function<typename Z::base, Z::safe_math>, Z>::value;
237
238//******************************************************************************
239// Common physics expressions.
240//******************************************************************************
241//------------------------------------------------------------------------------
253//------------------------------------------------------------------------------
254 template<jit::float_scalar T, bool SAFE_MATH=false>
255 static constexpr graph::shared_leaf<T, SAFE_MATH>
256 build_plasma_frequency(graph::shared_leaf<T, SAFE_MATH> n,
257 const T q,
258 const T m,
259 const T c,
260 const T epsilon0) {
261 return n*q*q/(epsilon0*m*c*c);
262 }
263
264//------------------------------------------------------------------------------
275//------------------------------------------------------------------------------
276 template<jit::float_scalar T, bool SAFE_MATH=false>
277 static constexpr graph::shared_leaf<T, SAFE_MATH>
278 build_cyclotron_frequency(const T q,
279 graph::shared_leaf<T, SAFE_MATH> b,
280 const T m,
281 const T c) {
282 return q*b/(m*c);
283 }
284
285//******************************************************************************
286// Dispersion function.
287//******************************************************************************
288//------------------------------------------------------------------------------
293//------------------------------------------------------------------------------
294 template<jit::float_scalar T, bool SAFE_MATH=false>
295 class dispersion_function {
296 public:
297//------------------------------------------------------------------------------
309//------------------------------------------------------------------------------
310 virtual graph::shared_leaf<T, SAFE_MATH>
311 D(graph::shared_leaf<T, SAFE_MATH> w,
312 graph::shared_leaf<T, SAFE_MATH> kx,
313 graph::shared_leaf<T, SAFE_MATH> ky,
314 graph::shared_leaf<T, SAFE_MATH> kz,
315 graph::shared_leaf<T, SAFE_MATH> x,
316 graph::shared_leaf<T, SAFE_MATH> y,
317 graph::shared_leaf<T, SAFE_MATH> z,
318 graph::shared_leaf<T, SAFE_MATH> t,
319 equilibrium::shared<T, SAFE_MATH> &eq) = 0;
320
322 typedef T base;
324 static constexpr bool safe_math = SAFE_MATH;
325 };
326
327//------------------------------------------------------------------------------
332//------------------------------------------------------------------------------
333 template<jit::float_scalar T, bool SAFE_MATH=false>
334 class stiff final : public dispersion_function<T, SAFE_MATH> {
335 public:
336//------------------------------------------------------------------------------
367//------------------------------------------------------------------------------
368 virtual graph::shared_leaf<T, SAFE_MATH>
380 };
381
382//------------------------------------------------------------------------------
387//------------------------------------------------------------------------------
388 template<jit::float_scalar T, bool SAFE_MATH=false>
389 class simple final : public dispersion_function<T, SAFE_MATH> {
390 public:
391//------------------------------------------------------------------------------
405//------------------------------------------------------------------------------
406 virtual graph::shared_leaf<T, SAFE_MATH>
407 D(graph::shared_leaf<T, SAFE_MATH> w,
408 graph::shared_leaf<T, SAFE_MATH> kx,
409 graph::shared_leaf<T, SAFE_MATH> ky,
410 graph::shared_leaf<T, SAFE_MATH> kz,
411 graph::shared_leaf<T, SAFE_MATH> x,
412 graph::shared_leaf<T, SAFE_MATH> y,
413 graph::shared_leaf<T, SAFE_MATH> z,
414 graph::shared_leaf<T, SAFE_MATH> t,
415 equilibrium::shared<T, SAFE_MATH> &eq) {
416 const T c = 1.0;
417
418 auto npar2 = kz*kz*c*c/(w*w);
419 auto nperp2 = (kx*kx + ky*ky)*c*c/(w*w);
420 return npar2 + nperp2 - c;
421 }
422 };
423
424//------------------------------------------------------------------------------
429//------------------------------------------------------------------------------
430 template<jit::float_scalar T, bool SAFE_MATH=false>
431 class physics : public dispersion_function<T, SAFE_MATH> {
432 protected:
433// Define some common constants.
435 const T epsilon0 = 8.8541878138E-12;
437 const T mu0 = M_PI*4.0E-7;
439 const T q = 1.602176634E-19;
441 const T me = 9.1093837015E-31;
443 const T c = static_cast<T> (1.0)/std::sqrt(epsilon0*mu0);
444 };
445
446//------------------------------------------------------------------------------
451//------------------------------------------------------------------------------
452 template<jit::float_scalar T, bool SAFE_MATH=false>
453 class bohm_gross final : public physics<T, SAFE_MATH> {
454 public:
455//------------------------------------------------------------------------------
471//------------------------------------------------------------------------------
472 virtual graph::shared_leaf<T, SAFE_MATH>
473 D(graph::shared_leaf<T, SAFE_MATH> w,
474 graph::shared_leaf<T, SAFE_MATH> kx,
475 graph::shared_leaf<T, SAFE_MATH> ky,
476 graph::shared_leaf<T, SAFE_MATH> kz,
477 graph::shared_leaf<T, SAFE_MATH> x,
478 graph::shared_leaf<T, SAFE_MATH> y,
479 graph::shared_leaf<T, SAFE_MATH> z,
480 graph::shared_leaf<T, SAFE_MATH> t,
481 equilibrium::shared<T, SAFE_MATH> &eq) {
482
483// Equilibrium quantities.
484 auto ne = eq->get_electron_density(x, y, z);
485 auto wpe2 = build_plasma_frequency(ne, physics<T, SAFE_MATH>::q,
486 physics<T, SAFE_MATH>::me,
487 physics<T, SAFE_MATH>::c,
488 physics<T, SAFE_MATH>::epsilon0);
489 auto te = eq->get_electron_temperature(x, y, z);
490// 2*1.602176634E-19 to convert eV to J.
491
492 auto vterm2 = static_cast<T> (2.0)*physics<T, SAFE_MATH>::q*te
493 / (physics<T, SAFE_MATH>::me *
494 physics<T, SAFE_MATH>::c *
495 physics<T, SAFE_MATH>::c);
496
497// Wave numbers should be parallel to B if there is a magnetic field. Otherwise
498// B should be zero.
499 auto b_vec = eq->get_magnetic_field(x, y, z);
500 auto k = kx*eq->get_esup1(x, y, z)
501 + ky*eq->get_esup2(x, y, z)
502 + kz*eq->get_esup3(x, y, z);
503 graph::shared_leaf<T, SAFE_MATH> kpara2;
504 if (b_vec->length()->is_match(graph::zero<T, SAFE_MATH> ())) {
505 kpara2 = k->dot(k);
506 } else {
507 auto b_hat = b_vec->unit();
508 auto kpara = b_hat->dot(k);
509 kpara2 = kpara*kpara;
510 }
511
512 return wpe2 + 3.0/2.0*kpara2*vterm2 - w*w;
513 }
514 };
515
516//------------------------------------------------------------------------------
521//------------------------------------------------------------------------------
522 template<jit::float_scalar T, bool SAFE_MATH=false>
523 class light_wave final : public physics<T, SAFE_MATH> {
524 public:
525//------------------------------------------------------------------------------
541//------------------------------------------------------------------------------
542 virtual graph::shared_leaf<T, SAFE_MATH>
543 D(graph::shared_leaf<T, SAFE_MATH> w,
544 graph::shared_leaf<T, SAFE_MATH> kx,
545 graph::shared_leaf<T, SAFE_MATH> ky,
546 graph::shared_leaf<T, SAFE_MATH> kz,
547 graph::shared_leaf<T, SAFE_MATH> x,
548 graph::shared_leaf<T, SAFE_MATH> y,
549 graph::shared_leaf<T, SAFE_MATH> z,
550 graph::shared_leaf<T, SAFE_MATH> t,
551 equilibrium::shared<T, SAFE_MATH> &eq) {
552
553// Equilibrium quantities.
554 auto ne = eq->get_electron_density(x, y, z);
555 auto wpe2 = build_plasma_frequency(ne,
556 physics<T, SAFE_MATH>::q,
557 physics<T, SAFE_MATH>::me,
558 physics<T, SAFE_MATH>::c,
559 physics<T, SAFE_MATH>::epsilon0);
560
561// Wave numbers should be parallel to B if there is a magnetic field. Otherwise
562// B should be zero.
563 assert(eq->get_magnetic_field(x, y, z)->length()->is_match(graph::zero<T, SAFE_MATH> ()) &&
564 "Expected equilibrium with no magnetic field.");
565
566 auto k = kx*eq->get_esup1(x, y, z)
567 + ky*eq->get_esup2(x, y, z)
568 + kz*eq->get_esup3(x, y, z);
569 auto k2 = k->dot(k);
570
571 return wpe2 + k2 - w*w;
572 }
573 };
574
575//------------------------------------------------------------------------------
580//------------------------------------------------------------------------------
581 template<jit::float_scalar T, bool SAFE_MATH=false>
582 class acoustic_wave final : public physics<T, SAFE_MATH> {
583 public:
584//------------------------------------------------------------------------------
600//------------------------------------------------------------------------------
601 virtual graph::shared_leaf<T, SAFE_MATH>
602 D(graph::shared_leaf<T, SAFE_MATH> w,
603 graph::shared_leaf<T, SAFE_MATH> kx,
604 graph::shared_leaf<T, SAFE_MATH> ky,
605 graph::shared_leaf<T, SAFE_MATH> kz,
606 graph::shared_leaf<T, SAFE_MATH> x,
607 graph::shared_leaf<T, SAFE_MATH> y,
608 graph::shared_leaf<T, SAFE_MATH> z,
609 graph::shared_leaf<T, SAFE_MATH> t,
610 equilibrium::shared<T, SAFE_MATH> &eq) {
611
612// Equilibrium quantities.
613 const T mi = eq->get_ion_mass(0);
614 auto te = eq->get_electron_temperature(x, y, z);
615 auto ti = eq->get_ion_temperature(0, x, y, z);
616 const T gamma = 3.0;
617 auto vs2 = (physics<T, SAFE_MATH>::q*te + gamma*physics<T, SAFE_MATH>::q*ti)
618 / (mi*physics<T, SAFE_MATH>::c*physics<T, SAFE_MATH>::c);
619
620// Wave numbers should be parallel to B if there is a magnetic field. Otherwise
621// B should be zero.
622 auto b_vec = eq->get_magnetic_field(x, y, z);
623 auto k = kx*eq->get_esup1(x, y, z)
624 + ky*eq->get_esup2(x, y, z)
625 + kz*eq->get_esup3(x, y, z);
626 graph::shared_leaf<T, SAFE_MATH> kpara2;
627 if (b_vec->length()->is_match(graph::zero<T, SAFE_MATH> ())) {
628 kpara2 = k->dot(k);
629 } else {
630 auto b_hat = b_vec->unit();
631 auto kpara = b_hat->dot(k);
632 kpara2 = kpara*kpara;
633 }
634
635 return kpara2*vs2 - w*w;
636 }
637 };
638
639//------------------------------------------------------------------------------
644//------------------------------------------------------------------------------
645 template<jit::float_scalar T, bool SAFE_MATH=false>
646 class gaussian_well final : public dispersion_function<T, SAFE_MATH> {
647 public:
648//------------------------------------------------------------------------------
662//------------------------------------------------------------------------------
663 virtual graph::shared_leaf<T, SAFE_MATH>
664 D(graph::shared_leaf<T, SAFE_MATH> w,
665 graph::shared_leaf<T, SAFE_MATH> kx,
666 graph::shared_leaf<T, SAFE_MATH> ky,
667 graph::shared_leaf<T, SAFE_MATH> kz,
668 graph::shared_leaf<T, SAFE_MATH> x,
669 graph::shared_leaf<T, SAFE_MATH> y,
670 graph::shared_leaf<T, SAFE_MATH> z,
671 graph::shared_leaf<T, SAFE_MATH> t,
672 equilibrium::shared<T, SAFE_MATH> &eq) {
673 const T c = 1.0;
674 auto well = c - 0.5*exp(-(x*x + y*y)/0.1);
675 auto npar2 = kz*kz*c*c/(w*w);
676 auto nperp2 = (kx*kx + ky*ky)*c*c/(w*w);
677 return npar2 + nperp2 - well;
678 }
679 };
680
681//------------------------------------------------------------------------------
686//------------------------------------------------------------------------------
687 template<jit::float_scalar T, bool SAFE_MATH=false>
688 class ion_cyclotron final : public physics<T, SAFE_MATH> {
689 public:
690//------------------------------------------------------------------------------
708//------------------------------------------------------------------------------
709 virtual graph::shared_leaf<T, SAFE_MATH>
710 D(graph::shared_leaf<T, SAFE_MATH> w,
711 graph::shared_leaf<T, SAFE_MATH> kx,
712 graph::shared_leaf<T, SAFE_MATH> ky,
713 graph::shared_leaf<T, SAFE_MATH> kz,
714 graph::shared_leaf<T, SAFE_MATH> x,
715 graph::shared_leaf<T, SAFE_MATH> y,
716 graph::shared_leaf<T, SAFE_MATH> z,
717 graph::shared_leaf<T, SAFE_MATH> t,
718 equilibrium::shared<T, SAFE_MATH> &eq) {
719// Equilibrium quantities.
720 const T mi = eq->get_ion_mass(0);
721
722 auto te = eq->get_electron_temperature(x, y, z);
723 auto ti = eq->get_ion_temperature(0, x, y, z);
724 const T gamma = 3.0;
725 auto vs2 = (physics<T, SAFE_MATH>::q*te +
726 gamma*physics<T, SAFE_MATH>::q*ti)
727 / (mi*physics<T, SAFE_MATH>::c *
728 physics<T, SAFE_MATH>::c);
729
730 auto b_vec = eq->get_magnetic_field(x, y, z);
731 auto wce = build_cyclotron_frequency(-physics<T, SAFE_MATH>::q,
732 b_vec->length(),
733 physics<T, SAFE_MATH>::me,
734 physics<T, SAFE_MATH>::c);
735
736// Wave numbers.
737 auto k = kx*eq->get_esup1(x, y, z)
738 + ky*eq->get_esup2(x, y, z)
739 + kz*eq->get_esup3(x, y, z);
740 auto b_hat = b_vec->unit();
741 auto kperp = b_hat->cross(k)->length();
742 auto kperp2 = kperp*kperp;
743
744 auto w2 = w*w;
745
746 return wce - kperp2*vs2 - w*w;
747 }
748 };
749
750//------------------------------------------------------------------------------
755//------------------------------------------------------------------------------
756 template<jit::float_scalar T, bool SAFE_MATH=false>
757 class ordinary_wave final : public physics<T, SAFE_MATH> {
758 public:
759//------------------------------------------------------------------------------
775//------------------------------------------------------------------------------
776 virtual graph::shared_leaf<T, SAFE_MATH>
777 D(graph::shared_leaf<T, SAFE_MATH> w,
778 graph::shared_leaf<T, SAFE_MATH> kx,
779 graph::shared_leaf<T, SAFE_MATH> ky,
780 graph::shared_leaf<T, SAFE_MATH> kz,
781 graph::shared_leaf<T, SAFE_MATH> x,
782 graph::shared_leaf<T, SAFE_MATH> y,
783 graph::shared_leaf<T, SAFE_MATH> z,
784 graph::shared_leaf<T, SAFE_MATH> t,
785 equilibrium::shared<T, SAFE_MATH> &eq) {
786// Equilibrium quantities.
787 auto ne = eq->get_electron_density(x, y, z);
788 auto wpe2 = build_plasma_frequency(ne,
789 physics<T, SAFE_MATH>::q,
790 physics<T, SAFE_MATH>::me,
791 physics<T, SAFE_MATH>::c,
792 physics<T, SAFE_MATH>::epsilon0);
793
794// Wave numbers.
795 auto n = (kx*eq->get_esup1(x, y, z) +
796 ky*eq->get_esup2(x, y, z) +
797 kz*eq->get_esup3(x, y, z))/w;
798 auto b_vec = eq->get_magnetic_field(x, y, z);
799 auto b_hat = b_vec->unit();
800 auto nperp = b_hat->cross(n);
801 auto nperp2 = nperp->dot(nperp);
802
803 auto w2 = w*w;
804
805 return 1.0 - wpe2/w2 - nperp2;
806 }
807 };
808
809//------------------------------------------------------------------------------
814//------------------------------------------------------------------------------
815 template<jit::float_scalar T, bool SAFE_MATH=false>
816 class extra_ordinary_wave final : public physics<T, SAFE_MATH> {
817 public:
818//------------------------------------------------------------------------------
839//------------------------------------------------------------------------------
840 virtual graph::shared_leaf<T, SAFE_MATH>
841 D(graph::shared_leaf<T, SAFE_MATH> w,
842 graph::shared_leaf<T, SAFE_MATH> kx,
843 graph::shared_leaf<T, SAFE_MATH> ky,
844 graph::shared_leaf<T, SAFE_MATH> kz,
845 graph::shared_leaf<T, SAFE_MATH> x,
846 graph::shared_leaf<T, SAFE_MATH> y,
847 graph::shared_leaf<T, SAFE_MATH> z,
848 graph::shared_leaf<T, SAFE_MATH> t,
849 equilibrium::shared<T, SAFE_MATH> &eq) {
850// Equilibrium quantities.
851 auto ne = eq->get_electron_density(x, y, z);
852 auto wpe2 = build_plasma_frequency(ne,
853 physics<T, SAFE_MATH>::q,
854 physics<T, SAFE_MATH>::me,
855 physics<T, SAFE_MATH>::c,
856 physics<T, SAFE_MATH>::epsilon0);
857
858 auto b_vec = eq->get_magnetic_field(x, y, z);
859 auto b_len = b_vec->length();
860 auto wec = build_cyclotron_frequency(-physics<T, SAFE_MATH>::q,
861 b_len,
862 physics<T, SAFE_MATH>::me,
863 physics<T, SAFE_MATH>::c);
864
865// Wave numbers.
866 auto n = (kx*eq->get_esup1(x, y, z) +
867 ky*eq->get_esup2(x, y, z) +
868 kz*eq->get_esup3(x, y, z))/w;
869 auto b_hat = b_vec->unit();
870 auto nperp = b_hat->cross(n);
871 auto nperp2 = nperp->dot(nperp);
872
873 auto wh = wpe2 + wec*wec;
874
875 auto w2 = w*w;
876
877 return 1.0 - wpe2/(w2)*(w2 - wpe2)/(w2 - wh) - nperp2;
878 }
879 };
880
881//------------------------------------------------------------------------------
886//------------------------------------------------------------------------------
887 template<jit::float_scalar T, bool SAFE_MATH=false>
888 class cold_plasma : public physics<T, SAFE_MATH> {
889 public:
890//------------------------------------------------------------------------------
926//------------------------------------------------------------------------------
927 virtual graph::shared_leaf<T, SAFE_MATH>
928 D(graph::shared_leaf<T, SAFE_MATH> w,
929 graph::shared_leaf<T, SAFE_MATH> kx,
930 graph::shared_leaf<T, SAFE_MATH> ky,
931 graph::shared_leaf<T, SAFE_MATH> kz,
932 graph::shared_leaf<T, SAFE_MATH> x,
933 graph::shared_leaf<T, SAFE_MATH> y,
934 graph::shared_leaf<T, SAFE_MATH> z,
935 graph::shared_leaf<T, SAFE_MATH> t,
936 equilibrium::shared<T, SAFE_MATH> &eq) {
937// Dielectric terms.
938// Frequencies
939 auto ne = eq->get_electron_density(x, y, z);
940 auto wpe2 = build_plasma_frequency(ne,
941 physics<T, SAFE_MATH>::q,
942 physics<T, SAFE_MATH>::me,
943 physics<T, SAFE_MATH>::c,
944 physics<T, SAFE_MATH>::epsilon0);
945 auto b_vec = eq->get_magnetic_field(x, y, z);
946 auto b_len = b_vec->length();
947 auto ec = build_cyclotron_frequency(-physics<T, SAFE_MATH>::q,
948 b_len,
949 physics<T, SAFE_MATH>::me,
950 physics<T, SAFE_MATH>::c);
951
952 auto w2 = w*w;
953 auto denome = 1.0 - ec*ec/w2;
954 auto e11 = 1.0 - (wpe2/w2)/denome;
955 auto e12 = ((ec/w)*(wpe2/w2))/denome;
956 auto e33 = wpe2;
957
958 for (size_t i = 0, ie = eq->get_num_ion_species(); i < ie; i++) {
959 const T mi = eq->get_ion_mass(i);
960 const T charge = static_cast<T> (eq->get_ion_charge(i))
961 * physics<T, SAFE_MATH>::q;
962
963 auto ni = eq->get_ion_density(i, x, y, z);
964 auto wpi2 = build_plasma_frequency(ni, charge, mi,
965 physics<T, SAFE_MATH>::c,
966 physics<T, SAFE_MATH>::epsilon0);
967 auto ic = build_cyclotron_frequency(charge, b_len, mi,
968 physics<T, SAFE_MATH>::c);
969
970 auto denomi = 1.0 - ic*ic/w2;
971 e11 = e11 - (wpi2/w2)/denomi;
972 e12 = e12 + ((ic/w)*(wpi2/w2))/denomi;
973 e33 = e33 + wpi2;
974 }
975
976 e12 = -1.0*e12;
977 e33 = 1.0 - e33/w2;
978
979// Wave numbers.
980 auto n = (kx*eq->get_esup1(x, y, z) +
981 ky*eq->get_esup2(x, y, z) +
982 kz*eq->get_esup3(x, y, z))/w;
983 auto b_hat = b_vec->unit();
984
985 auto npara = b_hat->dot(n);
986 auto npara2 = npara*npara;
987 auto nperp = b_hat->cross(n)->length();
988 auto nperp2 = nperp*nperp;
989
990// Determinate matrix elements
991 auto m11 = e11 - npara2;
992 auto m12 = e12;
993 auto m13 = npara*nperp;
994 auto m22 = e11 - npara2 - nperp2;
995 auto m33 = e33 - nperp2;
996
997 return (m11*m22 - m12*m12)*m33 - m22*(m13*m13);
998 }
999 };
1000
1001//------------------------------------------------------------------------------
1006//------------------------------------------------------------------------------
1007 template<jit::float_scalar T, bool SAFE_MATH=false>
1008 class cold_plasma_expansion : public physics<T, SAFE_MATH> {
1009 public:
1010//------------------------------------------------------------------------------
1033//------------------------------------------------------------------------------
1034 virtual graph::shared_leaf<T, SAFE_MATH>
1035 D(graph::shared_leaf<T, SAFE_MATH> w,
1036 graph::shared_leaf<T, SAFE_MATH> kx,
1037 graph::shared_leaf<T, SAFE_MATH> ky,
1038 graph::shared_leaf<T, SAFE_MATH> kz,
1039 graph::shared_leaf<T, SAFE_MATH> x,
1040 graph::shared_leaf<T, SAFE_MATH> y,
1041 graph::shared_leaf<T, SAFE_MATH> z,
1042 graph::shared_leaf<T, SAFE_MATH> t,
1043 equilibrium::shared<T, SAFE_MATH> &eq) {
1044// Setup plasma parameters.
1045 auto b_vec = eq->get_magnetic_field(x, y, z);
1046 auto b_len = b_vec->length();
1047 auto b_hat = b_vec/b_len;
1048 auto ne = eq->get_electron_density(x, y, z);
1049 auto te = eq->get_electron_temperature(x, y, z);
1050
1051 auto ve = graph::sqrt(static_cast<T> (2.0)*physics<T, SAFE_MATH>::q*te /
1052 physics<T, SAFE_MATH>::me)
1053 / physics<T, SAFE_MATH>::c;
1054
1055// Setup characteristic frequencies.
1056 auto ec = build_cyclotron_frequency(physics<T, SAFE_MATH>::q,
1057 b_len,
1058 physics<T, SAFE_MATH>::me,
1059 physics<T, SAFE_MATH>::c);
1060 auto wpe2 = build_plasma_frequency(ne, physics<T, SAFE_MATH>::q,
1061 physics<T, SAFE_MATH>::me,
1062 physics<T, SAFE_MATH>::c,
1063 physics<T, SAFE_MATH>::epsilon0);
1064
1065// Dispersion quantities.
1066 auto P = wpe2/(w*w);
1067 auto q = P/(2.0*(1.0 + ec/w));
1068
1069 auto n = (kx*eq->get_esup1(x, y, z) +
1070 ky*eq->get_esup2(x, y, z) +
1071 kz*eq->get_esup3(x, y, z))/w;
1072 auto n2 = n->dot(n);
1073 auto npara = n->dot(b_hat);
1074 auto npara2 = npara*npara;
1075 auto nperp = b_hat->cross(n);
1076 auto nperp2 = nperp->dot(nperp);
1077 auto n2nperp2 = n2*nperp2;
1078
1079 auto q_func = 1.0 - 2.0*q;
1080 auto n_func = n2 + npara2;
1081 auto p_func = 1.0 - P;
1082
1083 auto gamma1 = (1.0 - q)*n2nperp2
1084 + p_func*(n2*npara2 - (1.0 - q)*n_func)
1085 + q_func*(p_func - nperp2);
1086 auto gamma0 = nperp2*(n2 - 2.0*q_func) + p_func*(2.0*q_func - n_func);
1087
1088 return -P/2.0*(1.0 + ec/w)*gamma0 + (1.0 - ec*ec/(w*w))*gamma1;
1089 }
1090 };
1091
1092//------------------------------------------------------------------------------
1094//------------------------------------------------------------------------------
1095 template<jit::complex_scalar T, z_func Z, bool SAFE_MATH=false>
1096 class hot_plasma final : public physics<T, SAFE_MATH> {
1097 private:
1099 Z z;
1100
1101 public:
1102//------------------------------------------------------------------------------
1137//------------------------------------------------------------------------------
1138 virtual graph::shared_leaf<T, SAFE_MATH>
1139 D(graph::shared_leaf<T, SAFE_MATH> w,
1140 graph::shared_leaf<T, SAFE_MATH> kx,
1141 graph::shared_leaf<T, SAFE_MATH> ky,
1142 graph::shared_leaf<T, SAFE_MATH> kz,
1143 graph::shared_leaf<T, SAFE_MATH> x,
1144 graph::shared_leaf<T, SAFE_MATH> y,
1145 graph::shared_leaf<T, SAFE_MATH> z,
1146 graph::shared_leaf<T, SAFE_MATH> t,
1147 equilibrium::shared<T, SAFE_MATH> &eq) {
1148// Setup plasma parameters.
1149 auto b_vec = eq->get_magnetic_field(x, y, z);
1150 auto b_len = b_vec->length();
1151 auto b_hat = b_vec/b_len;
1152 auto ne = eq->get_electron_density(x, y, z);
1153 auto te = eq->get_electron_temperature(x, y, z);
1154
1155 auto ve = graph::sqrt(static_cast<T> (2.0)*physics<T, SAFE_MATH>::q*te /
1156 physics<T, SAFE_MATH>::me)
1157 / physics<T, SAFE_MATH>::c;
1158
1159// Setup characteristic frequencies.
1160 auto ec = build_cyclotron_frequency(physics<T, SAFE_MATH>::q,
1161 b_len,
1162 physics<T, SAFE_MATH>::me,
1163 physics<T, SAFE_MATH>::c);
1164 auto wpe2 = build_plasma_frequency(ne, physics<T, SAFE_MATH>::q,
1165 physics<T, SAFE_MATH>::me,
1166 physics<T, SAFE_MATH>::c,
1167 physics<T, SAFE_MATH>::epsilon0);
1168
1169// Dispersion quantities.
1170 auto P = wpe2/(w*w);
1171 auto q = P/(2.0*(1.0 + ec/w));
1172
1173 auto n = (kx*eq->get_esup1(x, y, z) +
1174 ky*eq->get_esup2(x, y, z) +
1175 kz*eq->get_esup3(x, y, z))/w;
1176 auto n2 = n->dot(n);
1177 auto npara = n->dot(b_hat);
1178 auto npara2 = npara*npara;
1179 auto nperp = b_hat->cross(n);
1180 auto nperp2 = nperp->dot(nperp);
1181
1182 auto zeta = (1.0 - ec/w)/(npara*ve);
1183 auto Z_func = this->z.Z(zeta);
1184 auto zeta_func = 1.0 + zeta*Z_func;
1185 auto F = ve*zeta*w/(2.0*npara*ec);
1186 auto isigma = P*Z_func/(2.0*npara*ve);
1187
1188 auto q_func = 1.0 - 2.0*q;
1189 auto n_func = n2 + npara2;
1190 auto p_func = 1.0 - P;
1191
1192 auto gamma5 = n2*npara2 - (1.0 - q)*n_func + q_func;
1193 auto gamma2 = (n2 - q_func)
1194 + P*w/(4.0*ec*npara2)*(n_func - 2.0*q_func);
1195 auto gamma1 = nperp2*((1.0 - q)*n2 - q_func)
1196 + p_func*(n2*npara2 - (1.0 - q)*n_func + q_func);
1197 auto gamma0 = nperp2*(n2 - 2.0*q_func) + p_func*(2.0*q_func - n_func);
1198
1199 return isigma*gamma0 + gamma1 + nperp2*P*w/ec*zeta_func*(gamma2 + gamma5*F);
1200 }
1201 };
1202
1203//------------------------------------------------------------------------------
1209//------------------------------------------------------------------------------
1210 template<jit::float_scalar T, class Z, bool SAFE_MATH=false>
1211 class hot_plasma_expansion final : public physics<T, SAFE_MATH> {
1212 private:
1214 Z z;
1215
1216 public:
1217//------------------------------------------------------------------------------
1250//------------------------------------------------------------------------------
1251 virtual graph::shared_leaf<T, SAFE_MATH>
1252 D(graph::shared_leaf<T, SAFE_MATH> w,
1253 graph::shared_leaf<T, SAFE_MATH> kx,
1254 graph::shared_leaf<T, SAFE_MATH> ky,
1255 graph::shared_leaf<T, SAFE_MATH> kz,
1256 graph::shared_leaf<T, SAFE_MATH> x,
1257 graph::shared_leaf<T, SAFE_MATH> y,
1258 graph::shared_leaf<T, SAFE_MATH> z,
1259 graph::shared_leaf<T, SAFE_MATH> t,
1260 equilibrium::shared<T, SAFE_MATH> &eq) {
1261// Setup plasma parameters.
1262 auto b_vec = eq->get_magnetic_field(x, y, z);
1263 auto b_hat = b_vec->unit();
1264 auto b_len = b_vec->length();
1265 auto ne = eq->get_electron_density(x, y, z);
1266 auto te = eq->get_electron_temperature(x, y, z);
1267
1268 auto ve = graph::sqrt(static_cast<T> (2.0)*physics<T, SAFE_MATH>::q*te /
1269 physics<T, SAFE_MATH>::me);
1270
1271// Setup characteristic frequencies.
1272 auto ec = build_cyclotron_frequency(physics<T, SAFE_MATH>::q, b_len,
1273 physics<T, SAFE_MATH>::me,
1274 physics<T, SAFE_MATH>::c);
1275 auto wpe2 = build_plasma_frequency(ne, physics<T, SAFE_MATH>::q,
1276 physics<T, SAFE_MATH>::me,
1277 physics<T, SAFE_MATH>::c,
1278 physics<T, SAFE_MATH>::epsilon0);
1279
1280// Dispersion quantities.
1281 auto P = wpe2/(w*w);
1282 auto q = P/(2.0*(1.0 + ec/w));
1283
1284 auto n = (kx*eq->get_esup1(x, y, z) +
1285 ky*eq->get_esup2(x, y, z) +
1286 kz*eq->get_esup3(x, y, z))/w;
1287 auto n2 = n->dot(n);
1288 auto npara = b_hat->dot(n);
1289 auto npara2 = npara*npara;
1290 auto nperp = b_hat->cross(n);
1291 auto nperp2 = nperp->dot(nperp);
1292
1293 auto vtnorm = ve/physics<T, SAFE_MATH>::c;
1294
1295 auto zeta = (1.0 - ec/w)/(npara*vtnorm);
1296 auto Z_func = this->z.Z(zeta);
1297
1298 auto q_func = 1.0 - 2.0*q;
1299 auto n_func = n2 + npara2;
1300 auto n2nperp2 = n2*nperp2;
1301 auto p_func = 1.0 - P;
1302
1303 auto gamma5 = P*(n2*npara2 - (1.0 - q)*n_func + q_func);
1304 auto gamma2 = P*w/ec*nperp2*(n2 - q_func)
1305 + P*P*w*w/(4.0*ec*ec)*(n_func - 2.0*q_func)*nperp2/npara2;
1306 auto gamma1 = (1.0 - q)*n2nperp2
1307 + p_func*(n2*npara2 - (1.0 - q)*n_func)
1308 + q_func*(p_func - nperp2);
1309
1310 return -(1.0 + ec/w)*npara*vtnorm *
1311 (gamma1 + gamma2 + nperp2/(2.0*npara)*(w*w/(ec*ec))*vtnorm*zeta*gamma5)*(1.0/Z_func + zeta);
1312 }
1313 };
1314
1315//******************************************************************************
1316// Dispersion interface.
1317//******************************************************************************
1319 template<class D>
1320 concept function = std::is_base_of<dispersion_function<typename D::base, D::safe_math>, D>::value;
1321
1322//------------------------------------------------------------------------------
1326//------------------------------------------------------------------------------
1327 template<function DISPERSION_FUNCTION>
1328 class dispersion_interface {
1329 protected:
1331 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1332 DISPERSION_FUNCTION::safe_math> D;
1333
1335 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1336 DISPERSION_FUNCTION::safe_math> dxdt;
1338 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1339 DISPERSION_FUNCTION::safe_math> dydt;
1341 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1342 DISPERSION_FUNCTION::safe_math> dzdt;
1344 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1345 DISPERSION_FUNCTION::safe_math> dkxdt;
1347 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1348 DISPERSION_FUNCTION::safe_math> dkydt;
1350 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1351 DISPERSION_FUNCTION::safe_math> dkzdt;
1353 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1354 DISPERSION_FUNCTION::safe_math> dsdt;
1355
1356 public:
1357//------------------------------------------------------------------------------
1369//------------------------------------------------------------------------------
1370 dispersion_interface(graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1371 DISPERSION_FUNCTION::safe_math> w,
1372 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1373 DISPERSION_FUNCTION::safe_math> kx,
1374 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1375 DISPERSION_FUNCTION::safe_math> ky,
1376 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1377 DISPERSION_FUNCTION::safe_math> kz,
1378 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1379 DISPERSION_FUNCTION::safe_math> x,
1380 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1381 DISPERSION_FUNCTION::safe_math> y,
1382 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1383 DISPERSION_FUNCTION::safe_math> z,
1384 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1385 DISPERSION_FUNCTION::safe_math> t,
1386 equilibrium::shared<typename DISPERSION_FUNCTION::base,
1387 DISPERSION_FUNCTION::safe_math> &eq) :
1388 D(DISPERSION_FUNCTION().D(w, kx, ky, kz, x, y, z, t, eq)) {
1389 auto dDdw = this->D->df(w);
1390 auto dDdkx = this->D->df(kx);
1391 auto dDdky = this->D->df(ky);
1392 auto dDdkz = this->D->df(kz);
1393 auto dDdx = this->D->df(x);
1394 auto dDdy = this->D->df(y);
1395 auto dDdz = this->D->df(z);
1396
1397 if (graph::pseudo_variable_cast(x).get()) {
1398 dDdw = dDdw->remove_pseudo();
1399 dDdkx = dDdkx->remove_pseudo();
1400 dDdky = dDdky->remove_pseudo();
1401 dDdkz = dDdkz->remove_pseudo();
1402 dDdx = dDdx->remove_pseudo();
1403 dDdy = dDdy->remove_pseudo();
1404 dDdz = dDdz->remove_pseudo();
1405 }
1406
1407 dxdt = -dDdkx/dDdw;
1408 dydt = -dDdky/dDdw;
1409 dzdt = -dDdkz/dDdw;
1410 dkxdt = dDdx/dDdw;
1411 dkydt = dDdy/dDdw;
1412 dkzdt = dDdz/dDdw;
1413
1414 dsdt = graph::vector(dxdt, dydt, dzdt)->length();
1415 }
1416
1417//------------------------------------------------------------------------------
1430//------------------------------------------------------------------------------
1431 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1432 DISPERSION_FUNCTION::safe_math>
1433 solve(graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1434 DISPERSION_FUNCTION::safe_math> x,
1435 graph::input_nodes<typename DISPERSION_FUNCTION::base,
1436 DISPERSION_FUNCTION::safe_math> inputs,
1437 const size_t index=0,
1438 const typename DISPERSION_FUNCTION::base tolerance = 1.0E-30,
1439 const size_t max_iterations = 1000) {
1440 auto x_var = graph::variable_cast(x);
1441
1442 workflow::manager<typename DISPERSION_FUNCTION::base,
1443 DISPERSION_FUNCTION::safe_math> work(index);
1444
1445 solver::newton(work, {x}, inputs, this->D,
1446 graph::shared_random_state<typename DISPERSION_FUNCTION::base,
1447 DISPERSION_FUNCTION::safe_math> (),
1448 tolerance, max_iterations);
1449
1450 work.compile();
1451 work.run();
1452
1453 work.copy_to_host(x, x_var->data());
1454
1455 return this->D*this->D;
1456 }
1457
1458//------------------------------------------------------------------------------
1462//------------------------------------------------------------------------------
1463 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1464 DISPERSION_FUNCTION::safe_math>
1465 get_residual() {
1466 return this->D*this->D;
1467 }
1468
1469//------------------------------------------------------------------------------
1473//------------------------------------------------------------------------------
1474 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1475 DISPERSION_FUNCTION::safe_math>
1476 get_d() {
1477 return this->D;
1478 }
1479
1480//------------------------------------------------------------------------------
1484//------------------------------------------------------------------------------
1485 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1486 DISPERSION_FUNCTION::safe_math>
1487 get_dsdt() {
1488 return this->dsdt;
1489 }
1490
1491//------------------------------------------------------------------------------
1495//------------------------------------------------------------------------------
1496 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1497 DISPERSION_FUNCTION::safe_math>
1498 get_dxdt() {
1499 return this->dxdt;
1500 }
1501
1502//------------------------------------------------------------------------------
1506//------------------------------------------------------------------------------
1507 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1508 DISPERSION_FUNCTION::safe_math>
1509 get_dydt() {
1510 return this->dydt;
1511 }
1512
1513//------------------------------------------------------------------------------
1517//------------------------------------------------------------------------------
1518 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1519 DISPERSION_FUNCTION::safe_math>
1520 get_dzdt() {
1521 return this->dzdt;
1522 }
1523
1524//------------------------------------------------------------------------------
1528//------------------------------------------------------------------------------
1529 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1530 DISPERSION_FUNCTION::safe_math>
1531 get_dkxdt() {
1532 return this->dkxdt;
1533 }
1534
1535//------------------------------------------------------------------------------
1539//------------------------------------------------------------------------------
1540 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1541 DISPERSION_FUNCTION::safe_math>
1542 get_dkydt() {
1543 return this->dkydt;
1544 }
1545
1546//------------------------------------------------------------------------------
1550//------------------------------------------------------------------------------
1551 graph::shared_leaf<typename DISPERSION_FUNCTION::base,
1552 DISPERSION_FUNCTION::safe_math>
1553 get_dkzdt() {
1554 return this->dkzdt;
1555 }
1556
1557//------------------------------------------------------------------------------
1559//------------------------------------------------------------------------------
1560 void print_dispersion() {
1561 D->to_latex();
1562 std::cout << std::endl;
1563 }
1564
1565//------------------------------------------------------------------------------
1567//------------------------------------------------------------------------------
1568 void print_dkxdt() {
1569 get_dkxdt()->to_latex();
1570 std::cout << std::endl;
1571 }
1572
1573//------------------------------------------------------------------------------
1575//------------------------------------------------------------------------------
1576 void print_dkydt() {
1577 get_dkydt()->to_latex();
1578 std::cout << std::endl;
1579 }
1580
1581//------------------------------------------------------------------------------
1583//------------------------------------------------------------------------------
1584 void print_dkzdt() {
1585 get_dkzdt()->to_latex();
1586 std::cout << std::endl;
1587 }
1588
1589//------------------------------------------------------------------------------
1591//------------------------------------------------------------------------------
1592 void print_dxdt() {
1593 get_dxdt()->to_latex();
1594 std::cout << std::endl;
1595 }
1596
1597//------------------------------------------------------------------------------
1599//------------------------------------------------------------------------------
1600 void print_dydt() {
1601 get_dydt()->to_latex();
1602 std::cout << std::endl;
1603 }
1604
1605//------------------------------------------------------------------------------
1607//------------------------------------------------------------------------------
1608 void print_dzdt() {
1609 get_dzdt()->to_latex();
1610 std::cout << std::endl;
1611 }
1612 };
1613}
1614
1615#endif /* dispersion_h */
dispersion::acoustic_wave
Ion wave dispersion function.
Definition dispersion.hpp:582
dispersion::acoustic_wave::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Ion acoustic wave function.
Definition dispersion.hpp:602
dispersion::bohm_gross
Bohm-Gross dispersion function.
Definition dispersion.hpp:453
dispersion::bohm_gross::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Bohm-Gross function.
Definition dispersion.hpp:473
dispersion::cold_plasma_expansion
Cold Plasma expansion dispersion function.
Definition dispersion.hpp:1008
dispersion::cold_plasma_expansion::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Cold Plasma expansion Dispersion function.
Definition dispersion.hpp:1035
dispersion::cold_plasma
Cold Plasma Dispersion function.
Definition dispersion.hpp:888
dispersion::cold_plasma::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Cold Plasma Dispersion function.
Definition dispersion.hpp:928
dispersion::dispersion_function
Interface for dispersion functions.
Definition dispersion.hpp:295
dispersion::dispersion_function::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)=0
Interface for a dispersion function.
dispersion::dispersion_function::base
T base
Type def to retrieve the backend base type.
Definition dispersion.hpp:322
dispersion::dispersion_function::safe_math
static constexpr bool safe_math
Retrieve template parameter of safe math.
Definition dispersion.hpp:324
dispersion::dispersion_interface
Class interface to build dispersion relation functions.
Definition dispersion.hpp:1328
dispersion::dispersion_interface::get_dydt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > get_dydt()
Provide right hand side for y update.
Definition dispersion.hpp:1509
dispersion::dispersion_interface::get_dkydt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > get_dkydt()
Provide right hand side for z update.
Definition dispersion.hpp:1542
dispersion::dispersion_interface::print_dxdt
void print_dxdt()
Print out the latex expression for the dxdt.
Definition dispersion.hpp:1592
dispersion::dispersion_interface::dkxdt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > dkxdt
Derivative with respect to kx.
Definition dispersion.hpp:1345
dispersion::dispersion_interface::dispersion_interface
dispersion_interface(graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > w, graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > kx, graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > ky, graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > kz, graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > x, graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > y, graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > z, graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > t, equilibrium::shared< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > &eq)
Construct a new dispersion_interface.
Definition dispersion.hpp:1370
dispersion::dispersion_interface::dxdt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > dxdt
Derivative with respect to kx.
Definition dispersion.hpp:1336
dispersion::dispersion_interface::print_dkxdt
void print_dkxdt()
Print out the latex expression for the dkxdt.
Definition dispersion.hpp:1568
dispersion::dispersion_interface::dydt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > dydt
Derivative with respect to ky.
Definition dispersion.hpp:1339
dispersion::dispersion_interface::get_dkxdt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > get_dkxdt()
Provide right hand side for z update.
Definition dispersion.hpp:1531
dispersion::dispersion_interface::print_dispersion
void print_dispersion()
Print out the latex expression for the dispersion relation.
Definition dispersion.hpp:1560
dispersion::dispersion_interface::dkydt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > dkydt
Derivative with respect to ky.
Definition dispersion.hpp:1348
dispersion::dispersion_interface::print_dkydt
void print_dkydt()
Print out the latex expression for the dkydt.
Definition dispersion.hpp:1576
dispersion::dispersion_interface::D
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > D
Dispersion function.
Definition dispersion.hpp:1332
dispersion::dispersion_interface::dzdt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > dzdt
Derivative with respect to kz.
Definition dispersion.hpp:1342
dispersion::dispersion_interface::get_dkzdt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > get_dkzdt()
Provide right hand side for z update.
Definition dispersion.hpp:1553
dispersion::dispersion_interface::get_dsdt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > get_dsdt()
Provide right hand side for s update.
Definition dispersion.hpp:1487
dispersion::dispersion_interface::print_dkzdt
void print_dkzdt()
Print out the latex expression for the dkzdt.
Definition dispersion.hpp:1584
dispersion::dispersion_interface::solve
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > solve(graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > x, graph::input_nodes< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > inputs, const size_t index=0, const typename DISPERSION_FUNCTION::base tolerance=1.0E-30, const size_t max_iterations=1000)
Solve the dispersion relation for x.
Definition dispersion.hpp:1433
dispersion::dispersion_interface::dkzdt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > dkzdt
Derivative with respect to kz.
Definition dispersion.hpp:1351
dispersion::dispersion_interface::dsdt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > dsdt
Derivative with respect to omega.
Definition dispersion.hpp:1354
dispersion::dispersion_interface::print_dydt
void print_dydt()
Print out the latex expression for the dydt.
Definition dispersion.hpp:1600
dispersion::dispersion_interface::get_d
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > get_d()
Get the dispersion function.
Definition dispersion.hpp:1476
dispersion::dispersion_interface::get_dzdt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > get_dzdt()
Provide right hand side for z update.
Definition dispersion.hpp:1520
dispersion::dispersion_interface::print_dzdt
void print_dzdt()
Print out the latex expression for the dzdt.
Definition dispersion.hpp:1608
dispersion::dispersion_interface::get_residual
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > get_residual()
Get the dispersion residual.
Definition dispersion.hpp:1465
dispersion::dispersion_interface::get_dxdt
graph::shared_leaf< typename DISPERSION_FUNCTION::base, DISPERSION_FUNCTION::safe_math > get_dxdt()
Provide right hand side for x update.
Definition dispersion.hpp:1498
dispersion::extra_ordinary_wave
Extra ordinary wave dispersion function.
Definition dispersion.hpp:816
dispersion::extra_ordinary_wave::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Dispersion relation for the X-Mode.
Definition dispersion.hpp:841
dispersion::gaussian_well
Gaussian Well dispersion function.
Definition dispersion.hpp:646
dispersion::gaussian_well::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Dispersion relation with a non uniform well.
Definition dispersion.hpp:664
dispersion::hot_plasma_expansion
Hot Plasma Expansion Dispersion function.
Definition dispersion.hpp:1211
dispersion::hot_plasma_expansion::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Hot plasma expansion dispersion function.
Definition dispersion.hpp:1252
dispersion::hot_plasma
Hot Plasma Dispersion function.
Definition dispersion.hpp:1096
dispersion::hot_plasma::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Hot Plasma Dispersion function.
Definition dispersion.hpp:1139
dispersion::ion_cyclotron
Electrostatic ion cyclotron wave dispersion function.
Definition dispersion.hpp:688
dispersion::ion_cyclotron::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Dispersion relation for the O mode.
Definition dispersion.hpp:710
dispersion::light_wave
Light Wave dispersion function.
Definition dispersion.hpp:523
dispersion::light_wave::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Light-wave function.
Definition dispersion.hpp:543
dispersion::ordinary_wave
Ordinary wave dispersion function.
Definition dispersion.hpp:757
dispersion::ordinary_wave::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Dispersion relation for the O mode.
Definition dispersion.hpp:777
dispersion::physics
Physics.
Definition dispersion.hpp:431
dispersion::physics::q
const T q
Fundamental charge.
Definition dispersion.hpp:439
dispersion::physics::me
const T me
Electron mass.
Definition dispersion.hpp:441
dispersion::physics::c
const T c
Speed of light.
Definition dispersion.hpp:443
dispersion::physics::mu0
const T mu0
Vacuum permeability.
Definition dispersion.hpp:437
dispersion::physics::epsilon0
const T epsilon0
Vacuum permittivity.
Definition dispersion.hpp:435
dispersion::simple
Simple dispersion function.
Definition dispersion.hpp:389
dispersion::simple::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Simple dispersion function.
Definition dispersion.hpp:407
dispersion::stiff
Stiff dispersion function.
Definition dispersion.hpp:334
dispersion::stiff::D
virtual graph::shared_leaf< T, SAFE_MATH > D(graph::shared_leaf< T, SAFE_MATH > w, graph::shared_leaf< T, SAFE_MATH > kx, graph::shared_leaf< T, SAFE_MATH > ky, graph::shared_leaf< T, SAFE_MATH > kz, graph::shared_leaf< T, SAFE_MATH > x, graph::shared_leaf< T, SAFE_MATH > y, graph::shared_leaf< T, SAFE_MATH > z, graph::shared_leaf< T, SAFE_MATH > t, equilibrium::shared< T, SAFE_MATH > &eq)
Stiff function.
Definition dispersion.hpp:369
dispersion::z_erfi
Class interface to build dispersion relation functions.
Definition dispersion.hpp:219
dispersion::z_erfi::Z
virtual graph::shared_leaf< T, SAFE_MATH > Z(graph::shared_leaf< T, SAFE_MATH > zeta)
Method to build the Z function.
Definition dispersion.hpp:228
dispersion::z_function
Class interface to build dispersion relation functions.
Definition dispersion.hpp:168
dispersion::z_function::safe_math
static constexpr bool safe_math
Retrieve template parameter of safe math.
Definition dispersion.hpp:182
dispersion::z_function::base
T base
Type def to retrieve the backend base type.
Definition dispersion.hpp:180
dispersion::z_function::Z
virtual graph::shared_leaf< T, SAFE_MATH > Z(graph::shared_leaf< T, SAFE_MATH > zeta)=0
Method to build the Z function.
dispersion::z_power_series
Class interface to build dispersion relation functions.
Definition dispersion.hpp:192
dispersion::z_power_series::Z
virtual graph::shared_leaf< T, SAFE_MATH > Z(graph::shared_leaf< T, SAFE_MATH > zeta)
Method to build the Z function.
Definition dispersion.hpp:201
workflow::manager
Class representing a workflow manager.
Definition workflow.hpp:171
dispersion::function
Dispersion concept.
Definition dispersion.hpp:1320
dispersion::z_func
Dispersion concept.
Definition dispersion.hpp:236
equilibrium.hpp
Class signature to implement plasma equilibrium.
assert
subroutine assert(test, message)
Assert check.
Definition f_binding_test.f90:38
dispersion
Name space for dispersion functions.
Definition dispersion.hpp:157
equilibrium::shared
std::shared_ptr< generic< T, SAFE_MATH > > shared
Convenience type alias for shared equilibria.
Definition equilibrium.hpp:472
graph::zero
constexpr shared_leaf< T, SAFE_MATH > zero()
Forward declare for zero.
Definition node.hpp:995
graph::vector
shared_vector< T, SAFE_MATH > vector(shared_leaf< T, SAFE_MATH > x, shared_leaf< T, SAFE_MATH > y, shared_leaf< T, SAFE_MATH > z)
Build a shared vector quantity.
Definition vector.hpp:133
graph::exp
shared_leaf< T, SAFE_MATH > exp(shared_leaf< T, SAFE_MATH > x)
Define exp convenience function.
Definition math.hpp:544
graph::shared_random_state
std::shared_ptr< random_state_node< T, SAFE_MATH > > shared_random_state
Convenience type alias for shared sqrt nodes.
Definition random.hpp:272
graph::input_nodes
std::vector< shared_variable< T, SAFE_MATH > > input_nodes
Convenience type alias for a vector of inputs.
Definition node.hpp:1731
graph::erfi
shared_leaf< T, SAFE_MATH > erfi(shared_leaf< T, SAFE_MATH > x)
Define erfi convenience function.
Definition math.hpp:1620
graph::variable_cast
shared_variable< T, SAFE_MATH > variable_cast(shared_leaf< T, SAFE_MATH > x)
Cast to a variable node.
Definition node.hpp:1747
graph::sqrt
shared_leaf< T, SAFE_MATH > sqrt(shared_leaf< T, SAFE_MATH > x)
Define sqrt convenience function.
Definition math.hpp:279
graph::shared_leaf
std::shared_ptr< leaf_node< T, SAFE_MATH > > shared_leaf
Convenience type alias for shared leaf nodes.
Definition node.hpp:674
graph::pseudo_variable_cast
shared_pseudo_variable< T, SAFE_MATH > pseudo_variable_cast(shared_leaf< T, SAFE_MATH > &x)
Cast to a pseudo variable node.
Definition node.hpp:1921
solver::newton
void newton(workflow::manager< T, SAFE_MATH > &work, graph::output_nodes< T, SAFE_MATH > vars, graph::input_nodes< T, SAFE_MATH > inputs, graph::shared_leaf< T, SAFE_MATH > func, graph::shared_random_state< T, SAFE_MATH > state, const T tolerance=1.0E-30, const size_t max_iterations=1000, const T step=1.0)
Determine the value of vars to minimize the loss function.
Definition newton.hpp:34
newton.hpp
Sets up the kernel for a newtons method.