Merge pull request #167 from ATFutures/proportions-only

add proportions_only parameter for #144

R/RcppExports.R

@@ -308,8 +308,8 @@ rcpp_get_paths <- function(graph, vert_map_in, fromi, toi_in, heap_type) {
 #' Implemented in parallal form only; no single-threaded version, and
 #' only for AStar (so graphs must be spatial).
 #' @noRd
-rcpp_get_sp_dists_proportional <- function(graph, vert_map_in, fromi, toi_in, heap_type) {
-    .Call(`_dodgr_rcpp_get_sp_dists_proportional`, graph, vert_map_in, fromi, toi_in, heap_type)
+rcpp_get_sp_dists_proportional <- function(graph, vert_map_in, fromi, toi_in, heap_type, proportions_only) {
+    .Call(`_dodgr_rcpp_get_sp_dists_proportional`, graph, vert_map_in, fromi, toi_in, heap_type, proportions_only)
 }
 
 #' rcpp_gen_hash
@@ -368,4 +368,3 @@ rcpp_points_index_par <- function(xy, pts) {
 rcpp_route_times <- function(graph, left_side, turn_penalty) {
     .Call(`_dodgr_rcpp_route_times`, graph, left_side, turn_penalty)
 }
-

R/dists-proportional.R

@@ -5,6 +5,10 @@
 #' graph which must have a column named "edge_type" which labels categories of
 #' edge types along which proportional distances are to be aggregated (see
 #' Note).
+#' @param proportions_only If `FALSE`, return full distance matrices for full
+#' distances and for each edge category; if `TRUE`, return single vector of
+#' proportional distances, like current `summary` function applied to full
+#' results. See Note.
 #' @return A list of distance matrices of equal dimensions (length(from),
 #' length(to)), the first of which ("distance") holds the final distances, while
 #' the rest are one matrix for each unique value of "edge_type", holding the
@@ -13,10 +17,17 @@
 #' @note The "edge_type" column in the graph can contain any kind of discrete or
 #' categorical values, although integer values of 0 are not permissible. `NA`
 #' values are ignored.
+#'
+#' @note Setting the `proportions_only` flag to `TRUE` may be advantageous for
+#' large jobs, because this avoids construction of the full matrices. Although
+#' this generally won't notably speed up calculations, it may make possible
+#' calculations which would otherwise require distance matrices too large to be
+#' directly stored.
 #' @export
 dodgr_dists_proportional <- function (graph,
                                       from = NULL,
                                       to = NULL,
+                                      proportions_only = FALSE,
                                       heap = "BHeap",
                                       quiet = TRUE) {
 
@@ -64,17 +75,29 @@ dodgr_dists_proportional <- function (graph,
                                          vert_map,
                                          from_index$index,
                                          to_index$index,
-                                         heap)
+                                         heap,
+                                         proportions_only)
+
+    n <- length (to)
+
+    if (!proportions_only) {
 
-    n <-length (to)
-    d0 <- list ("distances" = d [, seq (n)])
-    d <- lapply (seq_along (edge_type_table), function (i) {
-                     index <- i * n + seq (n) - 1
-                     d [, index]    })
-    names (d) <- names (edge_type_table)
+        d0 <- list ("distances" = d [, seq (n)])
+        d <- lapply (seq_along (edge_type_table), function (i) {
+                         index <- i * n + seq (n) - 1
+                         d [, index]    })
+        names (d) <- names (edge_type_table)
 
-    res <- c (d0, d)
-    class (res) <- append (class (res), "dodgr_dists_proportional")
+        res <- c (d0, d)
+        class (res) <- append (class (res), "dodgr_dists_proportional")
+
+    } else {
+
+        res <- apply (d, 2, sum)
+        res [2:length (res)] <- res [2:length (res)] / res [1]
+        res <- res [-1]
+        names (res) <- names (edge_type_table)
+    }
 
     return (res)
 }

src/RcppExports.cpp

@@ -243,8 +243,8 @@ BEGIN_RCPP
 END_RCPP
 }
 // rcpp_get_sp_dists_proportional
-Rcpp::NumericMatrix rcpp_get_sp_dists_proportional(const Rcpp::DataFrame graph, const Rcpp::DataFrame vert_map_in, Rcpp::IntegerVector fromi, Rcpp::IntegerVector toi_in, const std::string& heap_type);
-RcppExport SEXP _dodgr_rcpp_get_sp_dists_proportional(SEXP graphSEXP, SEXP vert_map_inSEXP, SEXP fromiSEXP, SEXP toi_inSEXP, SEXP heap_typeSEXP) {
+Rcpp::NumericMatrix rcpp_get_sp_dists_proportional(const Rcpp::DataFrame graph, const Rcpp::DataFrame vert_map_in, Rcpp::IntegerVector fromi, Rcpp::IntegerVector toi_in, const std::string& heap_type, const bool proportions_only);
+RcppExport SEXP _dodgr_rcpp_get_sp_dists_proportional(SEXP graphSEXP, SEXP vert_map_inSEXP, SEXP fromiSEXP, SEXP toi_inSEXP, SEXP heap_typeSEXP, SEXP proportions_onlySEXP) {
 BEGIN_RCPP
     Rcpp::RObject rcpp_result_gen;
     Rcpp::RNGScope rcpp_rngScope_gen;
@@ -253,7 +253,8 @@ BEGIN_RCPP
     Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type fromi(fromiSEXP);
     Rcpp::traits::input_parameter< Rcpp::IntegerVector >::type toi_in(toi_inSEXP);
     Rcpp::traits::input_parameter< const std::string& >::type heap_type(heap_typeSEXP);
-    rcpp_result_gen = Rcpp::wrap(rcpp_get_sp_dists_proportional(graph, vert_map_in, fromi, toi_in, heap_type));
+    Rcpp::traits::input_parameter< const bool >::type proportions_only(proportions_onlySEXP);
+    rcpp_result_gen = Rcpp::wrap(rcpp_get_sp_dists_proportional(graph, vert_map_in, fromi, toi_in, heap_type, proportions_only));
     return rcpp_result_gen;
 END_RCPP
 }
@@ -324,7 +325,7 @@ static const R_CallMethodDef CallEntries[] = {
     {"_dodgr_rcpp_get_iso", (DL_FUNC) &_dodgr_rcpp_get_iso, 5},
     {"_dodgr_rcpp_get_sp_dists", (DL_FUNC) &_dodgr_rcpp_get_sp_dists, 5},
     {"_dodgr_rcpp_get_paths", (DL_FUNC) &_dodgr_rcpp_get_paths, 5},
-    {"_dodgr_rcpp_get_sp_dists_proportional", (DL_FUNC) &_dodgr_rcpp_get_sp_dists_proportional, 5},
+    {"_dodgr_rcpp_get_sp_dists_proportional", (DL_FUNC) &_dodgr_rcpp_get_sp_dists_proportional, 6},
     {"_dodgr_rcpp_gen_hash", (DL_FUNC) &_dodgr_rcpp_gen_hash, 2},
     {"_dodgr_rcpp_sf_as_network", (DL_FUNC) &_dodgr_rcpp_sf_as_network, 2},
     {"_dodgr_rcpp_points_index_par", (DL_FUNC) &_dodgr_rcpp_points_index_par, 2},

src/run_sp_proportional.cpp

@@ -104,6 +104,91 @@ struct OneProportionalDist : public RcppParallel::Worker
 
 };
 
+// Modified version of OneProportionalDist to aggregate only a vector of
+// distances, one value for each `from`.
+struct OneProportion : public RcppParallel::Worker
+{
+    RcppParallel::RVector <int> dp_fromi;
+    const std::vector <size_t> toi;
+    const std::vector <size_t> edge_type;
+    const size_t nverts;
+    const std::vector <double> vx;
+    const std::vector <double> vy;
+    const std::shared_ptr <DGraph> g;
+    const std::string heap_type;
+    const size_t num_edge_types;
+
+    RcppParallel::RMatrix <double> dout;
+
+    // constructor
+    OneProportion (
+            const RcppParallel::RVector <int> fromi,
+            const std::vector <size_t> toi_in,
+            const std::vector <size_t> edge_type_in,
+            const size_t nverts_in,
+            const std::vector <double> vx_in,
+            const std::vector <double> vy_in,
+            const std::shared_ptr <DGraph> g_in,
+            const std::string & heap_type_in,
+            const size_t & num_edge_types_in,
+            RcppParallel::RMatrix <double> dout_in) :
+        dp_fromi (fromi), toi (toi_in), edge_type (edge_type_in),
+        nverts (nverts_in), vx (vx_in), vy (vy_in),
+        g (g_in), heap_type (heap_type_in), 
+        num_edge_types (num_edge_types_in),
+        dout (dout_in)
+    {
+    }
+
+    // Parallel function operator
+    void operator() (std::size_t begin, std::size_t end)
+    {
+        std::shared_ptr<PF::PathFinder> pathfinder =
+            std::make_shared <PF::PathFinder> (nverts,
+                    *run_sp::getHeapImpl (heap_type), g);
+        std::vector <double> w (nverts);
+        std::vector <double> d (nverts * (num_edge_types + 1));
+        std::vector <long int> prev (nverts);
+
+        std::vector <double> heuristic (nverts, 0.0);
+
+        const size_t nto = toi.size ();
+
+        for (std::size_t i = begin; i < end; i++)
+        {
+            size_t from_i = static_cast <size_t> (dp_fromi [i]);
+
+            // only implemented for spatial graphs
+            for (size_t j = 0; j < nverts; j++)
+            {
+                const double dx = vx [j] - vx [from_i],
+                    dy = vy [j] - vy [from_i];
+                heuristic [j] = sqrt (dx * dx + dy * dy);
+            }
+            pathfinder->AStarEdgeType (d, w, prev, heuristic, from_i, toi);
+
+            for (size_t j = 0; j < toi.size (); j++)
+            {
+                if (w [toi [j]] < INFINITE_DOUBLE)
+                {
+                    for (size_t k = 0; k < num_edge_types; k++)
+                    {
+                        const double dto = d [toi [j] + k * nverts];
+                        if (dto < INFINITE_DOUBLE) {
+                            //if (dout (i, k) == NA_REAL)
+                            if (Rcpp::NumericMatrix::is_na (dout (i, k)))
+                                dout (i, k) = dto;
+                            else
+                                dout (i, k) += dto;
+                        }
+                    }
+                }
+            }
+        }
+    }
+
+};
+
 size_t proportional::num_edge_types (const std::vector <size_t> &edge_type)
 {
     std::unordered_set <size_t> type_set;
@@ -222,7 +307,8 @@ Rcpp::NumericMatrix rcpp_get_sp_dists_proportional (const Rcpp::DataFrame graph,
         const Rcpp::DataFrame vert_map_in,
         Rcpp::IntegerVector fromi,
         Rcpp::IntegerVector toi_in,
-        const std::string& heap_type)
+        const std::string& heap_type,
+        const bool proportions_only)
 {
     std::vector <size_t> toi =
         Rcpp::as <std::vector <size_t> > ( toi_in);
@@ -253,19 +339,35 @@ Rcpp::NumericMatrix rcpp_get_sp_dists_proportional (const Rcpp::DataFrame graph,
     inst_graph (g, nedges, vert_map, from, to, edge_type, dist, wt);
 
     // One standard [nfrom, nto] matrix plus one extra for each value of edge_type:
-    Rcpp::NumericVector na_vec = Rcpp::NumericVector (nfrom * nto * (num_types + 1L),
+    size_t ncol;
+    if (proportions_only)
+        ncol = num_types + 1L;
+    else
+        ncol = nto * (num_types + 1L);
+
+    Rcpp::NumericVector na_vec = Rcpp::NumericVector (nfrom * ncol,
             Rcpp::NumericVector::get_na ());
     Rcpp::NumericMatrix dout (static_cast <int> (nfrom),
-            static_cast <int> (nto * (num_types + 1L)), na_vec.begin ());
+            static_cast <int> (ncol));
 
     // Create parallel worker
-    OneProportionalDist one_dist (RcppParallel::RVector <int> (fromi), toi,
-            edge_type, nverts, vx, vy,
-            g, heap_type, num_types,
-            RcppParallel::RMatrix <double> (dout));
-
     size_t chunk_size = run_sp::get_chunk_size (nfrom);
-    RcppParallel::parallelFor (0, nfrom, one_dist, chunk_size);
+    if (proportions_only)
+    {
+        OneProportion one_dist (RcppParallel::RVector <int> (fromi), toi,
+                edge_type, nverts, vx, vy,
+                g, heap_type, num_types,
+                RcppParallel::RMatrix <double> (dout));
+        RcppParallel::parallelFor (0, nfrom, one_dist, chunk_size);
+    } else
+    {
+        OneProportionalDist one_dist (RcppParallel::RVector <int> (fromi), toi,
+                edge_type, nverts, vx, vy,
+                g, heap_type, num_types,
+                RcppParallel::RMatrix <double> (dout));
+        RcppParallel::parallelFor (0, nfrom, one_dist, chunk_size);
+    }
+
 
     return (dout);
 }