Skip to content

Commit 97fecfa

Browse files
pdziekanpdziekan
authored
Merge pull request #393 from pdziekan/diag_cloud_time
[WIP] Diag cloud time
2 parents 80d4a13 + 726d0ad commit 97fecfa

21 files changed

+431
-58
lines changed

‎Readme.md‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -112,3 +112,4 @@ Some CMake hints:
112112
- the output of commands executed by "make test" can be viewed with:
113113
$ less Testing/Temporary/LastTest.log
114114

115+

‎bindings/python/lib.cpp‎

Lines changed: 9 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -23,13 +23,13 @@
2323
#define BP_ARR_FROM_BP_OBJ bp_array(bp::numpy::array(bp::object()))
2424
#endif
2525

26+
2627
BOOST_PYTHON_MODULE(libcloudphxx)
2728
{
2829
namespace bp = boost::python;
29-
using namespace libcloudphxx::python;
30-
3130
using real_t = double;
3231
using arr_t = blitz::Array<real_t, 1>;
32+
using namespace libcloudphxx::python;
3333

3434
#ifdef BPNUMERIC
3535
bp_array::set_module_and_type("numpy", "ndarray");
@@ -296,6 +296,7 @@ BOOST_PYTHON_MODULE(libcloudphxx)
296296
.def_readwrite("chem_rho", &lgr::opts_init_t<real_t>::chem_rho)
297297
.def_readwrite("RH_max", &lgr::opts_init_t<real_t>::RH_max)
298298
.def_readwrite("rng_seed", &lgr::opts_init_t<real_t>::rng_seed)
299+
.def_readwrite("diag_incloud_time", &lgr::opts_init_t<real_t>::diag_incloud_time)
299300
.add_property("w_LS", &lgrngn::get_w_LS<real_t>, &lgrngn::set_w_LS<real_t>)
300301
.add_property("SGS_mix_len", &lgrngn::get_SGS_mix_len<real_t>, &lgrngn::set_SGS_mix_len<real_t>)
301302
.add_property("kernel_parameters", &lgrngn::get_kp<real_t>, &lgrngn::set_kp<real_t>)
@@ -345,12 +346,16 @@ BOOST_PYTHON_MODULE(libcloudphxx)
345346
.def("diag_pressure",&lgr::particles_proto_t<real_t>::diag_pressure)
346347
.def("diag_temperature",&lgr::particles_proto_t<real_t>::diag_temperature)
347348
.def("diag_vel_div",&lgr::particles_proto_t<real_t>::diag_vel_div)
348-
.def("diag_dry_rng", &lgr::particles_proto_t<real_t>::diag_dry_rng)
349-
.def("diag_wet_rng", &lgr::particles_proto_t<real_t>::diag_wet_rng)
349+
.def("diag_dry_rng", &lgr::particles_proto_t<real_t>::diag_dry_rng)
350+
.def("diag_wet_rng", &lgr::particles_proto_t<real_t>::diag_wet_rng)
350351
.def("diag_kappa_rng", &lgr::particles_proto_t<real_t>::diag_kappa_rng)
352+
.def("diag_dry_rng_cons", &lgr::particles_proto_t<real_t>::diag_dry_rng_cons)
353+
.def("diag_wet_rng_cons", &lgr::particles_proto_t<real_t>::diag_wet_rng_cons)
354+
.def("diag_kappa_rng_cons", &lgr::particles_proto_t<real_t>::diag_kappa_rng_cons)
351355
.def("diag_dry_mom", &lgr::particles_proto_t<real_t>::diag_dry_mom)
352356
.def("diag_wet_mom", &lgr::particles_proto_t<real_t>::diag_wet_mom)
353357
.def("diag_kappa_mom", &lgr::particles_proto_t<real_t>::diag_kappa_mom)
358+
.def("diag_incloud_time_mom", &lgr::particles_proto_t<real_t>::diag_incloud_time_mom)
354359
.def("diag_wet_mass_dens", &lgr::particles_proto_t<real_t>::diag_wet_mass_dens)
355360
.def("diag_chem", &lgr::particles_proto_t<real_t>::diag_chem)
356361
.def("diag_precip_rate", &lgr::particles_proto_t<real_t>::diag_precip_rate)

‎include/libcloudph++/lgrngn/opts_init.hpp‎

Lines changed: 5 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -111,6 +111,9 @@ namespace libcloudphxx
111111
int sstp_chem;
112112
real_t chem_rho;
113113

114+
// do we want to track the time SDs spend inside clouds
115+
bool diag_incloud_time;
116+
114117
// RH threshold for calculating equilibrium condition at t=0
115118
real_t RH_max;
116119

@@ -165,7 +168,8 @@ namespace libcloudphxx
165168
n_sd_max(0),
166169
src_sd_conc(0),
167170
src_z1(0),
168-
rd_min(0.)
171+
rd_min(0.),
172+
diag_incloud_time(false)
169173
{}
170174

171175
// dtor (just to silence -Winline warnings)

‎include/libcloudph++/lgrngn/particles.hpp‎

Lines changed: 47 additions & 36 deletions
Original file line numberDiff line numberDiff line change
@@ -75,26 +75,38 @@ namespace libcloudphxx
7575
}
7676

