Sistemas de ecuaciones, matrices y determinantes

Utilizaremos como ejemplos los proporcionados por L. Merino y E. Santos en el primer bloque de su página de resolución de ejercicios tipo de su libro Álgebra lineal con métodos elementales.

Forma normal de Hermite de una matriz

Para calcular la forma normal de Hermite por filas, vamos a utilizar la librería sympy. Es interesante echar un vistazo a la documentación que ofrecen sobre matrices y álgebra lineal.

from sympy import Matrix

Repasemos cómo podemos hacer operaciones elementales con las filas y columnas de una matriz.

A=Matrix(3,4,range(1,13))
A

\(\displaystyle \left[\begin{matrix}1 & 2 & 3 & 4\\5 & 6 & 7 & 8\\9 & 10 & 11 & 12\end{matrix}\right]\)

Recordemos que los índices en python empiezan en 0. Así la primera fila es

A[0,:]

\(\displaystyle \left[\begin{matrix}1 & 2 & 3 & 4\end{matrix}\right]\)

Y la segunda columna

A[:,1]

\(\displaystyle \left[\begin{matrix}2\\6\\10\end{matrix}\right]\)

Podemos intercambiar columnas con col_swap.

A.col_swap(0,1)
A

\(\displaystyle \left[\begin{matrix}2 & 1 & 3 & 4\\6 & 5 & 7 & 8\\10 & 9 & 11 & 12\end{matrix}\right]\)

Y filas con row_swap.

A.row_swap(0,1)
A

\(\displaystyle \left[\begin{matrix}6 & 5 & 7 & 8\\2 & 1 & 3 & 4\\10 & 9 & 11 & 12\end{matrix}\right]\)

Para restarle a la primera fila tres veces la segunda, hacemos

A[0,:]=A[0,:]-3*A[1,:]
A

\(\displaystyle \left[\begin{matrix}0 & 2 & -2 & -4\\2 & 1 & 3 & 4\\10 & 9 & 11 & 12\end{matrix}\right]\)

Ejercicio

Calcular la forma de Hermite por filas y el rango de la matriz: \[ A= \left(\begin{array}{rrrr} -2 & -4 & 2 & 2 \\ 3 & 6 & -3 & -3 \\ 1 & 2 & 0 & 1 \end{array}\right). \]

Empezamos definiendo la matriz del enunciado.

A=Matrix([[-2,-4,2,2],[3,6,-3,-3],[1,2,0,1]])
A

\(\displaystyle \left[\begin{matrix}-2 & -4 & 2 & 2\\3 & 6 & -3 & -3\\1 & 2 & 0 & 1\end{matrix}\right]\)

Para calcular la forma normal de Hermite por filas, usamos rref (row reduced echelon form).

A.rref()
(Matrix([
 [1, 2, 0, 1],
 [0, 0, 1, 2],
 [0, 0, 0, 0]]),
 (0, 2))
_[0]

\(\displaystyle \left[\begin{matrix}1 & 2 & 0 & 1\\0 & 0 & 1 & 2\\0 & 0 & 0 & 0\end{matrix}\right]\)

Y el rango lo podemos calcular con rank.

A.rank()
2

Ejercicio

Calcular la forma de Hermite por columnas de la matriz \[ A= \left(\begin{array}{rrr} 1 & 2 & 1 \\ 2 & 1 & 0 \\ 4 & 5 & 2 \\ \end{array}\right). \]

A=Matrix([[1,2,1],[2,1,0],[4,5,2]])
A

\(\displaystyle \left[\begin{matrix}1 & 2 & 1\\2 & 1 & 0\\4 & 5 & 2\end{matrix}\right]\)

Para calcular la forma normal reducida por columnas, calculamos la transpuesta de la forma reducida por filas de la transpuesta de la matriz.

At=A.transpose()
rAt,piv =At.rref()
rAt

\(\displaystyle \left[\begin{matrix}1 & 0 & 2\\0 & 1 & 1\\0 & 0 & 0\end{matrix}\right]\)

rAt.transpose()

\(\displaystyle \left[\begin{matrix}1 & 0 & 0\\0 & 1 & 0\\2 & 1 & 0\end{matrix}\right]\)

Ejercicio

