Subversion Repositories wimsdev

Rev

Rev 10 | Rev 7847 | Go to most recent revision | Details | Compare with Previous | Last modification | View Log | RSS feed

Rev Author Line No. Line
10 reyssat 1
/*    Copyright (C) 1998-2003 XIAO, Gang of Universite de Nice - Sophia Antipolis
2
 *
3
 *  This program is free software; you can redistribute it and/or modify
4
 *  it under the terms of the GNU General Public License as published by
5
 *  the Free Software Foundation; either version 2 of the License, or
6
 *  (at your option) any later version.
7
 *
8
 *  This program is distributed in the hope that it will be useful,
9
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
10
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
11
 *  GNU General Public License for more details.
12
 *
13
 *  You should have received a copy of the GNU General Public License
14
 *  along with this program; if not, write to the Free Software
15
 *  Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
16
 */
17
 
18
/* log(-1) does not make sense in real */
19
#ifndef NAN
20
#define NAN log(-1)
21
#endif
22
 
7798 bpr 23
    /* Only two decimal points, less than 1 million.
24
     * No check of buffer length which should be at least 12.
25
     * returns the end of buffer. */
10 reyssat 26
char *moneyprint(char *p, double s)
27
{
28
    char *p1, *p2, buf[16];
29
    int t, t1, t2;
30
    if(s<0) {*p++='-'; s=-s;}
31
    if(s>999999) s=999999;
32
    t=floor(s*100+0.5); if(t>99999999) t=99999999; if(t<0) t=0;
33
    if(t==0) {*p++='0'; *p=0; return p;}
34
    t1=t/100; t2=t%100; p1=buf+10;
35
    for(*p1--=t1%10+'0',t1/=10;t1>0;*p1--=t1%10+'0',t1/=10);
36
    p2=buf+11;
37
    if(t2) {
7798 bpr 38
      *p2++='.';
39
      *p2++=t2/10+'0'; t2%=10;
40
      if(t2) *p2++=t2+'0';
10 reyssat 41
    }
42
    p1++; *p2=0; memmove(p,p1,p2-p1+1); p+=p2-p1;
43
    return p;
44
}
45
 
46
/* #define RAND_BUF_SIZE 4096
47
static char rand_buf[RAND_BUF_SIZE];
7798 bpr 48
*/
49
/* The trouble here is that httpd does not initialize
50
     * the variable RANDOM.
51
     * So I use time (microseconds) to get a quick solution. */
10 reyssat 52
void init_random(void)
53
{
54
    int r;
55
    struct timeval t;
56
/*    initstate(1,rand_buf,RAND_BUF_SIZE); */
57
    gettimeofday(&t,NULL);
58
    r=t.tv_usec+t.tv_sec*1000;
59
    if(r<0) r=-r; if(r==0) r=1;
60
    srandom(r);
61
}
62
 
7798 bpr 63
/* gives a double random number between 0 and m */
10 reyssat 64
double drand(double m)
65
{
66
    double r;
67
    r=((double) random()+(double) random()/(double) RAND_MAX);
68
    return (r/(double) RAND_MAX)*m;
69
}
70
 
7798 bpr 71
/* gives a random integer between 0 and n.
72
 * n maybe floating, but will be rounded */
10 reyssat 73
double irand(double n)
74
{
75
    int  end,r;
76
    if(n==0) return 0;
77
    if(n>0) end=n; else end=-n;
78
    r=(double) random()*end/RAND_MAX;
79
    if(r==n) r--;
7798 bpr 80
    if(n>0) return r; else return -r;
10 reyssat 81
}
82
 
7798 bpr 83
/* sign of d */
10 reyssat 84
double sign(double d)
85
{
86
    if(d==0) return 0;
87
    if(d<0) return -1;
88
    else return 1;
89
}
90
 
7798 bpr 91
/* rounding to integer: problem with half-way rounding */
10 reyssat 92
double myround(double d)
93
{
94
    long int t;
95
    if(d<0) t=d-0.5; else t=d+0.5;
96
    return t;
97
}
98
 
7798 bpr 99
/* log of base 2 */
10 reyssat 100
double mylog2(double d)
101
{
102
    return log(d)/log(2);
103
}
104
 
7798 bpr 105
/* sec */
10 reyssat 106
double sec(double d)
7798 bpr 107
{    return 1/cos(d);}
10 reyssat 108
 
7798 bpr 109
/* csc */
10 reyssat 110
double csc(double d)
7798 bpr 111
{    return 1/sin(d);}
10 reyssat 112
 
