Access to more information in aero_particle (+ execute rand_init befo…

…re each test, not just once per test session) (#286)

Co-authored-by: Sylwester Arabas <[email protected]>

src/aero_data.F90

@@ -185,4 +185,15 @@ subroutine f_aero_data_get_species_density(ptr_c, idx, val) bind(C)
 
   end subroutine
 
+  subroutine f_aero_data_n_source(ptr_c, n_source) bind(C)
+    type(aero_data_t), pointer :: ptr_f => null()
+    type(c_ptr), intent(in) :: ptr_c
+    integer(c_int) :: n_source
+
+    call c_f_pointer(ptr_c, ptr_f)
+
+    n_source = aero_data_n_source(ptr_f)
+
+  end subroutine
+
 end module

src/aero_data.hpp

@@ -14,6 +14,7 @@ extern "C" void f_aero_data_dtor(void *ptr) noexcept;
 extern "C" void f_aero_data_from_json(const void *ptr) noexcept;
 extern "C" void f_aero_data_spec_by_name(const void *ptr, int *value, const char *name_data, const int *name_size) noexcept;
 extern "C" void f_aero_data_len(const void *ptr, int *len) noexcept;
+extern "C" void f_aero_data_n_source(const void *ptr, int *len) noexcept;
 extern "C" void f_aero_data_set_frac_dim(void *ptr, const double*) noexcept;
 extern "C" void f_aero_data_get_frac_dim(const void *ptr, double*) noexcept;
 extern "C" void f_aero_data_set_vol_fill_factor(void *ptr, const double*) noexcept;
@@ -193,5 +194,16 @@ struct AeroData {
         return data;
     }
 
+    static std::size_t n_source(const AeroData &self) {
+        int len;
+        f_aero_data_n_source(
+            self.ptr.f_arg(),
+            &len
+        );
+        if (len == -1)
+            throw std::runtime_error("No sources defined.");
+        return len;
+    }
+
 };
 

src/aero_particle.F90

@@ -387,4 +387,141 @@ subroutine f_aero_particle_set_vols( &
     )
   end subroutine
 
+  subroutine f_aero_particle_absorb_cross_sect( &
+      aero_particle_ptr_c, &
+      absorb_cross_sect & 
+    ) bind(C)
+
+    type(aero_particle_t), pointer :: aero_particle_ptr_f => null()
+    type(c_ptr), intent(in) :: aero_particle_ptr_c
+    real(c_double), intent(out) :: absorb_cross_sect    
+
+    call c_f_pointer(aero_particle_ptr_c, aero_particle_ptr_f)
+
+    absorb_cross_sect = aero_particle_ptr_f%absorb_cross_sect
+
+  end subroutine
+
+  subroutine f_aero_particle_scatter_cross_sect( &
+      aero_particle_ptr_c, &
+      scatter_cross_sect &
+    ) bind(C)
+
+    type(aero_particle_t), pointer :: aero_particle_ptr_f => null()
+    type(c_ptr), intent(in) :: aero_particle_ptr_c
+    real(c_double), intent(out) :: scatter_cross_sect
+
+    call c_f_pointer(aero_particle_ptr_c, aero_particle_ptr_f)
+
+    scatter_cross_sect = aero_particle_ptr_f%scatter_cross_sect
+
+  end subroutine
+
+  subroutine f_aero_particle_asymmetry( &
+      aero_particle_ptr_c, &
+      asymmetry &
+    ) bind(C)
+
+    type(aero_particle_t), pointer :: aero_particle_ptr_f => null()
+    type(c_ptr), intent(in) :: aero_particle_ptr_c
+    real(c_double), intent(out) :: asymmetry
+
+    call c_f_pointer(aero_particle_ptr_c, aero_particle_ptr_f)
+
+    asymmetry = aero_particle_ptr_f%asymmetry
+
+  end subroutine
+
+  subroutine f_aero_particle_greatest_create_time( &
+      aero_particle_ptr_c, &
+      greatest_create_time &
+    ) bind(C)
+
+    type(aero_particle_t), pointer :: aero_particle_ptr_f => null()
+    type(c_ptr), intent(in) :: aero_particle_ptr_c
+    real(c_double), intent(out) :: greatest_create_time
+
+    call c_f_pointer(aero_particle_ptr_c, aero_particle_ptr_f)
+
+    greatest_create_time = aero_particle_ptr_f%greatest_create_time
+
+  end subroutine
+
+  subroutine f_aero_particle_least_create_time( &
+      aero_particle_ptr_c, &
+      least_create_time &
+    ) bind(C)
+
+    type(aero_particle_t), pointer :: aero_particle_ptr_f => null()
+    type(c_ptr), intent(in) :: aero_particle_ptr_c
+    real(c_double), intent(out) :: least_create_time
+
+    call c_f_pointer(aero_particle_ptr_c, aero_particle_ptr_f)
+
+    least_create_time = aero_particle_ptr_f%least_create_time
+
+  end subroutine
+
+  subroutine f_aero_particle_n_orig_part( &
+      aero_particle_ptr_c, &
+      n_orig_part, &
+      n_orig_part_size &
+    ) bind(C)
+
+    type(aero_particle_t), pointer :: aero_particle_ptr_f => null()
+    type(c_ptr), intent(in) :: aero_particle_ptr_c
+    integer(c_int), intent(in) :: n_orig_part_size
+    integer(c_int), dimension(n_orig_part_size), intent(out) :: n_orig_part
+
+    call c_f_pointer(aero_particle_ptr_c, aero_particle_ptr_f)
+
+    n_orig_part = aero_particle_ptr_f%n_orig_part
+
+  end subroutine
+
+  subroutine f_aero_particle_id( & 
+      aero_particle_ptr_c, &
+      id &
+    ) bind(C)
+
+    type(aero_particle_t), pointer :: aero_particle_ptr_f => null()
+    type(c_ptr), intent(in) :: aero_particle_ptr_c
+    integer(c_int), intent(out) :: id
+
+    call c_f_pointer(aero_particle_ptr_c, aero_particle_ptr_f)
+
+    id = aero_particle_ptr_f%id
+
+  end subroutine
+
+  subroutine f_aero_particle_refract_shell( &
+      aero_particle_ptr_c, &
+      refract_shell &
+    ) bind(C)
+
+    type(aero_particle_t), pointer :: aero_particle_ptr_f => null()
+    type(c_ptr), intent(in) :: aero_particle_ptr_c
+    complex(c_double_complex), intent(out) :: refract_shell
+
+    call c_f_pointer(aero_particle_ptr_c, aero_particle_ptr_f)
+
+    refract_shell = aero_particle_ptr_f%refract_shell
+
+  end subroutine
+
+  subroutine f_aero_particle_refract_core( &
+      aero_particle_ptr_c, &
+      refract_core &
+    ) bind(C)
+
+    type(aero_particle_t), pointer :: aero_particle_ptr_f => null()
+    type(c_ptr), intent(in) :: aero_particle_ptr_c
+    complex(c_double_complex), intent(out) :: refract_core
+
+    call c_f_pointer(aero_particle_ptr_c, aero_particle_ptr_f)
+
+    refract_core = aero_particle_ptr_f%refract_core
+
+  end subroutine
+
 end module