Calcular la forma de Hermite por filas \(H\) de la matriz \[ A= \left(\begin{array}{rrrr} -2 & -4 & 2 & 2 \\ 3 & 6 & -3 & -3 \\ 1 & 2 & 0 & 1 \end{array}\right) \] así como una matriz regular \(Q\) de forma que \(H=Q\cdot A\).

Empezamos definiendo la matriz del enunciado, y le añadimos la matriz identidad por columnas. En sympy, la identidad se puede crear con eye.

from sympy import eye
A=Matrix([[-2,-4,2,2],[3,6,-3,-3],[1,2,0,1]])
AI=A.row_join(eye(3))
AI

\(\displaystyle \left[\begin{matrix}-2 & -4 & 2 & 2 & 1 & 0 & 0\\3 & 6 & -3 & -3 & 0 & 1 & 0\\1 & 2 & 0 & 1 & 0 & 0 & 1\end{matrix}\right]\)

Que también podemos definir con Matrix.hstack (vstack para columnas).

Matrix.hstack(A,eye(3))

\(\displaystyle \left[\begin{matrix}-2 & -4 & 2 & 2 & 1 & 0 & 0\\3 & 6 & -3 & -3 & 0 & 1 & 0\\1 & 2 & 0 & 1 & 0 & 0 & 1\end{matrix}\right]\)

Ahora calculamos su forma normal reducida for filas.

hAI,piv=AI.rref()
hAI

\(\displaystyle \left[\begin{matrix}1 & 2 & 0 & 1 & 0 & 0 & 1\\0 & 0 & 1 & 2 & 0 & - \frac{1}{3} & 1\\0 & 0 & 0 & 0 & 1 & \frac{2}{3} & 0\end{matrix}\right]\)

Las primeras cuatro filas nos dan la forma normal reducida por filas.

hAI[:,0:4]

\(\displaystyle \left[\begin{matrix}1 & 2 & 0 & 1\\0 & 0 & 1 & 2\\0 & 0 & 0 & 0\end{matrix}\right]\)

Mientras que las restantes nos dan la matriz de paso.

Q=hAI[:,4:7]
Q

\(\displaystyle \left[\begin{matrix}0 & 0 & 1\\0 & - \frac{1}{3} & 1\\1 & \frac{2}{3} & 0\end{matrix}\right]\)

Comprobemos que efectivamente \(Q A\) es la forma reducida.

Q*A

\(\displaystyle \left[\begin{matrix}1 & 2 & 0 & 1\\0 & 0 & 1 & 2\\0 & 0 & 0 & 0\end{matrix}\right]\)

Ejercicio

Calcular la forma de Hermite por columnas \(C\) de la matriz \[ A= \left(\begin{array}{rrrr} -2 & -4 & 2 & 2 \\ 3 & 6 & -3 & -3 \\ 1 & 2 & 0 & 1 \\ \end{array}\right) \] así como una matriz regular \(Q\) de forma que \(C=A\cdot Q\).

Procedemos como en el ejercicio anterior, pero ahora añadiendo la identidad por columnas debajo de \(A\).

A=Matrix([[-2,-4,2,2],[3,6,-3,-3],[1,2,0,1]])
AI=A.col_join(eye(4))
AI

\(\displaystyle \left[\begin{matrix}-2 & -4 & 2 & 2\\3 & 6 & -3 & -3\\1 & 2 & 0 & 1\\1 & 0 & 0 & 0\\0 & 1 & 0 & 0\\0 & 0 & 1 & 0\\0 & 0 & 0 & 1\end{matrix}\right]\)

Calculamos la forma normal de hermite de AI por columnas (transponiendo, calculando la forma reducida por filas, y volviendo a transponer).

rtAI,piv=AI.transpose().rref()
rAI=rtAI.transpose()
rAI

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & 0\\- \frac{3}{2} & 0 & 0 & 0\\0 & 1 & 0 & 0\\0 & 0 & 1 & 0\\0 & 0 & 0 & 1\\\frac{1}{2} & -1 & 2 & 4\\0 & 1 & -1 & -2\end{matrix}\right]\)

Así la forma reducida por columnas es

H=rAI[0:3,:]
H

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & 0\\- \frac{3}{2} & 0 & 0 & 0\\0 & 1 & 0 & 0\end{matrix}\right]\)

Y la matriz \(Q\) es

Q=rAI[3:8,:]
Q

