[ Foro de Pascal ]

Programa sobre matrices

Llevo un tiempo (en mis escaso tiempo libre) creando un programa capaz de sumar, restar, multiplicar matrices, así como hallar la adjunta, inversa, traspuesta o el determinante de una matriz, aplicando para ello procedimientos, funciones, recursividad, etc...ya está casi terminado, faltando sólo un par de procedimientos. He intentado hacer un código fácil de leer, eficiente y limpio, pero hay partes del mismo que no soy capaz de hacer, o al menos, tan bien como me gustaría, aquí os muestro el código (es largo, casi 300 líneas) por si os resulta de ayuda, pues aunque faltan cosas, el programa compila y funciona sin problemas ni errores. También agradecería ayuda en todo lo que se refiera a mejorar el código, así como encontrar la mejor manera de escribir lo que falta:

-En la función del determinante, falta incorporar la recursividad, está ya creado el caso base, pero nada más, falta completarlo.

-Los procedimientos de matriz adjunta e inversa sólo devuelven una frase de prueba, no hacen nada.

Aquí tienen el código (desgraciadamente, se pierde mucha legibilidad por la eliminacion de la sangría, del editor de mensajes):

PROGRAM matrices;
Uses
crt;

Const
N=14; {Maxima dimension permitida}

type
Matriz = array [1..N,1..N] of real;
tpmatriz = ^Matriz;

Var
opcion2matrices,opcion1matriz: char;
opcion, dimension: integer;
a,b,prod: Matriz;

//---procedimiento del primer menú de opciones---//

Procedure menu;
Begin
writeln('Programa para el cálculo de matrices regulares');
Writeln('Esta versión del programa le permite multiplicar o sumar dos matrices');
Writeln('así como hallar la adjunta, inversa, traspuesta o el determinante de una matriz');
Writeln( 'Pulsa cualquier tecla para continuar..');
ReadKey;
ClrScr;
End; {menu}

//---procedimiento para leer la matriz---//

Procedure LeerMatriz(var mat:Matriz);
Var
fil, col: 1..N;
Begin
For fil:=1 to dimension do
For col:= 1 to dimension do
begin
Write(' Introduzca el componente a' ,fil,'.' , col, ': ');
Readln(mat[fil,col]);
end;
End; {Leer Matriz}

//---procedimiento para leer la dimensión de la matriz---//

Procedure leerdimension;
Var
correcto: Boolean;
Procedure Compruebadimension;
Begin
if dimension > N then
correcto := False;
End; {compruebadimension}
Begin
Repeat
correcto:=true;
writeln('Introduce la dimensión de la matriz (las matrices deben ser regulares)');
readln(dimension);
if dimension > N then
begin
writeln('La dimensión máxima permitida es 14x14 vuelve a intentarlo');
correcto:=false
end;
Compruebadimension;
Until
correcto;
End; {leerdimension}

//---procedimiento para sumar matrices---//

Procedure sumaMat (m1,m2:Matriz; var resul:Matriz);
Var i,j:1 .. N;
Begin
for i:=1 to dimension do
for j:=1 to dimension do
begin
resul[i,j]:=0;
resul[i,j]:= m1[i,j] + m2[i,j];
end;
End; {sumaMat}

//---procedimiento para restar matrices---//

Procedure restaMat (m1,m2:Matriz; var resul:Matriz);
Var i,j:1 .. N;
Begin
for i:=1 to dimension do
for j:=1 to dimension do
begin
resul[i,j]:=0;
resul[i,j]:= m1[i,j] - m2[i,j];
end;
End; {restaMat}

//---procedimiento para multiplicar matrices---//

Procedure MultiplicarMat(m1,m2: Matriz; var resul:Matriz);
Var i,j,k:1 .. N;
Begin
for i:=1 to dimension do
for j:=1 to dimension do
begin
resul[i,j]:=0;
for k:=1 to dimension do
resul[i,j]:=resul[i,j]+m1[i,k]*m2[k,j];
end;
End; {MultiplicarMat}

//---procedimiento para hallar la adjunta (sin hacer)---//

Procedure adjunta;
Begin
writeln('adjunta aún no disponible');
End;{adjunta}

//---procedimiento para hallar el determinante---//

