![]() ![]() | ![]() |
Passage de la branche 4.1 en branche stable.
1: /* 2: ================================================================================ 3: RPL/2 (R) version 4.1.0 4: Copyright (C) 1989-2011 Dr. BERTRAND Joël 5: 6: This file is part of RPL/2. 7: 8: RPL/2 is free software; you can redistribute it and/or modify it 9: under the terms of the CeCILL V2 License as published by the french 10: CEA, CNRS and INRIA. 11: 12: RPL/2 is distributed in the hope that it will be useful, but WITHOUT 13: ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or 14: FITNESS FOR A PARTICULAR PURPOSE. See the CeCILL V2 License 15: for more details. 16: 17: You should have received a copy of the CeCILL License 18: along with RPL/2. If not, write to info@cecill.info. 19: ================================================================================ 20: */ 21: 22: 23: #include "rpl-conv.h" 24: 25: 26: /* 27: ================================================================================ 28: Fonction 'qr' 29: ================================================================================ 30: Entrées : pointeur sur une structure struct_processus 31: -------------------------------------------------------------------------------- 32: Sorties : 33: -------------------------------------------------------------------------------- 34: Effets de bord : néant 35: ================================================================================ 36: */ 37: 38: void 39: instruction_qr(struct_processus *s_etat_processus) 40: { 41: complex16 registre; 42: complex16 *tau_complexe; 43: complex16 *vecteur_complexe; 44: 45: real8 *tau_reel; 46: real8 *vecteur_reel; 47: 48: struct_liste_chainee *registre_pile_last; 49: 50: struct_objet *s_copie_argument; 51: struct_objet *s_matrice_identite; 52: struct_objet *s_objet; 53: struct_objet *s_objet_argument; 54: struct_objet *s_objet_resultat; 55: 56: unsigned long i; 57: unsigned long j; 58: unsigned long k; 59: unsigned long nombre_reflecteurs_elementaires; 60: 61: void *tau; 62: 63: (*s_etat_processus).erreur_execution = d_ex; 64: 65: if ((*s_etat_processus).affichage_arguments == 'Y') 66: { 67: printf("\n QR "); 68: 69: if ((*s_etat_processus).langue == 'F') 70: { 71: printf("(décomposition QR)\n\n"); 72: } 73: else 74: { 75: printf("(QR décomposition)\n\n"); 76: } 77: 78: printf(" 1: %s, %s\n", d_MIN, d_MRL); 79: printf("-> 2: %s\n", d_MRL); 80: printf(" 1: %s\n\n", d_MRL); 81: 82: printf(" 1: %s\n", d_MCX); 83: printf("-> 2: %s\n", d_MCX); 84: printf(" 1: %s\n", d_MCX); 85: 86: return; 87: } 88: else if ((*s_etat_processus).test_instruction == 'Y') 89: { 90: (*s_etat_processus).nombre_arguments = -1; 91: return; 92: } 93: 94: if (test_cfsf(s_etat_processus, 31) == d_vrai) 95: { 96: if (empilement_pile_last(s_etat_processus, 1) == d_erreur) 97: { 98: return; 99: } 100: } 101: 102: if (depilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 103: &s_objet_argument) == d_erreur) 104: { 105: (*s_etat_processus).erreur_execution = d_ex_manque_argument; 106: return; 107: } 108: 109: if (((*s_objet_argument).type == MIN) || 110: ((*s_objet_argument).type == MRL)) 111: { 112: /* 113: * Matrice entière ou réelle 114: */ 115: 116: if ((s_copie_argument = copie_objet(s_etat_processus, 117: s_objet_argument, 'Q')) == NULL) 118: { 119: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 120: return; 121: } 122: 123: factorisation_qr(s_etat_processus, (*s_copie_argument).objet, &tau); 124: (*s_copie_argument).type = MRL; 125: 126: tau_reel = (real8 *) tau; 127: 128: if ((*s_etat_processus).erreur_systeme != d_es) 129: { 130: return; 131: } 132: 133: if (((*s_etat_processus).exception != d_ep) || 134: ((*s_etat_processus).erreur_execution != d_ex)) 135: { 136: free(tau); 137: liberation(s_etat_processus, s_objet_argument); 138: liberation(s_etat_processus, s_copie_argument); 139: return; 140: } 141: 142: if ((s_objet_resultat = copie_objet(s_etat_processus, 143: s_copie_argument, 'O')) == NULL) 144: { 145: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 146: return; 147: } 148: 149: // Matrice Q 150: 151: nombre_reflecteurs_elementaires = ((*((struct_matrice *) 152: (*s_copie_argument).objet)).nombre_colonnes < 153: (*((struct_matrice *) (*s_copie_argument).objet)) 154: .nombre_lignes) ? (*((struct_matrice *) 155: (*s_copie_argument).objet)).nombre_colonnes 156: : (*((struct_matrice *) (*s_copie_argument).objet)) 157: .nombre_lignes; 158: 159: registre_pile_last = NULL; 160: 161: if (test_cfsf(s_etat_processus, 31) == d_vrai) 162: { 163: registre_pile_last = (*s_etat_processus).l_base_pile_last; 164: (*s_etat_processus).l_base_pile_last = NULL; 165: } 166: 167: if ((s_objet = allocation(s_etat_processus, INT)) == NULL) 168: { 169: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 170: return; 171: } 172: 173: (*((integer8 *) (*s_objet).objet)) = (*((struct_matrice *) 174: (*s_copie_argument).objet)).nombre_lignes; 175: 176: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 177: s_objet) == d_erreur) 178: { 179: return; 180: } 181: 182: instruction_idn(s_etat_processus); 183: 184: if (((*s_etat_processus).erreur_systeme != d_es) || 185: ((*s_etat_processus).erreur_execution != d_ex) || 186: ((*s_etat_processus).exception != d_ep)) 187: { 188: liberation(s_etat_processus, s_copie_argument); 189: free(tau); 190: 191: if (empilement_pile_last(s_etat_processus, 0) == d_erreur) 192: { 193: return; 194: } 195: 196: (*s_etat_processus).l_base_pile_last = registre_pile_last; 197: return; 198: } 199: 200: if (depilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 201: &s_matrice_identite) == d_erreur) 202: { 203: (*s_etat_processus).erreur_execution = d_ex_manque_argument; 204: return; 205: } 206: 207: for(i = 0; i < nombre_reflecteurs_elementaires; i++) 208: { 209: // Calcul de H(i) = I - tau * v * v' 210: 211: if ((s_objet = copie_objet(s_etat_processus, s_matrice_identite, 212: 'P')) == NULL) 213: { 214: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 215: return; 216: } 217: 218: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 219: s_objet) == d_erreur) 220: { 221: return; 222: } 223: 224: if ((s_objet = allocation(s_etat_processus, REL)) == NULL) 225: { 226: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 227: return; 228: } 229: 230: (*((real8 *) (*s_objet).objet)) = tau_reel[i]; 231: 232: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 233: s_objet) == d_erreur) 234: { 235: return; 236: } 237: 238: if ((s_objet = allocation(s_etat_processus, MRL)) == NULL) 239: { 240: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 241: return; 242: } 243: 244: (*((struct_matrice *) (*s_objet).objet)).nombre_lignes = 245: (*((struct_matrice *) (*s_copie_argument).objet)) 246: .nombre_lignes; 247: (*((struct_matrice *) (*s_objet).objet)).nombre_colonnes = 248: (*((struct_matrice *) (*s_copie_argument).objet)) 249: .nombre_lignes; 250: 251: if ((vecteur_reel = malloc((*((struct_matrice *) (*s_objet).objet)) 252: .nombre_lignes * sizeof(real8))) == NULL) 253: { 254: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 255: return; 256: } 257: 258: for(j = 0; j < (*((struct_matrice *) (*s_objet).objet)) 259: .nombre_lignes; j++) 260: { 261: if (j < i) 262: { 263: vecteur_reel[j] = 0; 264: } 265: else if (j == i) 266: { 267: vecteur_reel[j] = 1; 268: } 269: else 270: { 271: vecteur_reel[j] = ((real8 **) (*((struct_matrice *) 272: (*s_copie_argument).objet)).tableau)[j][i]; 273: } 274: } 275: 276: if (((*((struct_matrice *) (*s_objet).objet)).tableau = 277: malloc((*((struct_matrice *) (*s_objet).objet)) 278: .nombre_lignes * sizeof(real8 *))) == NULL) 279: { 280: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 281: return; 282: } 283: 284: for(j = 0; j < (*((struct_matrice *) (*s_objet).objet)) 285: .nombre_lignes; j++) 286: { 287: if ((((real8 **) (*((struct_matrice *) (*s_objet).objet)) 288: .tableau)[j] = malloc((*((struct_matrice *) (*s_objet) 289: .objet)).nombre_lignes * sizeof(real8))) == NULL) 290: { 291: (*s_etat_processus).erreur_systeme = 292: d_es_allocation_memoire; 293: return; 294: } 295: 296: for(k = 0; k < (*((struct_matrice *) (*s_objet).objet)) 297: .nombre_colonnes; k++) 298: { 299: ((real8 **) (*((struct_matrice *) (*s_objet).objet)) 300: .tableau)[j][k] = vecteur_reel[j] * vecteur_reel[k]; 301: } 302: } 303: 304: free(vecteur_reel); 305: 306: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 307: s_objet) == d_erreur) 308: { 309: return; 310: } 311: 312: instruction_multiplication(s_etat_processus); 313: 314: if (((*s_etat_processus).erreur_systeme != d_es) || 315: ((*s_etat_processus).erreur_execution != d_ex) || 316: ((*s_etat_processus).exception != d_ep)) 317: { 318: liberation(s_etat_processus, s_copie_argument); 319: liberation(s_etat_processus, s_matrice_identite); 320: free(tau); 321: 322: if (empilement_pile_last(s_etat_processus, 0) == d_erreur) 323: { 324: return; 325: } 326: 327: (*s_etat_processus).l_base_pile_last = registre_pile_last; 328: return; 329: } 330: 331: instruction_moins(s_etat_processus); 332: 333: if (((*s_etat_processus).erreur_systeme != d_es) || 334: ((*s_etat_processus).erreur_execution != d_ex) || 335: ((*s_etat_processus).exception != d_ep)) 336: { 337: liberation(s_etat_processus, s_copie_argument); 338: liberation(s_etat_processus, s_matrice_identite); 339: free(tau); 340: 341: if (empilement_pile_last(s_etat_processus, 0) == d_erreur) 342: { 343: return; 344: } 345: 346: (*s_etat_processus).l_base_pile_last = registre_pile_last; 347: return; 348: } 349: 350: if (i > 0) 351: { 352: instruction_multiplication(s_etat_processus); 353: 354: if (((*s_etat_processus).erreur_systeme != d_es) || 355: ((*s_etat_processus).erreur_execution != d_ex) || 356: ((*s_etat_processus).exception != d_ep)) 357: { 358: liberation(s_etat_processus, s_copie_argument); 359: liberation(s_etat_processus, s_matrice_identite); 360: free(tau); 361: 362: if (empilement_pile_last(s_etat_processus, 0) == d_erreur) 363: { 364: return; 365: } 366: 367: (*s_etat_processus).l_base_pile_last = registre_pile_last; 368: return; 369: } 370: } 371: } 372: 373: if (test_cfsf(s_etat_processus, 31) == d_vrai) 374: { 375: if (empilement_pile_last(s_etat_processus, 0) == d_erreur) 376: { 377: return; 378: } 379: 380: (*s_etat_processus).l_base_pile_last = registre_pile_last; 381: } 382: 383: // Matrice R 384: 385: for(i = 0; i < (*((struct_matrice *) (*s_objet_resultat).objet)) 386: .nombre_lignes; i++) 387: { 388: for(j = 0; j < i; j++) 389: { 390: ((real8 **) (*((struct_matrice *) (*s_objet_resultat).objet)) 391: .tableau)[i][j] = 0; 392: } 393: } 394: 395: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 396: s_objet_resultat) == d_erreur) 397: { 398: return; 399: } 400: 401: liberation(s_etat_processus, s_matrice_identite); 402: liberation(s_etat_processus, s_copie_argument); 403: free(tau); 404: } 405: else if ((*s_objet_argument).type == MCX) 406: { 407: /* 408: * Matrice complexe 409: */ 410: 411: if ((s_copie_argument = copie_objet(s_etat_processus, 412: s_objet_argument, 'Q')) == NULL) 413: { 414: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 415: return; 416: } 417: 418: factorisation_qr(s_etat_processus, (*s_copie_argument).objet, &tau); 419: 420: tau_complexe = (complex16 *) tau; 421: 422: if ((*s_etat_processus).erreur_systeme != d_es) 423: { 424: return; 425: } 426: 427: if (((*s_etat_processus).exception != d_ep) || 428: ((*s_etat_processus).erreur_execution != d_ex)) 429: { 430: free(tau); 431: liberation(s_etat_processus, s_objet_argument); 432: liberation(s_etat_processus, s_copie_argument); 433: return; 434: } 435: 436: if ((s_objet_resultat = copie_objet(s_etat_processus, 437: s_copie_argument, 'O')) == NULL) 438: { 439: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 440: return; 441: } 442: 443: // Matrice Q 444: 445: nombre_reflecteurs_elementaires = ((*((struct_matrice *) 446: (*s_copie_argument).objet)).nombre_colonnes < 447: (*((struct_matrice *) (*s_copie_argument).objet)) 448: .nombre_lignes) ? (*((struct_matrice *) 449: (*s_copie_argument).objet)).nombre_colonnes 450: : (*((struct_matrice *) (*s_copie_argument).objet)) 451: .nombre_lignes; 452: 453: registre_pile_last = NULL; 454: 455: if (test_cfsf(s_etat_processus, 31) == d_vrai) 456: { 457: registre_pile_last = (*s_etat_processus).l_base_pile_last; 458: (*s_etat_processus).l_base_pile_last = NULL; 459: } 460: 461: if ((s_objet = allocation(s_etat_processus, INT)) == NULL) 462: { 463: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 464: return; 465: } 466: 467: (*((integer8 *) (*s_objet).objet)) = (*((struct_matrice *) 468: (*s_copie_argument).objet)).nombre_lignes; 469: 470: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 471: s_objet) == d_erreur) 472: { 473: return; 474: } 475: 476: instruction_idn(s_etat_processus); 477: 478: if (((*s_etat_processus).erreur_systeme != d_es) || 479: ((*s_etat_processus).erreur_execution != d_ex) || 480: ((*s_etat_processus).exception != d_ep)) 481: { 482: liberation(s_etat_processus, s_copie_argument); 483: free(tau); 484: 485: if (empilement_pile_last(s_etat_processus, 0) == d_erreur) 486: { 487: return; 488: } 489: 490: (*s_etat_processus).l_base_pile_last = registre_pile_last; 491: return; 492: } 493: 494: if (depilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 495: &s_matrice_identite) == d_erreur) 496: { 497: (*s_etat_processus).erreur_execution = d_ex_manque_argument; 498: return; 499: } 500: 501: for(i = 0; i < nombre_reflecteurs_elementaires; i++) 502: { 503: // Calcul de H(i) = I - tau * v * v' 504: 505: if ((s_objet = copie_objet(s_etat_processus, s_matrice_identite, 506: 'P')) == NULL) 507: { 508: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 509: return; 510: } 511: 512: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 513: s_objet) == d_erreur) 514: { 515: return; 516: } 517: 518: if ((s_objet = allocation(s_etat_processus, CPL)) == NULL) 519: { 520: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 521: return; 522: } 523: 524: (*((complex16 *) (*s_objet).objet)) = tau_complexe[i]; 525: 526: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 527: s_objet) == d_erreur) 528: { 529: return; 530: } 531: 532: if ((s_objet = allocation(s_etat_processus, MCX)) == NULL) 533: { 534: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 535: return; 536: } 537: 538: (*((struct_matrice *) (*s_objet).objet)).nombre_lignes = 539: (*((struct_matrice *) (*s_copie_argument).objet)) 540: .nombre_lignes; 541: (*((struct_matrice *) (*s_objet).objet)).nombre_colonnes = 542: (*((struct_matrice *) (*s_copie_argument).objet)) 543: .nombre_lignes; 544: 545: if ((vecteur_complexe = malloc((*((struct_matrice *) 546: (*s_objet).objet)).nombre_lignes * sizeof(complex16))) 547: == NULL) 548: { 549: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 550: return; 551: } 552: 553: for(j = 0; j < (*((struct_matrice *) (*s_objet).objet)) 554: .nombre_lignes; j++) 555: { 556: if (j < i) 557: { 558: vecteur_complexe[j].partie_reelle = 0; 559: vecteur_complexe[j].partie_imaginaire = 0; 560: } 561: else if (j == i) 562: { 563: vecteur_complexe[j].partie_reelle = 1; 564: vecteur_complexe[j].partie_imaginaire = 0; 565: } 566: else 567: { 568: vecteur_complexe[j].partie_reelle = ((complex16 **) 569: (*((struct_matrice *) (*s_copie_argument).objet)) 570: .tableau)[j][i].partie_reelle; 571: vecteur_complexe[j].partie_imaginaire = ((complex16 **) 572: (*((struct_matrice *) (*s_copie_argument).objet)) 573: .tableau)[j][i].partie_imaginaire; 574: } 575: } 576: 577: if (((*((struct_matrice *) (*s_objet).objet)).tableau = 578: malloc((*((struct_matrice *) (*s_objet).objet)) 579: .nombre_lignes * sizeof(complex16 *))) == NULL) 580: { 581: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire; 582: return; 583: } 584: 585: for(j = 0; j < (*((struct_matrice *) (*s_objet).objet)) 586: .nombre_lignes; j++) 587: { 588: if ((((complex16 **) (*((struct_matrice *) (*s_objet).objet)) 589: .tableau)[j] = malloc((*((struct_matrice *) (*s_objet) 590: .objet)).nombre_lignes * sizeof(complex16))) == NULL) 591: { 592: (*s_etat_processus).erreur_systeme = 593: d_es_allocation_memoire; 594: return; 595: } 596: 597: for(k = 0; k < (*((struct_matrice *) (*s_objet).objet)) 598: .nombre_colonnes; k++) 599: { 600: registre = vecteur_complexe[k]; 601: registre.partie_imaginaire = 602: -vecteur_complexe[k].partie_imaginaire; 603: 604: f77multiplicationcc_(&(vecteur_complexe[j]), 605: ®istre, &(((complex16 **) 606: (*((struct_matrice *) (*s_objet).objet)).tableau) 607: [j][k])); 608: } 609: } 610: 611: free(vecteur_complexe); 612: 613: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 614: s_objet) == d_erreur) 615: { 616: return; 617: } 618: 619: instruction_multiplication(s_etat_processus); 620: 621: if (((*s_etat_processus).erreur_systeme != d_es) || 622: ((*s_etat_processus).erreur_execution != d_ex) || 623: ((*s_etat_processus).exception != d_ep)) 624: { 625: liberation(s_etat_processus, s_copie_argument); 626: liberation(s_etat_processus, s_matrice_identite); 627: free(tau); 628: 629: if (empilement_pile_last(s_etat_processus, 0) == d_erreur) 630: { 631: return; 632: } 633: 634: (*s_etat_processus).l_base_pile_last = registre_pile_last; 635: return; 636: } 637: 638: instruction_moins(s_etat_processus); 639: 640: if (((*s_etat_processus).erreur_systeme != d_es) || 641: ((*s_etat_processus).erreur_execution != d_ex) || 642: ((*s_etat_processus).exception != d_ep)) 643: { 644: liberation(s_etat_processus, s_copie_argument); 645: liberation(s_etat_processus, s_matrice_identite); 646: free(tau); 647: 648: if (empilement_pile_last(s_etat_processus, 0) == d_erreur) 649: { 650: return; 651: } 652: 653: (*s_etat_processus).l_base_pile_last = registre_pile_last; 654: return; 655: } 656: 657: if (i > 0) 658: { 659: instruction_multiplication(s_etat_processus); 660: 661: if (((*s_etat_processus).erreur_systeme != d_es) || 662: ((*s_etat_processus).erreur_execution != d_ex) || 663: ((*s_etat_processus).exception != d_ep)) 664: { 665: liberation(s_etat_processus, s_copie_argument); 666: liberation(s_etat_processus, s_matrice_identite); 667: free(tau); 668: 669: if (empilement_pile_last(s_etat_processus, 0) == d_erreur) 670: { 671: return; 672: } 673: 674: (*s_etat_processus).l_base_pile_last = registre_pile_last; 675: return; 676: } 677: } 678: } 679: 680: if (test_cfsf(s_etat_processus, 31) == d_vrai) 681: { 682: if (empilement_pile_last(s_etat_processus, 0) == d_erreur) 683: { 684: return; 685: } 686: 687: (*s_etat_processus).l_base_pile_last = registre_pile_last; 688: } 689: 690: // Matrice R 691: 692: for(i = 0; i < (*((struct_matrice *) (*s_objet_resultat).objet)) 693: .nombre_lignes; i++) 694: { 695: for(j = 0; j < i; j++) 696: { 697: ((complex16 **) (*((struct_matrice *) (*s_objet_resultat) 698: .objet)).tableau)[i][j].partie_reelle = 0; 699: ((complex16 **) (*((struct_matrice *) (*s_objet_resultat) 700: .objet)).tableau)[i][j].partie_imaginaire = 0; 701: } 702: } 703: 704: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile), 705: s_objet_resultat) == d_erreur) 706: { 707: return; 708: } 709: 710: liberation(s_etat_processus, s_matrice_identite); 711: liberation(s_etat_processus, s_copie_argument); 712: free(tau); 713: } 714: 715: /* 716: * Type d'argument invalide 717: */ 718: 719: else 720: { 721: liberation(s_etat_processus, s_objet_argument); 722: 723: (*s_etat_processus).erreur_execution = d_ex_erreur_type_argument; 724: return; 725: } 726: 727: liberation(s_etat_processus, s_objet_argument); 728: 729: return; 730: } 731: 732: // vim: ts=4