src/aero_particle.hpp

@@ -10,6 +10,7 @@
 #include "aero_data.hpp"
 #include "env_state.hpp"
 #include "pybind11/stl.h"
+#include <complex>
 
 extern "C" void f_aero_particle_ctor(void *ptr) noexcept;
 extern "C" void f_aero_particle_dtor(void *ptr) noexcept;
@@ -35,6 +36,15 @@ extern "C" void f_aero_particle_crit_diameter(const void *aero_particle_ptr, con
 extern "C" void f_aero_particle_coagulate(const void *aero_particle_1_ptr, const void *aero_particle_2_ptr, void *new_particle_ptr) noexcept;
 extern "C" void f_aero_particle_zero(void *aero_particle_ptr, const void *aero_data_ptr) noexcept;
 extern "C" void f_aero_particle_set_vols(void *aero_particle_ptr, const int *vol_size, const void *volumes) noexcept;
+extern "C" void f_aero_particle_absorb_cross_sect(const void *aero_particle_ptr, double *val) noexcept;
+extern "C" void f_aero_particle_scatter_cross_sect(const void *aero_particle_ptr, double *val) noexcept;
+extern "C" void f_aero_particle_asymmetry(const void *aero_particle_ptr, double *val) noexcept;
+extern "C" void f_aero_particle_greatest_create_time(const void *aero_particle_ptr, double *val) noexcept;
+extern "C" void f_aero_particle_least_create_time(const void *aero_particle_ptr, double *val) noexcept;
+extern "C" void f_aero_particle_n_orig_part(const void *aero_particle_ptr, void *arr_data, const int *arr_size) noexcept;
+extern "C" void f_aero_particle_id(const void *aero_particle_ptr, int *val) noexcept;
+extern "C" void f_aero_particle_refract_shell(const void *aero_particle_ptr, std::complex<double> *val) noexcept;
+extern "C" void f_aero_particle_refract_core(const void *aero_particle_ptr, std::complex<double> *val) noexcept;
 
 namespace py = pybind11;
 struct AeroParticle {
@@ -272,4 +282,87 @@ struct AeroParticle {
         );
     }
 
-};
+    static auto scatter_cross_sect(const AeroParticle &self) {
+        double val;
+        f_aero_particle_scatter_cross_sect(
+            self.ptr.f_arg(),
+            &val
+        );
+        return val;
+    }
+
+    static auto absorb_cross_sect(const AeroParticle &self) {
+        double val;
+        f_aero_particle_absorb_cross_sect(
+            self.ptr.f_arg(),
+            &val
+        );
+        return val;
+    }
+
+    static auto asymmetry(const AeroParticle &self) {
+        double val;
+        f_aero_particle_asymmetry(
+            self.ptr.f_arg(),
+            &val
+        );
+        return val;
+    }
+
+    static auto n_orig_part(const AeroParticle &self) {
+        int len = AeroData::n_source(*self.aero_data);
+        std::valarray<int> data(len);
+
+        f_aero_particle_n_orig_part(
+            self.ptr.f_arg(),
+            begin(data),
+            &len
+        );
+        return data;
+    }
+
+    static auto least_create_time(const AeroParticle &self) {
+        double val;
+        f_aero_particle_least_create_time(
+            self.ptr.f_arg(),
+            &val
+        );
+        return val;
+    }
+
+    static auto greatest_create_time(const AeroParticle &self) {
+        double val;
+        f_aero_particle_greatest_create_time(
+            self.ptr.f_arg(),
+            &val
+        );
+        return val;
+    }
+
+    static auto id(const AeroParticle &self) {
+        int val;
+        f_aero_particle_id(
+            self.ptr.f_arg(),
+            &val
+        );
+        return val;
+    }
+
+    static auto refract_shell(const AeroParticle &self) {
+        std::complex<double> refract_shell;
+        f_aero_particle_refract_shell(
+            self.ptr.f_arg(),
+            &refract_shell
+        );
+        return refract_shell;
+    }
+
+    static auto refract_core(const AeroParticle &self) {
+        std::complex<double> refract_core;
+        f_aero_particle_refract_core(
+            self.ptr.f_arg(),
+            &refract_core
+        );
+        return refract_core;
+    }
+};