7798 bpr 113
/* cotangent function */
10 reyssat 114
double cotan(double d)
115
{
116
    return 1/tan(d);
117
}
118
 
7798 bpr 119
/* hyperbolic cotangent */
10 reyssat 120
double cotanh(double d)
121
{
122
    return 1/tanh(d);
123
}
124
 
7798 bpr 125
/* factorial of an integer */
10 reyssat 126
double factorial(double d)
127
{
128
    int i,n; double t;
129
    n=d;
130
    if(n<0 || n!=d) return NAN;
131
    if(n>1000) return HUGE_VAL;
132
    t=1; for(i=1;i<=n;i++) t=t*i;
133
    return t;
134
}
135
 
7798 bpr 136
/* binomial coefficient */
10 reyssat 137
double binomial(double d1,double d2)
138
{
139
    return factorial(d1)/(factorial(d2)*factorial(d1-d2));
140
}
141
 
7798 bpr 142
/* max and min */
10 reyssat 143
double max(double d1, double d2)
144
{
145
    if(!finite(d1) || !finite(d2)) return NAN;
146
    if(d1<d2) return d2; else return d1;
147
}
148
double min(double d1, double d2)
149
{
150
    if(!finite(d1) || !finite(d2)) return NAN;
151
    if(d1<d2) return d1; else return d2;
152
}
153
 
7798 bpr 154
/* gcd and lcm, not really checking errors. */
10 reyssat 155
double gcd(double n1, double n2)
156
{
157
    unsigned long long int l1, l2, ll;
158
    n1=abs(n1); n2=abs(n2);
159
    if(!finite(n1) || !finite(n2) || n1<0 || n2<0 ||
160
       n1>1E18 || n2>1E18) return NAN;
161
    l1=n1; l2=n2;
162
    if(l1<l2) {
7798 bpr 163
      ll=l1;l1=l2;l2=ll;
10 reyssat 164
    }
165
    if(l1==0) return NAN;
166
    while(l2>0) {
7798 bpr 167
      ll=l2;l2=l1%l2;l1=ll;
10 reyssat 168
    }
169
    return l1;
170
}
171
 
172
double lcm(double n1, double n2)
173
{
174
    return n1*n2/gcd(n1,n2);
175
}
176
 
