![]() ![]() | ![]() |
1.1 bertrand 1: /*
2: ================================================================================
1.4 ! bertrand 3: RPL/2 (R) version 4.0.12
1.1 bertrand 4: Copyright (C) 1989-2010 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 réalisation la factorisation LU d'une matrice quelconque
29: ================================================================================
30: Entrées : struct_matrice
31: --------------------------------------------------------------------------------
32: Sorties : décomposition A=PLU de la matrice d'entrée et drapeau d'erreur.
33: La matrice en entrée est écrasée. La matrice de sortie est réelle
34: si la matrice d'entrée est entière ou réelle, et complexe
35: si la matrice d'entrée est complexe.
36: La routine renvoie aussi une matrice d'entiers correspondant
37: à la matrice de permutation. Cette matrice est allouée par
38: la routine et vaut NULL sinon.
39: --------------------------------------------------------------------------------
40: Effets de bord : néant
41: ================================================================================
42: */
43:
44: void
45: factorisation_lu(struct_processus *s_etat_processus,
46: struct_matrice *s_matrice, struct_matrice **s_permutation)
47: {
48: integer4 dimension_vecteur_pivot;
49: integer4 erreur;
50: integer4 nombre_colonnes_a;
51: integer4 nombre_lignes_a;
52: integer4 *pivot;
53:
54: long n;
55:
56: unsigned long i;
57: unsigned long j;
58: unsigned long k;
59: unsigned long taille_matrice_f77;
60:
61: void *matrice_f77;
62: void *tampon;
63:
64: nombre_lignes_a = (integer4) (*s_matrice).nombre_lignes;
65: nombre_colonnes_a = (integer4) (*s_matrice).nombre_colonnes;
66: dimension_vecteur_pivot = (nombre_lignes_a < nombre_colonnes_a)
67: ? nombre_lignes_a : nombre_colonnes_a;
68: taille_matrice_f77 = nombre_lignes_a * nombre_colonnes_a;
69:
70: switch((*s_matrice).type)
71: {
72: case 'I' :
73: {
74: if ((matrice_f77 = (void *) malloc(taille_matrice_f77 *
75: sizeof(real8))) == NULL)
76: {
77: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
78: return;
79: }
80:
81: for(k = 0, i = 0; i < (*s_matrice).nombre_colonnes; i++)
82: {
83: for(j = 0; j < (*s_matrice).nombre_lignes; j++)
84: {
85: ((real8 *) matrice_f77)[k++] = ((integer8 **)
86: (*s_matrice).tableau)[j][i];
87: }
88: }
89:
90: for(i = 0; i < (*s_matrice).nombre_lignes; i++)
91: {
92: free(((integer8 **) (*s_matrice).tableau)[i]);
93: }
94:
95: free((integer8 **) (*s_matrice).tableau);
96:
97: if (((*s_matrice).tableau = (void **) malloc((*s_matrice)
98: .nombre_lignes * sizeof(real8 *))) == NULL)
99: {
100: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
101: return;
102: }
103:
104: for(i = 0; i < (*s_matrice).nombre_lignes; i++)
105: {
106: if ((((*s_matrice).tableau)[i] =
107: (real8 *) malloc((*s_matrice)
108: .nombre_colonnes * sizeof(real8))) == NULL)
109: {
110: (*s_etat_processus).erreur_systeme =
111: d_es_allocation_memoire;
112: return;
113: }
114: }
115:
116: (*s_matrice).type = 'R';
117:
118: if ((pivot = (integer4 *) malloc(dimension_vecteur_pivot *
119: sizeof(integer4))) == NULL)
120: {
121: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
122: return;
123: }
124:
125: dgetrf_(&nombre_lignes_a, &nombre_colonnes_a, matrice_f77,
126: &nombre_lignes_a, pivot, &erreur);
127:
128: if (erreur < 0)
129: {
130: (*s_etat_processus).erreur_execution =
131: d_ex_routines_mathematiques;
132:
133: free(pivot);
134: free(matrice_f77);
135: return;
136: }
137:
138: for(k = 0, i = 0; i < (*s_matrice).nombre_colonnes; i++)
139: {
140: for(j = 0; j < (*s_matrice).nombre_lignes; j++)
141: {
142: ((real8 **) (*s_matrice).tableau)[j][i] =
143: ((real8 *) matrice_f77)[k++];
144: }
145: }
146:
147: free(matrice_f77);
148:
149: (**s_permutation).type = 'I';
150: (**s_permutation).nombre_lignes = dimension_vecteur_pivot;
151: (**s_permutation).nombre_colonnes = dimension_vecteur_pivot;
152:
153: if (((**s_permutation).tableau = malloc((**s_permutation)
154: .nombre_lignes * sizeof(integer8 *))) == NULL)
155: {
156: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
157: return;
158: }
159:
160: for(i = 0; i < (**s_permutation).nombre_lignes; i++)
161: {
162: if ((((integer8 **) (**s_permutation).tableau)[i] =
163: (integer8 *) malloc((**s_permutation).nombre_colonnes *
164: sizeof(integer8))) == NULL)
165: {
166: (*s_etat_processus).erreur_systeme =
167: d_es_allocation_memoire;
168: return;
169: }
170:
171: for(j = 0; j < (**s_permutation).nombre_colonnes; j++)
172: {
173: ((integer8 **) (**s_permutation).tableau)[i][j] =
174: (j == i) ? 1 : 0;
175: }
176: }
177:
178: for(n = dimension_vecteur_pivot - 1; n >= 0; n--)
179: {
180: if ((pivot[n] - 1) != n)
181: {
182: tampon = ((integer8 **) (**s_permutation).tableau)[n];
183: ((integer8 **) (**s_permutation).tableau)[n] =
184: ((integer8 **) (**s_permutation).tableau)
185: [pivot[n] - 1];
186: ((integer8 **) (**s_permutation).tableau)[pivot[n] - 1] =
187: tampon;
188: }
189: }
190:
191: free(pivot);
192:
193: break;
194: }
195:
196: case 'R' :
197: {
198: if ((matrice_f77 = (void *) malloc(taille_matrice_f77 *
199: sizeof(real8))) == NULL)
200: {
201: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
202: return;
203: }
204:
205: for(k = 0, i = 0; i < (*s_matrice).nombre_colonnes; i++)
206: {
207: for(j = 0; j < (*s_matrice).nombre_lignes; j++)
208: {
209: ((real8 *) matrice_f77)[k++] = ((real8 **)
210: (*s_matrice).tableau)[j][i];
211: }
212: }
213:
214: if ((pivot = (integer4 *) malloc(dimension_vecteur_pivot *
215: sizeof(integer4))) == NULL)
216: {
217: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
218: return;
219: }
220:
221: dgetrf_(&nombre_lignes_a, &nombre_colonnes_a, matrice_f77,
222: &nombre_lignes_a, pivot, &erreur);
223:
224: if (erreur < 0)
225: {
226: (*s_etat_processus).erreur_execution =
227: d_ex_routines_mathematiques;
228:
229: free(pivot);
230: free(matrice_f77);
231: return;
232: }
233:
234: for(k = 0, i = 0; i < (*s_matrice).nombre_colonnes; i++)
235: {
236: for(j = 0; j < (*s_matrice).nombre_lignes; j++)
237: {
238: ((real8 **) (*s_matrice).tableau)[j][i] =
239: ((real8 *) matrice_f77)[k++];
240: }
241: }
242:
243: free(matrice_f77);
244:
245: (**s_permutation).type = 'I';
246: (**s_permutation).nombre_lignes = dimension_vecteur_pivot;
247: (**s_permutation).nombre_colonnes = dimension_vecteur_pivot;
248:
249: if (((**s_permutation).tableau = malloc((**s_permutation)
250: .nombre_lignes * sizeof(integer8 *))) == NULL)
251: {
252: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
253: return;
254: }
255:
256: for(i = 0; i < (**s_permutation).nombre_lignes; i++)
257: {
258: if ((((integer8 **) (**s_permutation).tableau)[i] =
259: (integer8 *) malloc((**s_permutation).nombre_colonnes *
260: sizeof(integer8))) == NULL)
261: {
262: (*s_etat_processus).erreur_systeme =
263: d_es_allocation_memoire;
264: return;
265: }
266:
267: for(j = 0; j < (**s_permutation).nombre_colonnes; j++)
268: {
269: ((integer8 **) (**s_permutation).tableau)[i][j] =
270: (j == i) ? 1 : 0;
271: }
272: }
273:
274: for(n = dimension_vecteur_pivot - 1; n >= 0; n--)
275: {
276: if ((pivot[n] - 1) != n)
277: {
278: tampon = ((integer8 **) (**s_permutation).tableau)[n];
279: ((integer8 **) (**s_permutation).tableau)[n] =
280: ((integer8 **) (**s_permutation).tableau)
281: [pivot[n] - 1];
282: ((integer8 **) (**s_permutation).tableau)[pivot[n] - 1] =
283: tampon;
284: }
285: }
286:
287: free(pivot);
288:
289: break;
290: }
291:
292: case 'C' :
293: {
294: if ((matrice_f77 = (void *) malloc(taille_matrice_f77 *
295: sizeof(complex16))) == NULL)
296: {
297: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
298: return;
299: }
300:
301: for(k = 0, i = 0; i < (*s_matrice).nombre_colonnes; i++)
302: {
303: for(j = 0; j < (*s_matrice).nombre_lignes; j++)
304: {
305: ((complex16 *) matrice_f77)[k++] = ((complex16 **)
306: (*s_matrice).tableau)[j][i];
307: }
308: }
309:
310: if ((pivot = (integer4 *) malloc(dimension_vecteur_pivot *
311: sizeof(integer4))) == NULL)
312: {
313: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
314: return;
315: }
316:
317: zgetrf_(&nombre_lignes_a, &nombre_colonnes_a, matrice_f77,
318: &nombre_lignes_a, pivot, &erreur);
319:
320: if (erreur < 0)
321: {
322: (*s_etat_processus).erreur_execution =
323: d_ex_routines_mathematiques;
324: free(pivot);
325: free(matrice_f77);
326: return;
327: }
328:
329: for(k = 0, i = 0; i < (*s_matrice).nombre_colonnes; i++)
330: {
331: for(j = 0; j < (*s_matrice).nombre_lignes; j++)
332: {
333: ((complex16 **) (*s_matrice).tableau)[j][i] =
334: ((complex16 *) matrice_f77)[k++];
335: }
336: }
337:
338: free(matrice_f77);
339:
340: (**s_permutation).type = 'I';
341: (**s_permutation).nombre_lignes = dimension_vecteur_pivot;
342: (**s_permutation).nombre_colonnes = dimension_vecteur_pivot;
343:
344: if (((**s_permutation).tableau = malloc((**s_permutation)
345: .nombre_lignes * sizeof(integer8 *))) == NULL)
346: {
347: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
348: return;
349: }
350:
351: for(i = 0; i < (**s_permutation).nombre_lignes; i++)
352: {
353: if ((((integer8 **) (**s_permutation).tableau)[i] =
354: (integer8 *) malloc((**s_permutation).nombre_colonnes *
355: sizeof(integer8))) == NULL)
356: {
357: (*s_etat_processus).erreur_systeme =
358: d_es_allocation_memoire;
359: return;
360: }
361:
362: for(j = 0; j < (**s_permutation).nombre_colonnes; j++)
363: {
364: ((integer8 **) (**s_permutation).tableau)[i][j] =
365: (j == i) ? 1 : 0;
366: }
367: }
368:
369: for(n = dimension_vecteur_pivot - 1; n >= 0; n--)
370: {
371: if ((pivot[n] - 1) != n)
372: {
373: tampon = ((integer8 **) (**s_permutation).tableau)[n];
374: ((integer8 **) (**s_permutation).tableau)[n] =
375: ((integer8 **) (**s_permutation).tableau)
376: [pivot[n] - 1];
377: ((integer8 **) (**s_permutation).tableau)[pivot[n] - 1] =
378: tampon;
379: }
380: }
381:
382: free(pivot);
383:
384: break;
385: }
386: }
387:
388: return;
389: }
390:
391:
392: /*
393: ================================================================================
394: Fonction réalisation la factorisation de Cholesky d'une matrice symétrique
395: définie positive
396: ================================================================================
397: Entrées : struct_matrice, mode (U ou L selon la décomposition demandée)
398: --------------------------------------------------------------------------------
399: Sorties : décomposition de Cholesky de la matrice d'entrée et drapeau
400: d'erreur. La matrice en entrée est écrasée.
401: --------------------------------------------------------------------------------
402: Effets de bord : néant
403: ================================================================================
404: */
405:
406: void
407: factorisation_cholesky(struct_processus *s_etat_processus,
408: struct_matrice *s_matrice, unsigned char mode)
409: {
410: integer4 erreur;
411: integer4 i;
412: integer4 j;
413: integer4 nombre_colonnes_a;
414: integer4 nombre_lignes_a;
415:
416: unsigned long taille_matrice_f77;
417:
418: void *matrice_f77;
419:
420: nombre_lignes_a = (integer4) (*s_matrice).nombre_lignes;
421: nombre_colonnes_a = (integer4) (*s_matrice).nombre_colonnes;
422:
423: if (nombre_lignes_a != nombre_colonnes_a)
424: {
425: (*s_etat_processus).erreur_execution = d_ex_dimensions_invalides;
426: return;
427: }
428:
429: taille_matrice_f77 = (nombre_lignes_a * (nombre_colonnes_a + 1)) / 2;
430:
431: switch((*s_matrice).type)
432: {
433: case 'I' :
434: {
435: /*
436: * Allocation du vecteur représentant la matrice triangulaire
437: */
438:
439: if ((matrice_f77 = (void *) malloc(taille_matrice_f77 *
440: sizeof(real8))) == NULL)
441: {
442: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
443: return;
444: }
445:
446: if (mode == 'U')
447: {
448: for(j = 1; j <= nombre_colonnes_a; j++)
449: {
450: for(i = 1; i <= j; i++)
451: {
452: ((real8 *) matrice_f77)[i + (((j - 1) * j) / 2) - 1] =
453: (real8) ((integer8 **) (*s_matrice).tableau)
454: [i - 1][j - 1];
455: }
456: }
457: }
458: else
459: {
460: for(j = 1; j <= nombre_colonnes_a; j++)
461: {
462: for(i = j; i <= nombre_colonnes_a; i++)
463: {
464: ((real8 *) matrice_f77)[(i - 1) + (((j - 1) *
465: ((nombre_colonnes_a * 2) - j)) / 2)] = (real8)
466: ((integer8 **) (*s_matrice).tableau)
467: [i - 1][j - 1];
468: }
469: }
470: }
471:
472: dpptrf_(&mode, &nombre_colonnes_a, matrice_f77, &erreur, 1);
473:
474: if (erreur != 0)
475: {
476: if (erreur > 0)
477: {
478: (*s_etat_processus).exception =
479: d_ep_matrice_non_definie_positive;
480: }
481: else
482: {
483: (*s_etat_processus).erreur_execution =
484: d_ex_routines_mathematiques;
485: }
486:
487: free(matrice_f77);
488: return;
489: }
490:
491: for(i = 0; i < (integer4) (*s_matrice).nombre_lignes; i++)
492: {
493: free(((integer8 **) (*s_matrice).tableau)[i]);
494: }
495:
496: free((integer8 **) (*s_matrice).tableau);
497:
498: if (((*s_matrice).tableau = (void **) malloc((*s_matrice)
499: .nombre_lignes * sizeof(real8 *))) == NULL)
500: {
501: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
502: return;
503: }
504:
505: for(i = 0; i < (integer4) (*s_matrice).nombre_lignes; i++)
506: {
507: if ((((*s_matrice).tableau)[i] =
508: (real8 *) malloc((*s_matrice)
509: .nombre_colonnes * sizeof(real8))) == NULL)
510: {
511: (*s_etat_processus).erreur_systeme =
512: d_es_allocation_memoire;
513: return;
514: }
515: }
516:
517: (*s_matrice).type = 'R';
518:
519: if (mode == 'U')
520: {
521: for(j = 1; j <= nombre_colonnes_a; j++)
522: {
523: for(i = 1; i <= j; i++)
524: {
525: ((real8 **) (*s_matrice).tableau)[i - 1][j - 1] =
526: ((real8 *) matrice_f77)
527: [i + (((j - 1) * j) / 2) - 1];
528: }
529: }
530:
531: for(j = 1; j <= nombre_colonnes_a; j++)
532: {
533: for(i = j + 1; i <= nombre_colonnes_a; i++)
534: {
535: ((real8 **) (*s_matrice).tableau)[i - 1][j - 1] = 0;
536: }
537: }
538: }
539: else
540: {
541: for(j = 1; j <= nombre_colonnes_a; j++)
542: {
543: for(i = j; i <= nombre_colonnes_a; i++)
544: {
545: ((real8 **) (*s_matrice).tableau)[i - 1][j - 1] =
546: ((real8 *) matrice_f77)[(i - 1) + (((j - 1) *
547: ((nombre_colonnes_a * 2) - j)) / 2)];
548: }
549: }
550:
551: for(j = 1; j <= nombre_colonnes_a; j++)
552: {
553: for(i = 1; i < j; i++)
554: {
555: ((real8 **) (*s_matrice).tableau)[i - 1][j - 1] = 0;
556: }
557: }
558: }
559:
560: free(matrice_f77);
561: break;
562: }
563:
564: case 'R' :
565: {
566: /*
567: * Allocation du vecteur représentant la matrice triangulaire
568: */
569:
570: if ((matrice_f77 = (void *) malloc(taille_matrice_f77 *
571: sizeof(real8))) == NULL)
572: {
573: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
574: return;
575: }
576:
577: if (mode == 'U')
578: {
579: for(j = 1; j <= nombre_colonnes_a; j++)
580: {
581: for(i = 1; i <= j; i++)
582: {
583: ((real8 *) matrice_f77)[i + (((j - 1) * j) / 2) - 1] =
584: ((real8 **) (*s_matrice).tableau)[i - 1][j - 1];
585: }
586: }
587: }
588: else
589: {
590: for(j = 1; j <= nombre_colonnes_a; j++)
591: {
592: for(i = j; i <= nombre_colonnes_a; i++)
593: {
594: ((real8 *) matrice_f77)[(i - 1) + (((j - 1) *
595: ((nombre_colonnes_a * 2) - j)) / 2)] = (real8)
596: ((real8 **) (*s_matrice).tableau)[i - 1][j - 1];
597: }
598: }
599: }
600:
601: dpptrf_(&mode, &nombre_colonnes_a, matrice_f77, &erreur, 1);
602:
603: if (erreur != 0)
604: {
605: if (erreur > 0)
606: {
607: (*s_etat_processus).exception =
608: d_ep_matrice_non_definie_positive;
609: }
610: else
611: {
612: (*s_etat_processus).erreur_execution =
613: d_ex_routines_mathematiques;
614: }
615:
616: free(matrice_f77);
617: return;
618: }
619:
620: if (mode == 'U')
621: {
622: for(j = 1; j <= nombre_colonnes_a; j++)
623: {
624: for(i = 1; i <= j; i++)
625: {
626: ((real8 **) (*s_matrice).tableau)[i - 1][j - 1] =
627: ((real8 *) matrice_f77)
628: [i + (((j - 1) * j) / 2) - 1];
629: }
630: }
631:
632: for(j = 1; j <= nombre_colonnes_a; j++)
633: {
634: for(i = j + 1; i <= nombre_colonnes_a; i++)
635: {
636: ((real8 **) (*s_matrice).tableau)[i - 1][j - 1] = 0;
637: }
638: }
639: }
640: else
641: {
642: for(j = 1; j <= nombre_colonnes_a; j++)
643: {
644: for(i = j; i <= nombre_colonnes_a; i++)
645: {
646: ((real8 **) (*s_matrice).tableau)[i - 1][j - 1] =
647: ((real8 *) matrice_f77)[(i - 1) + (((j - 1) *
648: ((nombre_colonnes_a * 2) - j)) / 2)];
649: }
650: }
651:
652: for(j = 1; j <= nombre_colonnes_a; j++)
653: {
654: for(i = 1; i < j; i++)
655: {
656: ((real8 **) (*s_matrice).tableau)[i - 1][j - 1] = 0;
657: }
658: }
659: }
660:
661: free(matrice_f77);
662: break;
663: }
664:
665: case 'C' :
666: {
667: /*
668: * Allocation du vecteur représentant la matrice triangulaire
669: */
670:
671: if ((matrice_f77 = (void *) malloc(taille_matrice_f77 *
672: sizeof(complex16))) == NULL)
673: {
674: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
675: return;
676: }
677:
678: if (mode == 'U')
679: {
680: for(j = 1; j <= nombre_colonnes_a; j++)
681: {
682: for(i = 1; i <= j; i++)
683: {
684: ((complex16 *) matrice_f77)[i + (((j - 1) * j) / 2)
685: - 1].partie_reelle = ((complex16 **)
686: (*s_matrice).tableau)[i - 1][j - 1]
687: .partie_reelle;
688: ((complex16 *) matrice_f77)[i + (((j - 1) * j) / 2)
689: - 1].partie_imaginaire = ((complex16 **)
690: (*s_matrice).tableau)[i - 1][j - 1]
691: .partie_imaginaire;
692: }
693: }
694: }
695: else
696: {
697: for(j = 1; j <= nombre_colonnes_a; j++)
698: {
699: for(i = j; i <= nombre_colonnes_a; i++)
700: {
701: ((complex16 *) matrice_f77)[(i - 1) + (((j - 1) *
702: ((nombre_colonnes_a * 2) - j)) / 2)]
703: .partie_reelle = ((complex16 **)
704: (*s_matrice).tableau)[i - 1][j - 1]
705: .partie_reelle;
706: ((complex16 *) matrice_f77)[(i - 1) + (((j - 1) *
707: ((nombre_colonnes_a * 2) - j)) / 2)]
708: .partie_imaginaire = ((complex16 **)
709: (*s_matrice).tableau)[i - 1][j - 1]
710: .partie_imaginaire;
711: }
712: }
713: }
714:
715: zpptrf_(&mode, &nombre_colonnes_a, matrice_f77, &erreur, 1);
716:
717: if (erreur != 0)
718: {
719: if (erreur > 0)
720: {
721: (*s_etat_processus).exception =
722: d_ep_matrice_non_definie_positive;
723: }
724: else
725: {
726: (*s_etat_processus).erreur_execution =
727: d_ex_routines_mathematiques;
728: }
729:
730: free(matrice_f77);
731: return;
732: }
733:
734: if (mode == 'U')
735: {
736: for(j = 1; j <= nombre_colonnes_a; j++)
737: {
738: for(i = 1; i <= j; i++)
739: {
740: ((complex16 **) (*s_matrice).tableau)[i - 1][j - 1]
741: .partie_reelle = ((complex16 *) matrice_f77)
742: [i + (((j - 1) * j) / 2) - 1].partie_reelle;
743: ((complex16 **) (*s_matrice).tableau)[i - 1][j - 1]
744: .partie_imaginaire = ((complex16 *) matrice_f77)
745: [i + (((j - 1) * j) / 2) - 1].partie_imaginaire;
746: }
747: }
748:
749: for(j = 1; j <= nombre_colonnes_a; j++)
750: {
751: for(i = j + 1; i <= nombre_colonnes_a; i++)
752: {
753: ((complex16 **) (*s_matrice).tableau)[i - 1][j - 1]
754: .partie_reelle = 0;
755: ((complex16 **) (*s_matrice).tableau)[i - 1][j - 1]
756: .partie_imaginaire = 0;
757: }
758: }
759: }
760: else
761: {
762: for(j = 1; j <= nombre_colonnes_a; j++)
763: {
764: for(i = j; i <= nombre_colonnes_a; i++)
765: {
766: ((complex16 **) (*s_matrice).tableau)[i - 1][j - 1]
767: .partie_reelle = ((complex16 *)
768: matrice_f77)[(i - 1) + (((j - 1) *
769: ((nombre_colonnes_a * 2) - j)) / 2)]
770: .partie_reelle;
771: ((complex16 **) (*s_matrice).tableau)[i - 1][j - 1]
772: .partie_imaginaire = ((complex16 *)
773: matrice_f77)[(i - 1) + (((j - 1) *
774: ((nombre_colonnes_a * 2) - j)) / 2)]
775: .partie_imaginaire;
776: }
777: }
778:
779: for(j = 1; j <= nombre_colonnes_a; j++)
780: {
781: for(i = 1; i < j; i++)
782: {
783: ((complex16 **) (*s_matrice).tableau)[i - 1][j - 1]
784: .partie_reelle = 0;
785: ((complex16 **) (*s_matrice).tableau)[i - 1][j - 1]
786: .partie_imaginaire = 0;
787: }
788: }
789: }
790:
791: free(matrice_f77);
792: break;
793: }
794: }
795:
796: return;
797: }
798:
799: // vim: ts=4