(%i3) eq1 : 3 * x - y - z + t = 2 $
eq2 : x + y - 2 * z - 5 * t = 1 $
sol : linsolve([ eq1, eq2],[ x, y, z, t]) ;
(sol) [ x = 3 %r2 + 4 %r1 + 3 4 , y = 5 %r2 + 16 %r1 + 1 4 , z = %r2 , t = %r1 ]

Expresemos las soluciones en forma matricial:
Un punto de la variedad lineal solución del sistema sería el dado con los parámetros igual a cero:

(%i4) P :[ ev( x, ev( sol[ 1], %r1 = 0, %r2 = 0)),
   ev( y, ev( sol[ 2], %r1 = 0, %r2 = 0)),
   ev( z, ev( sol[ 3], %r1 = 0, %r2 = 0)),
   ev( t, ev( sol[ 4], %r1 = 0, %r2 = 0))] ;
(P) [ 3 4 , 1 4 , 0 , 0 ]

Los vectores del espacio director de la variedad intersección son:

(%i6) v :[ ev( x, ev( sol[ 1], %r1 = 1, %r2 = 0)),
   ev( y, ev( sol[ 2], %r1 = 1, %r2 = 0)),
   ev( z, ev( sol[ 3], %r1 = 1, %r2 = 0)),
   ev( t, ev( sol[ 4], %r1 = 1, %r2 = 0))] - P ;
u :[ ev( x, ev( sol[ 1], %r1 = 0, %r2 = 1)),
   ev( y, ev( sol[ 2], %r1 = 0, %r2 = 1)),
   ev( z, ev( sol[ 3], %r1 = 0, %r2 = 1)),
   ev( t, ev( sol[ 4], %r1 = 0, %r2 = 1))] - P ;
(v) [ 1 , 4 , 0 , 1 ] (u) [ 3 4 , 5 4 , 1 , 0 ]

La forma matricial es:

(%i7) print( matrix([ x],[ y],[ z],[ t]), "=", transpose( matrix( P)), "+ λ",
   transpose( matrix( v)), "+ μ",
   transpose( matrix( u))) $ ;
[ x y z t ] = [ 3 4 1 4 0 0 ] + λ [ 1 4 0 1 ] + μ [ 3 4 5 4 1 0 ]

Created with wxMaxima.