\(\displaystyle \left[\begin{matrix}0 & 0 & 1 & 0\\0 & 0 & 0 & 1\\\frac{1}{2} & -1 & 2 & 4\\0 & 1 & -1 & -2\end{matrix}\right]\)

Comprobamos que el resultado es correcto.

A*Q==H
True

Sistemas de ecuaciones usando Gauss-Jordan

Ejercicio

Discutir y resolver, en su caso, el sistema de ecuaciones lineales: \[ \left\{ \begin{array}{rcl} 3x +2y +z&= & 1, \\ 2x+3y+z &= & 0,\\ 2x+ y+3z &= & 0. \end{array} \right. \]

Empezamos definiendo la matriz de coeficientes y la matriz ampliada.

A=Matrix([[3,2,1],[2,3,1],[2,1,3]])
b=Matrix([1,0,0])
A

\(\displaystyle \left[\begin{matrix}3 & 2 & 1\\2 & 3 & 1\\2 & 1 & 3\end{matrix}\right]\)

b

\(\displaystyle \left[\begin{matrix}1\\0\\0\end{matrix}\right]\)

Para resolverlo, usamos linsolve.

from sympy import linsolve
linsolve((A,b))

\(\displaystyle \left\{\left( \frac{2}{3}, \ - \frac{1}{3}, \ - \frac{1}{3}\right)\right\}\)

Por lo que la solución es \(x=2/3\), \(y=-1/3\) y \(z=-1/3\).

Usando la forma normal reducida de la matriz ampliada

Ab=A.row_join(b)
Ab

\(\displaystyle \left[\begin{matrix}3 & 2 & 1 & 1\\2 & 3 & 1 & 0\\2 & 1 & 3 & 0\end{matrix}\right]\)

rAb,piv=Ab.rref()
rAb

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & \frac{2}{3}\\0 & 1 & 0 & - \frac{1}{3}\\0 & 0 & 1 & - \frac{1}{3}\end{matrix}\right]\)

Obtenemos que la solución es única: \(x=2/3\), \(y=-1/3\) y \(z=-1/3\).

En este caso los rangos de la matriz de coeficientes y de la matriz ampliada coinciden y son máximos.

A.rank()
3
Ab.rank()
3

Ejercicio