7777
// method for accessing super-droplet statistics
78-
virtual void diag_sd_conc() { assert(false); }
79-
virtual void diag_pressure() { assert(false); }
80-
virtual void diag_temperature() { assert(false); }
81-
virtual void diag_RH() { assert(false); }
82-
virtual void diag_all() { assert(false); }
83-
virtual void diag_rw_ge_rc() { assert(false); }
84-
virtual void diag_RH_ge_Sc() { assert(false); }
85-
virtual void diag_dry_rng(const real_t&, const real_t&) { assert(false); }
86-
virtual void diag_wet_rng(const real_t&, const real_t&) { assert(false); }
87-
virtual void diag_dry_mom(const int&) { assert(false); }
88-
virtual void diag_wet_mom(const int&) { assert(false); }
89-
virtual void diag_wet_mass_dens(const real_t&, const real_t&) { assert(false); }
90-
virtual void diag_chem(const enum common::chem::chem_species_t&) { assert(false); }
91-
virtual void diag_precip_rate() { assert(false); }
92-
virtual void diag_kappa_mom(const int&) { assert(false); }
93-
virtual void diag_kappa_rng(const real_t&, const real_t&) { assert(false); }
94-
virtual void diag_max_rw() { assert(false); }
95-
virtual void diag_vel_div() { assert(false); }
96-
virtual std::map<libcloudphxx::common::output_t, real_t> diag_puddle() { assert(false); return std::map<libcloudphxx::common::output_t, real_t>(); }
97-
virtual real_t *outbuf() { assert(false); return NULL; }
78+
virtual void diag_sd_conc() { assert(false); }
79+
virtual void diag_pressure() { assert(false); }
80+
virtual void diag_temperature() { assert(false); }
81+
virtual void diag_RH() { assert(false); }
82+
virtual void diag_all() { assert(false); }
83+
virtual void diag_rw_ge_rc() { assert(false); }
84+
virtual void diag_RH_ge_Sc() { assert(false); }
85+
virtual void diag_dry_rng(const real_t&, const real_t&) { assert(false); }
86+
virtual void diag_wet_rng(const real_t&, const real_t&) { assert(false); }
87+
virtual void diag_kappa_rng(const real_t&, const real_t&) { assert(false); }
88+
// The 3 following functions are for consecutive selection of SDs.
89+
// It allows the user to select SDs based on multiple characteristics, e.g. wet radius (0.5, 1) and kappa (0.1, 0.2):
90+
// diag_wet_rng(0.5, 1); diag_kappa_rng_cons(0.1, 0.2);
91+
// NOTE: the call with "cons" needs to be right after the previous call to diag_X_rng!
92+
// otherwise some other function used in between can overwrite n_filtered used for selection of moments
93+
// NOTE2: We cannot implement this as an argument to diag_X_rng, because we would like it to default to null
94+
// and Boost Python does not work well with virtual member functions that have default arguments
95+
virtual void diag_dry_rng_cons(const real_t&, const real_t&) { assert(false); }
96+
virtual void diag_wet_rng_cons(const real_t&, const real_t&) { assert(false); }
97+
virtual void diag_kappa_rng_cons(const real_t&, const real_t&) { assert(false); }
98+
99+
virtual void diag_dry_mom(const int&) { assert(false); }
100+
virtual void diag_wet_mom(const int&) { assert(false); }
101+
virtual void diag_wet_mass_dens(const real_t&, const real_t&) { assert(false); }
102+
virtual void diag_chem(const enum common::chem::chem_species_t&) { assert(false); }
103+
virtual void diag_precip_rate() { assert(false); }
104+
virtual void diag_kappa_mom(const int&) { assert(false); }
105+
virtual void diag_incloud_time_mom(const int&) { assert(false); } // requires opts_init.diag_incloud_time==true
106+
virtual void diag_max_rw() { assert(false); }
107+
virtual void diag_vel_div() { assert(false); }
108+
virtual std::map<libcloudphxx::common::output_t, real_t> diag_puddle() { assert(false); return std::map<libcloudphxx::common::output_t, real_t>(); }
109+
virtual real_t *outbuf() { assert(false); return NULL; }
98110

