[Travelr-commits] r12 - pre-alpha pre-alpha/travelr pre-alpha/travelr/R pre-alpha/travelr/data pre-alpha/travelr/man pre-alpha/travelr/src pre-alpha/travelr/tests www
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Tue May 4 03:07:40 CEST 2010
Author: jrawbits
Date: 2010-05-04 03:07:40 +0200 (Tue, 04 May 2010)
New Revision: 12
Added:
pre-alpha/travelr/man/reaggregate.Rd
pre-alpha/travelr/tests/Test_05_SiouxFallsExample.R
pre-alpha/travelr_0.1-3.zip
Removed:
pre-alpha/travelr_0.1-2.zip
Modified:
pre-alpha/ReadMe.txt
pre-alpha/travelr/DESCRIPTION
pre-alpha/travelr/NAMESPACE
pre-alpha/travelr/R/Assignment_Set.R
pre-alpha/travelr/R/Highway_Assignment.R
pre-alpha/travelr/R/Iterative_Fitting.R
pre-alpha/travelr/R/zzz.R
pre-alpha/travelr/data/SiouxFalls.rda
pre-alpha/travelr/man/ipf.Rd
pre-alpha/travelr/man/skim.paths.Rd
pre-alpha/travelr/src/build_path_R.cc
pre-alpha/travelr/src/build_path_R.h
pre-alpha/travelr/src/short_path_R.cc
pre-alpha/travelr/tests/run.tests.bat
www/index.php
Log:
travelr package version 0.1-3
Complete assignment routines (lightly tested).
Some addition documentation, tests and examples.
Added a routine to reaggregate (aggregate/disaggregate) a matrix using correspondence tables.
Cleaned up the ipf (fratar) function so it works more smoothly.
Cleaned up the C++ code to conform a bit more cleanly with R's memory management and error handling.
Tidied the main web page a bit more.
Modified: pre-alpha/ReadMe.txt
===================================================================
--- pre-alpha/ReadMe.txt 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/ReadMe.txt 2010-05-04 01:07:40 UTC (rev 12)
@@ -7,7 +7,7 @@
for detailed information on what's in that package, where to get it, and how to
install it).
-Development (as of 4/29/2010) has been done with R 2.10.1 (even though 2.11.0 has
+Development (as of 5/3/2010) has been done with R 2.10.1 (even though 2.11.0 has
just been released) and a consistent Rtools package, and the package should
build with those tools.
Modified: pre-alpha/travelr/DESCRIPTION
===================================================================
--- pre-alpha/travelr/DESCRIPTION 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/DESCRIPTION 2010-05-04 01:07:40 UTC (rev 12)
@@ -1,5 +1,5 @@
Package: travelr
-Version: 0.1-2
+Version: 0.1-3
Date: 2009-04-26
Title: Functions for Travel Demand Modeling
Author: Jeremy Raw
Modified: pre-alpha/travelr/NAMESPACE
===================================================================
--- pre-alpha/travelr/NAMESPACE 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/NAMESPACE 2010-05-04 01:07:40 UTC (rev 12)
@@ -69,6 +69,7 @@
export(hwy.gamma.function)
export(rmse)
export(ipf)
+export(reaggregate.matrix)
S3method(print,iterative.fit)
# Low-level path functions
Modified: pre-alpha/travelr/R/Assignment_Set.R
===================================================================
--- pre-alpha/travelr/R/Assignment_Set.R 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/R/Assignment_Set.R 2010-05-04 01:07:40 UTC (rev 12)
@@ -40,7 +40,7 @@
build.paths<-function(aset,costs)
mapply(function(ac,cost,penalties) .shortest.paths(ac$network,cost,penalties),
- aset$classes,costs,MoreArgs=(penalties=aset$penalties),
+ aset$classes,costs,MoreArgs=list(penalties=aset$penalties),
USE.NAMES=TRUE, # Keep assignment class names on resulting list of path structures
SIMPLIFY=FALSE) # and leave results as a list
@@ -75,9 +75,13 @@
# skim paths returns a demand matrix by performing a function on a set of numeric values
# corresponding to the links on each path between origin and destination
skim.paths<-function( paths, costs, empty.val=0.0, FUN="sum", ... ) {
- if ( FUN=="sum" )
- return( .skim.paths(paths,costs,empty.val) )
- else
+ if ( FUN=="sum" ) {
+ (mapply(function(p,c) .skim.paths(p,c,empty.val),
+ paths, costs,
+ USE.NAMES=TRUE, # Keep assignment class names on resulting volume set list
+ SIMPLIFY=FALSE) # and leave the results as a list
+ )
+ } else
stop("Skimming with functions other than 'sum' is not yet implemented")
}
@@ -254,21 +258,23 @@
# Function to update an entire assignment class... Can't imagine why, but here it is!
hwy.update.class <- function( aset, name, aclass ) {
- if ( is.null(name) || is.null(class) )
+ if ( is.null(name) || is.null(aclass) )
stop("Cannot update assignment class unless name and new data are provided")
if ( class(aclass) != "highway.assignment.class" )
stop("Use make.assignment.class to construct the new class")
- if ( is.null(aset[[name]]) )
+ if ( is.null(aset$classes[[name]]) )
message("Adding assignment class named ",name)
aset$classes[[name]] <- aclass
}
# Function to update the demand for an assignment class once it's in the assignment set
hwy.update.demand <- function( aset, name, demand ) {
- if ( is.null(aset[[name]]) )
+ if ( is.null(name) )
+ stop("Cannot update assignment class demand unless name is provided")
+ if ( is.null(aset$classes[[name]]) )
warning("No class named ",name)
else if (!is.null(demand))
- aset[[name]]$demand<-demand
+ aset$classes[[name]]$demand<-demand
else
warning("No new demand: assignment set was unchanged")
return(aset)
Modified: pre-alpha/travelr/R/Highway_Assignment.R
===================================================================
--- pre-alpha/travelr/R/Highway_Assignment.R 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/R/Highway_Assignment.R 2010-05-04 01:07:40 UTC (rev 12)
@@ -75,7 +75,7 @@
)
parse.control<-function( control, control.element, default.value=NULL ) {
- cat("Parsing",control.element,"\n")
+# cat("Parsing",control.element,"\n")
element<-control[[control.element]]
if ( !is.null(element) ) {
result<-element
@@ -138,9 +138,9 @@
for ( i in test.min.names ) test.min <- c(test.min,parse.control(control,i))
test.max <- numeric(0)
for ( i in test.max.names ) test.max <- c(test.max,parse.control(control,i))
- cat("Convergence Test::\n")
- cat("Min:",paste(test.min.names,sep=","),"Max:",paste(test.max.names,sep=","),"\n")
- cat( "Min:",paste(test.min,collapse=","),"Max:",paste(test.max,collapse=","), "\n")
+# cat("Convergence Test::\n")
+# cat("Min:",paste(test.min.names,sep=","),"Max:",paste(test.max.names,sep=","),"\n")
+# cat( "Min:",paste(test.min,collapse=","),"Max:",paste(test.max,collapse=","), "\n")
function(results) {
any( unlist(results[min.names]) <= test.min ) ||
any( unlist(results[max.names]) >= test.max )
Modified: pre-alpha/travelr/R/Iterative_Fitting.R
===================================================================
--- pre-alpha/travelr/R/Iterative_Fitting.R 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/R/Iterative_Fitting.R 2010-05-04 01:07:40 UTC (rev 12)
@@ -1,6 +1,9 @@
####################################################################################################
# travelr\r\Iterative_Fitting.R by Jeremy Raw Copyright (C) 2010
-#
+#
+# Portions Copyright (C) 2002 Oregon Department of Transportation
+# and released under the GNU public license
+#
# This program is free software; you can redistribute it and/or modify
# it under the terms of the GNU General Public License as published by
# the Free Software Foundation; either version 2 of the License, or
@@ -24,29 +27,40 @@
# Gamma Function - For building gravity models with ipf
-hwy.gamma.function<-function(x,b,c) (x^b)*exp(x*c)
+hwy.gamma.function<-function(x,b,c) ifelse(x==0,0,(x^b)*exp(x*c))
# A function to perform iterative proportional fitting on an array with
-# 2 (or more?) dimensions -- how would that work for more than 2?
+# 2 (or more?) dimensions. The code has been tested with two dimensions
+# but in principle may have as many as you like.
ipf <- function(mtx,factors,method=NULL,max.rmse=1e-7,max.iter=50) {
- # a is an object with dimensions
- # Determine the dimensions of a, and return an error
+ # mtx is an object with dimensions
+ # Determine the dimensions of mtx, and return an error
# Note that default method (if method is not supplied) is "absolute"
- # Other options are "percent" or "fraction"
+ # Other options are "percent" or "fraction", in which case new factors
+ # are computed using the marginal totals of mtx
da<-dim(mtx)
la<-length(da)
dn<-dimnames(mtx)
- if ( is.null(da) || la!=2 )
- stop(sprintf(gettext("Parameter 'mtx' must have two dimensions, not %d"),length(la))) # for now...
+ if ( is.null(da) )
+ stop(sprintf("Parameter 'mtx' must have at least two dimensions"))
if ( is.list(factors) ) {
lf<-length(factors)
if ( ! all(sapply(factors,is.numeric)) )
stop("Growth factors must be numeric vectors")
- if ( lf==1 ) {
- factors[[2]]<-factors[[1]]
- lf=2
+ if ( lf<la ) {
+ warning(sprintf(gettext("Not enough growth factor vectors for %d dimensions: recycling"),la))
+ length.diff <- la-lf
+ copy.index <- 1
+ new.index <- lf+1
+ while ( length.diff>0 ) {
+ factors[[new.index]]<-factors[[copy.index]]
+ length.diff <- length.diff-1
+ new.index <- new.index+1
+ copy.index <- copy.index+1
+ }
}
+ lf <- length(factors)
} else if ( is.vector(factors) && is.numeric(factors) ) {
factors<-list(rows=factors,cols=factors)
lf<-2
@@ -56,7 +70,14 @@
if (lf!=la) {
stop('Must have as many growth factors as matrix dimensions')
}
- for ( i in seq_along(da) ) {
+ warn.zeroes<-0
+ dim.seq <- seq_along(da)
+ for ( i in dim.seq ) {
+ zero.factors<-which(factors[[i]]==0)
+ if ( any(zero.factors) ) {
+ factors[[i]][zero.factors]<-0.001
+ warn.zeroes<-warn.zeroes+1
+ }
if (da[i]!=length(factors[[i]])) {
warning(
sprintf(
@@ -65,10 +86,13 @@
)
)
}
- if ( any(factors[[i]]<0) ) {
- stop("All factors must be >= 0")
- }
+# if ( any(factors[[i]]<0) ) {
+# warning(paste(factors[[i]],collapse=", "),"\n","All factors must be >= 0")
+# factors[[i]][factors[[i]]<0]<-0
+# cat(factors[[i]])
+# }
}
+ if ( warn.zeroes>0 ) warning( warn.zeroes," zeroes in growth factors; adjusted to 0.001")
if ( !is.null(method) ) {
# For "growth factor" methods, the factors are specified relative to the existing margins
if ( method=="percent" ) {
@@ -82,29 +106,30 @@
} else {
method="absolute" # the core method factors the array so the margins match the factors
}
- if ( !isTRUE(all.equal(sum(factors[[1]]),sum(factors[[2]]))) ) {
- warning("Factors do not yield the same totals\nColumns will be adjusted to row totals")
- factors[[2]] <- factors[[2]]*(sum(factors[[1]])/sum(factors[[2]]))
- }
+# if ( !isTRUE(all.equal(sum(factors[[1]]),sum(factors[[2]]))) ) {
+# warning("Factors do not yield the same totals\nColumns will be adjusted to row totals")
+# factors[[2]] <- factors[[2]]*(sum(factors[[1]])/sum(factors[[2]]))
+# }
# Finally, we'll actually do the fitting...
iter<-0
- it.seq<-seq_along(factors)
- mtx<-t(mtx)
- new.rmse<-rmse(unlist(factors),c(colSums(mtx),rowSums(mtx))) # note that a is the transpose of the original matrix
+ new.rmse<-rmse(unlist(factors),c(rowSums(mtx),colSums(mtx)))
+ mtx.sums<-vector("list",0)
while( new.rmse>max.rmse && iter<max.iter ) {
- for (i in it.seq) mtx<-apply(mtx,1,function(v,rs)v*rs,rs=factors[[i]]/colSums(mtx))
- new.rmse<-rmse(unlist(factors),c(colSums(mtx),rowSums(mtx)))
+ for ( i in dim.seq ) {
+ mtx.sums[[i]]<-apply(mtx,i,sum)
+ mtx<-sweep(mtx,i,factors[[i]]/mtx.sums[[i]],"*")
+ }
+ new.rmse<-rmse(unlist(factors),unlist(mtx.sums))
iter<-iter+1
# cat("Iteration",iter,"RMSE:",new.rmse,"\n")
}
- final<-t(mtx)
- dimnames(final)<-dn
- attr(final,"Converged")<-iter<max.iter
- attr(final,"RMSE")<-new.rmse
- attr(final,"Iteration")<-iter
- class(final)<-c("iterative.fit",class(final)) # keep the original class identity (as matrix or whatever)
- return(final)
+ dimnames(mtx)<-dn
+ attr(mtx,"Converged")<-iter<max.iter
+ attr(mtx,"RMSE")<-new.rmse
+ attr(mtx,"Iteration")<-iter
+ class(mtx)<-c("iterative.fit",class(mtx)) # keep the original class identity (as matrix or whatever)
+ return(mtx)
}
print.iterative.fit <- function(x,...) {
@@ -116,3 +141,61 @@
cat("All",attr(x,"Iteration"),"Iterations Completed\n")
invisible(TRUE)
}
+
+###########################################################################
+# In keeping with the terms of the GNU Public License, source code for the
+# original ODOT code on which this version of ipf was based is here:
+#
+# # iterative proportional fitting function
+# ipf <- function(rowcontrol, colcontrol, seed, maxiter=50, closure=0.01){
+# # input data checks: sum of marginal totals equal and no zeros in marginal totals
+# #if(sum(rowcontrol) != sum(colcontrol)) stop("sum of rowcontrol must equal sum of colcontrol")
+# if(any(rowcontrol==0)){
+# numzero <- sum(rowcontrol==0)
+# rowcontrol[rowcontrol==0] <- 0.001
+# warning(paste(numzero, "zeros in rowcontrol argument replaced with 0.001", sep=" "))
+# }
+# if(any(colcontrol==0)){
+# numzero <- sum(colcontrol==0)
+# colcontrol[colcontrol==0] <- 0.001
+# warning(paste(numzero, "zeros in colcontrol argument replaced with 0.001", sep=" "))
+# }
+# # set initial values
+# result <- seed
+# rowcheck <- 1
+# colcheck <- 1
+# iter <- 0
+# # successively proportion rows and columns until closure or iteration criteria are met
+# while((rowcheck > closure) & (colcheck > closure) & (iter < maxiter))
+# {
+# rowtotal <- rowSums(result)
+# rowfactor <- rowcontrol/rowtotal
+# result <- sweep(result, 1, rowfactor, "*")
+# coltotal <- colSums(result)
+# colfactor <- colcontrol/coltotal
+# result <- sweep(result, 2, colfactor, "*")
+# rowcheck <- sum(abs(1-rowfactor))
+# colcheck <- sum(abs(1-colfactor))
+# iter <- iter + 1
+# }
+# result
+# }
+
+reaggregate.matrix <- function( mtx, eq, nrow=NULL, ncol=NULL ) {
+ # eq is a data.frame (or data matrix) that has columns i,j,o,p,fact
+ # to build such an equivalence table, the function expand.grid is useful
+ # example:
+ if ( is.null(nrow) ) nrow <- max(eq$o)
+ else if ( nrow < max(eq$o) ) eq <- eq[-which(eq$o > nrow),]
+ # else just use nrow as it is
+ if ( is.null(ncol) ) ncol<-max(eq$p)
+ else if ( ncol < max(eq$p) ) eq <- eq[-which(eq$p > ncol),]
+ # else just use ncol as it is
+
+ mtx.out <- matrix(0,nrow=nrow,ncol=ncol)
+ oval <- aggregate( data.frame(value=(mtx[data.matrix(eq[,c("i","j")])] * eq$fact)),
+ by=list(o=eq$o,p=eq$p),
+ sum)
+ mtx.out[data.matrix(oval[,c("o","p")])] <- oval$value
+ mtx.out
+}
Modified: pre-alpha/travelr/R/zzz.R
===================================================================
--- pre-alpha/travelr/R/zzz.R 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/R/zzz.R 2010-05-04 01:07:40 UTC (rev 12)
@@ -1,3 +1,3 @@
.onLoad <- function(libname,pkgname) {
- packageStartupMessage("travelr is SERIOUSLY in a 'pre-alpha' state.\n","Lots of Luck!")
+ packageStartupMessage("travelr is SERIOUSLY in a 'pre-alpha' state... Lots of Luck!")
}
\ No newline at end of file
Modified: pre-alpha/travelr/data/SiouxFalls.rda
===================================================================
(Binary files differ)
Modified: pre-alpha/travelr/man/ipf.Rd
===================================================================
--- pre-alpha/travelr/man/ipf.Rd 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/man/ipf.Rd 2010-05-04 01:07:40 UTC (rev 12)
@@ -13,13 +13,13 @@
}
\arguments{
\item{mtx}{
- A two-dimensinoal matrix that will be expanded to match the margins contained in \code{factors}
+ A matrix or array with at least two dimensions that will be expanded to match the margins contained in \code{factors}
}
\item{factors}{
- A vector, list, or \code{data.frame} of factors indicating how to grow the matrix
+ A vector, list, or \code{data.frame} of factors indicating how to grow the matri; See details
}
\item{method}{
- How to interpret the factors. Choices are \dQuote{absolute}, \dQuote{fraction} or \dQuote{percent}. See details.
+ How to interpret the factors. Choices are \dQuote{absolute}, \dQuote{fraction} or \dQuote{percent}; See details
}
\item{max.rmse}{
Convergence limit. Terminate fitting when the root mean square error between the expanded array margins and the target
@@ -44,21 +44,22 @@
\code{factors} must be a list.
The first element of \code{factors} corresponds to the desired row sums of the matrix, and the second element
-corresponds to the desired column sums of the matrix. If \code{factors} is a vector, the same vector will be used to
-establish the row and column sums.
+corresponds to the desired column sums of the matrix. Likewise if \code{mtx} has additional dimensions. If
+\code{factors} is a vector, the same vector will be used to establish the sums in all dimensions, but in that
+case all the dimensions must have the same lengths (so that is typically only useful for two-dimensional matrices).
-The \code{factors} must each be the same length as the corresponding dimension of the input matrix.
+The elements of \code{factors} must each be the same length as the corresponding dimension of the input matrix.
-The \code{method} parameter indicates how the \code{factors} are interpreted. The method may be \dQuote{absolute}, \dQuote{fraction},
-or \dQuote{percent}. If \code{method} is \dQuote{absolute}, the factors are taken to be the new margins, and the array is iteratively
-expanded until the margins match the new totals. If \code{method} is \dQuote{fraction} or \dQuote{percent} methods, the supplied
-factors are converted to absolute values by multiplying the corresponding row or column sums of the matrix by the
-corresponding factors.
+The \code{method} parameter indicates how the \code{factors} are interpreted. The method may be \dQuote{absolute},
+\dQuote{fraction}, or \dQuote{percent}. If \code{method} is \dQuote{absolute}, the factors are taken to be the new
+margins, and the array is iteratively expanded until the margins match the new totals. If \code{method} is
+\dQuote{fraction} or \dQuote{percent} methods, the supplied factors are converted to absolute values by multiplying the
+corresponding row or column sums of the matrix by the corresponding factors.
-Note that the \code{factors} must yield the same grand total for both rows and columns of \code{mtx}. After the
-absolute factors have been prepared according to the chosen \code{method}, if the total of the two sets of growth
-factors is not the same, the second growth factor (column) will itself be factored to match the total of the first
-growth factor and a warning will be issued.
+Note that the \code{factors} must yield the same grand total for all dimensions of \code{mtx}, which will be the same as
+\code{sum(mtx)} once the function is complete. After the absolute factors have been prepared according to the chosen
+\code{method}, if the sum of the any of the sets of growth factors is not the same, the second and subsequent growth
+factor s(column) will itself be factored to match the total of the first growth factor and a warning will be issued.
The print method shows computed array margins and convergence statistics.
}
Added: pre-alpha/travelr/man/reaggregate.Rd
===================================================================
--- pre-alpha/travelr/man/reaggregate.Rd (rev 0)
+++ pre-alpha/travelr/man/reaggregate.Rd 2010-05-04 01:07:40 UTC (rev 12)
@@ -0,0 +1,83 @@
+\name{reaggregate.matrix}
+\alias{reaggregate}
+\alias{reaggregate.matrix}
+\title{
+Aggregate or disaggregate at matrix
+}
+\description{
+Change the size of a matrix by adding or subdividing its cells.
+}
+\usage{
+reaggregate.matrix(mtx, eq, nrow=NULL, ncol=NULL)
+}
+\arguments{
+ \item{mtx}{
+ A matrix that will be refactored
+}
+ \item{eq}{
+ A \code{data.frame} which is a correspondence table describing how to compute output cells from input cells
+}
+ \item{nrow}{
+ An explict number of rows in the result matrix; the default is to use the maximum value of \code{eq$o}
+}
+ \item{ncol}{
+ An explict number of columns in the result matrix; the default is to use the maximum value of \code{eq$p}
+}
+}
+\details{
+
+\code{reaggregate.matrix} will compute a new matrix with revised dimensions from an input matrix, \code{mtx}, based on
+the correspondences described in the equivalence table, \code{eq}.
+
+The equivalence table must be a \code{data.frame} or \code{matrix} containing the following columns:
+\describe{
+ \item{i}{The row index of an input matrix cell}
+ \item{j}{The column index of an input matrix cell}
+ \item{o}{The row index of an output matrix cell}
+ \item{p}{The column index of an output matrix cell}
+ \item{fact}{The fraction of cell \code{mtx[i, j]} to be placed in cell \code{mtx.out[o, p]} }
+}
+
+Each row in the equivalence table describes how to move some or all of the value in the \code{[i,j]}
+cell in the input matrix to the corresponding $\code{[o,p]} cll in the output matrix. The value in
+the input cell is multiplied by \code{eq$fact} and all the values contributing to the value of cell
+$\code{[o,p]} in the output are added together to form the final result.
+
+If a particular \code{i,j} pair is not represented in the equivalence table then that pair will not contribute to the
+output matrix. Indices that are beyond the dimensions of the input matrix may introduce \code{NA}'s into an output
+cell, even if other correspondences introduce known values.
+
+With suitable correspondence factors, matrices may be expanded or reduced in size.
+}
+
+\value{
+A \code{matrix} with dimensions \code{[nrow,ncol]} computed by applying correspondence factors.
+}
+
+\author{ Jeremy Raw }
+
+\examples{
+ # The basic mechanics of aggregating a matrix...
+ m <- matrix(1:100,nrow=10,ncol=10)
+ eq <- expand.grid(i=1:10,j=1:10)
+ eq$fact <- 1
+ eq$o <- ceiling(eq$i/2)
+ eq$p <- ceiling(eq$j/2)
+ om <- reaggregate.matrix(m,eq)
+
+ # Disaggregating is trickier since factors are necessary
+ req <- eq[,c("i","j","fact","o","p")]
+ names(req) <- c("o","p","fact","i","j")
+ req$fact <- 0.25
+ rom <- reaggregate.matrix(om,req)
+
+ # Building an equivalence table using a correspondence table
+ # Used in travel modeling to convert zones to districts
+ ct <- data.frame( from=c( 1, 2, 3, 4, 5, 6, 7, 8, 9,10),
+ to=c( 2, 2, 1, 3, 3, 3, 4, 4, 5, 5) )
+ c.eq <- expand.grid(i=1:10,j=1:10)
+ c.eq$o <- ct$to[match(c.eq$i,ct$from)]
+ c.eq$p <- ct$to[match(c.eq$j,ct$from)]
+ c.eq$fact <- 1
+ c.om <- reaggregate.matrix(m,c.eq)
+}
Modified: pre-alpha/travelr/man/skim.paths.Rd
===================================================================
--- pre-alpha/travelr/man/skim.paths.Rd 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/man/skim.paths.Rd 2010-05-04 01:07:40 UTC (rev 12)
@@ -69,7 +69,6 @@
}
\examples{
\dontrun{
-# A "three-step" model (no mode choice) with trip distribution via gravity model
data(SiouxFalls)
# Trip Generation
@@ -77,19 +76,20 @@
attractions<-colSums(SiouxFalls.od)
# Highway Skims
-cost.function<-with(SiouxFalls.net$links,function(...)FFTime)
+cost.function<-with(SiouxFalls.net$Links,function(...)FFTime)
aclass <- make.assignment.class(SiouxFalls.net,"All",SiouxFalls.od)
aset <- new.assignment.set(SiouxFalls.net,list(All=aclass),cost.volume.type="vector",cost.function=cost.function)
paths <- build.paths(aset,aset$ff.cost)
-travel.times <- skim.paths(paths,aset$ff.cost)
+travel.times <- skim.paths(paths,aset$ff.cost)[["All"]] # only one purpose: "All trips"
# Trip Distribution (Gravity Model with gamma function)
base.distribution <- hwy.gamma.function(travel.times,-0.02,-0.123) # HBW coefficients from NCHRP 365
trip.table <- ipf(base.distribution,list(rows=productions, cols=attractions),method="absolute")
+print(round(trip.table,2))
aset <- hwy.update.demand(aset,"All",trip.table)
# Trip Assignment
-assignment.results <- highway.assign(aset,method="AON")
+assignment.results <- highway.assign(aset,method="Frank.Wolfe")
loaded.links <- assignment.results$volumes
}
}
Modified: pre-alpha/travelr/src/build_path_R.cc
===================================================================
--- pre-alpha/travelr/src/build_path_R.cc 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/src/build_path_R.cc 2010-05-04 01:07:40 UTC (rev 12)
@@ -125,7 +125,7 @@
static Candidate Null_Element;
public:
Heap(int capacity);
- ~Heap();
+// ~Heap();
void Clear() { size=0; }
int isEmpty() const { return size==0; }
int getSize() const { return size; }
@@ -137,16 +137,20 @@
Candidate PathBuilder::Heap::Null_Element; // Has lowest possible cost (0.0)
-PathBuilder::Heap::Heap( int capacity ) : capacity(capacity), size(0) {
- candidates = new Candidate[ capacity ];
- heap = new Candidate * [ capacity + 1 ];
+PathBuilder::Heap::Heap( int capacity ) : capacity(capacity), size(0), heap(0), candidates(0) {
+// Changed to use safer memory allocation (so R error trapping can occur safely)
+ candidates = (Candidate *) R_alloc(capacity, sizeof(Candidate));
+ heap = (Candidate **) R_alloc(capacity+1, sizeof(Candidate*));
+// candidates = new Candidate[ capacity ];
+// heap = new Candidate * [ capacity + 1 ];
heap[ 0 ] = &Null_Element; // Makes the math much easier!
}
-PathBuilder::Heap::~Heap() {
- delete [] candidates;
- delete [] heap;
-}
+// PathBuilder::Heap::~Heap() {
+// Changed for use with R_alloc -- memory will be freed when R .Call returns
+// delete [] candidates;
+// delete [] heap;
+// }
int PathBuilder::Heap::Push(Candidate * candidate) {
int i = ++size; // Set available position to new end of heap
Modified: pre-alpha/travelr/src/build_path_R.h
===================================================================
--- pre-alpha/travelr/src/build_path_R.h 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/src/build_path_R.h 2010-05-04 01:07:40 UTC (rev 12)
@@ -24,11 +24,11 @@
#undef DEBUG_TRAVELR
#ifdef DEBUG_TRAVELR
-#undef DEBUG_TRAVELR_PATH
-#define DEBUG_TRAVELR_LOAD
-#define DEBUG_TRAVELR_SKIM
-#define DEBUG_TRAVELR_WALK
-#define DEBUG_TRAVELR_ICPT
+#undef DEBUG_TRAVELR_PATH
+#undef DEBUG_TRAVELR_LOAD
+#undef DEBUG_TRAVELR_SKIM
+#undef DEBUG_TRAVELR_WALK
+#undef DEBUG_TRAVELR_ICPT
#endif
struct EdgeList {
Modified: pre-alpha/travelr/src/short_path_R.cc
===================================================================
--- pre-alpha/travelr/src/short_path_R.cc 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/src/short_path_R.cc 2010-05-04 01:07:40 UTC (rev 12)
@@ -299,9 +299,10 @@
}
extern "C"
- SEXP skim_paths( SEXP rSPF, SEXP rCostMatrix, SEXP rEmptyVal )
+ SEXP skim_paths( SEXP rSPF, SEXP rCostVector, SEXP rEmptyVal )
{
// Get number zones, number of nodes and number of links from tree attributes.
+
NetworkParameters np(rSPF);
double empty_val = REAL(rEmptyVal)[0];
@@ -312,7 +313,7 @@
dm.setAll(empty_val);
ShortestPathForest SPF(np.getNodes(),np.getZones(),np.getLinks(),INTEGER(rSPF));
- SkimVector<double> cv(REAL(rCostMatrix), np.getLinks());
+ SkimVector<double> cv(REAL(rCostVector), np.getLinks());
PathWalker<double,double> ( SPF, cv, dm ).Skim(empty_val); // accumulate cost vector into matrix cells along paths
Added: pre-alpha/travelr/tests/Test_05_SiouxFallsExample.R
===================================================================
--- pre-alpha/travelr/tests/Test_05_SiouxFallsExample.R (rev 0)
+++ pre-alpha/travelr/tests/Test_05_SiouxFallsExample.R 2010-05-04 01:07:40 UTC (rev 12)
@@ -0,0 +1,35 @@
+# A "three-step" model (no mode choice) with trip distribution via gravity model
+library(travelr)
+# load("data/SiouxFallsAset.Rdata")
+data(SiouxFalls)
+
+# Trip Generation
+str(SiouxFalls.net)
+productions<-rowSums(SiouxFalls.od)
+attractions<-colSums(SiouxFalls.od)
+print(data.frame(P=productions,A=attractions))
+cat("Total productions:",sum(productions),"\n")
+cat("Total attractions:",sum(attractions),"\n")
+
+# Highway Skims
+cost.function<-with(SiouxFalls.net$Links,function(...)FFTime)
+aclass <- make.assignment.class(SiouxFalls.net,"All",SiouxFalls.od)
+aset <- new.assignment.set(SiouxFalls.net,list(All=aclass),cost.volume.type="vector",cost.function=cost.function)
+paths <- build.paths(aset,aset$ff.cost)
+travel.times <- skim.paths(paths,aset$ff.cost)[["All"]] # only one purpose: "All trips"
+
+# Trip Distribution (Gravity Model with gamma function)
+options(width=180)
+base.distribution <- hwy.gamma.function(travel.times,-0.02,-0.123) # HBW coefficients from NCHRP 365
+trip.table <- ipf(base.distribution,list(rows=productions, cols=attractions),method="absolute")
+print(round(trip.table,2))
+
+# Note that either of the following is possible
+# The function version will give nicer error messages
+aset[[c("classes","All","demand")]]<-trip.table
+aset <- hwy.update.demand(aset,"All",trip.table)
+
+# Trip Assignment
+assignment.results <- highway.assign(aset,method="Frank.Wolfe")
+loaded.links <- assignment.results$volumes
+print(loaded.links)
Modified: pre-alpha/travelr/tests/run.tests.bat
===================================================================
--- pre-alpha/travelr/tests/run.tests.bat 2010-04-30 01:39:45 UTC (rev 11)
+++ pre-alpha/travelr/tests/run.tests.bat 2010-05-04 01:07:40 UTC (rev 12)
@@ -15,4 +15,6 @@
if not exist data\SiouxFallsAset.Rdata Rscript Test_02_Assignment_Classes.R
if exist data\SiouxFalls.rda if exist data\SiouxFallsAset.Rdata Rscript Test_03_Highway_Assignment.R
Rscript Test_04_ipf.R
-Rscript Test_Cost_Integrator.R
\ No newline at end of file
+Rscript Test_Cost_Integrator.R
+Rscript Test_05_SiouxFallsExample.R
+
Deleted: pre-alpha/travelr_0.1-2.zip
===================================================================
(Binary files differ)
Added: pre-alpha/travelr_0.1-3.zip
===================================================================
(Binary files differ)
Property changes on: pre-alpha/travelr_0.1-3.zip
___________________________________________________________________
Added: svn:mime-type
+ application/octet-stream
Modified: www/index.php
===================================================================
--- www/index.php 2010-04-30 01:39:45 UTC (rev 11)
+++ www/index.php 2010-05-04 01:07:40 UTC (rev 12)
@@ -46,13 +46,18 @@
analysis, and network skims.</p>
<h3>Project Status</h3>
-<p> The TravelR project's first package, not suprisingly named "<strong>travelr</strong>", is in "pre-alpha" state (the code is not quite
- complete, what is there sort of works, but what works and what doesn't is changing daily and sometimes hourly).
- A very preliminary version of the package is expected to be surreptitiously released some time in May 2010, and
- announced to a close circle of co-conspiratoRs. Subscribe to the travelr-announce list, or email the lead
- developer through the <a href="http://<?php echo $domain; ?>/projects/<?php echo $group_name; ?>/">project page</a>
- if you would like to join the inner circle.</p>
+ <p>
+ The definitive project status is maintained on the
+ <a href="http://<?php echo $domain; ?>/projects/<?php echo $group_name; ?>/">TravelR project page</a>.
+ </p>
+ <p> The TravelR project's first package, not suprisingly named "<strong>travelr</strong>", is in "pre-alpha" state
+(the code is not quite complete, what is there sort of works, but what works and what doesn't is changing daily and
+sometimes hourly). A very preliminary version of the package is expected to be surreptitiously released some time in
+May 2010, and announced to a close circle of co-conspiratoRs. Subscribe to the <a
+href="https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/travelr-announce">travelr-announce</a> e-mail list to
+see what's happening, or email the lead developer if you would like to join the inner circle.</p>
+
<h3>History, Motivation and Philosophy</h3>
<p> Some travel model developers and users have been using R very happily for years, but the ability to do a complete
model in R has eluded us since there hasn't been a good way to do
@@ -61,7 +66,9 @@
professional software developer).
Please be aware that TravelR is a "hobby" project for which no one is currently getting paid. Work on the project
- happens in Jeremy's (rather limited) spare time. Please attend to the
+ happens in Jeremy's (rather limited) spare time.</p>
+
+<p>Please attend to the
<a href="http://en.wikipedia.org/wiki/Yamas">Yamas</a> and
<a href="http://en.wikipedia.org/wiki/Niyamas">Niyamas</a>
as you contemplate this body of work.
@@ -79,29 +86,28 @@
Model developers who would like to build innovative and useful models in R for practical application.
</li>
</ul>
-
<h3>Design Goals</h3>
-<p>
- To support these two audiences, there are two design goals:
- <dl>
- <dt>Truly Open Source</dt>
- <dd>
- Source code is not open if it is incomprehensible. We want the <strong>travelr</strong> code to be clear,
- well-structured and well-documented; to use standard R approaches to data manipulation whenever possible;
- to be fast enough so small test problems can be solved almost instantaneously; and ultimately to be
- irresistible to researchers because you won't have to write any of the hard stuff all over again.<br/>(On the
- other hand, we would also like to be able to get into the Scorpion pose.)
- </dd>
- <dt>Industrial-strength features</dt>
- <dd>
- It is one thing to explore algorithms on simplified networks. It is something else again to perform real-world
- analyses. Because the primary users of R for travel demand modeling are what the industry calls "practitioners"
- (as opposed to "academics" or "researchers"), we find ourselves wanting to solve real problems with flexible,
- powerful tools. Certain indispensable features don't have open-source implementations at all, let alone ones
- that might actually be practical. TravelR aims to rectify that shortcoming.
- </dd>
- </dl>
-</p>
+<p>To support these two audiences, there are two design goals:</p>
+ <ul><li>Truly Open Source
+ <ul><li>
+ Source code is not open if it is incomprehensible. We want the <strong>travelr</strong> code to be clear,
+ well-structured and well-documented; to use standard R approaches to data manipulation whenever possible;
+ to be fast enough so small test problems can be solved almost instantaneously; and ultimately to be
+ irresistible to researchers because you won't have to write any of the hard stuff all over again.
+ </li>
+ </ul>
+ </li>
+ <li>Industrial-strength features
+ <ul><li>
[TRUNCATED]
To get the complete diff run:
svnlook diff /svnroot/travelr -r 12
More information about the Travelr-commits
mailing list