function determinante (m1:Matriz): real;
Begin
Case dimension of
1 : begin
determinante:= m1[1,1]{el caso más básico}
end;
2 : begin
determinante:= m1[1,1] * m1[2,2] - m1[1,2] * m1[2,1]{caso base para la recursividad}
end;
3 : begin
determinante:= (m1[1,1] * m1[2,2] * m1[3,3] + m1[1,3] * m1[2,1] * m1[3,2] + m1[1,2] * m1[2,3] * m1[3,1]) - (m1[1,3] * m1[2,2] * m1[3,1] + m1[1,2] * m1[2,1] * m1[3,3] + m1[1,1] * m1[2,3] * m1[3,2]){por sarrus}
end;
else
begin
determinante:= 0;
end;{else}
end;{Case-of}
End;{determinante}

//---procedimiento para hallar la inversa (sin hacer)---//

Procedure inversa;
Begin
writeln('inversa aún no disponible');
End;

//---procedimiento para hallar la traspuesta---//

Procedure traspuesta (m1:Matriz; var resul:Matriz);
Var i,j:1 .. N;
Begin
for i:=1 to dimension do
for j:=1 to dimension do
begin
resul[i,j]:=0;
resul[i,j]:= m1[j,i];
end;
End; {traspuesta}

//---procedimiento para escribir la matriz---//

Procedure EscribirMat(m:Matriz);
Var fil, col: 1..N;
Begin
For fil:=1 to dimension do
begin
write(' ( ');
For col:= 1 to dimension do
Write(m[fil,col]:4:0);
writeln(' ) ');
end;
End; {EscribirMat}

//-----------------------programa principal-----------------------//

Begin
menu;
repeat
writeln('Escribe 1 para cálculos con una matriz, 2 para cálculos con dos matrices');
writeln('y 0 para terminar');
readln(opcion);
case opcion of
1 : begin
leerdimension;
writeln('Introduce los componentes de la matriz A');
LeerMatriz(a);
repeat
writeln('Escribe "a" para la adjunta, "d" para el determinante, "i" para la inversa');
writeln('"t" para la traspuesta y "m" para volver al menú anterior. ');
readln(opcion1matriz);
Case opcion1matriz of
'a' : begin
adjunta;
end; {case 'a'}
'd' : begin
ClrScr;
writeln;
writeln(' Matriz A: ');
EscribirMat(a);
writeln;
writeln(' Determinante de A: ', determinante (a):5:2);
writeln(' Pulsa Enter para volver...');
writeln;
readln;
end; {case 'd'}
'i' : begin
inversa;
end; {case 'i'}
't' : begin
ClrScr;
traspuesta(a,prod);
writeln;
writeln(' Matriz A: ');
EscribirMat(a);
writeln;
writeln(' Traspuesta de A: ');
EscribirMat(prod);
writeln;
writeln(' Pulsa Enter para volver...');
writeln;
readln;
end; {case 't'}
end; {case opcion1matriz}
until opcion1matriz = 'm';
end; {case opcion = 1}
2 : begin
leerdimension;
writeln('Introduce los componentes de la matriz A');
LeerMatriz (a);
writeln('Introduce los componentes de la matriz B');
LeerMatriz (b);
repeat
writeln('Escribe "s" para sumar, "r" para restar, ');
writeln('"x" para multiplicar, y "m" para volver al menú anterior.');
readln(opcion2matrices);
Case opcion2matrices of
's' : begin
ClrScr;
SumaMat(a,b,prod);
writeln;
writeln(' Matriz A: ');
EscribirMat(a);
writeln;
writeln(' Matriz B: ');
EscribirMat(b);
writeln;
writeln(' Resultado A+B: ');
EscribirMat(prod);
writeln(' Pulsa Enter para volver...');
writeln;
readln;
end; {case 's'}
'r' : begin
ClrScr;
restaMat(a,b,prod);
writeln;
writeln(' Matriz A: ');
EscribirMat(a);
writeln;
writeln(' Matriz B: ');
EscribirMat(b);
writeln;
writeln(' Resultado A+B: ');
EscribirMat(prod);
writeln(' Pulsa Enter para volver...');
writeln;
readln;
end; {case 'r'}
'x' : begin
ClrScr;
MultiplicarMat(a,b,prod);
writeln;
writeln(' Matriz A: ');
EscribirMat(a);
writeln;
writeln(' Matriz B: ');
EscribirMat(b);
writeln;
writeln(' Matriz A * B: ');
EscribirMat(prod);
writeln;
writeln(' Pulsa Enter para volver...');
readln;
end; {case'm'}
end; {case opcion2matrices}
until opcion2matrices = 'm';
end; {case opcion = 2}
end {case "principal"}
until opcion = 0
End.

---

Yo lo veo bastante bien. Aun así, un par de comentarios:

a) No entiendo por qué haces:

resul[i,j]:=0;
resul[i,j]:= m1[j,i];

La primera operación es innecesaria, porque el valor que das se sobreescribe en la segunda.