Discutir y resolver, en su caso, el sistema de ecuaciones lineales: \[ \left\{ \begin{array}{rcl} 3x +2y +z&= & 1, \\ 2x+3y+z &= & 0,\\ x+ 4y+z &= & 0. \end{array} \right. \]

A=Matrix([[3,2,1],[2,3,1],[1,4,1]])
b=Matrix([1,0,0])
linsolve((A,b))

\(\displaystyle \emptyset\)

Que nos indica que el sistema no tiene solución.

Veamos cómo es la forma reducida por filas de la matriz ampliada.

Ab=A.row_join(b)
Ab.rref()[0]

\(\displaystyle \left[\begin{matrix}1 & 0 & \frac{1}{5} & 0\\0 & 1 & \frac{1}{5} & 0\\0 & 0 & 0 & 1\end{matrix}\right]\)

La incompatibilidad del sistema se aprecia en la última fila, que corresponde a una equación \(0=1\).

Podemos ver la incompatibilidad también estudiando los rangos de la matriz de coeficientes y de la ampliada.

A.rank()
2
Ab.rank()
3

Ejercicio

Discutir y resolver, en su caso, el sistema de ecuaciones lineales: \[ \left\{ \begin{array}{rcl} 3x +2y +z&= & 1, \\ 2x+3y+z &= & 2,\\ x+ 4y+z &= & 3. \end{array} \right. \]

A=Matrix([[3,2,1],[2,3,1],[1,4,1]])
b=Matrix([1,2,3])
linsolve((A,b))

\(\displaystyle \left\{\left( - \frac{\tau_{0}}{5} - \frac{1}{5}, \ \frac{4}{5} - \frac{\tau_{0}}{5}, \ \tau_{0}\right)\right\}\)

Por lo que el sistema tiene infinitas soluciones. Podemos tomar \(z\) como parámatro. Y las soluciones son \(x=-1/5-1/5z\), \(y=4/5-1/5z\).

Ab=A.row_join(b)

Veamos que en este caso los rangos de la matriz de coeficientes y de la matriz ampliada coinciden, pero no son máximos.

(A.rank(),Ab.rank())
(2, 2)

Matriz inversa usando operaciones elementales por filas

Ejercicio

Estudiar si la matriz \[ \left(\begin{array}{rrrr} 3 & 2 & 3 & 4 \\ 3 & 2 & 2 & 3 \\ 2 & 1 & 2 & 1 \\ 0 & 1 & 1 & 0 \end{array}\right) \] es regular y en tal caso, calcular su matriz inversa.

A=Matrix([[3,2,3,4],[3,2,2,3],[2,1,2,1],[0,1,1,0]])

Si \(A\) tiene inversa, la podemos calcular con inv o bien elevando a -1.

A.inv()

\(\displaystyle \left[\begin{matrix}- \frac{1}{2} & \frac{1}{2} & \frac{1}{2} & - \frac{1}{2}\\- \frac{1}{2} & \frac{5}{6} & - \frac{1}{2} & \frac{5}{6}\\\frac{1}{2} & - \frac{5}{6} & \frac{1}{2} & \frac{1}{6}\\\frac{1}{2} & - \frac{1}{6} & - \frac{1}{2} & - \frac{1}{6}\end{matrix}\right]\)

A**(-1)

\(\displaystyle \left[\begin{matrix}- \frac{1}{2} & \frac{1}{2} & \frac{1}{2} & - \frac{1}{2}\\- \frac{1}{2} & \frac{5}{6} & - \frac{1}{2} & \frac{5}{6}\\\frac{1}{2} & - \frac{5}{6} & \frac{1}{2} & \frac{1}{6}\\\frac{1}{2} & - \frac{1}{6} & - \frac{1}{2} & - \frac{1}{6}\end{matrix}\right]\)

Veamos cómo se puede calcular utilizando la forma normal de Hermite por filas.

AI=A.row_join(eye(4))
AI

\(\displaystyle \left[\begin{matrix}3 & 2 & 3 & 4 & 1 & 0 & 0 & 0\\3 & 2 & 2 & 3 & 0 & 1 & 0 & 0\\2 & 1 & 2 & 1 & 0 & 0 & 1 & 0\\0 & 1 & 1 & 0 & 0 & 0 & 0 & 1\end{matrix}\right]\)

rAI=AI.rref()[0]
rAI

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & 0 & - \frac{1}{2} & \frac{1}{2} & \frac{1}{2} & - \frac{1}{2}\\0 & 1 & 0 & 0 & - \frac{1}{2} & \frac{5}{6} & - \frac{1}{2} & \frac{5}{6}\\0 & 0 & 1 & 0 & \frac{1}{2} & - \frac{5}{6} & \frac{1}{2} & \frac{1}{6}\\0 & 0 & 0 & 1 & \frac{1}{2} & - \frac{1}{6} & - \frac{1}{2} & - \frac{1}{6}\end{matrix}\right]\)

Por tanto, la inversa es la submatriz formada por las últimas cuatro columnas.

rAI[:,4:9]

\(\displaystyle \left[\begin{matrix}- \frac{1}{2} & \frac{1}{2} & \frac{1}{2} & - \frac{1}{2}\\- \frac{1}{2} & \frac{5}{6} & - \frac{1}{2} & \frac{5}{6}\\\frac{1}{2} & - \frac{5}{6} & \frac{1}{2} & \frac{1}{6}\\\frac{1}{2} & - \frac{1}{6} & - \frac{1}{2} & - \frac{1}{6}\end{matrix}\right]\)

Ejercicio

Estudiar si la matriz \[ \left(\begin{array}{rrrr} 1 & 2 & 3 & 4 \\ 0 & 1 & 2 & 3 \\ 1 & 0 & 2 & 0 \\ 1 & 1 & 1 & 1 \end{array}\right) \] es regular y en tal caso, calcular su matriz inversa.

A=Matrix([(1,2,3,4),(0,1,2,3),(1,0,2,0),(1,1,1,1)])
A

\(\displaystyle \left[\begin{matrix}1 & 2 & 3 & 4\\0 & 1 & 2 & 3\\1 & 0 & 2 & 0\\1 & 1 & 1 & 1\end{matrix}\right]\)

A.rank()
3

Luego no tiene inversa.

Veamos qué aspecto tiene la forma normal reducida por filas cuando le añadimos la identidad.

