Add 3.2 challenge

1 files changed, 302 insertions, 0 deletions

diff --git a/3.2.py b/3.2.py
new file mode 100644
index 0000000..8006d21
--- /dev/null
+++ b/3.2.py
@@ -0,0 +1,302 @@
+#!/usr/bin/python2
+
+# Oppsieee https://github.com/mkutny/absorbing-markov-chains
+
+# Doomsday Fuel
+# =============
+
+# Making fuel for the LAMBCHOP's reactor core is a tricky process because of
+# the exotic matter involved. It starts as raw ore, then during processing,
+# begins randomly changing between forms, eventually reaching a stable form.
+# There may be multiple stable forms that a sample could ultimately reach, not
+# all of which are useful as fuel.
+
+# Commander Lambda has tasked you to help the scientists increase fuel creation
+# efficiency by predicting the end state of a given ore sample. You have
+# carefully studied the different structures that the ore can take and which
+# transitions it undergoes. It appears that, while random, the probability of
+# each structure transforming is fixed. That is, each time the ore is in 1
+# state, it has the same probabilities of entering the next state (which might
+# be the same state). You have recorded the observed transitions in a matrix.
+# The others in the lab have hypothesized more exotic forms that the ore can
+# become, but you haven't seen all of them.
+
+# Write a function solution(m) that takes an array of array of nonnegative ints
+# representing how many times that state has gone to the next state and return
+# an array of ints for each terminal state giving the exact probabilities of
+# each terminal state, represented as the numerator for each state, then the
+# denominator for all of them at the end and in simplest form. The matrix is at
+# most 10 by 10. It is guaranteed that no matter which state the ore is in,
+# there is a path from that state to a terminal state. That is, the processing
+# will always eventually end in a stable state. The ore starts in state 0. The
+# denominator will fit within a signed 32-bit integer during the calculation,
+# as long as the fraction is simplified regularly.
+
+# For example, consider the matrix m:
+# [
+#   [0,1,0,0,0,1],  # s0, the initial state, goes to s1 and s5 with equal probability
+#   [4,0,0,3,2,0],  # s1 can become s0, s3, or s4, but with different probabilities
+#   [0,0,0,0,0,0],  # s2 is terminal, and unreachable (never observed in practice)
+#   [0,0,0,0,0,0],  # s3 is terminal
+#   [0,0,0,0,0,0],  # s4 is terminal
+#   [0,0,0,0,0,0],  # s5 is terminal
+# ]
+# So, we can consider different paths to terminal states, such as:
+# s0 -> s1 -> s3
+# s0 -> s1 -> s0 -> s1 -> s0 -> s1 -> s4
+# s0 -> s1 -> s0 -> s5
+# Tracing the probabilities of each, we find that
+# s2 has probability 0
+# s3 has probability 3/14
+# s4 has probability 1/7
+# s5 has probability 9/14
+# So, putting that together, and making a common denominator, gives an answer in the form of
+# [s2.numerator, s3.numerator, s4.numerator, s5.numerator, denominator] which is
+# [0, 3, 2, 9, 14].
+
+
+from fractions import Fraction
+from fractions import gcd
+
+def num_of_transients(m):
+    if len(m) == 0:
+        raise Exception("Can't get transient states of empty matrix")
+
+    for r in range(len(m)):
+        for c in range(len(m[r])):
+            if m[r][c] != 0:
+                break
+        else:
+            return r
+    raise Exception("Not a valid AMC matrix: no absorbing (terminal) states")
+
+def decompose(m):
+    t = num_of_transients(m)
+    if t == 0:
+        raise Exception("No transient states. At least initial state is needed.")
+
+    Q = []
+    for r in range(t):
+        qRow = []
+        for c in range(t):
+            qRow.append(m[r][c])
+        Q.append(qRow)
+    if Q == []:
+        raise Exception("Not a valid AMC matrix: no transient states")
+
+    R = []
+    for r in range(t):
+        rRow = []
+        for c in range(t, len(m[r])):
+            rRow.append(m[r][c])
+        R.append(rRow)
+    if R == []:
+        raise Exception("Not a valid AMC matrix: missing absorbing states")
+    return Q, R
+
+def identity(t):
+    m = []
+    for i in range(t):
+        r = []
+        for j in range(t):
+            r.append(int(i == j))
+        m.append(r)
+    return m
+
+def isZero(m):
+    for r in range(len(m)):
+        for c in range(len(m[r])):
+            if m[r][c] != 0:
+                return False
+    return True
+
+def swap(m, i, j):
+    n = []
+    s = len(m)
+
+    if s != len(m[0]):
+        raise Exception("Cannot swap non-square matrix")
+
+    if i == j:
+        return m
+
+    for r in range(s):
+        nRow = []
+        tmpRow = m[r]
+        if r == i:
+            tmpRow = m[j]
+        if r == j:
+            tmpRow = m[i]
+        for c in range(s):
+            tmpEl = tmpRow[c]
+            if c == i:
+                tmpEl = tmpRow[j]
+            if c == j:
+                tmpEl = tmpRow[i]
+            nRow.append(tmpEl)
+        n.append(nRow)
+    return n
+
+def sort(m):
+    size = len(m)
+
+    zero_row = -1
+    for r in range(size):
+        sum = 0
+        for c in range(size):
+            sum += m[r][c]
+        if sum == 0:
+            zero_row = r
+        if sum != 0 and zero_row > -1:
+            n = swap(m, r, zero_row)
+            return sort(n)
+    return m
+
+def normalize(m, use_fractions=False ):
+    n = []
+    for r in range(len(m)):
+        sum = 0
+        cols = len(m[r])
+        for c in range(cols):
+            sum += m[r][c]
+
+        nRow = []
+
+        if sum == 0:
+            nRow = m[r]
+        else:
+            for c in range(cols):
+                if use_fractions:
+                    nRow.append(Fraction(m[r][c], sum))
+                else:
+                    nRow.append(float(m[r][c])/sum)
+        n.append(nRow)
+    return n
+
+def subtract(i, q):
+    if len(i) != len(i[0]) or len(q) != len(q[0]):
+        raise Exception("non-square matrices")
+
+    if len(i) != len(q) or len(i[0]) != len(q[0]):
+        raise Exception("Cannot subtract matrices of different sizes")
+
+    s = []
+    for r in range(len(i)):
+        sRow = []
+        for c in range(len(i[r])):
+            sRow.append(i[r][c] - q[r][c])
+        s.append(sRow)
+    return s
+
+def multiply(a, b):
+    if a == [] or b == []:
+        raise Exception("Cannot multiply empty matrices")
+
+    if len(a[0]) != len(b):
+        raise Exception("Cannot multiply matrices of incompatible sizes")
+
+    m = []
+    rows = len(a)
+    cols = len(b[0])
+    iters = len(a[0])
+
+    for r in range(rows):
+        mRow = []
+        for c in range(cols):
+            sum = 0
+            for i in range(iters):
+                sum += a[r][i]*b[i][c]
+            mRow.append(sum)
+        m.append(mRow)
+    return m
+
+def transposeMatrix(m):
+    t = []
+    for r in range(len(m)):
+        tRow = []
+        for c in range(len(m[r])):
+            if c == r:
+                tRow.append(m[r][c])
+            else:
+                tRow.append(m[c][r])
+        t.append(tRow)
+    return t
+
+def getMatrixMinor(m,i,j):
+    return [row[:j] + row[j+1:] for row in (m[:i]+m[i+1:])]
+
+def getMatrixDeternminant(m):
+    if len(m) == 2:
+        return m[0][0]*m[1][1]-m[0][1]*m[1][0]
+
+    d = 0
+    for c in range(len(m)):
+        d += ((-1)**c)*m[0][c]*getMatrixDeternminant(getMatrixMinor(m,0,c))
+
+    return d
+
+def getMatrixInverse(m):
+    d = getMatrixDeternminant(m)
+
+    if d == 0:
+        raise Exception("Cannot get inverse of matrix with zero determinant")
+
+    if len(m) == 2:
+        return [[m[1][1]/d, -1*m[0][1]/d],
+                [-1*m[1][0]/d, m[0][0]/d]]
+
+    cofactors = []
+    for r in range(len(m)):
+        cofactorRow = []
+        for c in range(len(m)):
+            minor = getMatrixMinor(m,r,c)
+            cofactorRow.append(((-1)**(r+c)) * getMatrixDeternminant(minor))
+        cofactors.append(cofactorRow)
+    cofactors = transposeMatrix(cofactors)
+    for r in range(len(cofactors)):
+        for c in range(len(cofactors)):
+            cofactors[r][c] = cofactors[r][c]/d
+    return cofactors
+
+
+def find_lcm(a):
+    lcm = a[0]
+    for i in a[1:]:
+        lcm = lcm*i//gcd(lcm, i)
+    return lcm
+    
+def simplify(b):
+    a = [i.denominator for i in b]
+    lcm = find_lcm(a)
+
+    c = []
+    for i in range(len(b)):
+        c.append(lcm * b[i].numerator // b[i].denominator)
+    c.append(lcm)
+    return c
+
+def solution(m):
+    m = sort(m)
+    n = normalize(m,use_fractions=True)
+    (q, r) = decompose(n)
+    i = identity(len(q))
+    s = subtract(i, q)
+    v = getMatrixInverse(s)
+    b = multiply(v, r)
+
+    return simplify(b[0])
+
+print(solution([[0, 2, 1, 0, 0],
+                [0, 0, 0, 3, 4],
+                [0, 0, 0, 0, 0],
+                [0, 0, 0, 0, 0],
+                [0, 0, 0, 0, 0]]))
+#   [7, 6, 8, 21]
+
+print(solution([[0, 1, 0, 0, 0, 1],
+                [4, 0, 0, 3, 2, 0],
+                [0, 0, 0, 0, 0, 0],
+                [0, 0, 0, 0, 0, 0],
+                [0, 0, 0, 0, 0, 0],
+                [0, 0, 0, 0, 0, 0]]))
+#    [0, 3, 2, 9, 14]