time-varying aero specs in Scenario (#207)

Co-authored-by: Jeff Curtis <[email protected]>

Author: slayoo

src/aero_dist.hpp

@@ -54,6 +54,10 @@ struct AeroDist {
         gimmick_ptr().reset();
     }
 
+    AeroDist() :
+        ptr(f_aero_dist_ctor, f_aero_dist_dtor)
+    {}
+
     static auto get_n_mode(const AeroDist &self) {
         int n_mode;
         f_aero_dist_n_mode(self.ptr.f_arg(), &n_mode);

src/fake_spec_file.F90

@@ -82,6 +82,11 @@ subroutine c_spec_file_read_real_named_array_data( &
             integer, intent(in) :: row, vals_size, names_size
         end subroutine
 
+        subroutine c_spec_file_read_line_data_size(sz) bind(C)
+            import c_int
+            integer(kind=c_int), intent(out) :: sz
+        end subroutine
+
         subroutine c_spec_file_read_line( &
             name_data, name_size, &
             data0_data, data0_size, &
@@ -100,7 +105,11 @@ subroutine spec_file_check_line_name(file, line, name)
         type(spec_file_t), intent(in) :: file
         type(spec_line_t), intent(in) :: line
         character(len=*), intent(in) :: name
-        if (line%name /= name) then
+        integer :: substr_len
+
+        substr_len = len(name)
+
+        if (line%name(1:substr_len) /= name(1:substr_len)) then
             print*, 'line must begin with: ' // trim(name) // ' not: ' // trim(line%name)
             call pmc_stop(462932478)
         end if
@@ -115,10 +124,14 @@ subroutine spec_file_read_real_named_array(file, max_lines, names, vals)
         integer :: row, col, n_rows, n_cols, name_size;
 
         call c_spec_file_read_real_named_array_size(n_rows, n_cols)
+
+        if (max_lines .ne. 0 .and. max_lines .lt. n_rows) then
+            n_rows = max_lines
+        end if
+
         allocate(names(n_rows))
         allocate(vals(n_rows, n_cols))
         allocate(vals_row(n_cols))
-        ! TODO #112: handle max_lines
         do row = 1, n_rows
             name_size = len(names(row))
             call c_spec_file_read_real_named_array_data( &
@@ -146,7 +159,9 @@ subroutine spec_file_check_line_length(file, line, length)
         type(spec_line_t), intent(in) :: line
         integer, intent(in) :: length
         if (size(line%data) /= length) then
-            print*, 'expected ' // trim(integer_to_string(length)) // ' data items on line'
+            print*, 'expected ' // trim(integer_to_string(length)) &
+                                // ' data items on line, not ' &
+                                // trim(integer_to_string(size(line%data)))
             call pmc_stop(189339129)
         end if
     end subroutine
@@ -197,17 +212,21 @@ subroutine spec_file_read_line(file, line, eof)
         type(spec_line_t), intent(inout) :: line
         logical, intent(out) :: eof
 
-        integer :: name_size, data0_size
+        integer :: name_size, data0_size, data_size, i
         logical(c_bool) :: c_eof
 
-        allocate(line%data(1))
+        call c_spec_file_read_line_data_size(data_size)
+
+        allocate(line%data(data_size))
         name_size = len(line%name)
         data0_size = len(line%data(1))
         call c_spec_file_read_line(line%name, name_size, line%data(1), data0_size, c_eof)
         eof = c_eof
         if (.not. eof) then
             line%name = line%name(1:name_size)
-            line%data(1) = line%data(1)(1:data0_size)
+            do i = 1, data_size
+                line%data(i) = line%data(1)(1:data0_size)
+            end do
         else
             deallocate(line%data)
         end if

src/fake_spec_file.cpp

@@ -205,6 +205,12 @@ void c_spec_file_read_real_named_array_data(
     );
 }
 
+/*********************************************************************************/
+extern "C"
+void c_spec_file_read_line_data_size(int *size) noexcept {
+    *size = gimmick_ptr()->read_line_data_size_and_start_enumerating();
+}
+
 /*********************************************************************************/
 extern "C"
 void c_spec_file_read_line(

src/gimmicks.hpp

@@ -26,6 +26,8 @@ struct Gimmick {
     }
 
   protected:
+    size_t index = 0;
+
     Gimmick() {}
 
     Gimmick(const nlohmann::json &json) {
@@ -73,9 +75,15 @@ struct Gimmick {
         auto it = this->json->is_array() 
             ? this->json->at(this->json->size()-1).find(sub)
             : this->json->find(sub);
-        // TODO #112: handle errors
+
+        assert(*it != NULL);
         this->json_parent.push(this->json);
-        this->set_current_json_ptr(&(*it));
+
+        if (sub == "dist")  // TODO #112: do better than hardcoding "dist"...
+            this->set_current_json_ptr(&(it->at(this->n_elements(sub) - this->index--)));
+        else
+            this->set_current_json_ptr(&(*it));
+
     }
 
     void zoom_out() noexcept {
@@ -84,7 +92,7 @@ struct Gimmick {
         this->json_parent.pop();
     }
 
-    auto zoom_level() noexcept {
+    auto zoom_level() const noexcept {
         auto level = this->json_parent.size();
         return level;
     }
@@ -108,7 +116,7 @@ struct Gimmick {
         return name;
     }
 
-    auto n_elements(const bpstd::string_view &name) noexcept {
+    std::size_t n_elements(const bpstd::string_view &name) const noexcept {
         std::size_t n_elem = 0;
         for (auto i=0u; i<this->json->size(); ++i) {
             for (auto &entry : this->json->at(i).items()) {
@@ -203,6 +211,8 @@ struct Gimmick {
     virtual std::string str() const = 0;
 
     virtual bool read_line(std::string &name, std::string &data) = 0;
+    
+    virtual int read_line_data_size_and_start_enumerating() = 0;
 };
 
 struct InputGimmick: Gimmick {
@@ -218,7 +228,8 @@ struct InputGimmick: Gimmick {
         const std::string key_name = "",
         const std::size_t max_zoom_level = 3
     ) : Gimmick(json), key_cond(key_cond), key_name(key_name), max_zoom_level(max_zoom_level)
-    {}
+    {
+    }
 
     std::string str() const {
         throw std::logic_error("str() called on InputGimmick!");
@@ -231,7 +242,10 @@ struct InputGimmick: Gimmick {
 
             auto key = this->next_dict_key(this->last_read_line_key, this->key_cond);
             if (key == "") {
-                this->last_read_line_key = "";
+                if (this->index == 0)
+                    this->last_read_line_key = "";
+                else
+                    this->last_read_line_key = this->key_cond;
                 return true;
             }
             else
@@ -250,6 +264,14 @@ struct InputGimmick: Gimmick {
 
         return false;
     }
+
+    int read_line_data_size_and_start_enumerating() {
+        if (this->zoom_level() == this->max_zoom_level - 2) {
+            this->index = this->n_elements(this->key_cond);
+            return this->index;
+        }
+        return 1;
+    }
 };
 
 struct OutputGimmick: Gimmick {
@@ -266,6 +288,10 @@ struct OutputGimmick: Gimmick {
     bool read_line(std::string &, std::string &) {
         throw std::logic_error("read_line() called on OutputGimmick!");
     }
+
+    int read_line_data_size_and_start_enumerating() {
+        throw std::logic_error("read_line_data_size_and_start_enumerating() called on OutputGimmick!");
+    }
 };
 
 std::unique_ptr<Gimmick> &gimmick_ptr();

src/pypartmc.cpp

@@ -281,6 +281,19 @@ PYBIND11_MODULE(_PyPartMC, m) {
             "returns a string with JSON representation of the object")
         .def("init_env_state", Scenario::init_env_state,
             "initializes the EnvState")
+        .def("aero_emissions", Scenario::get_dist,
+            "returns aero_emissions AeroDists at a given index")
+        .def_property_readonly("aero_emissions_n_times", Scenario::get_emissions_n_times,
+            "returns the number of times specified for emissions")
+        .def_property_readonly("aero_emissions_rate_scale", Scenario::emission_rate_scale)
+        .def_property_readonly("aero_emissions_time", Scenario::emission_time)
+        .def("aero_background", Scenario::get_aero_background_dist,
+            "returns aero_background AeroDists at a given index")
+        .def_property_readonly("aero_dilution_n_times", Scenario::get_aero_dilution_n_times,
+            "returns the number of times specified for dilution")
+        .def_property_readonly("aero_dilution_rate", Scenario::aero_dilution_rate)
+        .def_property_readonly("aero_dilution_time", Scenario::aero_dilution_time)
+
     ;
 
     py::class_<GasState>(m,

src/scenario.F90

@@ -109,4 +109,112 @@ subroutine f_scenario_init_env_state(scenario_ptr_c, env_state_ptr_c, &
 
   end subroutine
 
+  subroutine f_scenario_aero_dist_emission(scenario_ptr_c, aero_dist_ptr_c, index) bind(C)
+
+    type(c_ptr), intent(in) :: scenario_ptr_c, aero_dist_ptr_c
+    integer(c_int), intent(in) :: index
+    type(scenario_t), pointer :: scenario_ptr_f => null()
+    type(aero_dist_t), pointer :: aero_dist_ptr_f => null()
+
+    call c_f_pointer(scenario_ptr_c, scenario_ptr_f)
+    call c_f_pointer(aero_dist_ptr_c, aero_dist_ptr_f)
+
+    aero_dist_ptr_f = scenario_ptr_f%aero_emission(index + 1)
+
+  end subroutine
+
+  subroutine f_scenario_aero_emission_n_times(scenario_ptr_c, n_times) bind(C)
+
+    type(c_ptr), intent(in) :: scenario_ptr_c
+    integer(c_int), intent(out) :: n_times
+    type(scenario_t), pointer :: scenario_ptr_f => null()
+
+    call c_f_pointer(scenario_ptr_c, scenario_ptr_f)
+
+    n_times = size(scenario_ptr_f%aero_emission_time)
+
+  end subroutine
+
+  subroutine f_scenario_emission_rates(scenario_ptr_c, emission_rates, &
+    n_times) bind(C)
+
+    type(c_ptr), intent(in) :: scenario_ptr_c
+    type(scenario_t), pointer :: scenario_ptr_f => null()
+    integer(c_int) :: n_times
+    real(c_double) :: emission_rates(n_times)
+
+    call c_f_pointer(scenario_ptr_c, scenario_ptr_f)
+
+    emission_rates = scenario_ptr_f%aero_emission_rate_scale
+
+  end subroutine
+
+  subroutine f_scenario_emission_time(scenario_ptr_c, emission_time, &
+    n_times) bind(C)
+
+    type(c_ptr), intent(in) :: scenario_ptr_c
+    type(scenario_t), pointer :: scenario_ptr_f => null()
+    integer(c_int) :: n_times
+    real(c_double) :: emission_time(n_times)
+
+    call c_f_pointer(scenario_ptr_c, scenario_ptr_f)
+
+    emission_time = scenario_ptr_f%aero_emission_time
+
+  end subroutine
+
+  subroutine f_scenario_aero_dist_background(scenario_ptr_c, aero_dist_ptr_c, index) bind(C)
+
+    type(c_ptr), intent(in) :: scenario_ptr_c, aero_dist_ptr_c
+    integer(c_int), intent(in) :: index
+    type(scenario_t), pointer :: scenario_ptr_f => null()
+    type(aero_dist_t), pointer :: aero_dist_ptr_f => null()
+
+    call c_f_pointer(scenario_ptr_c, scenario_ptr_f)
+    call c_f_pointer(aero_dist_ptr_c, aero_dist_ptr_f)
+
+    aero_dist_ptr_f = scenario_ptr_f%aero_background(index + 1)
+
+  end subroutine
+
+  subroutine f_scenario_aero_background_n_times(scenario_ptr_c, n_times) bind(C)
+
+    type(c_ptr), intent(in) :: scenario_ptr_c
+    integer(c_int), intent(out) :: n_times
+    type(scenario_t), pointer :: scenario_ptr_f => null()
+
+    call c_f_pointer(scenario_ptr_c, scenario_ptr_f)
+
+    n_times = size(scenario_ptr_f%aero_dilution_time)
+
+  end subroutine
+
+  subroutine f_scenario_aero_background_rate_scale(scenario_ptr_c, &
+    aero_background_rate_scale, n_times) bind(C)
+
+    type(c_ptr), intent(in) :: scenario_ptr_c
+    type(scenario_t), pointer :: scenario_ptr_f => null()
+    integer(c_int) :: n_times
+    real(c_double) :: aero_background_rate_scale(n_times)
+
+    call c_f_pointer(scenario_ptr_c, scenario_ptr_f)
+
+    aero_background_rate_scale = scenario_ptr_f%aero_dilution_rate
+
+  end subroutine
+
+  subroutine f_scenario_aero_background_time(scenario_ptr_c, aero_background_time, &
+    n_times) bind(C)
+
+    type(c_ptr), intent(in) :: scenario_ptr_c
+    type(scenario_t), pointer :: scenario_ptr_f => null()
+    integer(c_int) :: n_times
+    real(c_double) :: aero_background_time(n_times)
+
+    call c_f_pointer(scenario_ptr_c, scenario_ptr_f)
+
+    aero_background_time = scenario_ptr_f%aero_dilution_time
+
+  end subroutine
+
 end module