Annotation of rpl/src/calcul_integral.c, revision 1.1
1.1 ! bertrand 1: /*
! 2: ================================================================================
! 3: RPL/2 (R) version 4.0.9
! 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 'integrale_romberg'
! 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: integrale_romberg(struct_processus *s_etat_processus,
! 40: struct_objet *s_expression, unsigned char *variable,
! 41: real8 a, real8 b, real8 precision)
! 42: {
! 43: real8 erreur;
! 44: real8 h;
! 45: real8 **t;
! 46: real8 **t_tampon;
! 47: real8 *vecteur;
! 48: real8 *vecteur_precedent;
! 49: real8 x;
! 50:
! 51: logical1 validite;
! 52: logical1 erreur_memoire;
! 53:
! 54: struct_objet *s_objet;
! 55:
! 56: struct_variable s_variable;
! 57:
! 58: unsigned long i;
! 59: unsigned long j;
! 60: unsigned long k;
! 61: unsigned long n;
! 62: unsigned long p;
! 63: unsigned long nombre_elements;
! 64: unsigned long nombre_termes;
! 65: unsigned long taille_vecteur_precedent;
! 66:
! 67: /*
! 68: * Création d'une variable locale représentant la variable d'intégration
! 69: */
! 70:
! 71: (*s_etat_processus).niveau_courant++;
! 72: s_variable.niveau = (*s_etat_processus).niveau_courant;
! 73:
! 74: if ((s_variable.nom = (unsigned char *) malloc((strlen(variable) + 1) *
! 75: sizeof(unsigned char))) == NULL)
! 76: {
! 77: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 78: return;
! 79: }
! 80:
! 81: strcpy(s_variable.nom, variable);
! 82:
! 83: if ((s_variable.objet = allocation(s_etat_processus, REL)) == NULL)
! 84: {
! 85: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 86: return;
! 87: }
! 88:
! 89: if (creation_variable(s_etat_processus, &s_variable, 'V', 'P') == d_erreur)
! 90: {
! 91: return;
! 92: }
! 93:
! 94: /*
! 95: * Algorithme de Romberg
! 96: */
! 97:
! 98: if ((vecteur = (real8 *) malloc(2 * sizeof(real8))) == NULL)
! 99: {
! 100: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 101: return;
! 102: }
! 103:
! 104: evaluation_romberg(s_etat_processus, s_expression, variable, &a,
! 105: &(vecteur[nombre_termes = 0]), &validite);
! 106: nombre_termes += (validite == d_vrai) ? 1 : 0;
! 107:
! 108: evaluation_romberg(s_etat_processus, s_expression, variable, &b,
! 109: &(vecteur[nombre_termes]), &validite);
! 110: nombre_termes += (validite == d_vrai) ? 1 : 0;
! 111:
! 112: if ((t = (real8 **) malloc(sizeof(real8 *))) == NULL)
! 113: {
! 114: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 115: return;
! 116: }
! 117:
! 118: if ((t[0] = (real8 *) malloc(sizeof(real8))) == NULL)
! 119: {
! 120: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 121: return;
! 122: }
! 123:
! 124: t[0][n = 0] = ((b - a) / 2) * sommation_vecteur_reel(vecteur,
! 125: &nombre_termes, &erreur_memoire);
! 126:
! 127: if (erreur_memoire == d_vrai)
! 128: {
! 129: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 130: return;
! 131: }
! 132:
! 133: free(vecteur);
! 134:
! 135: vecteur_precedent = NULL;
! 136: taille_vecteur_precedent = 0;
! 137:
! 138: do
! 139: {
! 140: h = (b - a) / (nombre_elements = (1ULL << (++n)));
! 141: x = a;
! 142:
! 143: /*
! 144: * Nouvelle allocation de t
! 145: */
! 146:
! 147: t_tampon = t;
! 148:
! 149: if ((t = (real8 **) malloc((n + 1) * sizeof(real8 *))) == NULL)
! 150: {
! 151: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 152: return;
! 153: }
! 154:
! 155: for(i = 0; i <= n; i++)
! 156: {
! 157: if ((t[i] = (real8 *) malloc((n + 1) * sizeof(real8))) == NULL)
! 158: {
! 159: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 160: return;
! 161: }
! 162:
! 163: if (i != n)
! 164: {
! 165: for(t[i][n] = 0, j = 0; j < n; j++)
! 166: {
! 167: t[i][j] = t_tampon[i][j];
! 168: }
! 169: }
! 170: else
! 171: {
! 172: for(j = 0; j <= n; t[i][j++] = 0);
! 173: }
! 174: }
! 175:
! 176: for(i = 0; i < n; free(t_tampon[i++]));
! 177: free(t_tampon);
! 178:
! 179: /*
! 180: * Boucle principale
! 181: */
! 182:
! 183: if ((vecteur = (real8 *) malloc((nombre_elements + 1) *
! 184: sizeof(real8))) == NULL)
! 185: {
! 186: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 187: return;
! 188: }
! 189:
! 190:
! 191: if (vecteur_precedent == NULL)
! 192: {
! 193: /*
! 194: * Première itération
! 195: */
! 196:
! 197: evaluation_romberg(s_etat_processus, s_expression, variable, &a,
! 198: &(vecteur[nombre_termes = 0]), &validite);
! 199: nombre_termes += (validite == d_vrai) ? 1 : 0;
! 200:
! 201: evaluation_romberg(s_etat_processus, s_expression, variable, &b,
! 202: &(vecteur[nombre_termes]), &validite);
! 203: nombre_termes += (validite == d_vrai) ? 1 : 0;
! 204:
! 205: if (nombre_termes == 2)
! 206: {
! 207: vecteur[0] += vecteur[1];
! 208: vecteur[0] /= 2;
! 209: nombre_termes--;
! 210: }
! 211:
! 212: for(i = 1; i < nombre_elements; i++)
! 213: {
! 214: x += h;
! 215:
! 216: evaluation_romberg(s_etat_processus, s_expression, variable, &x,
! 217: &(vecteur[nombre_termes]), &validite);
! 218: nombre_termes += (validite == d_vrai) ? 1 : 0;
! 219: }
! 220: }
! 221: else
! 222: {
! 223: /*
! 224: * Pour les itérations suivantes, la moitié des points
! 225: * a déjà été calculée.
! 226: */
! 227:
! 228: for(i = 0; i < taille_vecteur_precedent; i++)
! 229: {
! 230: vecteur[i] = vecteur_precedent[i];
! 231: }
! 232:
! 233: free(vecteur_precedent);
! 234:
! 235: for(nombre_termes = taille_vecteur_precedent, i = 1;
! 236: i < nombre_elements; x += h, i += 2)
! 237: {
! 238: x += h;
! 239:
! 240: evaluation_romberg(s_etat_processus, s_expression, variable, &x,
! 241: &(vecteur[nombre_termes]), &validite);
! 242: nombre_termes += (validite == d_vrai) ? 1 : 0;
! 243: }
! 244: }
! 245:
! 246: t[0][n] = h * sommation_vecteur_reel(vecteur, &nombre_termes,
! 247: &erreur_memoire);
! 248:
! 249: if (erreur_memoire == d_vrai)
! 250: {
! 251: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 252: return;
! 253: }
! 254:
! 255: vecteur_precedent = vecteur;
! 256: taille_vecteur_precedent = nombre_termes;
! 257:
! 258: for(p = 1; p <= n; p++)
! 259: {
! 260: k = n - p;
! 261: t[p][k] = ((pow(4, p) * t[p - 1][k + 1]) - t[p - 1][k]) /
! 262: ((real8) (pow(4, p) - 1));
! 263: }
! 264: } while(((erreur = fabs(t[n][0] - t[n - 1][0])) > precision) &&
! 265: ((*s_etat_processus).var_volatile_requete_arret == 0));
! 266:
! 267: free(vecteur_precedent);
! 268:
! 269: /*
! 270: * Empilement du résultat et liberation de t
! 271: */
! 272:
! 273: if ((s_objet = allocation(s_etat_processus, REL)) == NULL)
! 274: {
! 275: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 276: return;
! 277: }
! 278:
! 279: (*((real8 *) (*s_objet).objet)) = t[n][0];
! 280:
! 281: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile),
! 282: s_objet) == d_erreur)
! 283: {
! 284: return;
! 285: }
! 286:
! 287: if ((s_objet = allocation(s_etat_processus, REL)) == NULL)
! 288: {
! 289: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 290: return;
! 291: }
! 292:
! 293: (*((real8 *) (*s_objet).objet)) = erreur;
! 294:
! 295: if (empilement(s_etat_processus, &((*s_etat_processus).l_base_pile),
! 296: s_objet) == d_erreur)
! 297: {
! 298: return;
! 299: }
! 300:
! 301: for(i = 0; i <= n; free(t[i++]));
! 302: free(t);
! 303:
! 304: /*
! 305: * Libération de la variable locale
! 306: */
! 307:
! 308: (*s_etat_processus).niveau_courant--;
! 309:
! 310: if (retrait_variable(s_etat_processus, s_variable.nom, 'L')
! 311: == d_erreur)
! 312: {
! 313: if ((*s_etat_processus).erreur_systeme != d_es)
! 314: {
! 315: if ((*s_etat_processus).erreur_systeme ==
! 316: d_es_variable_introuvable)
! 317: {
! 318: (*s_etat_processus).erreur_systeme = d_es;
! 319: }
! 320: else
! 321: {
! 322: /*
! 323: * Erreur système
! 324: */
! 325:
! 326: return;
! 327: }
! 328: }
! 329:
! 330: (*s_etat_processus).erreur_execution = d_ex_variable_non_definie;
! 331: return;
! 332: }
! 333:
! 334: return;
! 335: }
! 336:
! 337:
! 338: void
! 339: evaluation_romberg(struct_processus *s_etat_processus,
! 340: struct_objet *s_expression, unsigned char *variable, real8 *point,
! 341: real8 *valeur, logical1 *validite)
! 342: {
! 343: long hauteur_pile;
! 344:
! 345: struct_liste_pile_systeme *l_position_normale;
! 346:
! 347: struct_objet *s_objet;
! 348:
! 349: /*
! 350: * Vérification du type de la variable
! 351: */
! 352:
! 353: (*validite) = d_faux;
! 354: (*valeur) = 0;
! 355:
! 356: if (recherche_variable(s_etat_processus, variable) == d_vrai)
! 357: {
! 358: if (((*s_etat_processus).s_liste_variables[(*s_etat_processus)
! 359: .position_variable_courante]).variable_verrouillee == d_vrai)
! 360: {
! 361: (*s_etat_processus).erreur_execution = d_ex_variable_verrouillee;
! 362: return;
! 363: }
! 364:
! 365: if ((*((*s_etat_processus).s_liste_variables
! 366: [(*s_etat_processus).position_variable_courante])
! 367: .objet).type != REL)
! 368: {
! 369: liberation(s_etat_processus, ((*s_etat_processus).s_liste_variables
! 370: [(*s_etat_processus).position_variable_courante]).objet);
! 371:
! 372: if ((((*s_etat_processus).s_liste_variables
! 373: [(*s_etat_processus).position_variable_courante]).objet =
! 374: allocation(s_etat_processus, REL)) == NULL)
! 375: {
! 376: (*s_etat_processus).erreur_systeme = d_es_allocation_memoire;
! 377: return;
! 378: }
! 379: }
! 380:
! 381: (*((real8 *) (*((*s_etat_processus).s_liste_variables
! 382: [(*s_etat_processus).position_variable_courante])
! 383: .objet).objet)) = (*point);
! 384: }
! 385: else
! 386: {
! 387: /*
! 388: * La variable d'intégration étant locale, l'utilisateur ne peut
! 389: * pas l'effacer. On ne doit donc jamais passer par ici. Si
! 390: * c'est néanmoins le cas, une erreur système est générée
! 391: * provoquant l'arrêt du programme même dans une
! 392: * structure IFERR.
! 393: */
! 394:
! 395: return;
! 396: }
! 397:
! 398: hauteur_pile = (*s_etat_processus).hauteur_pile_operationnelle;
! 399: l_position_normale = (*s_etat_processus).l_base_pile_systeme;
! 400:
! 401: (*s_etat_processus).erreur_execution = d_ex;
! 402: (*s_etat_processus).exception = d_ep;
! 403:
! 404: if (evaluation(s_etat_processus, s_expression, 'N') == d_erreur)
! 405: {
! 406: if ((*s_etat_processus).erreur_systeme != d_es)
! 407: {
! 408: /*
! 409: * Erreur système
! 410: */
! 411:
! 412: return;
! 413: }
! 414: else
! 415: {
! 416: /*
! 417: * Restauration de la pile initiale
! 418: */
! 419:
! 420: while((*s_etat_processus).hauteur_pile_operationnelle >
! 421: (unsigned long) hauteur_pile)
! 422: {
! 423: if (depilement(s_etat_processus, &((*s_etat_processus)
! 424: .l_base_pile), &s_objet) == d_erreur)
! 425: {
! 426: (*s_etat_processus).erreur_execution =
! 427: d_ex_manque_argument;
! 428: break;
! 429: }
! 430:
! 431: liberation(s_etat_processus, s_objet);
! 432: }
! 433:
! 434: /*
! 435: * Restauration de la pile système
! 436: */
! 437:
! 438: while((*s_etat_processus).l_base_pile_systeme !=
! 439: l_position_normale)
! 440: {
! 441: depilement_pile_systeme(s_etat_processus);
! 442:
! 443: if ((*s_etat_processus).erreur_systeme != d_es)
! 444: {
! 445: /*
! 446: * Une pile vide provoque une erreur système
! 447: */
! 448:
! 449: return;
! 450: }
! 451: }
! 452: }
! 453:
! 454: (*s_etat_processus).erreur_execution = d_ex;
! 455: (*s_etat_processus).exception = d_ep;
! 456: }
! 457: else
! 458: {
! 459: /*
! 460: * Retrait de la valeur de la pile
! 461: */
! 462:
! 463: if (depilement(s_etat_processus, &((*s_etat_processus).l_base_pile),
! 464: &s_objet) == d_erreur)
! 465: {
! 466: (*s_etat_processus).erreur_execution = d_ex_manque_argument;
! 467: return;
! 468: }
! 469:
! 470: if ((*s_objet).type == INT)
! 471: {
! 472: (*valeur) = (*((integer8 *) (*s_objet).objet));
! 473: (*validite) = d_vrai;
! 474: }
! 475: else if ((*s_objet).type == REL)
! 476: {
! 477: (*valeur) = (*((real8 *) (*s_objet).objet));
! 478: (*validite) = d_vrai;
! 479: }
! 480:
! 481: liberation(s_etat_processus, s_objet);
! 482: }
! 483:
! 484: return;
! 485: }
! 486:
! 487: // vim: ts=4
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to calcul_integral.c CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / src