From b8b7acdb2de02fe13ed4559b73327297da44b560 Mon Sep 17 00:00:00 2001 From: Andy Eschbacher Date: Wed, 29 Mar 2017 18:55:49 -0400 Subject: [PATCH] refactoring to limit the number of attributes --- .../crankshaft/optimization/optim.py | 137 +++++++++--------- 1 file changed, 66 insertions(+), 71 deletions(-) diff --git a/src/py/crankshaft/crankshaft/optimization/optim.py b/src/py/crankshaft/crankshaft/optimization/optim.py index 01124a9..3b3ec35 100644 --- a/src/py/crankshaft/crankshaft/optimization/optim.py +++ b/src/py/crankshaft/crankshaft/optimization/optim.py @@ -6,7 +6,13 @@ from cvxopt.glpk import ilp from crankshaft.analysis_data_provider import AnalysisDataProvider class Optim(object): - """Linear optimization class for logistics cost minimization""" + """Linear optimization class for logistics cost minimization + Optimization for logistics + based on models: + - waste_per_person * (1 - recycle_rate) * population + - waste_in_area * (marginal_cost + transport_cost * distance) + That is, `cost ~ population * distance` + """ def __init__(self, drain_table, source_table, capacity_column, production_column, marginal_column, **kwargs): @@ -14,24 +20,24 @@ class Optim(object): # set data provider (defaults to SQL database access self.data_provider = kwargs.get('data_provider', AnalysisDataProvider()) - - # optional params - self.waste_per_person = kwargs.get('waste_per_person', 0.01) - self.recycle_rate = kwargs.get('recycle_rate', 0.0) - self.dist_cost = kwargs.get('dist_cost', 0.15) - - # data sources - self.drain_table = drain_table - self.source_table = source_table - + # model parameters + self.model_params = { + 'waste_per_person': kwargs.get('waste_per_person', 0.01), + 'dist_cost': kwargs.get('dist_cost', 0.15), + 'recycle_rate': kwargs.get('recycle_rate', 0.0) + } # model data - self.plant_capacity = self.data_provider.get_column(drain_table, - capacity_column) - self.waste_in_area = (0.01 * (1. - self.recycle_rate) * + self.model_data = { + 'plant_capacity': self.data_provider.get_column(drain_table, + capacity_column), + 'waste_in_area': (self.model_params['waste_per_person'] * + (1. - self.model_params['recycle_rate']) * self.data_provider.get_column(source_table, - production_column)) - self.marginal_cost = self.data_provider.get_column(drain_table, + production_column)), + 'marginal_cost': self.data_provider.get_column(drain_table, marginal_column) + } + # database ids self.drain_ids = self.data_provider.get_column(drain_table, 'cartodb_id', @@ -40,44 +46,38 @@ class Optim(object): 'cartodb_id', dtype=int) # derivative data - self.distances = self.data_provider.get_pairwise_distances(source_table, - drain_table) - self.n_areas = len(self.waste_in_area) - self.n_plants = len(self.distances) - self.cost = self.calc_cost() + self.n_sources = len(self.source_ids) + self.n_drains = len(self.drain_ids) + self.cost = self.calc_cost(source_table, + drain_table) + + def _check_constraints(self): + """Check if inputs are within constraints""" + if self.model_data['waste_in_area'].sum() > self.model_data['plant_capacity'].sum(): + plpy.error("Solution not possible. Drain capacity is smaller " + "than total source production.") def output(self): """...""" # n_drains x n_sources matrix (row, column) assignments = self.optim() - # - plpy.notice("self.cost.shape: {}".format(str(self.cost.shape))) - plpy.notice("self.cost: {}".format(self.cost)) - # - plpy.notice(assignments) - plpy.notice("assignments.shape: {}".format(str(assignments.shape))) - + # crosswalks for matrix index -> cartodb_id drain_id_crosswalk = {} for idx, cid in enumerate(self.drain_ids): # matrix index -> cartodb_id drain_id_crosswalk[idx] = cid - # plpy.notice(drain_id_crosswalk) - + source_id_crosswalk = {} for idx, cid in enumerate(self.source_ids): # matrix index -> cartodb_id source_id_crosswalk[idx] = cid - # plpy.notice(source_id_crosswalk) - + # find non-zero entries nonzeros = np.nonzero(assignments) - plpy.notice("nonzeros: {}".format(str(nonzeros))) source_index, drain_index = nonzeros[0], nonzeros[1] - # - plpy.notice(source_index) - plpy.notice(drain_index) + # assigned_costs = [(drain_id_crosswalk[drain_index[i]], source_id_crosswalk[source_index[i]], self.cost[drain_index[i], @@ -85,42 +85,38 @@ class Optim(object): for i in range(len(source_index))] return assigned_costs - def test(self): - """ - just plpy.notice the stored information - """ - - plpy.notice(self.source_table) - plpy.notice(self.drain_table) - plpy.notice(self.distances) - plpy.notice(self.plant_capacity) - plpy.notice(self.waste_in_area) - return None - def cost_func(self, distance, waste, marginal): """ cost equation - """ - return waste * (marginal + self.dist_cost * distance) - def calc_cost(self): + :param distance: distance (in km) + :type distance: float + :param waste: number of tons of waste. This was previously calculated + as self.model_params['waste_per_person'] * number of people minus the recycle_rate + :type waste: numeric + :param marginal: intrinsic cost per ton of a plant + :type marginal: numeric + :returns: cost + :rtype: numeric + + Note: dist_cost is the cost per ton (0.15 GBP/ton) + """ + return waste * (marginal + self.model_params['dist_cost'] * distance) + + def calc_cost(self, source_table, drain_table): """ Populate an d x s matrix according to the cost equation :returns: d x s matrix of costs from area i to plant j :rtype: NumPy matrix """ - plpy.notice('self.waste_in_area: {}'.format(str(self.waste_in_area.shape))) - plpy.notice(self.waste_in_area) - plpy.notice('self.marginal_cost: {}'.format(str(self.marginal_cost.shape))) - plpy.notice(self.marginal_cost) - plpy.notice('self.distances: {}'.format(str(self.distances.shape))) - plpy.notice(self.distances) + distances = self.data_provider.get_pairwise_distances(source_table, + drain_table) costs = np.array([self.cost_func(distance, - self.waste_in_area[pair[1]], - self.marginal_cost[pair[0]]) - for pair, distance in np.ndenumerate(self.distances)]) - return costs.reshape(self.distances.shape) + self.model_data['waste_in_area'][pair[1]], + self.model_data['marginal_cost'][pair[0]]) + for pair, distance in np.ndenumerate(distances)]) + return costs.reshape(distances.shape) def optim(self): """solve linear optimization problem @@ -139,22 +135,21 @@ class Optim(object): # equality constraint variables # each area is serviced once A = cvxopt.spmatrix(1., - [i // self.n_plants - for i in range(self.n_plants * self.n_areas)], - range(self.n_plants * self.n_areas)) - b = cvxopt.matrix(np.ones((self.n_areas, 1)), tc='d') + [i // self.n_drains + for i in range(self.n_drains * self.n_sources)], + range(self.n_drains * self.n_sources)) + b = cvxopt.matrix(np.ones((self.n_sources, 1)), tc='d') # inequality constraint variables # each plant never goes over capacity - h = cvxopt.matrix(self.plant_capacity, tc='d') - G = cvxopt.spmatrix(np.repeat(self.waste_in_area, self.n_plants), - [i % self.n_plants - for i in range(self.n_plants * self.n_areas)], - range(self.n_plants * self.n_areas)) + h = cvxopt.matrix(self.model_data['plant_capacity'], tc='d') + G = cvxopt.spmatrix(np.repeat(self.model_data['waste_in_area'], self.n_drains), + [i % self.n_drains + for i in range(self.n_drains * self.n_sources)], + range(self.n_drains * self.n_sources)) binary_entries = set(range(len(c))) # solve (sol, x) = ilp(c=c, G=G, h=h, A=A, b=b, B=binary_entries) - # assignment = np.array(x).reshape((self.n_areas, self.n_plants)) if sol != 'optimal': raise Exception("No solution possible: {}".format(sol)) x_shape = (self.cost.shape[1], self.cost.shape[0])