99111
// storing a pointer to opts_init (e.g. for interrogatin about
100112
// dimensions in Python bindings)
@@ -160,26 +172,23 @@ namespace libcloudphxx
160172
void diag_pressure();
161173
void diag_temperature();
162174
void diag_RH();
163-
void diag_dry_rng(
164-
const real_t &r_mi, const real_t &r_mx
165-
);
166-
void diag_wet_rng(
167-
const real_t &r_mi, const real_t &r_mx
168-
);
169-
void diag_kappa_rng(
170-
const real_t &r_mi, const real_t &r_mx
171-
);
175+
void diag_dry_rng(const real_t &r_mi, const real_t &r_mx);
176+
void diag_wet_rng(const real_t &r_mi, const real_t &r_mx);
177+
void diag_kappa_rng(const real_t &r_mi, const real_t &r_mx);
178+
void diag_dry_rng_cons(const real_t &r_mi, const real_t &r_mx);
179+
void diag_wet_rng_cons(const real_t &r_mi, const real_t &r_mx);
180+
void diag_kappa_rng_cons(const real_t &r_mi, const real_t &r_mx);
172181
void diag_dry_mom(const int &k);
173182
void diag_wet_mom(const int &k);
174183
void diag_kappa_mom(const int &k);
184+
void diag_incloud_time_mom(const int &k);
175185
void diag_wet_mass_dens(const real_t&, const real_t&);
176186

177187
void diag_chem(const enum common::chem::chem_species_t&);
178188
void diag_rw_ge_rc();
179189
void diag_RH_ge_Sc();
180190
void diag_all();
181191
void diag_precip_rate();
182-
void diag_kappa(const int&);
183192
void diag_max_rw();
184193
void diag_vel_div();
185194
std::map<libcloudphxx::common::output_t, real_t> diag_puddle();
@@ -256,14 +265,16 @@ namespace libcloudphxx
256265
void diag_pressure();
257266
void diag_temperature();
258267
void diag_RH();
259-
void diag_dry_rng(
260-
const real_t &r_mi, const real_t &r_mx
261-
);
262-
void diag_wet_rng(
263-
const real_t &r_mi, const real_t &r_mx
264-
);
268+
void diag_dry_rng(const real_t &r_mi, const real_t &r_mx);
269+
void diag_wet_rng(const real_t &r_mi, const real_t &r_mx);
270+
void diag_kappa_rng(const real_t &r_mi, const real_t &r_mx);
271+
void diag_dry_rng_cons(const real_t &r_mi, const real_t &r_mx);
272+
void diag_wet_rng_cons(const real_t &r_mi, const real_t &r_mx);
273+
void diag_kappa_rng_cons(const real_t &r_mi, const real_t &r_mx);
265274
void diag_dry_mom(const int &k);
266275
void diag_wet_mom(const int &k);
276+
void diag_kappa_mom(const int&);
277+
void diag_incloud_time_mom(const int&);
267278
void diag_wet_mass_dens(const real_t&, const real_t&);
268279
real_t *outbuf();
269280

‎src/impl/particles_impl.ipp‎

Lines changed: 7 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -85,10 +85,11 @@ namespace libcloudphxx
8585
sstp_tmp_rv, // either rv_old or advection-caused change in water vapour mixing ratio
8686
sstp_tmp_th, // ditto for theta
8787
sstp_tmp_rh, // ditto for rho
88-
sstp_tmp_p; // ditto for pressure
88+
sstp_tmp_p, // ditto for pressure
89+
incloud_time; // time this SD has been within a cloud
8990

9091
const int no_of_n_vctrs_copied = 1;
91-
const int no_of_real_vctrs_copied = 15;
92+
const int no_of_real_vctrs_copied = opts_init.diag_incloud_time ? 16 : 15;
9293

