1: /*
2: ================================================================================
3: RPL/2 (R) version 4.0.10
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 le tri d'un vecteur de réels du plus petit au plus
29: grand en valeur absolue (algorithme de tri dit de Shell-Metzner)
30: ================================================================================
31: Entrées : pointeur sur une structure struct_processus
32: --------------------------------------------------------------------------------
33: Sorties :
34: --------------------------------------------------------------------------------
35: Effets de bord : néant
36: ================================================================================
37: */
38:
39: void
40: tri_vecteur(real8 *vecteur, unsigned long taille)
41: {
42: logical1 terminaison_boucle_1;
43: logical1 terminaison_boucle_2;
44: logical1 terminaison_boucle_3;
45:
46: signed long indice_i;
47: signed long indice_j;
48: signed long indice_k;
49: signed long indice_l;
50:
51: unsigned long ecartement;
52:
53: ecartement = taille;
54: terminaison_boucle_1 = d_faux;
55:
56: do
57: {
58: ecartement = ecartement / 2;
59:
60: if (ecartement >= 1)
61: {
62: indice_j = 0;
63: indice_k = taille - ecartement;
64: terminaison_boucle_2 = d_faux;
65:
66: do
67: {
68: indice_i = indice_j;
69: terminaison_boucle_3 = d_faux;
70:
71: do
72: {
73: indice_l = indice_i + ecartement;
74:
75: if ((indice_i > 0) && (indice_l > 0))
76: {
77: if (fabs(vecteur[indice_i - 1]) >
78: fabs(vecteur[indice_l - 1]))
79: {
80: swap((void *) &(vecteur[indice_i - 1]),
81: (void *) &(vecteur[indice_l - 1]),
82: sizeof(real8));
83:
84: indice_i -= ecartement;
85:
86: if (indice_i < 1)
87: {
88: terminaison_boucle_3 = d_vrai;
89: }
90: }
91: else
92: {
93: terminaison_boucle_3 = d_vrai;
94: }
95: }
96: else
97: {
98: terminaison_boucle_3 = d_vrai;
99: }
100: } while(terminaison_boucle_3 == d_faux);
101:
102: indice_j++;
103:
104: if (indice_j > indice_k)
105: {
106: terminaison_boucle_2 = d_vrai;
107: }
108: } while(terminaison_boucle_2 == d_faux);
109: }
110: else
111: {
112: terminaison_boucle_1 = d_vrai;
113: }
114: } while(terminaison_boucle_1 == d_faux);
115: }
116:
117:
118: /*
119: ================================================================================
120: Fonction réalisation la sommation d'un vecteur de réel sans perte
121: de précision
122: ================================================================================
123: Entrées : pointeur sur une structure struct_processus
124: --------------------------------------------------------------------------------
125: Sorties :
126: --------------------------------------------------------------------------------
127: Effets de bord : néant
128: ================================================================================
129: */
130:
131: real8
132: sommation_vecteur_reel(real8 *vecteur, unsigned long *taille,
133: logical1 *erreur_memoire)
134: {
135: #if 0
136: unsigned long nombre_elements;
137: unsigned long pointeur;
138:
139: /*
140: * Sommation des termes en commençant par le plus petit.
141: * Algorithme optimal mais NP-complet...
142: */
143:
144: nombre_elements = (*taille);
145: (*erreur_memoire) = d_faux;
146:
147: while(nombre_elements != 1)
148: {
149: pointeur = (*taille) - nombre_elements;
150: tri_vecteur(&(vecteur[pointeur]), nombre_elements);
151: vecteur[pointeur + 1] += vecteur[pointeur];
152: nombre_elements--;
153: }
154:
155: return(vecteur[(*taille) - 1]);
156: #else
157: real8 erreur;
158: real8 somme;
159: real8 registre;
160: real8 tampon;
161:
162: unsigned long i;
163:
164: somme = 0;
165: erreur = 0;
166:
167: (*erreur_memoire) = d_faux;
168:
169: for(i = 0; i < (*taille); i++)
170: {
171: tampon = somme;
172: registre = vecteur[i] + erreur;
173: somme = tampon + registre;
174: erreur = (tampon - somme) + registre;
175: }
176:
177: return(somme);
178: #endif
179: }
180:
181:
182: /*
183: ================================================================================
184: Fonction réalisation la sommation d'un vecteur de complexes sans perte
185: de précision
186: ================================================================================
187: Entrées : pointeur sur une structure struct_processus
188: --------------------------------------------------------------------------------
189: Sorties :
190: --------------------------------------------------------------------------------
191: Effets de bord : néant
192: ================================================================================
193: */
194:
195: complex16
196: sommation_vecteur_complexe(complex16 *vecteur, unsigned long *taille,
197: logical1 *erreur_memoire)
198: {
199: complex16 cumul;
200:
201: real8 *tampon;
202:
203: unsigned long i;
204: unsigned long nombre_elements;
205:
206: if ((tampon = malloc((*taille) * sizeof(real8))) == NULL)
207: {
208: (*erreur_memoire) = d_vrai;
209:
210: cumul.partie_reelle = 0;
211: cumul.partie_imaginaire = 0;
212:
213: return(cumul);
214: }
215:
216: (*erreur_memoire) = d_faux;
217:
218: /*
219: * Sommation des termes en commençant par le plus petit
220: */
221:
222: for(i = 0, nombre_elements = (*taille); i < nombre_elements;
223: tampon[i] = vecteur[i].partie_reelle, i++);
224:
225: cumul.partie_reelle = sommation_vecteur_reel(tampon, taille,
226: erreur_memoire);
227:
228: /*
229: * Même traitement, mais sur la partie imaginaire
230: */
231:
232: for(i = 0, nombre_elements = (*taille); i < nombre_elements;
233: tampon[i] = vecteur[i].partie_imaginaire, i++);
234:
235: cumul.partie_imaginaire = sommation_vecteur_reel(tampon, taille,
236: erreur_memoire);
237:
238: free(tampon);
239:
240: return(cumul);
241: }
242:
243: // vim: ts=4
CVSweb interface <joel.bertrand@systella.fr>