(%i3) A : matrix([ 1, - 2, 1, - 1],[ 2, - 3, 2, - 3],[ 3, - 5, 3, 4],[ - 1, 1, - 1, 2]) ;
b : transpose( matrix([ 4, - 1, 3, 5])) $
Ab : addcol( A, b) ;
(A) ( 1 2 1 1 2 3 2 3 3 5 3 4 1 1 1 2 ) (Ab) ( 1 2 1 1 4 2 3 2 3 1 3 5 3 4 3 1 1 1 2 5 )
(%i4) [ rank( A), rank( Ab)] ;
(%o4) [ 3 , 3 ]

Como los rangos coinciden, tiene solución y esta dependerá de las filas linealmente independientes. Para ver cuáles son, escalonaremos la matriz ampliada:

(%i6) for i : 2 while i< 5 do(
   Ab : rowop( Ab, i, 1, Ab[ i, 1] / A[ 1, 1])
) $
Ab ;
0 errores, 0 advertencias (%o6) ( 1 2 1 1 4 0 1 0 1 9 0 1 0 7 9 0 1 0 1 9 )
(%i8) for i : 3 while i< 5 do(
   Ab : rowop( Ab, i, 2, Ab[ i, 2] / Ab[ 2, 2])
) $
Ab ;
(%o8) ( 1 2 1 1 4 0 1 0 1 9 0 0 0 8 0 0 0 0 0 0 )

Consideramos z como parámetro, z=λ

(%i10) A2 : submatrix( 4, Ab, 3, 5) ;
Ab2 : submatrix( 4, Ab, 3, 2, 1, 4) - submatrix( 4, Ab, 1, 2, 4, 5) *λ ;
(A2) ( 1 2 1 0 1 1 0 0 8 ) (Ab2) ( 4 λ 9 0 )
(%i11) A3 : invert( A2) ;
(A3) ( 1 2 3 8 0 1 1 8 0 0 1 8 )
(%i12) rat( A3. Ab2) ;
(%o12)/R/ ( λ 14 9 0 )

Las ecuaciones paramétricas serán:

(%i13) print( transpose( matrix([ x, y, t, z])), "=", addrow( %,[λ])) $
( x y t z ) = ( λ 14 9 0 λ )

Created with wxMaxima.