Merge pull request #605 from HSLdevcom/feat/agent-based-park-and-ride-aggregated

Park and ride for agent-based model

Author: hsl-petrhaj

Scripts/datatypes/purpose.py

@@ -7,6 +7,8 @@
 import openmatrix as omx
 from datahandling.zonedata import ZoneData
 
+if TYPE_CHECKING:
+    from datatypes.person import Person
 from models.park_and_ride_model import ParkAndRideModel, ParkAndRidePseudoPurpose
 import parameters.zone as param
 from parameters.destination_choice import secondary_destination_threshold, destination_choice
@@ -195,6 +197,38 @@ def calc_basic_prob(self, impedance):
         self.model.calc_basic_prob(impedance)
         self.dist = impedance["car"]["dist"]
 
+    def calc_pnr_demand(self, population, estimation_mode=False, log=None):
+        """Calculate pnr demand matrices for agent based model.
+              
+        Returns
+        -------
+        dict
+            Mode (car/transit/bike) : dict
+                Demand matrix for whole day : Demand
+        """
+        #Initiate OD matrix
+        gen_zones = numpy.zeros_like(self.zone_numbers)
+        attr_zones = numpy.zeros_like(self.zone_data.zone_numbers)
+        od_matrix = numpy.zeros((gen_zones.size,attr_zones.size))
+        #Aggregate tour to matrix
+        for person in population:
+            for tour in person.tours:
+                if tour.purpose_name == self.name and tour.mode == "park_and_ride":
+                    od_matrix[self.zone_data.zone_index(tour.orig),
+                              self.zone_data.zone_index(tour.dest)] += 1
+        log.info(f"Matrix contains {od_matrix.sum()} park and ride tours")
+        demand = {}
+        # if estimation_mode:
+        #     omx_file = omx.open_file(f"{self.resultdata.path}/estimation/demand_{self.name}.omx","w")
+        #     omx_file.create_mapping("zone_number",self.zone_data.all_zone_numbers)
+
+        car_demand, transit_demand = self.park_and_ride_model.distribute_demand(od_matrix)
+        pnr_purpose = ParkAndRidePseudoPurpose(self)
+        demand["pnr_car"] = Demand(pnr_purpose, "car", car_demand)
+        demand["pnr_transit"] = Demand(pnr_purpose, "transit", transit_demand)
+
+        return demand
+    
     def calc_demand(self, estimation_mode=False, add_sec_dest: bool = True):
         """Calculate purpose specific demand matrices.
               

Scripts/modelsystem.py

@@ -731,6 +731,22 @@ def _add_internal_demand(self, previous_iter_impedance, is_last_iteration, estim
                 thread.join()
         for purpose in self.dm.tour_purposes:
             purpose.print_data()
+            #Park and ride
+            if purpose.park_and_ride_model is not None:
+                # Apply penalty for overcrowded park and ride facilities.
+                MAX_PNR_ITERATIONS = 5 # Maximum number of iterations. Set to 0 for no penalty
+                for i in range(MAX_PNR_ITERATIONS):
+                    modified = purpose.park_and_ride_model.apply_crowding_penalty()
+                    purpose.calc_basic_prob(saved_pnr_impedance[purpose.name])
+                    demand = purpose.calc_pnr_demand(self.dm.population,estimation_mode=estimation_mode, log=log)
+                    log.debug(f"Park and ride crowding penalty iteration {i+1} modified {modified} facilities.")
+                    if modified < 1:                                
+                        break
+                log.debug("Park and ride demand calculation completed.")
+        #Add park and ride to time periods
+        self.dtm.add_demand(demand["pnr_car"])
+        self.dtm.add_demand(demand["pnr_transit"])
+        #Add all tours to time periods
         for person in self.dm.population:
             for tour in person.tours:
                 self.dtm.add_demand(tour)