b) En principio, parece más razonable que en cada operación (por ejemplo, la suma) incluyas (cuatro o cinco) parámetros en vez de 3: la primera matriz, la segunda, la de resultado, pero también el tamaño de la matriz (el ancho, si es cuadrada, o ancho y alto si es rectangular). Eso permite que puedas reservar espacio para matrices de 100x100 (por ejemplo) pero que aun asi las operaciones no se realicen sobre 10.000 elementos, sino sólo sobre los que realmente estás usando (por ejemplo, los 9 elementos de una matriz 3x3), y así se gana en rapidez.

---

a)
Imagino que cuando lo hice, fué simplemente para darle un valor inicial y así evitar usarlo antes de tener un valor, aunque ahora que lo veo, si puede ser innecesario.

b)
mmmmmm, no estoy seguro ¿se refiere a hacer algo asi?

Variables globales -> f,c: integer; {f=filas, c=columnas}

Procedure sumaMat (m1,m2:Matriz; f,c:Integer; var resul:Matriz);
Var i,j:1 .. N;
Begin
for i:=1 to f do
for j:=1 to c do
begin
resul[i,j]:=0;
resul[i,j]:= m1[i,j] + m2[i,j];
end;
End; {sumaMat}

---

He modificado el procedimiento del determinante, ya puede utilizarse cualquier dimensión, éste era el viejo:

function determinante (m1:Matriz): real;
Begin
Case dimension of
1 : begin
determinante:= m1[1,1]{el caso más básico}
end;
2 : begin
determinante:= m1[1,1] * m1[2,2] - m1[1,2] * m1[2,1]{caso base para la recursividad}
end;
3 : begin
determinante:= (m1[1,1] * m1[2,2] * m1[3,3] + m1[1,3] * m1[2,1] * m1[3,2] + m1[1,2] * m1[2,3] * m1[3,1]) - (m1[1,3] * m1[2,2] * m1[3,1] + m1[1,2] * m1[2,1] * m1[3,3] + m1[1,1] * m1[2,3] * m1[3,2]){por sarrus}
end;
else
// begin
determinante:= 0;
// end;{else}
end;{Case-of}
End;{determinante}

Y ésta la nueva función:

function determinante(m1:Matriz;dimension:integer):real;
var
n,i,j:integer;
det:real;
B:Matriz;
begin
if dimension=2 then (*caso base para la recursividad*)
determinante:= m1[1,1] * m1[2,2] - m1[1,2] * m1[2,1]
else
begin
det:= 0;
for n:= 1 to dimension do
begin
for i:= 2 to dimension do (*carga la matriz B con los valores*)
begin (*que quedan eliminando la fila y columna*)
for j:= 1 to n-1 do (*correspondiente de la matriz m1*)
B[i-1,j]:= m1[i,j];
for j:= n+1 to dimension do
B[i-1,j-1]:= m1[i,j];
end;
if (1+n) mod 2= 0 then i:=1 (*Signo del complemento algebraico*)
else i:= -1;
det:= det + i * m1[1,n] * determinante(B,dimension-1); (*Llamada recursiva*)
end;
determinante:= det;
end;
end; (*determinante*)

La llamada a la función se hace con -> determinante (a, dimension)

Y con esto, ya solo quedan los procedimientos de la inversa y la adjunta...¿Alguien me echa una mano? :P

P.D. ¿Que diferencia hay, a efectos prácticos, en escribir un comentario entre {} y hacerlo entre (* *)?

---

mmmmm no me ha gustado la nueva "versión" de la función determinante, (entre otras cosas, me parece menos "limpio", y si buscas el determinante de rango 1, te da siempre valor 0), creo que como mejor funciona es "mezclando" ambas funciones, creo que esta podría ser la versión definitiva de la función (con la sangría en condiciones, se ve mucho mejor):

function determinante (m1:Matriz;dimension:integer): real;

Var

Begin
End; (*Determinante*)

---

Esa última versión me gusta más: lo de separar los casos especiales (que además son los más habituales a nivel "académico"), hace que sea más legible, y más eficiente en muchos casos.

En cuanto a la adjunta y la inversa... si ya está correcto la función que calcula determinantes, tienes la mitad del trabajo hecho...  ;-)

---

El procedimiento de la adjunta creo que podré sacarlo solito en cuanto le dedique un poco más tiempo (o al menos, eso quiero pensar :-P ) y sí, evidemente, con todos los demás procedimientos hechos (calculo del determinante, adjunta y traspuesta), hallar la inversa sería muy sencillo aplicando la ecuación:

A^-1 = 1/Det(A) A^*t