src/pypartmc.cpp

@@ -7,6 +7,7 @@
 #include "pybind11/pybind11.h"
 #include "nlohmann/json.hpp"
 #include "pybind11_json/pybind11_json.hpp"
+#include "pybind11/complex.h"
 
 #include "util.hpp"
 #include "rand.hpp"
@@ -76,6 +77,7 @@ PYBIND11_MODULE(_PyPartMC, m) {
         .def(py::init<const nlohmann::json&>())
         .def("spec_by_name", AeroData::spec_by_name)
         .def("__len__", AeroData::__len__)
+        .def_property_readonly("n_source", AeroData::n_source)
         .def_property("frac_dim", &AeroData::get_frac_dim, &AeroData::set_frac_dim)
         .def_property("vol_fill_factor", &AeroData::get_vol_fill_factor, &AeroData::set_vol_fill_factor)
         .def_property("prime_radius", &AeroData::get_prime_radius, &AeroData::set_prime_radius)
@@ -128,6 +130,23 @@ PYBIND11_MODULE(_PyPartMC, m) {
             "Returns the average of the solute kappas (1).")
         .def_property_readonly("moles", AeroParticle::moles,
             "Total moles in the particle (1).")
+        .def_property_readonly("absorb_cross_sect", AeroParticle::absorb_cross_sect,
+            "Absorption cross-section (m^-2).")
+        .def_property_readonly("scatter_cross_sect", AeroParticle::scatter_cross_sect,
+            "Scattering cross-section (m^-2).")
+        .def_property_readonly("asymmetry", AeroParticle::asymmetry,
+            "Asymmetry parameter (1).")
+        .def_property_readonly("refract_shell", AeroParticle::refract_shell,
+            "Refractive index of the shell (1).")
+        .def_property_readonly("refract_core", AeroParticle::refract_core,
+            "Refractive index of the core (1).")
+        .def_property_readonly("n_orig_part", AeroParticle::n_orig_part,
+            "Number of original particles from each source that coagulated to form particle.")
+        .def_property_readonly("least_create_time", AeroParticle::least_create_time,
+            "First time a constituent was created (s).")
+        .def_property_readonly("greatest_create_time", AeroParticle::greatest_create_time,
+            "Last time a constituent was created (s).")
+        .def_property_readonly("id", AeroParticle::id, "Unique ID number.")
         .def("mobility_diameter", AeroParticle::mobility_diameter,
             "Mobility diameter of the particle (m).")
         .def_property_readonly("density", AeroParticle::density,

tests/conftest.py

@@ -1,3 +1,8 @@
+import pytest
+
 import PyPartMC as ppmc
 
-ppmc.rand_init(44)
+
+@pytest.fixture(autouse=True)
+def rand_init():
+    ppmc.rand_init(44)

tests/test_aero_data.py

@@ -362,3 +362,15 @@ def test_ctor_error_on_nonunique_keys():
 
         # assert
         assert str(exc_info.value) == "Species names must be unique"
+
+    @staticmethod
+    def test_n_source_uninitialized():
+        # arrange
+        sut = ppmc.AeroData(AERO_DATA_CTOR_ARG_FULL)
+
+        # act
+        with pytest.raises(Exception) as exc_info:
+            _ = sut.n_source
+
+        # assert
+        assert str(exc_info.value) == "No sources defined."