A.row_join(eye(4)).rref()[0]

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & - \frac{4}{3} & 0 & - \frac{2}{3} & \frac{1}{3} & \frac{2}{3}\\0 & 1 & 0 & \frac{5}{3} & 0 & \frac{1}{3} & - \frac{2}{3} & \frac{2}{3}\\0 & 0 & 1 & \frac{2}{3} & 0 & \frac{1}{3} & \frac{1}{3} & - \frac{1}{3}\\0 & 0 & 0 & 0 & 1 & -1 & 0 & -1\end{matrix}\right]\)

La última fila muestra que la matriz no puede tener inversa.

Matriz de paso entre dos matrices equivalentes por filas

Ejercicio

Estudiar si las matrices: \[ A= \left(\begin{array}{rrrr} 1 & 1 & 1 & 1 \\ 0 & 1 & 1 & 1 \\ 0 & 0 & 1 & 1 \\ \end{array}\right) ; \;\; B= \left(\begin{array}{rrrr} 1 & 0 & 0 & 0 \\ 1 & 1 & 0 & 0 \\ 1 & 1 & 1 & 1 \\ \end{array}\right) \] son equivalentes por filas, y en tal caso determinar una matriz regular \(Q\) de forma que \(A=Q\cdot B\).

A=Matrix(3,4, [1,1,1,1,0,1,1,1,0,0,1,1])
A

\(\displaystyle \left[\begin{matrix}1 & 1 & 1 & 1\\0 & 1 & 1 & 1\\0 & 0 & 1 & 1\end{matrix}\right]\)

B=Matrix([(1,0,0,0),(1,1,0,0),(1,1,1,1)])
B

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & 0\\1 & 1 & 0 & 0\\1 & 1 & 1 & 1\end{matrix}\right]\)

Veamos cómo son sus formas normales de Hermite por filas.

A.rref()[0]

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & 0\\0 & 1 & 0 & 0\\0 & 0 & 1 & 1\end{matrix}\right]\)

B.rref(pivots=False)

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & 0\\0 & 1 & 0 & 0\\0 & 0 & 1 & 1\end{matrix}\right]\)

Ahora calculemos la matriz \(Q\) del enunciado. Para ello buscamos las matrices \(Q_A\) y \(Q_B\) tales que \(Q_A A\) y \(Q_B B\) estén en forma normal reducida por filas.

rAI=A.row_join(eye(3)).rref(pivots=False)
rAI

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & 0 & 1 & -1 & 0\\0 & 1 & 0 & 0 & 0 & 1 & -1\\0 & 0 & 1 & 1 & 0 & 0 & 1\end{matrix}\right]\)

QA=rAI[:,4:]
QA

\(\displaystyle \left[\begin{matrix}1 & -1 & 0\\0 & 1 & -1\\0 & 0 & 1\end{matrix}\right]\)

QB=(B.row_join(eye(3)).rref(pivots=False))[:,4:]
QB

\(\displaystyle \left[\begin{matrix}1 & 0 & 0\\-1 & 1 & 0\\0 & -1 & 1\end{matrix}\right]\)

Sabemos que \(Q_A A = Q_B B\), por lo que la matriz \(Q\) que buscamos es \(Q_A^{-1}Q_B\).

Q=QA**(-1)*QB
Q

\(\displaystyle \left[\begin{matrix}0 & 0 & 1\\-1 & 0 & 1\\0 & -1 & 1\end{matrix}\right]\)

Finalmente, comprobamos que \(QB\) es \(A\).

Q*B==A
True

Determinantes y operaciones elementales

Ejercicio

Calcular el determinante \[ \left| \begin{array}{rrrr} 1 & 1/2 & 1/3 & 1/5\\ -1 & 1/2 & -1/3 & -1/5\\ 1 & 1/2 & 4/3 & 1/5\\ 2 & 1 & 2/3 & 11/5\\ \end{array} \right|. \]

A=Matrix([(1,1/2,1/3,1/5),(-1,1/2,-1/3,-1/5),(1,1/2,4/3,1/5),(2,1,2/3,11/5)])
A

\(\displaystyle \left[\begin{matrix}1 & 0.5 & 0.333333333333333 & 0.2\\-1 & 0.5 & -0.333333333333333 & -0.2\\1 & 0.5 & 1.33333333333333 & 0.2\\2 & 1 & 0.666666666666667 & 2.2\end{matrix}\right]\)