Ya que tenemos todos los procedimientos necesarios. Peeero...no era eso lo que tenía en mente. Quería utilizar otro método para su cálculo (y así incluir otra caracteristica al programa, la del cálculo de la forma normal de Hermite). Me explico:

consideremos la matriz A:


pues si buscamos su forma normal de Hermite "incrustando" en su parte derecha la identidad (A | I ) tendríamos:

( 1 2 0 1 | 1 0 0)
( 1 1 2 0 | 0 1 0)
( 1 0 2 1 | 0 0 1)

Que al operar, quedaría en

( 1 0 0 3 | 1 -2 2 )
( 0 1 0 -1 | 0 1 -1 )
(0 0 1 -1 | -1/2 1 -1/2 )

Teniendo la inversa en la parte derecha, y la resolución de Hermite en la parte izquierda. Este procedimiento podría aprovecharse así no sólo para calcular la inversa, sino también para la resolucion de la matriz por Hermite. Y ahí es donde tengo el problema.

---

¿Pero cual de las cosas que comentas es "el problema"? ;-) Me explico:

- Lo de ampliar la matriz para convertirla en una rectangular, que contenga la matriz identidad en su parte derecha, es trivial: Ya estás guardando un "ancho" y un "alto" de la matriz, sólo falta que añadas un dato adicional, el "anchoAmpliado", o que tengas en cuenta que si el ancho original era "n", ahora tienes además otras "n-1" columnas que representan lo que inicialmente era la matriz identidad.

- Y si te refieres a resolver la matriz... no recuerdo cual era el método de Hermite, pero por Gauss es sencillo: tomas un pivote (típicamente el número de abajo a la izquierda), que te permita conseguir que todo lo que hay en la primera columna por debajo de la diagonal principal sea 0, luego haces lo mismo en la segunda columna, y así sucesivamente, hasta tener una matriz diagonal superior. Cuando ya tienes una matriz diagonal superior, podrías aplicar sustitución si sólo quisieras saber las soluciones, o bien, si quieres seguir hasta tener la inversa en la parte derecha, aplicar el mismo método por encima de la diagonal principal, para hacer ceros también en esa zona.

---

Me refiero simplemente a que no encuentro la manera correcta de realizar el algoritmo del método de Hermite......el método es "sencillamente" aplicar transformaciones elementales (Cambiar dos filas de lugar/Multiplicar una fila por una constante distinta de cero/Sumar a una fila otra multiplicada por una constante) por filas a una matriz dada, de forma que dé como resultado una matriz escalonada reducida, esto es aquélla que verifica las siguientes propiedades:
1. Todas las filas llenas de ceros se encuentran en la parte inferior de la matriz.
2. De izquierda a derecha, el primer número distinto de cero (llamado pivote) es siempre igual a 1.
3. El pivote de cada fila está más a la derecha que el pivote de la fila de encima.
4. Todos los elementos de una columna en que hay un pivote son cero, excepto él mismo.

En el ejemplo anterior, pasar de:

a:


Entonces quiero aprovechar esto, para hacer tanto la inversa, como la resolucion de las matrices, en vez de realizar la inversa con la conocida ecuación.

---

Los pasos que yo daría para aplicar ese algoritmo son:

- Tomas como pivote el primer elemento de la primera fila, si no es cero (si lo fuera, tendrías que intercambiar esa fila por la siguiente en que no lo fuera).

Ej:

3 6 8
6 2 3
0 4 5

- Divides todos los elementos de esa fila entre ese valor, para que ese pivote pase a ser uno.

1 2 8/3
6 2 3
0 4 5

- Para cada una de las filas inferiores, si el primer elemento no es cero sino n, deberás restarle a toda la fila la primera multiplicada por -n:

[6 2 3] - 6 [1 2 8/3] = [6 2 3] - [6 12 16] = [0 -10 -13]

luego

1 2 8/3
0 -10 -13
0 4 5

- Ya están los 0 por debajo del primer elemento de la diagonal. Ahora conservas la primera fila tal cual, ya repites los pasos con el elemento (2,2):

Primero preparas un nuevo pivote:

1 2 8/3
0 1 13/10
0 4 5

Y luego haces ceros por debajo de él

[0 4 5] - 4 [0 1 13/10] = [0 4 5] -  [0 4 26/5] = [0 0 -1/5]

obteniendo

1 2 8/3
0 1 13/10
0 0 -1/5

Y así sucesivamente.


Así obtienes una matriz diagonal superior. Luego aplicas el mismo sistema para hacer ceros por encima de la diagonal, y así consigues finalmente una matriz diagonal.

---

---