############   capacitated P-median.py #######################
## Goal:    Computes capacitated p-median on CMS dataset
## Inputs:  demand nodes & actual demand
##          schools & capacities
## Outputs: CPLEX file
###########################################################

import string, os,  math, time, sys
import arcpy
from arcpy import env
arcpy.OverWriteOutput = True

def distance(x1, y1, x2, y2):
    dx = x2 - x1
    dy = y2 - y1
    d = math.sqrt(dx**2 + dy**2)
    return d

#########################################################
#### STEP 1 ######### INPUTS  ###########################
#########################################################

demand = "demand.shp"
schools = "schools.shp"
p = 16
maxCapacity = .75 #85% of its capacity
outputFile = open('capacitatedPmedian.txt', 'w')

# reading demand information
originFID=[]
popInfo=[]
originX=[]
originY=[]
for ori in arcpy.SearchCursor(demand):
    originFID.append(ori.getValue("FID"))
    popInfo.append(ori.getValue("Weight"))

for row in arcpy.da.SearchCursor(demand, ["SHAPE@XY"]):
    originX.append(row[0][0])
    originY.append(row[0][1])

print("...origin information read....")

# reading destination information
destinationFID = []
destinationX=[]
destinationY=[]
destinationCapacity=[]
for dest in arcpy.SearchCursor(schools):
    destinationFID.append(dest.getValue("FID"))
    destinationCapacity.append(dest.getValue("Capacity"))

for row in arcpy.da.SearchCursor(schools, ["SHAPE@XY"]):
    destinationX.append(row[0][0])
    destinationY.append(row[0][1])

print("...destination information read....")


OD = []
i=0
while i<len(originX):
    j=0
    while j<len(destinationX):
        OD.append(str(distance(originX[i], originY[i], destinationX[j], destinationY[j])))
        j=j+1
    i=i+1

arcpy.AddMessage("euclidean OD matrix built......")


#  WRITING L-P file
#  Model headings

# Inputs for objective function
# CPLEX Heading Requirements
# Objective function headings
print("writing model output...")
outputFile.write("Minimize\n")

i=0
while i< len(originX):
    j=0
    while j < len(destinationX):
        b = float(OD[(i*len(destinationX))+j])
        outputFile.write("+")
        outputFile.write(str(b))
        outputFile.write("X")
        outputFile.write(str(int(originFID[i])) + "_"+ str(int(destinationFID[j])))
        j=j+1
    i=i+1

#  CONSTRAINTS  #

# demand node assigned if school is open.
outputFile.write("\nSubject to")
outputFile.write("\n")
i=0
while i< len(originFID):
    j=0
    while j < len(destinationFID):
        outputFile.write("+")
        outputFile.write("X" + str(int(originFID[i])) + "_"+ str(int(destinationFID[j])))
        outputFile.write("-Y" + str(int(destinationFID[j])))
        outputFile.write("<=0\n")
        j=j+1
    i=i+1
print("constraint a demand node is assigned only if a facility is open...done")


# Constraint each demand node is assigned -
outputFile.write("\n")
i=0
while i< len(originFID):
    j=0
    while j < len(destinationFID):
        outputFile.write("+")
        outputFile.write("X" + str(int(originFID[i])) + "_"+ str(int(destinationFID[j])))
        j=j+1
    outputFile.write("=1\n")
    i=i+1
print("constraint that each demand node is assigned ...done")

###  MAXIMAL CAPACITY !Maximal Capacity Constraint!
outputFile.write("\n")
##
CapacityField = destinationCapacity
j=0
while j<len(destinationFID):
    outputFile.write(str(maxCapacity * float(CapacityField[j]))+ " Y"+str(destinationFID[j]))
    i=0
    while i<len(originFID):
        a = popInfo[i]
        if a>0:
            outputFile.write(" -")
            outputFile.write(str(int(a)))
            outputFile.write(" X")
            outputFile.write(str(int(originFID[i]))+"_"+str(destinationFID[j]))
        i=i+1
    outputFile.write(">=0\n")
    j=j+1
outputFile.write("\n\n")
print("capacity constraint....done")

###  p-facilities

j=0
while j < len(destinationFID):
    outputFile.write("+")
    outputFile.write(" Y")
    outputFile.write(str(int(destinationFID[j])))
    j=j+1
outputFile.write("<=" +str(p))
print("keep track of how many schools")
outputFile.write("\n\n")



### INTEGER CONSTRAINTS
outputFile.write("\n") ##REMOVE LATER
outputFile.write("Binary\n")
i=0
while i<len(originFID):
    j=0
    while j<len(destinationFID):
        outputFile.write("X")
        outputFile.write(str(int(originFID[i]))+"_"+str(int(destinationFID[j]))+ "\n")
        j=j+1
    i=i+1
outputFile.write("\n")
print("X integers all set")

j=0
while j<len(destinationFID):
    outputFile.write("Y")
    outputFile.write(str(int(destinationFID[j]))+"\n")
    j=j+1
print("Y integers all set")


outputFile.write("END\n")
outputFile.close()

ModuleNotFoundError: No module named 'arcpy'