from sympy import Rational,S,Integer
1/2
0.5
Integer(1)/2

\(\displaystyle \frac{1}{2}\)

S(1)/2

\(\displaystyle \frac{1}{2}\)

Rational(1/2)

\(\displaystyle \frac{1}{2}\)

Rational("1/2")

\(\displaystyle \frac{1}{2}\)

Rational(1,2)

\(\displaystyle \frac{1}{2}\)

En sympy podemos utilizar det para calcular el determinante de una matriz.

A.det()

\(\displaystyle 1.8\)

u=Integer(1)
A=Matrix([(u,u/2,u/3,u/5),(-1,u/2,-u/3,-u/5),(u,u/2,u*4/3,u/5),(2,1,u*2/3,u*11/5)])
A

\(\displaystyle \left[\begin{matrix}1 & \frac{1}{2} & \frac{1}{3} & \frac{1}{5}\\-1 & \frac{1}{2} & - \frac{1}{3} & - \frac{1}{5}\\1 & \frac{1}{2} & \frac{4}{3} & \frac{1}{5}\\2 & 1 & \frac{2}{3} & \frac{11}{5}\end{matrix}\right]\)

A.det()

\(\displaystyle \frac{9}{5}\)

Hagamos operaciones elementales por filas que no alteren el determinante de forma que lleguemos a una matriz triangular.

A

\(\displaystyle \left[\begin{matrix}1 & \frac{1}{2} & \frac{1}{3} & \frac{1}{5}\\-1 & \frac{1}{2} & - \frac{1}{3} & - \frac{1}{5}\\1 & \frac{1}{2} & \frac{4}{3} & \frac{1}{5}\\2 & 1 & \frac{2}{3} & \frac{11}{5}\end{matrix}\right]\)

Vamos a conseguir ceros debajo de la primera posición de la primera fila. Para ello empezamos cambiando la segunda fila por la suma de las dos primeras files (lo que no altera el valor del determinante). Hay que tener cuidado con los índices, pues en python empiezan a contar desde cero.

El método elementary_row_operation nos permite hacer esto.

A.elementary_row_op(op="n->n+km",k=1,row=1,row2=0)

\(\displaystyle \left[\begin{matrix}1 & \frac{1}{2} & \frac{1}{3} & \frac{1}{5}\\0 & 1 & 0 & 0\\1 & \frac{1}{2} & \frac{4}{3} & \frac{1}{5}\\2 & 1 & \frac{2}{3} & \frac{11}{5}\end{matrix}\right]\)

Pero si no queremos conservar la matriz original, podemos simplemente sumar las filas correspondientes.

A[1,:]=A[1,:]+A[0,:]
A

\(\displaystyle \left[\begin{matrix}1 & \frac{1}{2} & \frac{1}{3} & \frac{1}{5}\\0 & 1 & 0 & 0\\1 & \frac{1}{2} & \frac{4}{3} & \frac{1}{5}\\2 & 1 & \frac{2}{3} & \frac{11}{5}\end{matrix}\right]\)

Ahora cambiamos la tercera por la diferencia de la tercera con la primera.

A[2,:]=A[2,:]-A[0,:]
A

\(\displaystyle \left[\begin{matrix}1 & \frac{1}{2} & \frac{1}{3} & \frac{1}{5}\\0 & 1 & 0 & 0\\0 & 0 & 1 & 0\\2 & 1 & \frac{2}{3} & \frac{11}{5}\end{matrix}\right]\)

Y la cuarta por la diferencia de la cuarta por dos veces la primera.

A[3,:]=A[3,:]-2*A[0,:]
A

\(\displaystyle \left[\begin{matrix}1 & \frac{1}{2} & \frac{1}{3} & \frac{1}{5}\\0 & 1 & 0 & 0\\0 & 0 & 1 & 0\\0 & 0 & 0 & \frac{9}{5}\end{matrix}\right]\)

La matriz ya es triangular, con lo que el determinante es el producto de los elementos de la diagonal.

A.diagonal()

\(\displaystyle \left[\begin{matrix}1 & 1 & 1 & \frac{9}{5}\end{matrix}\right]\)

A.det()

\(\displaystyle \frac{9}{5}\)

from functools import reduce
reduce(lambda x,y:x*y,list(A.diagonal()))

\(\displaystyle \frac{9}{5}\)

Ejercicio