9394
// dry radii distribution characteristics
9495
real_t log_rd_min, // logarithm of the lower bound of the distr
@@ -382,6 +383,7 @@ namespace libcloudphxx
382383
void init_ijk();
383384
void init_xyz();
384385
void init_kappa(const real_t &);
386+
void init_incloud_time();
385387
void init_count_num_sd_conc(const real_t & = 1);
386388
void init_count_num_const_multi(const common::unary_function<real_t> &);
387389
void init_count_num_const_multi(const common::unary_function<real_t> &, const thrust_size_t &);
@@ -432,7 +434,8 @@ namespace libcloudphxx
432434
);
433435
void moms_rng(
434436
const real_t &min, const real_t &max,
435-
const typename thrust_device::vector<real_t>::iterator &vec_bgn
437+
const typename thrust_device::vector<real_t>::iterator &vec_bgn,
438+
const bool cons
436439
);
437440
void moms_calc(
438441
const typename thrust_device::vector<real_t>::iterator &vec_bgn,
@@ -469,6 +472,7 @@ namespace libcloudphxx
469472
void update_th_rv(thrust_device::vector<real_t> &);
470473
void update_state(thrust_device::vector<real_t> &, thrust_device::vector<real_t> &);
471474
void update_pstate(thrust_device::vector<real_t> &, thrust_device::vector<real_t> &);
475+
void update_incloud_time();
472476

473477
void coal(const real_t &dt, const bool &turb_coal);
474478

‎src/impl/particles_impl_coal.ipp‎

Lines changed: 37 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -14,6 +14,22 @@ namespace libcloudphxx
1414
{
1515
enum{na_ge_nb = -2, nb_gt_na = -1};
1616

17+
struct selector // keep the max value of some parameter in the SD that represents droplets that collided
18+
{
19+
template<class tpl_t>
20+
BOOST_GPU_ENABLED
21+
void operator()(tpl_t tpl)
22+
{
23+
#if !defined(__NVCC__)
24+
using std::max;
25+
#endif
26+
if(thrust::get<2>(tpl) <= 0) return; // do nothing if no collisions or first one passed was a SD with an uneven number in the cell
27+
thrust::get<3>(tpl) == na_ge_nb ? // does the first SD of the pair have greater multiplicity?
28+
thrust::get<1>(tpl) = max(thrust::get<0>(tpl), thrust::get<1>(tpl)): // set val = max(val_SD_A, val_SD_B) in the one with smaller multiplicity
29+
thrust::get<0>(tpl) = max(thrust::get<0>(tpl), thrust::get<1>(tpl)); // set val = max(val_SD_A, val_SD_B) in the one with smaller multiplicity
30+
}
31+
};
32+
1733
struct summator
1834
{
1935
template<class tpl_t>
@@ -466,6 +482,27 @@ namespace libcloudphxx
466482
);
467483
nancheck(kpa, "kpa - post coalescence");
468484
}
485+
486+
// update incloud time
487+
if(opts_init.diag_incloud_time)
488+
{
489+
thrust::for_each(
490+
thrust::make_zip_iterator(thrust::make_tuple(
491+
thrust::make_permutation_iterator(incloud_time.begin(), sorted_id.begin()),
492+
thrust::make_permutation_iterator(incloud_time.begin(), sorted_id.begin())+1,
493+
col.begin(),
494+
col.begin()+1
495+
)),
496+
thrust::make_zip_iterator(thrust::make_tuple(
497+
thrust::make_permutation_iterator(incloud_time.begin(), sorted_id.begin()),
498+
thrust::make_permutation_iterator(incloud_time.begin(), sorted_id.begin())+1,
499+
col.begin(),
500+
col.begin()+1
501+
)) + n_part -1,
502+
detail::selector()
503+
);
504+
nancheck(incloud_time, "incloud_time - post coalescence");
505+
}
469506
}
470507
};
471508
};

‎src/impl/particles_impl_hskpng_remove.ipp‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -56,6 +56,9 @@ namespace libcloudphxx
5656
real_t_vctrs.insert(&sstp_tmp_p);
5757
}
5858

59+
if(opts_init.diag_incloud_time)
60+
real_t_vctrs.insert(&incloud_time);
61+
5962
namespace arg = thrust::placeholders;
6063

6164
// remove chemical stuff

‎src/impl/particles_impl_hskpng_resize.ipp‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -69,6 +69,9 @@ namespace libcloudphxx
6969
sstp_tmp_p.resize(n_part);
7070
}
7171
}
72+
73+
if(opts_init.diag_incloud_time)
74+
incloud_time.resize(n_part);
7275
}
7376
};
7477
};

‎src/impl/particles_impl_init_hskpng_npart.ipp‎

Lines changed: 3 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -52,6 +52,9 @@ namespace libcloudphxx
5252
n.reserve(opts_init.n_sd_max);
5353
kpa.reserve(opts_init.n_sd_max);
5454

55+
if(opts_init.diag_incloud_time)
56+
incloud_time.reserve(opts_init.n_sd_max);
57+
5558
if(opts_init.chem_switch || opts_init.sstp_cond > 1 || n_dims >= 2)
5659
{
5760
tmp_device_real_part1.reserve(opts_init.n_sd_max);

‎src/impl/particles_impl_init_incloud_time.ipp‎

Lines changed: 20 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,20 @@
1+
// vim:filetype=cpp
2+
/** @file
3+
* @copyright University of Warsaw
4+
* @section LICENSE
5+
* GPLv3+ (see the COPYING file or http://www.gnu.org/licenses/)
6+
* @brief initialisation routine for super droplets
7+
*/
8+
9+
namespace libcloudphxx
10+
{
11+
namespace lgrngn
12+
{
13+
template <typename real_t, backend_t device>
14+
void particles_t<real_t, device>::impl::init_incloud_time()
15+
{
16+
thrust::fill(incloud_time.begin() + n_part_old, incloud_time.end(), 0);
17+
}
18+
};
19+
};
20+

0 commit comments

Comments
 (0)