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