Calcular el determinante dependiendo de los valores de los parámetros \(a\) y \(b\). Decidir para qué valores es cero. \[ \left| \begin{array}{rrr} 1-b & 1-a & a-1\\ -1 & 2-a-b & a\\ -1 & 1-a & a+1-b\\ \end{array} \right|. \]

Empezamos definiendo dos símbolos, a y b

from sympy import symbols
a,b =symbols("a,b")
a,b
(a, b)
A=Matrix([[1-b,1-a,a-1],[-1,2-a-b,a],[-1,1-a,a+1-b]])
A

\(\displaystyle \left[\begin{matrix}1 - b & 1 - a & a - 1\\-1 & - a - b + 2 & a\\-1 & 1 - a & a - b + 1\end{matrix}\right]\)

p=A.det()
p

\(\displaystyle - b^{3} + 4 b^{2} - 5 b + 2\)

Usamos factor para factorizar.

p.factor()

\(\displaystyle - \left(b - 2\right) \left(b - 1\right)^{2}\)

Por lo que sabemos que el determinante se anula para los valores de \(b=2\) y \(b=1\), pero no depende de \(a\).

Podríamos usar el comando solve para encontrar las raíces del determinante de la matriz.

from sympy import solve
solve(p)
[1, 2]

Veamos ahora cómo podemos hacer operaciones elementales por filas y columnas con A tal y como lo hemos definido. Vamos a crear una nueva matriz para no perder la matriz original. No podemos hacer B=A, pues cualquier cambio en B se haría en A (al ser básicamente listas); usamos el método copy.

B=A.copy()
B

\(\displaystyle \left[\begin{matrix}1 - b & 1 - a & a - 1\\-1 & - a - b + 2 & a\\-1 & 1 - a & a - b + 1\end{matrix}\right]\)

Le restamos a la tercera fila la segunda.

B[2,:]=B[2,:]-B[1,:]
B

\(\displaystyle \left[\begin{matrix}1 - b & 1 - a & a - 1\\-1 & - a - b + 2 & a\\0 & b - 1 & 1 - b\end{matrix}\right]\)

Vemos que A queda intacta.

A

\(\displaystyle \left[\begin{matrix}1 - b & 1 - a & a - 1\\-1 & - a - b + 2 & a\\-1 & 1 - a & a - b + 1\end{matrix}\right]\)

Ahora podemos cambiar la segunda columna por su suma con la tercera.

B[:,1]=B[:,1]+B[:,2]
B

\(\displaystyle \left[\begin{matrix}1 - b & 0 & a - 1\\-1 & 2 - b & a\\0 & 0 & 1 - b\end{matrix}\right]\)

Ahora el determinante es much más sencillo de calcular, pues se puede desarrollar por la segunda columna.

Ejercicio

Discutir y resolver en función del parámetro \(a\) el sistema de ecuaciones lineales: \[ \left\{ \begin{array}{rcl} x+ay+az&= & 1, \\ x+2ay+(a+1)z&= & 1, \\ 2x+ay+az & = & 2. \\ \end{array} \right. \]

A=Matrix([[1,a,a],[1,2*a,a+1],[2,a,a]])
A

\(\displaystyle \left[\begin{matrix}1 & a & a\\1 & 2 a & a + 1\\2 & a & a\end{matrix}\right]\)

b=Matrix([1,1,2])
Ab=A.row_join(b)
Ab

\(\displaystyle \left[\begin{matrix}1 & a & a & 1\\1 & 2 a & a + 1 & 1\\2 & a & a & 2\end{matrix}\right]\)

Hagamos una copia de la matriz extendida y reduzcamos por filas.

B=Ab.copy()
B

\(\displaystyle \left[\begin{matrix}1 & a & a & 1\\1 & 2 a & a + 1 & 1\\2 & a & a & 2\end{matrix}\right]\)

B[1,:]=B[1,:]-B[0,:]
B

\(\displaystyle \left[\begin{matrix}1 & a & a & 1\\0 & a & 1 & 0\\2 & a & a & 2\end{matrix}\right]\)

B[2,:]=B[2,:]-2*B[0,:]
B

\(\displaystyle \left[\begin{matrix}1 & a & a & 1\\0 & a & 1 & 0\\0 & - a & - a & 0\end{matrix}\right]\)

B[2,:]=B[2,:]+B[1,:]
B

\(\displaystyle \left[\begin{matrix}1 & a & a & 1\\0 & a & 1 & 0\\0 & 0 & 1 - a & 0\end{matrix}\right]\)

Los valores de \(a\) que tenemos que distinguir son 0 y 1, que son precisamente los ceros del determinante de A como polinomio en \(a\).

p=A.det()
p

\(\displaystyle - a^{2} + a\)

p.factor()

\(\displaystyle - a \left(a - 1\right)\)

solve(p)
[0, 1]

Veamos qué ocurre para \(a=0\). Usamos subs (que toma como argumento un diccionario con las substituciones).

Ab0=Ab.subs({a:0})
Ab0.rref(pivots=False)

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & 1\\0 & 0 & 1 & 0\\0 & 0 & 0 & 0\end{matrix}\right]\)

