1 #ifndef COL_SPECIES_HPP
2 #define COL_SPECIES_HPP
5 #include <Kokkos_Core.hpp>
6 #include <Kokkos_DualView.hpp>
38 - sqrt(std::abs(1e-5 + pow(log(
T_avg)-2,2)/16));
40 double ti_ev,massi,densi,chargei,te_ev,masse,dense;
59 if(ti_ev*masse/massi<te_ev){
60 if(te_ev<1e1*(chargei*chargei)){
61 lambda = 2.3e1 - log(chargei*sqrt(dense*1e-6/(te_ev*te_ev*te_ev)));
63 lambda = 2.4e1 - log(sqrt(dense*1e-6)/te_ev);
66 lambda = 3e1-log((chargei*chargei)*sqrt(densi*1e-6/(ti_ev*ti_ev*ti_ev))/massi);
80 template<
class Device>
85 Kokkos::DualView<CollisionSpeciesScalars**,Device>
s;
87 Kokkos::DualView<double****,Kokkos::LayoutRight,Device>
pdf_n;
88 Kokkos::DualView<double****,Kokkos::LayoutRight,Device>
pdf_np1;
98 int col_species_cnt=0;
99 for (
int isp = 0; isp<all_species.size(); isp++){
100 if (!all_species[isp].is_adiabatic){
101 for (
int mesh_ind = 0; mesh_ind<mb_n_nodes; mesh_ind++){
102 s.h_view(col_species_cnt, mesh_ind).gid = isp;
103 s.h_view(col_species_cnt, mesh_ind).is_electron = all_species[isp].is_electron;
104 s.h_view(col_species_cnt, mesh_ind).grid_index = all_species[isp].collision_grid_index;
116 s.h_view(isp,mesh_ind).mass = species.
mass;
117 s.h_view(isp,mesh_ind).mass_au =
s.h_view(isp,mesh_ind).mass/
PROTON_MASS;
123 double conv_fac_temp = species.
f0.fg_temp_ev_h(local_node_ind);
124 s.h_view(isp,mesh_ind).conv_factor = 1.0/sqrt(conv_fac_temp*(conv_inv_mass*conv_inv_mass*conv_inv_mass));
128 for (
int imu = 0; imu <
pdf_n.h_view.extent(3); imu++){
129 for (
int ivp = 0; ivp<
pdf_n.h_view.extent(2); ivp++){
135 for (
int imu = 0; imu <
pdf_n.h_view.extent(3); imu++){
136 double vnorm_to_SI =
s.h_view(isp,mesh_ind).conv_factor/
f.
get_smu_n(imu);
137 for (
int ivp = 0; ivp<
pdf_n.h_view.extent(2); ivp++){
138 pdf_n.h_view(mesh_ind,isp,ivp,imu) *= vnorm_to_SI;
143 for (
int iz = 0; iz <
pdf_np1.h_view.extent(2); iz++){
144 for (
int ir = 0; ir <
pdf_np1.h_view.extent(3); ir++){
145 pdf_np1.h_view(mesh_ind,isp,iz,ir) =
pdf_n.h_view(mesh_ind,isp,iz,ir);
152 for (
int mesh_ind = 0; mesh_ind<
s.extent(1); mesh_ind++){
153 for (
int isp = 0; isp <
s.extent(0); isp++){
154 setup_one(isp, mesh_ind, all_species[
s.h_view(isp,mesh_ind).gid], moments, mesh_nodes(mesh_batch_ind,mesh_ind));
160 Kokkos::deep_copy(
s.d_view,
s.h_view);
161 Kokkos::deep_copy(
pdf_n.d_view,
pdf_n.h_view);
165 void lambda_gamma_pair(
const View<double*,CLayout,HostType>& dt, Kokkos::View<double***,Device>& gammac){
168 auto gammac_h = Kokkos::create_mirror_view(gammac);
170 for (
int mesh_ind = 0; mesh_ind<
s.extent(1); mesh_ind++){
172 for (
int isp = 0; isp <
s.extent(0); isp++){
173 for (
int jsp = 0; jsp <
s.extent(0); jsp++){
174 gammac_h(mesh_ind,jsp,isp) = dt(mesh_ind)*
s.h_view(isp,mesh_ind).lambda_gamma_pair(
s.h_view(jsp,mesh_ind));
178 Kokkos::deep_copy(gammac, gammac_h);
181 inline int n()
const{
KOKKOS_INLINE_FUNCTION int n_vr() const
Definition: vgrid_distribution.hpp:155
constexpr double PROTON_MASS
Definition: constants.hpp:7
Distribution< Device > f0
Species distribution in velocity space on local mesh nodes.
Definition: species.hpp:127
View< double **, CLayout, HostType > temp_local
temperature, , local to compute node [Nsp,Nlocalnode]
Definition: moments.hpp:17
Kokkos::DualView< double ****, Kokkos::LayoutRight, Device > pdf_np1
Definition: col_species.hpp:88
double mass_au
Definition: col_species.hpp:17
double get_smu_n(int imu) const
Definition: vgrid_distribution.hpp:172
double den_lambda_gamma
Definition: col_species.hpp:19
bool is_electron
Definition: col_species.hpp:14
Definition: moments.hpp:13
int gid
Definition: col_species.hpp:12
void setup_one(int isp, int mesh_ind, const Species< Device > &species, const Moments &moments, int local_node_ind)
Definition: col_species.hpp:115
double conv_factor
Definition: col_species.hpp:20
int nonadiabatic_idx
Index of species skipping adiabatic species (for compatibility with fortran arrays) ...
Definition: species.hpp:81
Kokkos::DualView< CollisionSpeciesScalars **, Device > s
Definition: col_species.hpp:85
Definition: col_species.hpp:11
KOKKOS_INLINE_FUNCTION int n_species() const
Definition: vgrid_distribution.hpp:151
double mass
Particle mass.
Definition: species.hpp:83
double numeric_t
Definition: col_species.hpp:22
KOKKOS_INLINE_FUNCTION int n_vz() const
Definition: vgrid_distribution.hpp:163
Kokkos::DualView< double ****, Kokkos::LayoutRight, Device > pdf_n
Definition: col_species.hpp:87
void lambda_gamma_pair(const View< double *, CLayout, HostType > &dt, Kokkos::View< double ***, Device > &gammac)
Definition: col_species.hpp:165
double mass
Definition: col_species.hpp:16
double charge
Particle charge.
Definition: species.hpp:84
View< double **, CLayout, HostType > den_local
density, local to compute node [Nsp,Nlocalnode]
Definition: moments.hpp:14
double charge_eu
Definition: col_species.hpp:18
void setup_all(const std::vector< Species< Device >> &all_species, const Moments &moments, Kokkos::View< int **, HostType > &mesh_nodes, int mesh_batch_ind)
Definition: col_species.hpp:150
constexpr double UNIT_CHARGE
Charge of an electron (C)
Definition: constants.hpp:4
double numeric_vth2
Definition: col_species.hpp:21
int grid_index
Definition: col_species.hpp:13
int n() const
Definition: col_species.hpp:181
double T_avg
Definition: col_species.hpp:15
double dens
Definition: col_species.hpp:23
double mom
Definition: col_species.hpp:25
const VGridDistribution< HostType > f
Definition: col_species.hpp:83
Definition: species.hpp:75
Definition: col_species.hpp:81
constexpr double PI
Definition: constants.hpp:8
constexpr double TWOPI
Definition: constants.hpp:9
double ens
Definition: col_species.hpp:24
KOKKOS_INLINE_FUNCTION double lambda_gamma_pair(const CollisionSpeciesScalars &sp_b) const
Definition: col_species.hpp:28