177
struct {
178
    char *name;
179
    int type;
180
    double val;
181
    double (*f1) (double parm);
182
    double (*f2) (double parm1, double parm2);
183
} evalname[]={
184
      {"Argch", 1,      0,      acosh,  NULL},
185
      {"Argsh", 1,      0,      asinh,  NULL},
186
      {"Argth", 1,      0,      atanh,  NULL},
187
      {"E",     0,      M_E,    NULL,   NULL},
188
      {"EULER", 0,      0.57721566490153286,    NULL,   NULL},
189
      {EV_S,    0,      0,      NULL,   NULL},
190
      {EV_T,    0,      0,      NULL,   NULL},
191
      {EV_X,    0,      0,      NULL,   NULL},
192
      {EV_Y,    0,      0,      NULL,   NULL},
193
      {"Euler", 0,      0.57721566490153286,    NULL,   NULL},
194
      {"Inf",   0,      1,      log,    NULL},
195
      {"NaN",   0,      0,      log,    NULL},
196
      {"PI",    0,      M_PI,   NULL,   NULL},
197
      {"Pi",    0,      M_PI,   NULL,   NULL},
198
      {"abs",   1,      0,      fabs,   NULL},
199
      {"acos",  1,      0,      acos,   NULL},
200
      {"acosh", 1,      0,      acosh,  NULL},
201
      {"arccos",1,      0,      acos,   NULL},
202
      {"arcsin",1,      0,      asin,   NULL},
203
      {"arctan",1,      0,      atan,   NULL},
204
      {"arctg", 1,      0,      atan,   NULL},
205
      {"argch", 1,      0,      acosh,  NULL},
206
      {"argsh", 1,      0,      asinh,  NULL},
207
      {"argth", 1,      0,      atanh,  NULL},
208
      {"asin",  1,      0,      asin,   NULL},
209
      {"asinh", 1,      0,      asinh,  NULL},
210
      {"atan",  1,      0,      atan,   NULL},
211
      {"atanh", 1,      0,      atanh,  NULL},
212
      {"binomial",2,    0,      NULL,   binomial},
213
      {"ceil",  1,      0,      ceil,   NULL}, /* round-up integer */
214
      {"ch",    1,      0,      cosh,   NULL},
215
      {"cos",   1,      0,      cos,    NULL},
216
      {"cosh",  1,      0,      cosh,   NULL},
217
      {"cot",   1,      0,      cotan,  NULL},
218
      {"cotan", 1,      0,      cotan,  NULL},
219
      {"cotanh",1,      0,      cotanh, NULL},
220
      {"coth",  1,      0,      cotanh, NULL},
221
      {"csc",   1,      0,      csc,    NULL},
222
      {"ctg",   1,      0,      cotan,  NULL},
223
      {"cth",   1,      0,      cotanh, NULL},
224
      {"drand", 1,      0,      drand,  NULL},
225
      {"e",     0,      M_E,    NULL,   NULL},
226
      {"erf",   1,      0,      erf,    NULL},
227
      {"erfc",  1,      0,      erfc,   NULL},
228
      {"euler", 0,      0.57721566490153286,    NULL,   NULL},
229
      {"exp",   1,      0,      exp,    NULL},
230
      {"factorial",1,   0,      factorial,      NULL},
231
      {"floor", 1,      0,      floor,  NULL},
232
      {"gcd",   2,      0,      NULL,   gcd},
233
      {"irand", 1,      0,      irand,  NULL},
234
/*      {"j0",  1,      0,      j0,     NULL}, */ /* Bessel functions */
235
/*      {"j1",  1,      0,      j1,     NULL}, */
236
      {"lcm",   2,      0,      NULL,   lcm},
237
      {"lg",    1,      0,      log10,  NULL},
238
      {"lgamma",1,      0,      lgamma, NULL}, /* log of Gamma function */
239
      {"ln",    1,      0,      log,    NULL},
240
      {"log",   1,      0,      log,    NULL},
241
      {"log10", 1,      0,      log10,  NULL},
242
      {"log2",  1,      0,      mylog2, NULL},
243
      {"max",   2,      0,      NULL,   max},
244
      {"min",   2,      0,      NULL,   min},
245
      {"pi",    0,      M_PI,   NULL,   NULL},
246
      {"pow",   2,      0,      NULL,   pow},
247
      {"rand",  1,      0,      drand,  NULL},
248
      {"randdouble",1,  0,      drand,  NULL},
249
      {"randfloat",1,   0,      drand,  NULL},
250
      {"randint",1,     0,      irand,  NULL},
251
      {"random",1,      0,      drand,  NULL},
252
      {"randreal",1,    0,      drand,  NULL},
253
      {"rint",  1,      0,      myround,        NULL}, /* closest integer */
254
      {"round", 1,      0,      myround,        NULL}, /* closest integer */
255
      {"sec",   1,      0,      sec,    NULL},
256
      {"sgn",   1,      0,      sign,   NULL}, /* sign of the value */
257
      {"sh",    1,      0,      sinh,   NULL},
258
      {"sign",  1,      0,      sign,   NULL}, /* sign of the value */
259
      {"sin",   1,      0,       sin,   NULL},
260
      {"sinh",  1,      0,      sinh,   NULL},
261
      {"sqrt",  1,      0,      sqrt,   NULL},
262
      {"tan",   1,      0,      tan,    NULL},
263
      {"tanh",  1,      0,      tanh,   NULL},
264
      {"tg",    1,      0,      tan,    NULL},
265
      {"th",    1,      0,      tanh,   NULL},
266
/*      {"y0",  1,      0,      y0,     NULL}, */
267
/*      {"y1",  1,      0,      y1,     NULL},  */
268
};
269
#define evalname_no (sizeof(evalname)/sizeof(evalname[0]))
270
 
271
int get_evalcnt(void) {return evalname_no;}
272
char *get_evalname(int i) {return evalname[i].name;}
273
int get_evaltype(int i) {return evalname[i].type;}
274
int evaltab_verify(void) {return verify_order(evalname,evalname_no,sizeof(evalname[0]));}
275
int search_evaltab(char *p) {
276
    return search_list(evalname,evalname_no,sizeof(evalname[0]),p);
277
}
278
 
279
static char *evalue_pt;
280
int evalue_error;
281
 
282
int get_evalue_error(void) { return evalue_error; }
283
void set_evalue_error(int e) {evalue_error=e; return;}
284
 
285
/* prepare pointer for evaluation */
286
void set_evalue_pointer(char *p)
287
{
288
    evalue_pt=p;
289
}
290
 