Por lo que en este caso es un sistema compatible indeterminado. Las soluciones son \(x=1\), \(z=0\) e \(y\) puede tomar cualquier valor. Esto también lo podemos conseguir con linsolve.

linsolve((A.subs({a:0}),b))

\(\displaystyle \left\{\left( 1, \ \tau_{0}, \ 0\right)\right\}\)

Veamos ahora qué ocurre con \(a=1\).

Ab1=Ab.subs({a:1})
Ab1

\(\displaystyle \left[\begin{matrix}1 & 1 & 1 & 1\\1 & 2 & 2 & 1\\2 & 1 & 1 & 2\end{matrix}\right]\)

Ab1.rref(pivots=False)

\(\displaystyle \left[\begin{matrix}1 & 0 & 0 & 1\\0 & 1 & 1 & 0\\0 & 0 & 0 & 0\end{matrix}\right]\)

linsolve((A.subs({a:1}),b))

\(\displaystyle \left\{\left( 1, \ - \tau_{0}, \ \tau_{0}\right)\right\}\)

En este caso el sistema también es compatible indeterminado: \(x=1\), \(y=-z\) y \(z\) puede tomar cualquier valor.

Para el resto de valores de \(a\), el sistema es compatible determinado, pues el rango de \(A\) y el de la matriz ampliada es el mismo y es máximo.

Por lo que la solución es \(x=1\), \(y=z=0\).

linsolve((A,b))

\(\displaystyle \left\{\left( 1, \ 0, \ 0\right)\right\}\)

Ejercicio

Discutir, en función de los valores del parámetro \(a\), el sistema \[ \left\{ \begin{array}{rrrl} x&+y& +z&=a,\\ x&+ay& + a^2z&=a,\\ x&+a^2y& +az&=a,\\ x&+y& +az&=0.\\ \end{array} \right. \]

A=Matrix([[1,1,1],[1,a,a**2],[1,a**2,a],[1,1,a]])
A

\(\displaystyle \left[\begin{matrix}1 & 1 & 1\\1 & a & a^{2}\\1 & a^{2} & a\\1 & 1 & a\end{matrix}\right]\)

b=Matrix([a,a,a,0])
Ab=A.row_join(b)
Ab

\(\displaystyle \left[\begin{matrix}1 & 1 & 1 & a\\1 & a & a^{2} & a\\1 & a^{2} & a & a\\1 & 1 & a & 0\end{matrix}\right]\)

Observemos que la matriz ampliada tiene cuatro columnas y cuatro filas, mientras que A tiene sólo tres columnas. Esto quiere decir que como mucho el rango de A es tres, mientras que si el determinante de la matriz ampliada es no nulo, entonces su rango sería cuatro, dando lugar a un sistema incompatible. Veamos pues para qué valores de \(a\) se anula el determinante de Ab.

p=Ab.det()
p.factor()

\(\displaystyle a^{2} \left(a - 1\right)^{2} \left(a + 2\right)\)

Por tanto si \(a\not\in \{-2,0,1\}\) el sistema será incompatible. Veamos qué ocurre para el resto de los casos.

Para \(a=-2\), obtenemos un sistema compatible determinado, y \(x=y=z=-2/3\).

Para \(a=0\), también (\(x=y=z=0\)).

Por último, para \(a=1\), obtenemos un sistema incompatible.

[linsolve((A.subs({a:i}),b.subs({a:i}))) for i in [-2,0,1]]
[{(-2/3, -2/3, -2/3)}, {(0, 0, 0)}, EmptySet]