Advertisement
brandblox

numpy lab(22/04/2024)(simplex)

Apr 21st, 2024 (edited)
627
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
Python 4.71 KB | Source Code | 0 0
  1. # -*- coding: utf-8 -*-
  2. """
  3. Created on Mon Apr 22 09:45:34 2024
  4.  
  5. @author: lab
  6. """
  7. import numpy as np
  8. from fractions import Fraction # so that numbers are not displayed in decimal.
  9.  
  10. print("\n****Simplex Algorithm ****\n\n")
  11.  
  12. # ‘A’ will contain the coefficients of the constraints
  13. A = np.array([[1, 1, 0, 1], [2, 1, 1, 0]])
  14. # b will contain the amount of resources
  15. b = np.array([8, 10])
  16. # c will contain coefficients of objective function Z
  17. c = np.array([1, 1, 0, 0])
  18. # ‘B’ will contain the basic variables that make identity matrix
  19. cb = np.array(c[3])
  20. B = np.array([[3], [2]])
  21. # cb contains their corresponding coefficients in Z
  22. cb = np.vstack((cb, c[2]))
  23. xb = np.transpose([b])
  24.  
  25. # combine matrices B and cb
  26. table = np.hstack((B, cb))
  27. table = np.hstack((table, xb))
  28. # combine matrices B, cb and xb
  29. # finally combine matrix A to form the complete simplex table
  30. table = np.hstack((table, A))
  31. # change the type of table to float
  32. table = np.array(table, dtype ='float')
  33. MIN = 0
  34.  
  35. print("Table at itr = 0")
  36. print("B \tCB \tXB \ty1 \ty2 \ty3 \ty4")
  37. for row in table:
  38.     for el in row:
  39.         print(Fraction(str(el)).limit_denominator(100), end ='\t')
  40.     print()
  41. print()
  42.  
  43. print("Simplex Working....")
  44. # when optimality reached it will be made 1
  45. reached = 0
  46. itr = 1
  47. unbounded = 0
  48. alternate = 0
  49. while reached == 0:
  50.     print("Iteration: ", end =' ')
  51.     print(itr)
  52.     print("B \tCB \tXB \ty1 \ty2 \ty3 \ty4")
  53.     for row in table:
  54.         for el in row:
  55.             print(Fraction(str(el)).limit_denominator(100), end ='\t')
  56.         print()
  57.    
  58.     # calculate Relative profits-> cj - zj for non-basics
  59.     i = 0
  60.     rel_prof = []
  61.     while i<len(A[0]):
  62.         rel_prof.append(c[i] - np.sum(table[:, 1]*table[:, 3 + i]))
  63.         i = i + 1
  64.     print("rel profit: ", end =" ")
  65.     for profit in rel_prof:
  66.         print(Fraction(str(profit)).limit_denominator(100), end =", ")
  67.     print()
  68.    
  69.     i = 0
  70.     b_var = table[:, 0]
  71.     # checking for alternate solution
  72.     while i<len(A[0]):
  73.         j = 0
  74.         present = 0
  75.         while j<len(b_var):
  76.             if int(b_var[j]) == i:
  77.                 present = 1
  78.                 break
  79.             j += 1
  80.         if present == 0:
  81.             if rel_prof[i] == 0:
  82.                 alternate = 1
  83.                 print("Case of Alternate found")
  84.         i += 1
  85.     print()
  86.    
  87.     flag = 0
  88.     for profit in rel_prof:
  89.         if profit>0:
  90.             flag = 1
  91.             break
  92.     # if all relative profits <= 0
  93.     if flag == 0:
  94.         print("All profits are <= 0, optimality reached")
  95.         reached = 1
  96.         break
  97.    
  98.     # kth var will enter the basis
  99.     k = rel_prof.index(max(rel_prof))
  100.     min_val = 99999
  101.     i = 0
  102.     r = -1
  103.     # min ratio test (only positive values)
  104.     while i<len(table):
  105.         if (table[:, 2][i]>0 and table[:, 3 + k][i]>0):
  106.             val = table[:, 2][i] / table[:, 3 + k][i]
  107.             if val < min_val:
  108.                 min_val = val
  109.                 r = i # leaving variable
  110.         i += 1
  111.     # if no min ratio test was performed
  112.     if r == -1:
  113.         unbounded = 1
  114.         print("Case of Unbounded")
  115.         break
  116.    
  117.     print("pivot element index:", end =' ')
  118.     print(np.array([r, 3 + k]))
  119.     pivot = table[r][3 + k]
  120.     print("pivot element: ", end =" ")
  121.     print(Fraction(pivot).limit_denominator(100))
  122.    
  123.     # perform row operations
  124.     # divide the pivot row with the pivot element
  125.     table[r, 2:len(table[0])] = table[r, 2:len(table[0])] / pivot
  126.     # do row operation on other rows
  127.     i = 0
  128.     while i<len(table):
  129.         if i != r:
  130.             table[i, 2:len(table[0])] = table[i, 2:len(table[0])] - table[i][3 + k] * table[r, 2:len(table[0])]
  131.         i += 1
  132.     # assign the new basic variable
  133.     table[r][0] = k
  134.     table[r][1] = c[k]
  135.     print()
  136.     print()
  137.     itr += 1
  138.     print()
  139.     print("***************************************************************")
  140.  
  141. if unbounded == 1:
  142.     print("UNBOUNDED LPP")
  143.     exit()
  144. if alternate == 1:
  145.     print("ALTERNATE Solution")
  146.  
  147. print("optimal table:")
  148. print("B \tCB \tXB \ty1 \ty2 \ty3 \ty4")
  149. for row in table:
  150.     for el in row:
  151.         print(Fraction(str(el)).limit_denominator(100), end ='\t')
  152.     print()
  153.  
  154. print()
  155. print("value of Z at optimality: ", end =" ")
  156. basis = []
  157. i = 0
  158. sum_val = 0
  159. while i < len(table):
  160.     sum_val += c[int(table[i][0])] * table[i][2]
  161.     temp = "x" + str(int(table[i][0]) + 1)
  162.     basis.append(temp)
  163.     i += 1
  164. # if MIN problem make z negative
  165. if MIN == 1:
  166.     print(-Fraction(str(sum_val)).limit_denominator(100))
  167. else:
  168.     print(Fraction(str(sum_val)).limit_denominator(100))
  169. print("Final Basis: ", end =" ")
  170. print(basis)
  171. print("Simplex Finished...")
  172. print()
  173.  
  174.  
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement