+/* ec_lineal.h
+ ***********************
+ * Análisis Numérico I *
+ ***********************
+ */
+
+// Tipos básicos
+typedef unsigned int Indice;
+typedef float Numero;
+
+// Cantidad máxima de elementos
+const Indice N = 5000;
+
+// Tipos compuestos
+typedef Indice VectorPermutaciones[N];
+typedef Numero Vector[N];
+typedef Numero Matriz[N][N];
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+
+// Declaración de funciones
+void generar_matriz( Matriz& );
+void generar_vector( Vector& );
+void generar_vector_permutaciones( VectorPermutaciones& );
+void imprimir_matriz( Matriz& );
+void imprimir_vector( Vector& );
+void imprimir_matriz_permutada( Matriz&, VectorPermutaciones& );
+void imprimir_vector_permutado( Vector&, VectorPermutaciones& );
+void imprimir_matriz_L( Matriz&, VectorPermutaciones& );
+void imprimir_matriz_U( Matriz&, VectorPermutaciones& );
+void sustitucion_directa( Matriz&, Vector&, VectorPermutaciones& );
+void sustitucion_inversa( Matriz&, Vector&, VectorPermutaciones& );
+
+#include <time.h>
+
+// Constantes de error
+const int MATRIZ_SINGULAR = 1;
+
+int main( int argc, char *argv[] ) {
+
+ Matriz A;
+ Vector b;
+ VectorPermutaciones p;
+
+ // Reinicio los numeros aleatorios
+ //srand( time( NULL ) );
+
+ // Genero matriz al azar
+ generar_matriz( A );
+ // Genero vector al azar
+ generar_vector( b );
+ // Genero vector de permutaciones
+ generar_vector_permutaciones( p );
+
+ // La imprimo
+ printf( "Matriz A:\n=========\n" );
+ imprimir_matriz( A );
+ printf( "\n" );
+ // La imprimo
+ printf( "Vector b:\n=========\n" );
+ imprimir_vector( b );
+ printf( "\n\n" );
+
+ // Recorro por columna (o diagonal).
+ for ( Indice j = 0; j < N; j++ ) {
+ // Busco máximo en la columna para usar de pivote
+ Indice maxi = j;
+ Numero max = fabs( A[j][j] );
+ for ( Indice i = j + 1; i < N; i++ )
+ if ( fabs( A[i][j] ) > max ) {
+ maxi = i;
+ max = fabs( A[i][j] );
+ }
+
+ // Nos fijamos si la matriz es singular
+ if ( max == 0 )
+ return MATRIZ_SINGULAR;
+
+ // Nos fijamos si hay permutación
+ if ( maxi != j ) {
+ // Intercambio vector de permutaciones
+ Indice aux = p[maxi];
+ p[maxi] = p[j];
+ p[j] = aux;
+ }
+
+ // Reduzco las filas
+ for ( Indice i = j + 1; i < N; i++ ) {
+ // Calculo el m, almacenandolo en A
+ A[p[i]][j] /= A[p[j]][j];
+ // Recorro columna por columna, operando
+ for ( Indice k = j + 1; i < N; i++ )
+ A[p[i]][j] -= A[p[i]][j] * A[p[j]][k]; // Aij = Aij - m * Ajk
+ }
+ }
+
+ // Imprimo la matriz modificada
+ printf( "Matriz A modificada:\n====================\n" );
+ imprimir_matriz_permutada( A, p );
+ printf( "\n" );
+ // Imprimo la matriz L
+ printf( "Matriz L:\n========\n" );
+ imprimir_matriz_L( A, p );
+ printf( "\n" );
+ // Imprimo la matriz U
+ printf( "Matriz U:\n========\n" );
+ imprimir_matriz_U( A, p );
+ printf( "\n" );
+
+ // Calcula solucion de Ly=b
+ sustitucion_directa( A, b, p );
+
+ // Imprimo el vector solucion y
+ printf( "y (Ly=b):\n=========\n" );
+ imprimir_vector_permutado( b, p );
+ printf( "\n" );
+
+ // Calcula solucion de Ux=y
+ sustitucion_inversa( A, b, p );
+
+ // Imprimo el vector solucion x
+ printf( "x (Ux=y):\n=========\n" );
+ imprimir_vector_permutado( b, p );
+ printf( "\n" );
+
+ return EXIT_SUCCESS;
+
+}
+
+void generar_matriz( Matriz& m ) {
+ for ( Indice i = 0; i < N; i++ )
+ for ( Indice j = 0; j < N; j++ )
+ m[i][j] = Numero( rand() ) * 3.0 / RAND_MAX;
+}
+
+void generar_vector( Vector& v ) {
+ for ( Indice i = 0; i < N; i++ )
+ v[i] = Numero( rand() ) * 3.0 / RAND_MAX;
+}
+
+void generar_vector_permutaciones( VectorPermutaciones& v ) {
+ for ( Indice i = 0; i < N; i++ )
+ v[i] = i;
+}
+
+void imprimir_matriz( Matriz& m ) {
+ for ( Indice i = 0; i < N; i++ ) {
+ for ( Indice j = 0; j < N; j++ )
+ printf( "%.02f ", m[i][j] );
+ printf( "\n" );
+ }
+}
+
+void imprimir_vector( Vector& v ) {
+ for ( Indice i = 0; i < N; i++ )
+ printf( "%.02f ", v[i] );
+ printf( "\n" );
+}
+
+void imprimir_matriz_permutada( Matriz& m, VectorPermutaciones& p ) {
+ for ( Indice i = 0; i < N; i++ ) {
+ for ( Indice j = 0; j < N; j++ )
+ printf( "%.02f ", m[p[i]][j] );
+ printf( "\n" );
+ }
+}
+
+void imprimir_vector_permutado( Vector& v, VectorPermutaciones& p ) {
+ for ( Indice i = 0; i < N; i++ )
+ printf( "%.02f ", v[p[i]] );
+ printf( "\n" );
+}
+
+void imprimir_matriz_L( Matriz& m, VectorPermutaciones& p ) {
+ for ( Indice i = 0; i < N; i++ ) {
+ for ( Indice j = 0; j < N; j++ )
+ if ( i < j )
+ printf( "0.00 " );
+ else if ( i == j )
+ printf( "1.00 " );
+ else
+ printf( "%.02f ", m[p[i]][j] );
+ printf( "\n" );
+ }
+}
+
+void imprimir_matriz_U( Matriz& m, VectorPermutaciones& p ) {
+ for ( Indice i = 0; i < N; i++ ) {
+ for ( Indice j = 0; j < N; j++ )
+ if ( i > j )
+ printf( "0.00 " );
+ else
+ printf( "%.02f ", m[p[i]][j] );
+ printf( "\n" );
+ }
+}
+
+void sustitucion_directa( Matriz& A, Vector& b, VectorPermutaciones& p ) {
+ for ( Indice i = 1; i < N; i++ ) {
+ Numero sum = 0;
+ for ( Indice j = 0; j < i - 1; j++ )
+ sum += A[p[i]][j] * b[p[i]];
+ b[p[i]] -= sum;
+ }
+}
+
+void sustitucion_inversa( Matriz& A, Vector& b, VectorPermutaciones& p ) {
+ for ( Indice i = N; i > 0; i-- ) {
+ Numero sum = 0;
+ for ( Indice j = i; j < N; j++ )
+ sum += A[p[i-1]][j] * b[p[j]];
+ b[p[i-1]] = ( b[p[i-1]] - sum ) / A[p[i-1]][i];
+ }
+}
projects / z.facultad / 75.12 / ejercicios.git / commitdiff