7798 bpr 291
/* get position of name in nametable */
10 reyssat 292
int eval_getpos(char *name)
293
{
294
    return search_list(evalname,evalname_no,sizeof(evalname[0]),name);
295
}
296
 
7798 bpr 297
/* set value to name */
10 reyssat 298
void eval_setval(int pos, double v)
299
{
300
    if(pos>=0 && pos<evalname_no) evalname[pos].val=v;
301
}
302
 
7798 bpr 303
/* get string pointer (after evaluation) */
10 reyssat 304
char *get_evalue_pointer(void)
305
{
306
    return evalue_pt;
307
}
308
 
309
double _evalue(int ord)
310
{
311
    double d,dd;
312
    int i,k;
313
    char buf[32];
314
 
7798 bpr 315
 
10 reyssat 316
    if(evalue_error) return NAN;
317
    d=0;
318
    while(*evalue_pt=='+') evalue_pt++;
319
    if(*evalue_pt==0) return 0; /* empty string */
320
    switch(*evalue_pt) {
7798 bpr 321
      case '(':
322
       evalue_pt++; d=_evalue(')');goto vld;
323
      case '|':
324
       if(ord=='|') {
325
        evalue_pt++; return 0;
326
       }
327
       evalue_pt++; d=fabs(_evalue('|'));goto vld;
328
      case '-':
329
       evalue_pt++; d=-_evalue(6);goto vld;
10 reyssat 330
    }
7798 bpr 331
    if((128&*evalue_pt)!=0) {/* special character */
332
      k=(*evalue_pt)&255; evalue_pt++;
333
      if(k>=130 && k<140) {
334
          i=(k-130)*200; k=(*evalue_pt)&255; evalue_pt++;
335
          if(k<33 || k>=233) goto badeval;
336
          i+=k-33; if(i<0 || i>=evalname_no) goto badeval;
337
          goto ename;
338
      }
339
      if(k>=140 && k<150) {
340
          i=(k-140)*200; k=(*evalue_pt)&255; evalue_pt++;
341
          if(k<33 || k>=233) goto badeval;
342
          if(ev_var==NULL || ev_varcnt==NULL) goto badeval;
343
          i+=k-33; if(i<0 || i>=*ev_varcnt) goto badeval;
344
          goto vname;
345
      }
346
      evalue_pt++; goto badeval;
10 reyssat 347
    }
348
    if(*evalue_pt=='.' || myisdigit(*evalue_pt))
7798 bpr 349
      {d=strtod(evalue_pt,&evalue_pt);goto binary;}
10 reyssat 350
    for(i=0;myisalnum(*(evalue_pt+i)) && i<16; i++)
351
      buf[i]=*(evalue_pt+i);
352
    buf[i]=0; evalue_pt+=i;
353
    if(i==0) goto badeval;
354
    if(ev_varcnt!=NULL && ev_var!=NULL && *ev_varcnt>0)
355
      for(i=0;i<*ev_varcnt;i++) {
7798 bpr 356
        if(strcmp(buf,ev_var[i].name)==0) {
357
            vname: d=ev_var[i].value; goto vld;
358
        }
10 reyssat 359
      }
360
    i=search_list(evalname,evalname_no,sizeof(evalname[0]),buf);
7798 bpr 361
    ename: if(i>=0) switch(evalname[i].type) { /* evaluation of expressions */
362
      case 0: {
363
          d=evalname[i].val;
364
          if(evalname[i].f1!=NULL) {
365
            if(d==0) d=NAN; if(d==1) d=HUGE_VAL;
366
          }
367
          break;
368
      }
369
      case 1: {
370
          if(*evalue_pt!='(') return NAN;
371
          evalue_pt++;
372
          d=evalname[i].f1(_evalue(')')); break;/* evaluation of function */
373
      }
374
      case 2: {
375
          double parm1,parm2;
376
          if(*evalue_pt!='(') return NAN;
377
          evalue_pt++;
378
          parm1=_evalue(',');parm2=_evalue(')');
379
          d=evalname[i].f2(parm1,parm2); break;
380
      }
381
      default: {      /* This is impossible. */
382
          return NAN;
383
      }
10 reyssat 384
    }
385
    else {
7798 bpr 386
      badeval: evalue_error=-1; return NAN;
10 reyssat 387
    }
388
  vld:
389
    if(evalue_error) return NAN;
7798 bpr 390
  binary: /*evaluation des expressions */
10 reyssat 391
    if(*evalue_pt=='!') {
7798 bpr 392
      evalue_pt++; d=factorial(d);
10 reyssat 393
    }
7798 bpr 394
    if(*evalue_pt==ord) {evalue_pt++;goto ok;}/* */
395
    if(*evalue_pt==0 || /* chaine de caractere finie*/
10 reyssat 396
       (ord<10 && (*evalue_pt==',' || *evalue_pt==';' || *evalue_pt==')'
7798 bpr 397
               || *evalue_pt=='|')))
10 reyssat 398
       goto ok;
399
    switch(*evalue_pt) {
7798 bpr 400
      case '+':
401
        if(ord<=8) break;
402
        evalue_pt++; d+=_evalue(8);goto vld;
403
      case '-':
404
        if(ord<=8) break;
405
        evalue_pt++; d-=_evalue(8);goto vld;
406
      case '*':
407
        if(ord<=6) break;
408
        evalue_pt++; d*=_evalue(6);goto vld;
409
      case '/':
410
        if(ord<=6) break;
411
        evalue_pt++; dd=_evalue(6);
412
        if(dd==0) {evalue_error=10;return NAN;}
413
            d/=dd;goto vld;
414
      case '%': {
415
        int di, ddi;
416
        if(ord<=6) break;
417
        evalue_pt++; dd=_evalue(6);
418
        if(dd==0) {evalue_error=10;return NAN;}
419
          di=d; ddi=dd; d=di%ddi;goto vld;
420
      }
421
      case '^': {
422
       if(ord<5) break;
423
       evalue_pt++; d=pow(d,_evalue(5));goto vld;
424
      }
425
      default : {
426
          return NAN;
427
      }
10 reyssat 428
    }
429
    ok: return d;
430
}
431
 
