Tema: Transposicion de una matriz MxN
Ver Mensaje Individual
  #5 (permalink)  
Antiguo 05/06/2012, 13:25
mikelromero
 
Fecha de Ingreso: junio-2012
Mensajes: 7
Antigüedad: 12 años, 8 meses
Puntos: 0
Respuesta: Transposicion de una matriz MxN
Mas que darmela al revés, no me la acaba de completar, es decir, me rellena las primeras 4 columnas y luego ya no funciona mas. Os voy a pasar todo el codigo (es larguisimo, pero solo hay que compilarlo). El metodo de la transpuesta esta en la linea 166, y lo llamo en la línea 344.
Codigo C:
Código c: 
Ver original
  1. #include <stdio.h>
  2. #include <QtCore/QCoreApplication>
  3. #include <iostream>
  4. #include <cmath>
  5. using namespace std;
  6.  
  7.  
  8. void affichage_matrice(double**matrice, int n);
  9. double* algo_solution (double **U, double**L, double*b, int n);
  10. double** cholesky (double **matrice, int taille);
  11. double* fcourbe(double *x,double *y,double *v,int longueur);
  12. double* gnewton(double *x,double *y, double *v, int longueur, int longueur1);
  13. double** transposee(double **chol,int col, int fil);
  14.  
  15.  
  16. int main(int argc, char *argv[]){
  17.     QCoreApplication a(argc, argv);
  18.  
  19.     int longueur;
  20.     int longueur1;
  21.     double *x;
  22.     double *y;
  23.     double *v;
  24.  
  25.  
  26.     //Allocation de memoire dynamique
  27.     x =(double*)malloc(sizeof(double)*longueur);
  28.     y =(double*)malloc(sizeof(double)*longueur);
  29.     v =(double*)malloc(sizeof(double)*longueur1);
  30.  
  31.  
  32.     longueur=10;
  33.     longueur1=4;
  34.  
  35.     x[0]=0.1;
  36.     x[1]=0.3;
  37.     x[2]=0.7;
  38.     x[3]=1.2;
  39.     x[4]=1.6;
  40.     x[5]=2.2;
  41.     x[6]=2.7;
  42.     x[7]=3.1;
  43.     x[8]=3.5;
  44.     x[9]=3.9;
  45.  
  46.     y[0]=0.558;
  47.     y[1]=0.569;
  48.     y[2]=0.176;
  49.     y[3]=-0.207;
  50.     y[4]=-0.133;
  51.     y[5]=0.132;
  52.     y[6]=0.055;
  53.     y[7]=-0.09;
  54.     y[8]=-0.069;
  55.     y[9]=0.027;
  56.  
  57.     v[0]=1;
  58.     v[1]=2;
  59.     v[2]=2;
  60.     v[3]=1;
  61.  
  62.     double *model;
  63.     model=gnewton(x,y,v,longueur,longueur1);
  64.  
  65.     int z;
  66.     printf("Le vecteur v optimise est:\n");
  67.     for(z=0;z<longueur1;z++){
  68.         printf("%lf " ,model[z]);
  69.         printf("\n");
  70.     }
  71.  
  72.  
  73.     return (int)a.exec;
  74.  
  75.  
  76.  
  77. }
  78.  
  79. double** cholesky (double **matrice, int taille){
  80.  
  81.     int i, j, m;
  82.     double b, a, b1, a1;
  83.     double **chol;
  84.  
  85.     chol =(double**)malloc(sizeof(double*)*taille);
  86.  
  87.     for(i=0;i<taille;i++){
  88.  
  89.         chol[i]=(double*)calloc(taille,sizeof(double));
  90.  
  91.  
  92. }
  93.     int k;
  94.  
  95.     double *matest;
  96.     double *Test;
  97.     double Tf;
  98.     int cont=1;
  99.     int cont1=1;
  100.     matest =(double*)malloc(sizeof(double)*taille);
  101.     Test =(double*)malloc(sizeof(double)*taille);
  102.  
  103. //Initialisation des tableaux
  104.     int z=0;
  105.     for(z=0;z<taille;z++){
  106.         matest[z]=1;
  107. }
  108.  
  109.  
  110. //On teste si la matrice entrée est symétrique définie positive
  111.  
  112.     for (i=0;i<taille;i++){
  113.         for(j=0;j<taille;j++)
  114.     {
  115.         if (abs(matrice[i][j]-matrice[j][i])>0.01 )
  116.         cont=0;
  117.  
  118.     }
  119.  
  120. }
  121.     if (!cont)
  122.          printf ("\n La matrice n'est pas symetrique!");
  123.  
  124.     for(int i=0;i<taille;i++){
  125.         Test[i]=0;
  126.         for(int j=0;j<taille;j++){
  127.  
  128.             Test[i]=Test[i]+(matest[i]*matrice[i][j]);
  129.             Tf=Tf+(Test[i]*matest[i]);
  130.             if(Tf<0)
  131.                 cont1=0;
  132.         }
  133.  
  134.     }
  135.  
  136.  
  137.     if (!cont1)
  138.      printf ("\n La matrice n'est pas definie positive!\n");
  139.  
  140.   if(cont && cont1){
  141.  
  142.     //On calcule les valeurs des cases de la matrice decomposée L
  143.     for(i=0;i<taille;i++){
  144.         b=0;
  145.         for(j=0;j<=i-1;j++){
  146.             b1=chol[i][j]*chol[i][j];
  147.             b=b+b1;
  148.     }
  149.         chol[i][i]=sqrt(matrice[i][i]-b);
  150.  
  151.  
  152.     for(k=i+1;k<taille;k++){
  153.             a=0;
  154.         for(m=0;m<=k-1;m++){
  155.                 a1=chol[i][m]*chol[k][m];
  156.                 a=a+a1;
  157.             }
  158.             chol[k][i]=(matrice[k][i]-a)/chol[i][i];
  159.         }
  160.     }
  161. }
  162.     return chol;
  163. }
  164.  
  165. //ALGORITMO DE CALCULO DE LA TRANSPUESTA
  166. double** transposee(double **chol,int col, int fil){
  167.  
  168.     double **cholin;
  169.     cholin = (double**)malloc(fil * sizeof(*cholin));
  170.     for(int m = 0; m < fil; m++)
  171.     {
  172.         cholin[m] = (double*)malloc(col * sizeof(**cholin));
  173.     }
  174.  
  175.     /* Initialisation */
  176.     for (int i=0; i < fil; i++){
  177.         for (int j=0; j < col; j++){
  178.             cholin[i][j] = chol[j][i];
  179.         }
  180.     }
  181.  
  182.     return cholin;
  183.  
  184. }
  185.  
  186.  
  187. void affichage_matrice(double**matrice, int n)
  188. {
  189.     for (int i=0;i<n;i++)
  190.     {
  191.         for (int j=0;j<n;j++)
  192.         {
  193.             cout<<matrice[i][j]<<"   ";
  194.         }
  195.         cout<<endl;
  196.     }
  197.     cout<<endl;
  198. }
  199.  
  200.  
  201.  
  202.  
  203. double* algo_remontee (double**U, double*b, int n)
  204. {
  205.     double*x;
  206.     x=new double[n];
  207.     double somme=0;
  208.  
  209.     x[n-1]=b[n-1]/U[n-1][n-1];
  210.  
  211.     for(int i=n-2;i>=0;i--)
  212.     {
  213.         somme=0;
  214.         for(int j=i+1;j<n;j++)
  215.         {
  216.             somme+=U[i][j]*x[j];
  217.         }
  218.  
  219.         x[i]=(b[i]-somme)/U[i][i];
  220.     }
  221.     return x;
  222. }
  223.  
  224.  
  225.  
  226.  
  227. double* algo_descente (double**L, double*b, int n)
  228. {
  229.     double*x;
  230.     x=new double[n];
  231.     double somme=0;
  232.  
  233.     x[0]=b[0]/L[0][0];
  234.  
  235.     for(int i=1;i<n;i++)
  236.     {
  237.         somme=0;
  238.         for(int j=0;j<i;j++)
  239.         {
  240.             somme+=L[i][j]*x[j];
  241.         }
  242.  
  243.         x[i]=(b[i]-somme)/L[i][i];
  244.     }
  245.     return x;
  246. }
  247.  
  248.  
  249.  
  250.  
  251. double* algo_solution(double **U, double**L, double*b, int n)
  252. {
  253.  
  254.     // on a A, on décompose en LLt, L=L et Lt=U
  255.    double*y=new double[n];
  256.  
  257.     y=algo_descente(L,b,n);
  258.     return algo_remontee(U,y,n);
  259. }
  260.  
  261.  
  262.  
  263. double* fcourbe(double *x,double *y,double *v,int longueur){
  264.     double *Fc=new double[longueur];
  265.     for(int i=0; i<longueur;i++){
  266.  
  267.         Fc[i]=y[i]-(v[0]*exp(-v[1]*x[i])*sin((v[2]*x[i])+v[3]));
  268.  
  269.  
  270.     }
  271.        return Fc;
  272. }
  273.  
  274.  
  275. double** jacobien_courbe(double *x,double *v,int longueur, int longueur1){
  276.  
  277.     double **J;
  278.  
  279.  
  280.         J = (double**)malloc(longueur * sizeof(*J));
  281.         for(int m = 0; m < longueur; m++)
  282.         {
  283.             J[m] = (double*)malloc(longueur1 * sizeof(**J));
  284.         }
  285.  
  286.     for(int i=0; i<longueur; i++){
  287.         J[i][0]= -exp(-v[1]*x[i])*sin((v[2]*x[i])+v[3]);
  288.         J[i][1]= x[i]*v[0]*exp(-v[1]*x[i])*sin((v[2]*x[i])+v[3]);
  289.         J[i][2]= -x[i]*v[0]*exp(-v[1]*x[i])*cos((v[2]*x[i])+v[3]);
  290.         J[i][3]= -v[0]*exp(-v[1]*x[i])*cos((v[2]*x[i])+v[3]);
  291.     }
  292.     return J;
  293.     int e; int d;
  294.     printf("La matrice J est:");
  295.  
  296.     for(e=0;e<longueur;e++){
  297.         for(d=0;d<longueur1;d++){
  298.             printf("%lf " ,J [e][d]);
  299.  
  300.         }
  301.         printf("\n");
  302.     }
  303.  
  304. }
  305.  
  306.  
  307. double* gnewton(double *x,double *y, double *v, int longueur, int longueur1){
  308.     double *F;
  309.     double **J;
  310.     double **A;
  311.  
  312.     A = (double**)malloc(longueur1 * sizeof(*A));
  313.     for(int m = 0; m < longueur1; m++){
  314.         A[m] = (double*)malloc(longueur1 * sizeof(**A));
  315.     }
  316.     double *B=new double[longueur1];
  317.     double *C=new double[longueur1];
  318.     double **L;
  319.     double **Lt;
  320.     double *S;
  321.     double **Jt;
  322.     int p;
  323.  
  324.  
  325.  
  326.  
  327.  
  328.         do{
  329.  
  330.  
  331.         F=fcourbe(x,y,v,longueur);
  332.  
  333.         J=jacobien_courbe(x,v,longueur,longueur1);
  334.          int e; int d;
  335.          printf("La matrice J est:\n");
  336.  
  337.         for(e=0;e<longueur;e++){
  338.             for(d=0;d<longueur1;d++){
  339.                 printf("%lf " ,J [e][d]);
  340.  
  341.             }
  342.             printf("\n");
  343.         }
  344.         Jt=transposee(J,longueur,longueur1);
  345.  
  346.         printf("La matrice Jt est:\n");
  347.  
  348.         for(e=0;e<longueur;e++){
  349.             for(d=0;d<longueur1;d++){
  350.                 printf("%lf " ,Jt [e][d]);
  351.  
  352.             }
  353.             printf("\n");
  354.         }
  355.  
  356.         for(int i=0;i<longueur1;i++){
  357.             for(int j=0;j<longueur1;j++){
  358.                 A[i][j]=0;
  359.                 for(p=0;p<longueur;p++){
  360.  
  361.                     A[i][j]=A[i][j]+(Jt[i][p]*J[p][j]);
  362.  
  363.                 }
  364.             }
  365.  
  366.         }
  367.  
  368.         for(int i=0;i<longueur1;i++){
  369.             B[i]=0;
  370.             for(int j=0;j<longueur;j++){
  371.                 B[i]=B[i]+(Jt[i][j]*(F[j]));
  372.  
  373.             }
  374.             C[i]=-B[i];
  375.         }
  376.  
  377.         printf("La matrice symetrique J*Jt est:");printf("\n");
  378.  
  379.  
  380.         for(e=0;e<longueur1;e++){
  381.             for(d=0;d<longueur1;d++){
  382.                 printf("%lf " ,A [e][d]);
  383.  
  384.             }
  385.             printf("\n");
  386.         }
  387.  
  388.  
  389.             L=cholesky(A,longueur1);
  390.  
  391.             Lt=transposee(L,longueur1,longueur1);
  392.  
  393.  
  394.             S=algo_solution(L,Lt,C,longueur1);
  395.  
  396.  
  397.             printf("S est:");printf("\n");
  398.             for(e=0;e<longueur1;e++){
  399.  
  400.                     printf("%lf " ,S [e]);
  401.             }
  402.             printf("\n");
  403.  
  404.             for(int k=0; k<longueur1;k++){
  405.  
  406.  
  407.                 v[k]=v[k]+S[k];
  408.             }
  409.             printf("v est:");
  410.             printf("\n");
  411.  
  412.  
  413.             for(d=0;d<longueur1;d++){
  414.  
  415.                     printf("%lf " ,v [d]);
  416.             }
  417.             printf("\n");
  418.  
  419.  
  420.     }while(sqrt(B[0]*B[0]+B[1]*B[1]+B[2]*B[2]+B[3]*B[3])>(pow(10,-6)));
  421.  
  422.  
  423.  
  424.  
  425.     return v;
  426. }
Última edición por mikelromero; 05/06/2012 a las 13:35