A simple randomized python model on the Free Water Scheme announced by the AAP in Delhi. Back calculated the mean usage as 16 units and Standard deviation as 4 units to reach a cost impact of 166 Cr between the schemes.

free_water.py

import numpy
import matplotlib.pylab as plt

__author__ = 'sujeet'

# Define Constants

METERED_CONNECTIONS = 900000
TOTAL_CONNECTIONS = 1800000
TOTAL_HOUSEHOLDS = 5000000

CURRENT_PRICING = [
	{'start':0,'end':10,'rate':2.42, 'service_charge':60.50},
	{'start':10,'end':20,'rate':3.63, 'service_charge':121.00},
	{'start':20,'end':30,'rate':18.15, 'service_charge':181.50},
	{'start':30,'end':999999,'rate':30.25, 'service_charge':242.00},
]
SEWERAGE_CESS = 0.6

MEAN_USAGE_UNITS = 16
STDDEV_USAGE_UNITS = 4


# Create a normalized Usage sample
usage_sample = numpy.random.normal(MEAN_USAGE_UNITS,STDDEV_USAGE_UNITS,METERED_CONNECTIONS)

# Function to calculate cost of water currently
def calculate_price(x):
	volumetric_charge = sum(map(lambda p: (min(x,p['end'])-p['start'])*p['rate'] if x>p['start'] else 0,
	                            CURRENT_PRICING))
	fixed_charge = sum(map(lambda p: p['service_charge'] if x>p['start'] and x<=p['end'] else 0, CURRENT_PRICING))
	return volumetric_charge * (1+SEWERAGE_CESS) + fixed_charge

# Function to calculate new cost of water
def calculate_new_price(x):
	if x>20:
		return calculate_price(x)
	else:
		return 0

# Vectorize functions
v_calculate_price = numpy.vectorize(calculate_price)
v_calculate_new_price = numpy.vectorize(calculate_new_price)

# Apply on usage
current_price = v_calculate_price(usage_sample)
new_price = v_calculate_new_price(usage_sample)

# Display differences
print 'Current: '+str(numpy.sum(current_price)*12/10000000)
print 'New: '+str(numpy.sum(new_price)*12/10000000)
print 'Difference: '+str((numpy.sum(new_price)-numpy.sum(current_price))*12/10000000)

# Plot distributions
plt.subplot(311)
plt.hist(usage_sample,50)
plt.title('Mean Usage Units = '+str(MEAN_USAGE_UNITS)+' Std Deviation = '+str(STDDEV_USAGE_UNITS))
plt.subplot(312)
plt.hist(current_price,30)
plt.title('Total Cost = '+str(numpy.sum(current_price)*12/10000000)+' Crores')
plt.subplot(313)
plt.hist(new_price,30)
plt.title('Total Cost = '+str(numpy.sum(new_price)*12/10000000)+' Crores')
plt.subplots_adjust(bottom = 0.1)