7798 bpr 432
/* substitute variable names by their environment strings
433
 * The buffer pointed to by p must have enough space
434
 * (defined by MAX_LINELEN). */
10 reyssat 435
char *_substit(char *p)
436
{
437
    return p;
438
}
439
 
440
char *(*substitute) (char *p)=_substit;
441
 
7798 bpr 442
/* evalue a string to double */
10 reyssat 443
double strevalue(char *p)
444
{
445
    char buf[MAX_LINELEN+1];
7798 bpr 446
 
10 reyssat 447
    if(p==NULL) return 0;
448
    mystrncpy(buf,p,sizeof(buf));
449
    substitute(buf); nospace(buf);
450
    if(check_parentheses(buf,0)) return NAN;
451
    set_evalue_error(0);
452
    set_evalue_pointer(buf);
453
    return _evalue(10);
454
}
455
 
7798 bpr 456
/* compile an expression for faster evaluation
457
 * returns -1 if cannot be compiled.
458
 * else returns the number of compilations.
459
 */
10 reyssat 460
int evalue_compile(char *p)
461
{
462
    char *p1, *p2, *pe, name[256], buf[8];
463
    int i,k;
7798 bpr 464
 
10 reyssat 465
    k=0;
466
    for(p1=p; *p1; p1++) if((128&*p1)!=0) return -1;
467
    nospace(p);
468
    for(p1=find_mathvar_start(p); *p1; p1=find_mathvar_start(pe)) {
7798 bpr 469
      pe=find_mathvar_end(p1);
470
      if(!myisalpha(*p1)) continue;
471
      p2=pe; if(p2-p1>16) continue;
472
      memmove(name,p1,p2-p1); name[p2-p1]=0;
473
/* replace the variables by a number
474
 * at most 2000 variables on two characters :
475
 * variable: 140 <= integer <150,  number between 33 and 233
476
 * function: 130 <= integer < 140, number between 33 and 233
477
 */
478
    if(ev_varcnt!=NULL && ev_var!=NULL && *ev_varcnt>0) {
479
        for(i=0;i<*ev_varcnt && strcmp(name,ev_var[i].name)!=0;i++);
480
        if(i<*ev_varcnt && i<2000) {
481
        buf[0]=i/200+140; buf[1]=i%200+33; buf[2]=0;
482
        string_modify(p,p1,p2,"%s",buf);
483
        pe=p1+2; k++; continue;
484
        }
10 reyssat 485
    }
7798 bpr 486
    i=search_list(evalname,evalname_no,sizeof(evalname[0]),name);
487
    if(i>=0 && i<2000) {
488
        buf[0]=i/200+130; buf[1]=i%200+33; buf[2]=0;
489
        string_modify(p,p1,p2,"%s",buf);
490
        pe=p1+2; k++; continue;
491
    }
492
  }
493
  return k;
10 reyssat 494
}