1 #ifndef COL_SPECIES_HPP
2 #define COL_SPECIES_HPP
5 #include <Kokkos_Core.hpp>
6 #include <Kokkos_DualView.hpp>
40 - sqrt(std::abs(1e-5 + pow(log(
T_avg)-2,2)/16));
42 double ti_ev,massi,densi,chargei,te_ev,masse,dense;
61 if(ti_ev*masse/massi<te_ev){
62 if(te_ev<1e1*(chargei*chargei)){
63 lambda = 2.3e1 - log(chargei*sqrt(dense*1e-6/(te_ev*te_ev*te_ev)));
65 lambda = 2.4e1 - log(sqrt(dense*1e-6)/te_ev);
68 lambda = 3e1-log((chargei*chargei)*sqrt(densi*1e-6/(ti_ev*ti_ev*ti_ev))/massi);
82 template<
class Device>
87 Kokkos::DualView<CollisionSpeciesScalars**,Device>
s;
89 Kokkos::DualView<double****,Kokkos::LayoutRight,Device>
pdf_n;
90 Kokkos::DualView<double****,Kokkos::LayoutRight,Device>
pdf_np1;
100 int col_species_cnt=0;
101 for (
int isp = 0; isp<all_species.size(); isp++){
102 if (!all_species[isp].is_adiabatic){
103 for (
int mesh_ind = 0; mesh_ind<mb_n_nodes; mesh_ind++){
104 s.h_view(col_species_cnt, mesh_ind).gid = isp;
105 s.h_view(col_species_cnt, mesh_ind).is_electron = all_species[isp].is_electron;
106 s.h_view(col_species_cnt, mesh_ind).grid_index = all_species[isp].collision_grid_index;
118 s.h_view(isp,mesh_ind).mass = species.
mass;
119 s.h_view(isp,mesh_ind).mass_au =
s.h_view(isp,mesh_ind).mass/
PROTON_MASS;
123 s.h_view(isp,mesh_ind).den_lambda_gamma =
get_den_at_node(local_node_ind + first_node,
s.h_view(isp,mesh_ind).gid);
124 s.h_view(isp,mesh_ind).T_avg =
get_t_ev_at_node(local_node_ind + first_node,
s.h_view(isp,mesh_ind).gid);
125 s.h_view(isp,mesh_ind).conv_factor = 1.0/sqrt(
s.h_view(isp,mesh_ind).T_avg*(conv_inv_mass*conv_inv_mass*conv_inv_mass));
129 for (
int imu = 0; imu <
pdf_n.h_view.extent(3); imu++){
131 for (
int invp = 0; invp<
pdf_n.h_view.extent(2); invp++){
132 pdf_n.h_view(mesh_ind,isp,invp,imu) = std::max(
f(species.
nonadiabatic_idx, imu,local_node_ind,invp),0.0)*
s.h_view(isp,mesh_ind).conv_factor/smu_n;
137 for (
int iz = 0; iz <
pdf_np1.h_view.extent(2); iz++){
138 for (
int ir = 0; ir <
pdf_np1.h_view.extent(3); ir++){
139 pdf_np1.h_view(mesh_ind,isp,iz,ir) =
pdf_n.h_view(mesh_ind,isp,iz,ir);
144 void setup_all(
const std::vector<
Species<Device>>& all_species, Kokkos::View<int**,HostType>& mesh_nodes,
int mesh_batch_ind,
int first_node){
146 for (
int mesh_ind = 0; mesh_ind<
s.extent(1); mesh_ind++){
147 for (
int isp = 0; isp <
s.extent(0); isp++){
148 setup_one(isp, mesh_ind, all_species[
s.h_view(isp,mesh_ind).gid], mesh_nodes(mesh_batch_ind,mesh_ind), first_node);
154 Kokkos::deep_copy(
s.d_view,
s.h_view);
155 Kokkos::deep_copy(
pdf_n.d_view,
pdf_n.h_view);
162 auto gammac_h = Kokkos::create_mirror_view(gammac);
164 for (
int mesh_ind = 0; mesh_ind<
s.extent(1); mesh_ind++){
166 for (
int isp = 0; isp <
s.extent(0); isp++){
167 for (
int jsp = 0; jsp <
s.extent(0); jsp++){
168 gammac_h(mesh_ind,jsp,isp) = dt*
s.h_view(isp,mesh_ind).lambda_gamma_pair(
s.h_view(jsp,mesh_ind));
172 Kokkos::deep_copy(gammac, gammac_h);
175 inline int n()
const{
void setup_all(const std::vector< Species< Device >> &all_species, Kokkos::View< int **, HostType > &mesh_nodes, int mesh_batch_ind, int first_node)
Definition: col_species.hpp:144
Kokkos::DualView< double ****, Kokkos::LayoutRight, Device > pdf_np1
Definition: col_species.hpp:90
double mass_au
Definition: col_species.hpp:19
double get_smu_n(int imu) const
Definition: vgrid_distribution.hpp:63
double den_lambda_gamma
Definition: col_species.hpp:21
bool is_electron
Definition: col_species.hpp:16
int gid
Definition: col_species.hpp:14
double conv_factor
Definition: col_species.hpp:22
int nonadiabatic_idx
Index of species skipping adiabatic species (for compatibility with fortran arrays) ...
Definition: species.hpp:77
Kokkos::DualView< CollisionSpeciesScalars **, Device > s
Definition: col_species.hpp:87
Definition: col_species.hpp:13
void setup_one(int isp, int mesh_ind, const Species< Device > &species, int local_node_ind, int first_node)
Definition: col_species.hpp:117
double mass
Particle mass.
Definition: species.hpp:79
double numeric_t
Definition: col_species.hpp:24
Kokkos::DualView< double ****, Kokkos::LayoutRight, Device > pdf_n
Definition: col_species.hpp:89
double mass
Definition: col_species.hpp:18
double charge
Particle charge.
Definition: species.hpp:80
constexpr double PROTON_MASS
Definition: globals.hpp:121
double charge_eu
Definition: col_species.hpp:20
constexpr double TWOPI
Definition: globals.hpp:123
int n_vz() const
Definition: vgrid_distribution.hpp:58
double numeric_vth2
Definition: col_species.hpp:23
int grid_index
Definition: col_species.hpp:15
int n() const
Definition: col_species.hpp:175
double T_avg
Definition: col_species.hpp:17
double get_den_at_node(int node, int species)
double dens
Definition: col_species.hpp:25
constexpr double PI
Definition: globals.hpp:122
double mom
Definition: col_species.hpp:27
const VGridDistribution< HostType > f
Definition: col_species.hpp:85
int n_species() const
Definition: vgrid_distribution.hpp:46
Definition: species.hpp:71
Definition: col_species.hpp:83
void lambda_gamma_pair(double dt, Kokkos::View< double ***, Device > &gammac)
Definition: col_species.hpp:159
double get_t_ev_at_node(int node, int species)
int n_vr() const
Definition: vgrid_distribution.hpp:50
constexpr double UNIT_CHARGE
Charge of an electron (C)
Definition: globals.hpp:118
double ens
Definition: col_species.hpp:26
KOKKOS_INLINE_FUNCTION double lambda_gamma_pair(const CollisionSpeciesScalars &sp_b) const
Definition: col_species.hpp:30