SageWeek08.sagews

Q1.

y = var('y')
solve([10 * x + 2 * y == 0, 2 * x - y == 154], [x, y])

[[x == 22, y == -110]]

Q2.

b = vector([0,154])
b

(0, 154)

A = matrix([[10, 2], [2, -1]])
A

[10 2] [ 2 -1]

y = var('y')
X = vector([x, y])
A * X

(10*x + 2*y, 2*x - y)

X = vector([22, -110])
A * X  # This checks that above is a solution

(0, 154)

Q3.

A = matrix([[4, -2, 3], [-1, -5, -8], [1, 1, 1]])
b = vector([10, 9, 1])
A \ b

(53/19, -16/19, -18/19)

Q4.

A = matrix([[1,2],[3,4]])
B = matrix([[7,8],[9,10]])
A * B

[25 28] [57 64]

A = matrix([[-1,1],[-1,-1]])
B = matrix([[-1,4],[1,1]])
A * B

[ 2 -3] [ 0 -5]

A = matrix([[0,144],[-2,1]])
B = matrix([[1,0],[0,1]])
A * B

[ 0 144] [ -2 1]

A = matrix([[1,0],[0,1]])
B = matrix([[-2,0],[-1,-17]])
A * B

[ -2 0] [ -1 -17]

A = matrix([[0,-1],[3,1]])
B = matrix([[1/3,1/3],[-1,0]])
A * B

[1 0] [0 1]

Q5.

A = matrix([[10, 2], [2, -1]])
b = vector([0,154])
A.inverse() * b

(22, -110)

A = matrix([[4, -2, 3], [-1, -5, -8], [1, 1, 1]])
b = vector([10, 9, 1])
A.inverse() * b

(53/19, -16/19, -18/19)

Q6.

def attemptinv(matrix):
    """
    A function to try to calculate the inverse of a matrix

    Arguments: a matrix

    Outputs: the inverse of the matrix if it exists, False otherwise
    """
    try:
        return matrix.inverse()
    except:
        return False

matrixlist = [matrix([[1/2, 0 , 0, -1, 1], [-1, -1, 1, -1/2, 2], [0, -1, 0, -2, 0] ,[0, 0, 1/2, -1, 0], [-1, 0, -2, 2, 0]]),
              matrix([[-1, -1, 0, 0, -1], [2, 1, 0, 1, 1], [-2, 0, 1, 2, 2], [-1/2, 0, -1/2, 0, 1/2], [0, 0, 0, 1/2, -1]]),
              matrix([[-1/2, -1/2], [-2, -1]]),
              matrix([[2, -2, 1], [6, -1, 1], [12, -2, 2]]),
              matrix([[1, 2], [2, 0]])]
[[attemptinv(M), M.det()] for M in matrixlist]  # Use list comprehensions to obtain inverse and determinant.

[[[ 8/19 -4/19 4/19 -36/19 -11/19] [ 8/19 -4/19 -15/19 40/19 8/19] [ -8/19 4/19 -4/19 -2/19 -8/19] [ -4/19 2/19 -2/19 -20/19 -4/19] [ 11/19 4/19 -4/19 -2/19 3/38], -19/4], [[ 7/13 7/13 -2/13 -4/13 -6/13] [-23/13 -10/13 1/13 2/13 16/13] [ -4/13 -4/13 3/13 -20/13 -4/13] [ 6/13 6/13 2/13 4/13 6/13] [ 3/13 3/13 1/13 2/13 -10/13], -13/4], [[ 2 -1] [-4 1], -1/2], [False, 0], [[ 0 1/2] [ 1/2 -1/4], -4]]

Q7.

matrixsamples = [random_matrix(QQ, 5) for i in range(10000)]  # Creating a 10000 random matrices

matriceswithnoinverse = [m for m in matrixsamples if not attemptinv(m)]  # From sample find matrices with no inverse
dets = [m.det() for m in matriceswithnoinverse]  # Calculate determinant of these matrices
print max(dets) == min(dets)  # Returns true if all dets are equal
dets[0]  # Shows that all matrices with no inverse in our sample have determinant 0

True 0

matriceswith0det = [m for m in matrixsamples if m.det() == 0]  # From sample find matrices with 0 determinant
inverses = [attemptinv(m) for m in matriceswith0det]  # Calculate inverse
any(inverses)  # `any` would return True if not all our elements where False. It returns False so all matrices with 0 determinant have no inverse

False

Q8.

import csv

f = open('W08_D01.csv', 'r')
data = csv.reader(f)
data = [row for row in data]
f.close()

rawdata = [[int(r) for r in row] for row in data[1:]]

class eqnsys():
    """
    A class for a system of equations

    Attributes:
        A: The matrix
        b: The right hand side vector
        x: the solution
        y: the solution
        norm: sqrt(x^2+y^2)

    Methods:
        init
        solve: Solves the system and calculates norm
    """
    def __init__(self, a, b, c, d, f, g):
        self.A = matrix([[a, b], [d, f]])
        self.b = vector([c, g])
    def solve(self):
        """
        A function to solve the system.

        Arguments: self

        Outputs: NA
        """
        try:
            self.x, self.y = self.A.inverse() * self.b
            self.norm = sqrt(self.x ^ 2 + self.y ^ 2)
        except:
            self.x, self.y, self.norm = False, False, False

data = [eqnsys(*row) for row in rawdata]  # Converts all data to instances of eqnsys. Note that this uses some Python kung-fu the following is what I expect students to do and is completely fine:
data = [eqnsys(row[0], row[1], row[2], row[3], row[4], row[5]) for row in rawdata]

outputdata = []
for m in data:
    m.solve()
    outputdata.append([m.x, m.y, m.norm])

nbr = len(data)
nbrwithsol = len([r for r in outputdata if r[0]])
print "%s have a solution" % (float(nbrwithsol) / nbr)  # Need to force real division here (len returns python objects  and not sage objects)

0.9845 have a solution

Q9.

import csv

f = open('W08_D02.csv', 'r')
data = csv.reader(f)
data = [row for row in data]
f.close()

data = [[eval(r) for r in row] for row in data]
mat = matrix(data)

mat.dimensions()

(501, 501)

mat.plot()