Code cleanup, moved main simulation struct to simulation.cpp

generate.sh

@@ -13,10 +13,25 @@ function run_sim {
 }
 
 function Two_Real {
+        REAL_LANES=4
+        VIRTUAL_LANES=0
+        LANE_CHANGE=1
+        DIR=Two_Real
+        mkdir $DIR
+        cd $DIR
+
+        for car_ratio in `seq 0 0.05 1`;
+        do
+                run_sim
+        done
+        cd ..
+}
+
+function Two_Real_NoLaneChange {
 	REAL_LANES=4
 	VIRTUAL_LANES=0
 	LANE_CHANGE=0
-	DIR=Two_Real
+	DIR=Two_Real_NoLaneChange
 	mkdir $DIR
 	cd $DIR
 
@@ -25,7 +40,6 @@ function Two_Real {
                 run_sim
         done
         cd ..
-
 }
 
 

hdf_save.h

@@ -1,6 +1,10 @@
 #ifndef _HDF_SAVE_COMPRESS_H
 #define _HDF_SAVE_COMPRESS_H
 
+#if (_MSC_VER >= 1400)
+	#define _CRT_SECURE_NO_WARNINGS
+#endif
+
 #ifdef _WIN32
 	#include "cpp/H5Cpp.h"
 #elif __linux

makefile

@@ -7,7 +7,7 @@ LFLAGS = -Wall -g -fopenmp -Wl,-rpath='$$ORIGIN' -L/usr/lib
 TARGET = traffic_simulation
 
 LIBS = -lgsl -lgslcblas -lhdf5 -lhdf5_cpp
-OBJS = traffic_simulation.o vehicles.o parameters.o hdf_save.o
+OBJS = traffic_simulation.o vehicles.o parameters.o hdf_save.o simulation.o
 
 
 $(TARGET): $(OBJS)
@@ -22,7 +22,10 @@ vehicles.o: parameters.o vehicles.h vehicles.cpp
 hdf_save.o: hdf_save.h hdf_save.cpp
 	$(CPP) -std=c++11 $(CFLAGS) hdf_save.cpp $(LIBS)
 
-traffic_simulation.o: hdf_save.o vehicles.o traffic_simulation.cpp
+simulation.o: simulation.h simulation.cpp
+	$(CPP) -std=c++11 $(CFLAGS) simulation.cpp $(LIBS)
+
+traffic_simulation.o: hdf_save.o vehicles.o simulation.o parameters.o traffic_simulation.cpp
 	$(CPP) -std=c++11 $(CFLAGS) traffic_simulation.cpp $(LIBS)
 
 clean:

parameters.cpp

@@ -6,7 +6,7 @@ int VIRTUAL_LANES = 0; //-v
 int LANES = REAL_LANES + VIRTUAL_LANES;
 int V_MAX = 5;
 int TIMESTEPS = 10; //-t
-int TRIALS = 10; //-T
+int TRIALS = 50; //-T
 double SLOWDOWN = 0.3;
 double LANE_CHANGE_PROB = 0.8;
 double car_ratio = 0;

plot_throughput.py

@@ -6,10 +6,8 @@
 @author: Damian
 """
 import numpy as np
-#import cPickle
 import matplotlib
 matplotlib.use("Agg")
-import matplotlib.lines as lines
 import matplotlib.pyplot as plt
 import glob
 import re
@@ -18,31 +16,36 @@
 import h5py
 
 
-REAL_LANES = 4
+#REAL_LANES = 4
 ROADLENGTH = 50
 TRIALS = 50
 
-AREA = 1*(REAL_LANES)*ROADLENGTH
+#AREA = 1 * (REAL_LANES) * ROADLENGTH
 SPEED = 0
 SIZE = -1
 LAST = -1
 
 
-def throughput(_trial):
-    '''Computes the throughput of one trial of the simulation.'''
-    return np.sum([vehicles[SPEED]/ROADLENGTH for vehicles in _trial[LAST]])
+def throughput(dist_trial):
+    """Computes the throughput of one trial of the simulation."""
+    return np.sum([vehicles[SPEED]/ROADLENGTH for vehicles in dist_trial[LAST]])
+
+
+def velocities(velo_trial):
+    """Computes the average velocity of one trial of the simulation"""
+    return np.average([vehicles[SPEED] for vehicles in velo_trial[LAST]])
 
 
 def load(_ratio, _density):
-    '''Loads data from the hdf5 dataset.'''
+    """Loads data from the hdf5 dataset."""
     vehicledata = np.array([], dtype=np.int8)
     filename = "CarRatio.%.2f_Density.%.2f.h5" % (_ratio, _density)
-    call(['bunzip2', filename+'.bz2'])
+    call(['bunzip2', filename + '.bz2'])
     fid = h5py.File(filename, 'r')
     for n in xrange(TRIALS):
         group = "CarRatio::%.2f/Density::%.2f/" % (_ratio, _density)
-        _trial = "Trial::%04d" % (n+1)
-        dset = fid[group+_trial]
+        _trial = "Trial::%04d" % (n + 1)
+        dset = fid[group + _trial]
         vehicledata = np.append(vehicledata, dset)
     vehicledata = np.reshape(vehicledata, (TRIALS, dset.shape[0],
                                            dset.shape[1],
@@ -51,40 +54,58 @@ def load(_ratio, _density):
     call(["bzip2", "-6", filename])
     return vehicledata
 
+
 if __name__ == "__main__":
+    __plot_ratios__ = [0, 0.25, 0.5, 0.75, 1]
     FILES = glob.glob("*.bz2")
     DIRNAME = os.path.split(os.getcwd())[1]
     DENSITIES = np.arange(0.05, 1, 0.05)
     RATIOS = np.array([])
     for i in FILES:
         car_ratio = float(re.findall(r"CarRatio.(\d.\d+)", i)[0])
-        if car_ratio not in RATIOS:
+        if car_ratio not in RATIOS and car_ratio in __plot_ratios__:
             RATIOS = np.append(RATIOS, car_ratio)
     RATIOS.sort()
-    THROUGHPUT = np.zeros([len(RATIOS), len(DENSITIES)])
-    STDEV = np.zeros([len(RATIOS), len(DENSITIES)])
+    THROUGHPUT = [[None for i in range(len(DENSITIES))]
+                  for j in range(len(RATIOS))]
+    VELOCITIES = [[None for i in range(len(DENSITIES))]
+                  for j in range(len(RATIOS))]
     for x, ratio in enumerate(RATIOS):
         for y, density in enumerate(DENSITIES):
             data = load(ratio, density)
-            data = np.cumsum(data[:, :, :, :], axis=1)
+            datasum = np.cumsum(data[:, :, :, :], axis=1)
             flux = [throughput(trial)
-                    for trial in data]
-            THROUGHPUT[x, y] = np.mean(flux)
-            STDEV[x, y] = np.std(flux)
-
-    __markers__ = list(lines.Line2D.markers.keys())
-    __markers__.sort()
-    N = len(__markers__)
-    __labels__ = RATIOS
+                    for trial in datasum]
+            THROUGHPUT[x][y] = flux
+            VELOCITIES[x][y] = np.average([velocities(trial) for trial in data])
+    MEDIANS = [[np.median(_density) for _density in _car_ratio] for _car_ratio in THROUGHPUT]
     fig = plt.figure(1)
     ax = fig.add_subplot(111)
-    for ydata, label, i in zip(THROUGHPUT, __labels__, range(len(__labels__))):
-        marker = __markers__[i % N]
-        ax.errorbar(DENSITIES, ydata,
-                    yerr=STDEV[i, :], label=label, marker=marker)
-    lgd = ax.legend(loc="upper left", bbox_to_anchor=(1.02, 1))
+    for ydata, median, label, i in zip(THROUGHPUT, MEDIANS, RATIOS, range(len(RATIOS))):
+        color = "%s" % (i*0.15)
+        ax.plot(DENSITIES, median)
+        bp = ax.boxplot(ydata, positions=DENSITIES, widths=0.02)
+        plt.setp(bp['boxes'], color=color)
+        plt.setp(bp['whiskers'], color=color)
+        plt.setp(bp['fliers'], color=color)
+        plt.setp(bp['medians'], color=color)
     ax.set_xlabel('road density')
     ax.set_ylabel('number of exiting vehicles')
     ax.set_title(DIRNAME)
-    fig.savefig('throughput_%s.png'%DIRNAME,
-                bbox_extra_artists=(lgd,),bbox_inches='tight')
+    ax.set_xlim(0, 1)
+    ax.set_xticks(DENSITIES[1::2])
+    fig.savefig('throughput_%s.png' % DIRNAME, bbox_inches='tight')
+    ax.cla()
+
+    fig2 = plt.figure(2)
+    ax2 = fig2.add_subplot(111)
+    for ydata, median, label, i in zip(VELOCITIES, MEDIANS, RATIOS, range(len(RATIOS))):
+        color = "%s" % (i*0.15)
+        ax2.plot(DENSITIES, ydata,color=color)
+    ax2.set_xlabel('road density')
+    ax2.set_ylabel('average velocity')
+    ax2.set_title(DIRNAME)
+    ax2.set_xlim(0, 1)
+    ax2.set_xticks(DENSITIES[1::2])
+    fig2.savefig('velocities_%s.png' % DIRNAME, bbox_inches='tight')
+    ax2.cla()

traffic_simulation.cpp

@@ -1,7 +1,3 @@
-#if (_MSC_VER >= 1400)
-	#define _CRT_SECURE_NO_WARNINGS
-#endif
-
 #include <iostream>
 #include <stdio.h>
 #include <gsl/gsl_rng.h>
@@ -14,116 +10,11 @@
 #include <bzlib.h>
 
 #include "hdf_save.h"
-#include "parameters.h"
-#include "vehicles.h"
+#include "simulation.h"
 
 using namespace std;
 
-int DATAPOINTS = 1000;
-
-struct Simulation{
-public:
-	road_arr road;
-	vector<vehicle> vehicle_array;
-	vector<int> throughput;
-	vector<vector<vector<int> > > vehicle_data;
-	int number_vehicles;
-
-	Simulation(void) {}
-	~Simulation(void) {}
-
-	void evolve(float density, float car_ratio){
-		int passed_vehicles;
-		initialize(density, car_ratio);
-		vector<int> permutation(vehicle_array.size());
-		throughput.resize(TIMESTEPS);
-		for (unsigned i = 0; i < permutation.size(); i++) permutation[i] = i;
-		random_shuffle(permutation.begin(), permutation.end());
-		for (int t = 0; t < TIMESTEPS; t++){
-			vector<vector<int> > vehicle_stats;
-			passed_vehicles = 0;
-			for (int i : permutation){
-				vehicle_array[i].accelerate();
-				if (LANE_CHANGE){
-					vehicle_array[i].change_lane(road);
-					vehicle_array[i].decelerate(road);
-					if (!vehicle_array[i].changed_lane) vehicle_array[i].random_slow();
-				}
-				else{
-					vehicle_array[i].decelerate(road);
-					vehicle_array[i].random_slow();
-				}
-				vehicle_array[i].move(road);
-				if (vehicle_array[i].exit_road == true) passed_vehicles += 1;
-			}
-			/* Eliminate the transient 2000 steps */
-			if (t >= TIMESTEPS-DATAPOINTS){
-				for (vehicle vehicle : vehicle_array) vehicle_stats.push_back(vehicle.stats());
-				vehicle_data.push_back(vehicle_stats);
-				vector<vector<int> >().swap(vehicle_stats);
-			}
-			random_shuffle(permutation.begin(), permutation.end());
-		}
-	}
-
-private:
-	void initialize(float density, float car_ratio){
-		int pos, lane, counter;
-		float motor_ratio = 1 - car_ratio;
-		number_vehicles = density*ROADLENGTH*(REAL_LANES) /
-			(car().size*car_ratio + motorcycle().size*motor_ratio);
-		int number_car = car_ratio*number_vehicles;
-		int number_motorcycle;
-		for (int i=0; i < LANES; i++) this->road.push_back(vector<int>(ROADLENGTH,0));
-		vector<vehicle> car_array;
-		/* Initializes Cars */
-		if (car_ratio>0){
-			vector<int> lane_choice(LANES / car().width);
-			for (int i = 0; i < lane_choice.size(); i++) lane_choice[i] = i * 2;
-			for (counter = 0; counter < number_car; counter++){
-				int iterations = 0;
-				pos = gsl_rng_uniform_int(generator, ROADLENGTH / 2) * 2 + 1;
-				lane = lane_choice[gsl_rng_uniform_int(generator, 2)];
-				while (!place_check(pos, lane, car().length, car().width, road, ROADLENGTH)){
-					pos = gsl_rng_uniform_int(generator, ROADLENGTH / 2) * 2 + 1;
-					lane = lane_choice[gsl_rng_uniform_int(generator, lane_choice.size())];
-					if (iterations > 500) break;
-					iterations += 1;
-				}
-				car_array.push_back(car());
-				car_array[counter].pos = pos;
-				car_array[counter].lane = lane;
-				car_array[counter].vel = gsl_rng_uniform_int(generator, V_MAX + 1);
-				car_array[counter].place(road);
-			}
-		}
-		/* Initializes Motorcycles */
-		if (motor_ratio > 0){
-			number_motorcycle = number_vehicles - car_array.size();
-			vector<vehicle> moto_array(number_motorcycle, motorcycle());
-			for (counter = 0; counter < number_motorcycle; counter++){
-				pos = 0;
-				lane = 0;
-				while (!place_check(pos, lane, motorcycle().length, motorcycle().width, road, ROADLENGTH)){
-					pos = gsl_rng_uniform_int(generator, ROADLENGTH);
-					lane = gsl_rng_uniform_int(generator, REAL_LANES);
-				}
-				moto_array[counter].pos = pos;
-				moto_array[counter].lane = lane;
-				moto_array[counter].vel = gsl_rng_uniform_int(generator, V_MAX + 1);
-				moto_array[counter].place(road);
-			}
-			car_array.insert(car_array.end(), moto_array.begin(), moto_array.end());
-		}
-		else vector<vehicle> moto_array(0);
-		/* Cleanup Operations */
-
-		vehicle_array.swap(car_array);
-		//vector<vehicle>().swap(car_array);
-		//vector<vehicle>().swap(moto_array);
-	}
 
-};
 
 void BZIP(char* _filename){
 	char message[100];
@@ -286,10 +177,9 @@ int main(int argc, char* argv[]){
 		densities[i] = first;
 		first += 0.05;
 	}
-	int j = 4;
 		omp_set_num_threads(2);
 		#pragma omp parallel for
-		for (int i = 0; i < densities.size(); i++){
+		for (unsigned i = 0; i < densities.size(); i++){
 			if (LOAD_SEED == false) seed = time(NULL) * 123456789;
 			runmsg[i] = start(densities[i], car_ratio, seed);
 			printstat(runmsg, densities, car_ratio);