transportation_problem.py

# -*- coding: utf-8 -*-
"""transportation_problem.ipynb

Automatically generated by Colaboratory.

Original file is located at
    https://colab.research.google.com/drive/12T0rOZw_O_gmA0LM83Da-KfPL48982dV
"""

import numpy as np

def find_north_west_corner(temp_source, temp_dest):
    for i, s in enumerate(temp_source):
        if s != 0:
            for j, d in enumerate(temp_dest):
                if d != 0:
                    return i, j
    return False

def north_west_corner(sources, dest, costs):
    temp_source = sources.copy()
    temp_dest = dest.copy()
    n_source, n_dests = costs.shape
    assigned = np.zeros(costs.shape)

    while find_north_west_corner(temp_source, temp_dest):
        x, y = find_north_west_corner(temp_source, temp_dest)
        min_cost = min(temp_source[x], temp_dest[y])
        temp_source[x] -= min_cost
        temp_dest[y] -= min_cost
        assigned[x][y] = min_cost
    return assigned

def find_uv_values(sources, dest, costs, assigned):
    U = [None for _ in range(len(sources))]
    V = [None for _ in range(len(dest))]
    U[0] = 0
    filled = False

    while not filled:
        for i, row in enumerate(assigned):
            for j, row_el in enumerate(row):
                if U[i] is not None:
                    if row_el > 0:
                        if V[j] is None:
                            V[j] = costs[i][j] - U[i]

        assigned_T = np.transpose(np.array(assigned))
        for i, col in enumerate(assigned_T):
            for j, col_el in enumerate(col):
                if V[i] is not None:
                    if col_el > 0:
                        if U[j] is None:
                            U[j] = costs[i][j] - V[i]

        filled_u = all(u is not None for u in U)
        filled_v = all(v is not None for v in V)
        filled = filled_u and filled_v
    return U, V

def get_penalties(U, V, cost, assigned):
    penalty = np.zeros(cost.shape)
    for i, u in enumerate(U):
        for j, v in enumerate(V):
            if assigned[i][j] == 0:
                penalty[i][j] = u + v - cost[i][j]
    return penalty

def check_optimality(penalty):
    return not any(el > 0 for row in penalty for el in row)

def get_entering_variable_position(penalty):
    max_value = np.max(penalty)
    position = np.where(penalty == max_value)
    return int(position[0]), int(position[1])

def get_next_cells(matrix, position):
    x, y = position
    next_cells = [(x - 1, y), (x + 1, y), (x, y - 1), (x, y + 1)]
    next_cells = [
        cell
        for cell in next_cells
        if 0 <= cell[0] < matrix.shape[0]
        and 0 <= cell[1] < matrix.shape[1]
        and matrix[cell[0], cell[1]] != 0
    ]
    return next_cells

def find_closed_loop(assigned, start, visited, path, check=0):
    i, j = start
    visited[i, j] = 1

    next_cells = get_next_cells(assigned, start)
    for cell in next_cells:
        if visited[cell[0], cell[1]] == 0:
            if cell in path:
                return path
            else:
                visited[cell[0], cell[1]] = 1
                result = find_closed_loop(assigned, cell, visited, path + [cell], check + 1)
                if result is not None:
                    return result
        elif check >= 2:
            return path

def get_closed_loop(assigned, start):
    visited = np.zeros(assigned.shape)
    path = [start]
    final_path = find_closed_loop(assigned, start, visited, path)
    return final_path

def pivoting(cost, assigned, closed_loop):
    pivot_cell = closed_loop[1]
    assigned_copy = np.array(assigned)
    pivot_cost = assigned[pivot_cell[0], pivot_cell[1]]
    for i in range(len(closed_loop)):
        cell = closed_loop[i]
        if i % 2 == 0:
            assigned_copy[cell[0], cell[1]] += pivot_cost
        else:
            assigned_copy[cell[0], cell[1]] -= pivot_cost
    return assigned_copy

def display_results(cost_matrix, allocations, u_values, v_values):
    num_sources, num_destinations = cost_matrix.shape
    print("Cost Matrix:")
    for i in range(num_sources):
        for j in range(num_destinations):
            print("{:<5}".format(cost_matrix[i, j]), end=" ")
        print()

    print("\nAssigned Values:")
    for i in range(num_sources):
        for j in range(num_destinations):
            print("{:<5}".format(allocations[i, j]), end=" ")
        print()

    print("\nU Values:")
    for i in range(num_sources):
        print("U[{}]: {}".format(i, u_values[i]))

    print("\nV Values:")
    for j in range(num_destinations):
        print("V[{}]: {}".format(j, v_values[j]))

def Modi_method(sources, dest, costs):
    assigned = north_west_corner(sources, dest, costs)
    print("Initial Assigned Values:\n", assigned)
    print("Initial Cost:", get_cost(costs, assigned))

    U, V = find_uv_values(sources, dest, costs, assigned)
    print("Initial U Values:", U)
    print("Initial V Values:", V)

    penalty = get_penalties(U, V, costs, assigned)
    print("Initial Penalty Matrix:\n", penalty)

    while not check_optimality(penalty):
        start = get_entering_variable_position(penalty)
        loop = get_closed_loop(assigned, start)
        assigned = pivoting(costs, assigned, loop)
        U, V = find_uv_values(sources, dest, costs, assigned)
        penalty = get_penalties(U, V, costs, assigned)

        display_results(costs, assigned, U, V)
        print("Penalty Matrix:\n", penalty)

    print("\nFinal Assigned Values:\n", assigned)
    return assigned

def main():
    num_sources = int(input("Enter the number of sources: "))
    num_destinations = int(input("Enter the number of destinations: "))

    sources = np.array(list(map(int, input("Enter the supply for each source: ").split())))
    dest = np.array(list(map(int, input("Enter the demand for each destination: ").split())))

    costs = []
    print("Enter the transportation costs:")
    for _ in range(num_sources):
        row = list(map(int, input().split()))
        costs.append(row)
    costs = np.array(costs)

    assigned = Modi_method(sources, dest, costs)
    print("\nFinal Cost:", get_cost(costs, assigned))

if __name__ == "__main__":
    main()