1/* BEGIN LICENSE BLOCK 2 * Version: CMPL 1.1 3 * 4 * The contents of this file are subject to the Cisco-style Mozilla Public 5 * License Version 1.1 (the "License"); you may not use this file except 6 * in compliance with the License. You may obtain a copy of the License 7 * at www.eclipse-clp.org/license. 8 * 9 * Software distributed under the License is distributed on an "AS IS" 10 * basis, WITHOUT WARRANTY OF ANY KIND, either express or implied. See 11 * the License for the specific language governing rights and limitations 12 * under the License. 13 * 14 * The Original Code is The ECLiPSe Constraint Logic Programming System. 15 * The Initial Developer of the Original Code is Cisco Systems, Inc. 16 * Portions created by the Initial Developer are 17 * Copyright (C) 1989-2006 Cisco Systems, Inc. All Rights Reserved. 18 * 19 * Contributor(s): 20 * 21 * END LICENSE BLOCK */ 22 23/* 24 * VERSION $Id: bip_arith.c,v 1.24 2015/04/02 03:35:08 jschimpf Exp $ 25 */ 26 27/* 28 * Overview of the arithmetic system 29 * ---------------------------------- 30 * 31 * The externals corresponding to the arithmetic built-in predicates are 32 * either here or (for the more critical ones) in the emulator. E.g. 33 * BIAdd in the emulator is the implementation of +/3, and p_sin() here 34 * is the implementation of sin/1. 35 * 36 * Although the external functions sometimes handle TINT and TDBL arguments 37 * directly, the general mechanism is to go through tag-indexed function 38 * tables that handle type coercion, i.e. 39 * 40 * tag_desc[<type>].coerce_to[<othertype>] 41 * functions in this table are called _<type>_<othertype> 42 * 43 * or dispatch to the routine that implements a particular arithmetic 44 * operation for a particular type: 45 * 46 * tag_desc[<type>].arith_op[<operation>] 47 * functions in this table are called _<type>_<operation> 48 * 49 * Other table functions that must be implemented for each numeric type are 50 * 51 * tag_desc[<type>].arith_compare 52 * functions in this table are called _arith_compare_<type> 53 * tag_desc[<type>].compare 54 * functions in this table are called _compare_<type> 55 * tag_desc[<type>].equal 56 * functions in this table are called _equal_<type> 57 * 58 * All table functions for TINT and TDBL are here, the functions for TBIG 59 * and TRAT are in bigrat.c, and the ones for breals (TIVL) in intervals.c. 60 * In principle, it should be possible to build eclipse without bigrat.c 61 * and intervals.c (the corresponding functions then remain initialised 62 * to dummies). 63 * 64 * Sufficiently regular arithmetic externals share code by calling 65 * unary_arith_op() or binary_arith_op(). These do the following: 66 * 67 * 1. check for instantiation fault/delay 68 * 2. make argument types equal, and/or lift argument type(s) if necessary 69 * 3. invoke the actual operation via arith_op table on that type 70 * 4. unify the result 71 * 72 */ 73 74/* 75 * To make delaying and error reporting consistent in all combinations 76 * of inline-expanded, interpreted and coroutining arithmetic, the 77 * policy is as follows. Note that we differ from the usual rule 78 * of checking for type errors first. 79 * 80 * 1. evaluate sub-expressions 81 * This is needed for compatibility with inline expansion 82 * 2. check delay conditions 83 * In non-coroutine mode we simply report an instantiation 84 * fault instead of delaying. Thus the same code can be used 85 * for delaying and non-delaying builtins. 86 * 3. check for type errors 87 * 4. compute the result 88 */ 89 90#include <math.h> 91 92#include "config.h" 93#include "sepia.h" 94#include "types.h" 95#include "embed.h" 96#include "error.h" 97#include "mem.h" 98#include "dict.h" 99#include "emu_export.h" 100#include "rounding_control.h" /* for init_rounding_modes() */ 101 102/* range checks for functions */ 103 104#define OneMOne(x) ( (x) >= -1.0 && (x) <= 1.0 ) 105#define Positive(x) ( (x) > 0.0 ) 106#define NonNegative(x) ( (x) >= 0.0 ) 107 108#define NonIntNum(t) (IsDouble(t) || IsRational(t) || IsInterval(t)) 109 110#define BITS_PER_WORD (8*SIZEOF_WORD) 111 112static int 113 _reverse_times(word x, word y, value zval, type ztag); 114 115 116#if defined(i386) && defined(__GNUC__) 117double (*pow_ptr_to_avoid_buggy_inlining)(double,double) = pow; 118#endif 119 120 121/*------------------------------------------------------------------------ 122 * The multi-directional builtins succ/2, plus/3 and times/3 123 * Some don't work with bignums yet! Maybe write them in Prolog? 124 *-----------------------------------------------------------------------*/ 125 126static int 127p_succ(value x, type tx, value y, type ty) 128{ 129 pword result; 130 131 if (!(IsRef(tx) || IsNumber(tx))) { Bip_Error(ARITH_TYPE_ERROR) } 132 if (NonIntNum(tx)) { Bip_Error(TYPE_ERROR) } 133 if (!(IsRef(ty) || IsNumber(ty))) { Bip_Error(ARITH_TYPE_ERROR) } 134 if (NonIntNum(ty)) { Bip_Error(TYPE_ERROR) } 135 136 result.tag.kernel = TINT; 137 if (IsInteger(tx)) 138 { 139 if (x.nint == MAX_S_WORD) { 140 Bip_Error(INTEGER_OVERFLOW); 141 } 142 if (x.nint < 0) { 143 Fail_; 144 } 145 result.val.nint = x.nint + 1; 146 Return_Numeric(y, ty, result); 147 } 148 else if (IsRef(tx)) 149 { 150 if (IsInteger(ty)) { 151 if (y.nint <= 0) { 152 Fail_ 153 } 154 result.val.nint = y.nint - 1; 155 Return_Numeric(x, tx, result); 156 } else if (IsRef(ty)) { 157 return PDELAY_1_2; 158 } 159 return unary_arith_op(y, ty, x, tx, ARITH_PREV, TINT); 160 } 161 else 162 return unary_arith_op(x, tx, y, ty, ARITH_NEXT, TINT); 163} 164 165 166static int 167p_plus(value x, type tx, value y, type ty, value z, type tz) 168{ 169 if (IsRef(tx)) 170 { 171 if (IsRef(ty)) 172 { 173 return PDELAY_1_2; 174 } 175 else if (IsRef(tz)) 176 { 177 return PDELAY_1_3; 178 } 179 else if (IsInteger(ty) && IsInteger(tz)) 180 { 181 Kill_DE; 182 Return_Unify_Integer(x, tx, z.nint - y.nint); 183 } 184 } 185 else if (IsRef(ty)) 186 { 187 if (IsRef(tz)) 188 { 189 return PDELAY_2_3; 190 } 191 else if (IsInteger(tx) && IsInteger(tz)) 192 { 193 Kill_DE; 194 Return_Unify_Integer(y, ty, z.nint - x.nint); 195 } 196 } 197 else if (IsInteger(tx) && IsInteger(ty) && (IsRef(tz) || IsInteger(tz))) 198 { 199 Kill_DE; 200 Return_Unify_Integer(z, tz, x.nint + y.nint); 201 } 202 203 if (NonIntNum(tx) || NonIntNum(ty) || NonIntNum(tz)) 204 {Bip_Error(TYPE_ERROR);} 205 else if (!IsNumber(tx) || !IsNumber(ty) || !IsNumber(tz)) 206 {Bip_Error(ARITH_TYPE_ERROR);} 207 else if (IsBignum(tx) || !IsBignum(ty) || !IsBignum(tz)) 208 {Bip_Error(RANGE_ERROR);} 209 else 210 {Bip_Error(TYPE_ERROR);} 211} 212 213 214static int 215p_times(value x, type tx, value y, type ty, value z, type tz) 216{ 217 if ((IsRef(tx) || IsInteger(tx)) && 218 (IsRef(ty) || IsInteger(ty)) && 219 (IsRef(tz) || IsInteger(tz))) 220 { 221 if (IsRef(tx)) 222 { 223 if (IsRef(ty)) 224 { 225 if (x.ptr == y.ptr && IsInteger(tz) && z.nint == 0) { 226 Kill_DE; 227 Return_Unify_Integer(x, tx, 0) 228 } else { 229 return PDELAY_1_2; 230 } 231 } 232 else if (y.nint == 0) 233 { 234 Kill_DE; 235 Return_Unify_Integer(z, tz, 0); 236 } 237 else if (y.nint == 1) 238 { 239 Kill_DE; 240 Return_Unify_Pw(z, tz, x, tx) 241 } 242 else if (IsRef(tz)) 243 { 244 return PDELAY_1_3; 245 } 246 else 247 return _reverse_times(z.nint, y.nint, x, tx); 248 } 249 else if (x.nint == 0) 250 { 251 Kill_DE; 252 Return_Unify_Integer(z, tz, 0); 253 } 254 else if (x.nint == 1) 255 { 256 Kill_DE; 257 Return_Unify_Pw(z, tz, y, ty) 258 } 259 else if (IsRef(ty)) 260 { 261 if (IsRef(tz)) 262 { 263 if (y.ptr == z.ptr) { 264 if (x.nint != 1) { 265 Kill_DE; 266 Return_Unify_Integer(z, tz, 0) 267 } 268 } 269 return PDELAY_2_3; 270 } 271 else 272 return _reverse_times(z.nint, x.nint, y, ty); 273 } 274 else 275 { 276 Kill_DE; 277 Return_Unify_Integer(z, tz, x.nint * y.nint); 278 } 279 } 280 if (NonIntNum(tx) || NonIntNum(ty) || NonIntNum(tz)) 281 {Bip_Error(TYPE_ERROR);} 282 else if (!IsNumber(tx) || !IsNumber(ty) || !IsNumber(tz)) 283 {Bip_Error(ARITH_TYPE_ERROR);} 284 else if (IsBignum(tx) || !IsBignum(ty) || !IsBignum(tz)) 285 {Bip_Error(RANGE_ERROR);} 286 else 287 {Bip_Error(TYPE_ERROR);} 288} 289 290 291/* 292 * _reverse_times is an auxiliary function for times 293 * used when times(INT,VAR,INT) or times(VAR,INT,INT) are called. 294 * receives two integers x,y 295 * returns_unifies z= x/y if this is an integer otherwise fails 296 * for times(X, 0, 0) we delay since this may be true (X=0) or false 297 */ 298 299static int 300_reverse_times(word x, word y, value zval, type ztag) 301{ 302 if (y == 0) 303 if (x == 0) 304 { 305 Push_var_delay(zval.ptr, ztag.all); 306 return PDELAY; 307 } 308 else 309 { 310 Fail_ 311 } 312 else if (x % y != 0) 313 { 314 Fail_; 315 } 316 else 317 { 318 Kill_DE; 319 Return_Unify_Integer(zval,ztag, x/y); 320 } 321} 322 323 324/*------------------------------------------------------------------------ 325 * Other arithmetic-related built-ins 326 *------------------------------------------------------------------------ */ 327 328/* 329 * between(Min, Max, Step, Index) 330 */ 331static int 332p_between(value vmi, type tmi, value vma, type tma, value vs, type ts, value vi, type ti) 333{ 334 value v; 335 336 Check_Integer(tmi) 337 Check_Integer(tma) 338 Check_Integer(ts) 339 if (vs.nint == 0) { 340 Bip_Error(RANGE_ERROR) 341 } 342 if (IsInteger(ti)) { 343 Cut_External 344 if (vs.nint > 0) { 345 Succeed_If( 346 vi.nint >= vmi.nint 347 && vi.nint <= vma.nint 348 && (vi.nint - vmi.nint) % vs.nint == 0 349 ) 350 } else { 351 Succeed_If( 352 vi.nint <= vmi.nint 353 && vi.nint >= vma.nint 354 && (vmi.nint - vi.nint) % -vs.nint == 0 355 ) 356 } 357 } else { 358 Check_Output_Integer(ti) 359 } 360 if (vs.nint > 0) { 361 if (vmi.nint >= vma.nint) { 362 Cut_External 363 if (vmi.nint > vma.nint) { 364 Fail_ 365 } 366 } 367 } else { 368 if (vmi.nint <= vma.nint) { 369 Cut_External 370 if (vmi.nint < vma.nint) { 371 Fail_ 372 } 373 } 374 } 375 v.nint = vmi.nint + vs.nint; 376 Remember(1, v, tmi) 377 Return_Unify_Integer(vi, ti, vmi.nint) 378} 379 380 381/* 382 * is_zero/1 and collect/3 383 * support externals for lib(r) - largely obsolete 384 */ 385 386static int 387p_is_zero(value v, type t) 388{ 389 pword result; 390 Succeed_If(tag_desc[TagType(t)].arith_op[ARITH_SGN](v, &result) == PSUCCEED 391 && result.val.nint == 0); 392} 393 394 395/* 396 * collect(+LinNormExpr, -LinNornSimplified, -ZeroVars) 397 * This was for library(r) and is pretty much obsolete. 398 */ 399 400#define OFF_C 0 401#define OFF_V 1 402 403static int 404p_collect(value vin, type tin, value vout, type tout, value vzero, type tzero) 405{ 406 register pword *curr_var, *curr_tail, *zero_tail, *new_tail; 407 register pword *pcoeff, *pvar, *pw; 408 pword in_list, out_list, zero_list, new_coeff; 409 int err; 410 Prepare_Requests; 411 412 Check_List(tin); 413 Check_Output_List(tout); 414 Check_Output_List(tzero); 415 416 in_list.val = vin; 417 in_list.tag = tin; 418 curr_tail = &in_list; 419 new_tail = &out_list; 420 zero_tail = &zero_list; 421 422 new_coeff.tag.kernel = TINT; 423 new_coeff.val.nint = 0; 424 curr_var = 0; /* 0 stands for constant sequence */ 425 426 while (IsList(curr_tail->tag)) 427 { 428 pw = curr_tail->val.ptr; 429 curr_tail = &pw[1]; 430 Dereference_(curr_tail); 431 Dereference_(pw); 432 pw = pw->val.ptr; 433 pcoeff = &pw[OFF_C]; 434 Dereference_(pcoeff); 435 pvar = &pw[OFF_V]; 436 Dereference_(pvar); 437 if (IsRef(pvar->tag)) /* mono(..., Var) */ 438 { 439 if (pvar == curr_var) /* inside a sequence */ 440 { 441 err = bin_arith_op(pcoeff->val, pcoeff->tag, 442 new_coeff.val, new_coeff.tag, &new_coeff, ARITH_ADD); 443 if (err != PSUCCEED) goto _error_; 444 } 445 else /* end of a sequence */ 446 { 447 if (p_is_zero(new_coeff.val, new_coeff.tag) == PSUCCEED) 448 { 449 if (curr_var) 450 { 451 pw = TG; 452 TG += 4; 453 Check_Gc; 454 Make_List(zero_tail, pw); 455 Make_List(&pw[0], &pw[2]); 456 zero_tail = &pw[1]; 457 pw = pw + 2; 458 pw[OFF_C].val.nint = 0; 459 pw[OFF_C].tag.kernel = TINT; 460 pw[OFF_V].val.ptr = curr_var; 461 pw[OFF_V].tag.kernel = TREF; 462 } 463 } 464 else 465 { 466 pw = TG; 467 TG += 4; 468 Check_Gc; 469 Make_List(new_tail, pw); 470 Make_List(&pw[0], &pw[2]); 471 new_tail = &pw[1]; 472 pw = pw + 2; 473 pw[OFF_C] = new_coeff; 474 if (curr_var) 475 { 476 pw[OFF_V].val.ptr = curr_var; 477 pw[OFF_V].tag.kernel = TREF; 478 } 479 else /* build new constant mono */ 480 { 481 Make_Integer(&pw[OFF_V], 1); 482 } 483 } 484 curr_var = pvar; 485 new_coeff = *pcoeff; 486 } 487 } 488 else /* in the constant part */ 489 { 490 pword product; 491 err = bin_arith_op(pcoeff->val, pcoeff->tag, 492 pvar->val, pvar->tag, &product, ARITH_MUL); 493 if (err != PSUCCEED) goto _error_; 494 err = bin_arith_op(product.val, product.tag, 495 new_coeff.val, new_coeff.tag, &new_coeff, ARITH_ADD); 496 if (err != PSUCCEED) goto _error_; 497 } 498 } 499 500 /* end of last sequence */ 501 if (p_is_zero(new_coeff.val, new_coeff.tag) == PSUCCEED) 502 { 503 if (curr_var) 504 { 505 pw = TG; 506 TG += 4; 507 Check_Gc; 508 Make_List(zero_tail, pw); 509 Make_List(&pw[0], &pw[2]); 510 zero_tail = &pw[1]; 511 pw = pw + 2; 512 pw[OFF_C].val.nint = 0; 513 pw[OFF_C].tag.kernel = TINT; 514 pw[OFF_V].val.ptr = curr_var; 515 pw[OFF_V].tag.kernel = TREF; 516 } 517 } 518 else 519 { 520 pw = TG; 521 TG += 4; 522 Check_Gc; 523 Make_List(new_tail, pw); 524 Make_List(&pw[0], &pw[2]); 525 new_tail = &pw[1]; 526 pw = pw + 2; 527 pw[OFF_C] = new_coeff; 528 if (curr_var) 529 { 530 pw[OFF_V].val.ptr = curr_var; 531 pw[OFF_V].tag.kernel = TREF; 532 } 533 else /* build new constant mono */ 534 { 535 Make_Integer(&pw[OFF_V], 1); 536 } 537 } 538 539 if (IsNil(curr_tail->tag)) 540 { 541 Make_Nil(new_tail); 542 Make_Nil(zero_tail); 543 Request_Unify_Pw(vout, tout, out_list.val, out_list.tag); 544 Request_Unify_Pw(vzero, tzero, zero_list.val, zero_list.tag); 545 Return_Unify 546 } 547 548 err = TYPE_ERROR; 549_error_: 550 /* may have to pop incomplete junk on the stack */ 551 Bip_Error(err); 552} 553 554 555/* 556 * collapse_linear(+LinDenormalSorted, -LinNormalised) 557 * 558 * The linear normal form is [C0*1,C1*X1,...,Cn*Xn] with numbers Ci and 559 * variables Xi. The Ci are nonzero for i>0, and the Xi are distinct. 560 * The C0*1 term is always present. 561 * 562 * Our input is expected to be denormalised due to variable instantiation 563 * and var-var bindings. But it is expected to have been sorted using 564 * sort(2, >=, LinDenormal, LinDenormalSorted). I.e. there may be several 565 * constant terms in the front (but not necessarily starting with C0*1), 566 * and sequences of terms with identical variables. 567 * 568 * This is a bit elaborate because we try to reuse parts of the input: 569 * - the whole tail of the input list that remains unchanged 570 * - the singleton monomial terms Ci*Xi 571 * This should reduce garbage collection time for large applications. 572 */ 573 574#define OFF_CONST 1 575#define OFF_VAR 2 576#define SIZE_MONO 3 577#define SIZE_LIST 2 578 579static int 580p_collapse_linear(value vin, type tin, value vout, type tout) 581{ 582 pword *in_tail, *out_tail; 583 pword *seq_mono, *seq_var, seq_coeff; 584 pword *in_reuse_tail, *out_reuse_tail, *reuse_tg = TG; 585 pword in_list, out_list, unit_var; 586 pword *old_tg = TG; 587 int err, const_seen = 0; 588 589 Check_List(tin); 590 in_list.val = vin; 591 in_list.tag = tin; 592 in_tail = in_reuse_tail = &in_list; 593 out_tail = out_reuse_tail = &out_list; 594 Make_Integer(&seq_coeff, 0); 595 Make_Integer(&unit_var, 1); 596 seq_var = &unit_var; 597 seq_mono = 0; 598 599 for(;;) 600 { 601 pword *cur_mono, *cur_var, *cur_coeff; 602 pword *cur_tail = in_tail; 603 604 if (IsList(in_tail->tag)) 605 { 606 /* read next input list element */ 607 pword *pw = in_tail->val.ptr; 608 in_tail = &pw[1]; 609 Dereference_(in_tail); 610 611 /* read the monomial Coeff*VarOrConst */ 612 Dereference_(pw); 613 if (!IsStructure(pw->tag)) 614 { 615 err= IsRef(pw->tag) ? INSTANTIATION_FAULT : TYPE_ERROR; 616 goto _error_; 617 } 618 cur_mono = pw->val.ptr; 619 if (cur_mono->val.did != d_.times) 620 { 621 err=TYPE_ERROR; 622 goto _error_; 623 } 624 cur_coeff = &cur_mono[OFF_CONST]; 625 Dereference_(cur_coeff); 626 cur_var = &cur_mono[OFF_VAR]; 627 Dereference_(cur_var); 628 629 if (!IsRef(cur_var->tag)) /* still part of the constant sequence */ 630 { 631 if (seq_var != &unit_var) 632 { 633 err = RANGE_ERROR; 634 goto _error_; /* malformed: constant after var */ 635 } 636 if (!const_seen && IsInteger(cur_var->tag) && cur_var->val.nint==1) 637 { 638 /* first (maybe only) C0*1 term, potentially reusable */ 639 seq_mono = cur_mono; 640 seq_coeff = *cur_coeff; 641 const_seen = 1; 642 continue; 643 } 644 else /* general constant monomial */ 645 { 646 pword product; 647 err = bin_arith_op(cur_coeff->val, cur_coeff->tag, 648 cur_var->val, cur_var->tag, &product, ARITH_MUL); 649 if (err != PSUCCEED) goto _error_; 650 err = bin_arith_op(product.val, product.tag, 651 seq_coeff.val, seq_coeff.tag, &seq_coeff, ARITH_ADD); 652 if (err != PSUCCEED) goto _error_; 653 seq_mono = 0; /* not reusable */ 654 const_seen = 1; 655 continue; 656 } 657 } 658 else if (cur_var == seq_var) /* still part of a variable sequence */ 659 { 660 err = bin_arith_op(cur_coeff->val, cur_coeff->tag, 661 seq_coeff.val, seq_coeff.tag, &seq_coeff, ARITH_ADD); 662 if (err != PSUCCEED) goto _error_; 663 seq_mono = 0; /* not reusable */ 664 continue; 665 } 666 /* else start of a new sequence */ 667 } 668 else if (IsNil(in_tail->tag)) 669 { 670 cur_mono = 0; 671 } 672 else 673 { 674 err = IsRef(in_tail->tag) ? INSTANTIATION_FAULT : TYPE_ERROR; 675 goto _error_; 676 } 677 678 /* emit monomial for finished sequence (unless it is 0*Var) */ 679 if (seq_var != &unit_var && p_is_zero(seq_coeff.val, seq_coeff.tag) == PSUCCEED) 680 { 681 /* an element was dropped: can't reuse earlier input list */ 682 in_reuse_tail = cur_tail; 683 out_reuse_tail = out_tail; 684 reuse_tg = TG; 685 } 686 else /* do emit the monomial, reusing the input one if possible */ 687 { 688 pword *pw = TG; 689 TG += SIZE_LIST; 690 Make_List(out_tail, pw); 691 out_tail = &pw[1]; 692 if (seq_mono) 693 { 694 Make_Struct(&pw[0], seq_mono); /* reuse the old mono */ 695 } 696 else 697 { 698 Make_Struct(&pw[0], TG); /* make new Coeff*VarOr1 */ 699 pw = TG; 700 TG += SIZE_MONO; 701 pw->val.did = d_.times; 702 pw->tag.kernel = TDICT; 703 pw[OFF_CONST] = seq_coeff; 704 pw[OFF_VAR] = *seq_var; 705 706 /* elements were collapsed: can't reuse earlier input list */ 707 in_reuse_tail = cur_tail; 708 out_reuse_tail = out_tail; 709 reuse_tg = TG; 710 } 711 Check_Gc; 712 } 713 if (!cur_mono) /* finished! */ 714 break; 715 716 /* init next sequence */ 717 seq_mono = cur_mono; 718 seq_var = cur_var; 719 seq_coeff = *cur_coeff; 720 } 721 722 /* Reuse reusable tail of input list, discard already constructed copy */ 723 *out_reuse_tail = *in_reuse_tail; 724 TG = reuse_tg; 725 726 Return_Unify_Pw(vout, tout, out_list.val, out_list.tag); 727 728_error_: 729 TG = old_tg; 730 Bip_Error(err); 731} 732 733 734/*------------------------------------------------------------------------ 735 * Standard arithmetic bips 736 * More frequently used ones are in the emulator 737 *-----------------------------------------------------------------------*/ 738 739int 740p_sgn(value v1, type t1, value v, type t) 741{ 742 int err; 743 pword result; 744 if (IsRef(t1)) { Bip_Error(PDELAY_1) } 745 err = tag_desc[TagType(t1)].arith_op[ARITH_SGN](v1, &result); 746 if (err != PSUCCEED) return err; 747 Return_Numeric(v, t, result) 748} 749 750static int 751p_min(value v1, type t1, value v2, type t2, value v, type t) 752{ 753 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_MIN); 754} 755 756static int 757p_max(value v1, type t1, value v2, type t2, value v, type t) 758{ 759 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_MAX); 760} 761 762static int 763p_gcd(value v1, type t1, value v2, type t2, value v, type t) 764{ 765 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_GCD); 766} 767 768static int 769p_gcd_ext(value v1, type t1, value v2, type t2, value s, type ts, value t, type tt, value g, type tg) 770{ 771 pword res1,res2,res3; 772 int err; 773 Prepare_Requests; 774 if (IsRef(t1)) { Bip_Error(PDELAY_1) } 775 if (IsRef(t2)) { Bip_Error(PDELAY_2) } 776 Check_Integer_Or_Bignum(t1) 777 Check_Integer_Or_Bignum(t2) 778 Check_Output_Integer_Or_Bignum(ts) 779 Check_Output_Integer_Or_Bignum(tt) 780 Check_Output_Integer_Or_Bignum(tg) 781 /* we don't have a TINT implementation, always compute via bignums */ 782 err = tag_desc[TagType(t1)].coerce_to[TBIG](v1, &v1); 783 if (err != PSUCCEED) return(err); 784 err = tag_desc[TagType(t2)].coerce_to[TBIG](v2, &v2); 785 if (err != PSUCCEED) return(err); 786 err = tag_desc[TBIG].arith_op[ARITH_GCD_EXT](v1, v2, &res1, &res2, &res3); 787 if (err != PSUCCEED) return err; 788 Kill_DE; /* in case it's a demon */ 789 Request_Unify_Pw(s, ts, res1.val, res1.tag); 790 Request_Unify_Pw(t, tt, res2.val, res2.tag); 791 Request_Unify_Pw(g, tg, res3.val, res3.tag); 792 Return_Unify 793} 794 795static int 796p_lcm(value v1, type t1, value v2, type t2, value v, type t) 797{ 798 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_LCM); 799} 800 801static int 802p_uplus(value v1, type t1, value v, type t) 803{ 804 return unary_arith_op(v1, t1, v, t, ARITH_PLUS, TINT); 805} 806 807static int 808p_abs(value v1, type t1, value v, type t) 809{ 810 return unary_arith_op(v1, t1, v, t, ARITH_ABS, TINT); 811} 812 813static int 814p_sin(value v1, type t1, value v, type t) 815{ 816 return unary_arith_op(v1, t1, v, t, ARITH_SIN, TDBL); 817} 818 819static int 820p_cos(value v1, type t1, value v, type t) 821{ 822 return unary_arith_op(v1, t1, v, t, ARITH_COS, TDBL); 823} 824 825static int 826p_tan(value v1, type t1, value v, type t) 827{ 828 return unary_arith_op(v1, t1, v, t, ARITH_TAN, TDBL); 829} 830 831static int 832p_asin(value v1, type t1, value v, type t) 833{ 834 return unary_arith_op(v1, t1, v, t, ARITH_ASIN, TDBL); 835} 836 837static int 838p_acos(value v1, type t1, value v, type t) 839{ 840 return unary_arith_op(v1, t1, v, t, ARITH_ACOS, TDBL); 841} 842 843static int 844p_atan(value v1, type t1, value v, type t) 845{ 846 return unary_arith_op(v1, t1, v, t, ARITH_ATAN, TDBL); 847} 848 849static int 850p_atan2(value v1, type t1, value v2, type t2, value v, type t) 851{ 852 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_ATAN2); 853} 854 855static int 856p_exp(value v1, type t1, value v, type t) 857{ 858 return unary_arith_op(v1, t1, v, t, ARITH_EXP, TDBL); 859} 860 861static int 862p_ln(value v1, type t1, value v, type t) 863{ 864 return unary_arith_op(v1, t1, v, t, ARITH_LN, TDBL); 865} 866 867static int 868p_sqrt(value v1, type t1, value v, type t) 869{ 870 return unary_arith_op(v1, t1, v, t, ARITH_SQRT, TDBL); 871} 872 873static int 874p_round(value v1, type t1, value v, type t) 875{ 876 return unary_arith_op(v1, t1, v, t, ARITH_ROUND, TINT); 877} 878 879static int 880p_floor(value v1, type t1, value v, type t) 881{ 882 return unary_arith_op(v1, t1, v, t, ARITH_FLOOR, TINT); 883} 884 885static int 886p_ceiling(value v1, type t1, value v, type t) 887{ 888 return unary_arith_op(v1, t1, v, t, ARITH_CEIL, TINT); 889} 890 891static int 892p_truncate(value v1, type t1, value v, type t) 893{ 894 return unary_arith_op(v1, t1, v, t, ARITH_TRUNCATE, TINT); 895} 896 897static int 898p_fix(value v1, type t1, value v, type t) 899{ 900 return unary_arith_op(v1, t1, v, t, ARITH_FIX, TINT); 901} 902 903static int 904p_integer2(value v1, type t1, value v, type t) 905{ 906 return unary_arith_op(v1, t1, v, t, ARITH_INT, TINT); 907} 908 909static int 910p_numerator(value v1, type t1, value v, type t) 911{ 912 return unary_arith_op(v1, t1, v, t, ARITH_NUM, TRAT); 913} 914 915static int 916p_denominator(value v1, type t1, value v, type t) 917{ 918 return unary_arith_op(v1, t1, v, t, ARITH_DEN, TRAT); 919} 920 921static int 922p_pi(value v, type t) 923{ 924 pword result; 925 Make_Float(&result, 3.1415926535897932) 926 Return_Numeric(v, t, result) 927} 928 929static int 930p_e(value v, type t) 931{ 932 pword result; 933 Make_Float(&result, 2.7182818284590455) 934 Return_Numeric(v, t, result) 935} 936 937/*------------------------------------------------------------------------ 938 * Irregular arithmetic Bips 939 * They have special rules concerning their argument and result types 940 *-----------------------------------------------------------------------*/ 941 942static int 943p_rational2(value v1, type t1, value v, type t) 944{ 945 int err; 946 pword result; 947 if (IsRef(t1)) { Bip_Error(PDELAY_1) } 948 result.tag.kernel = TRAT; 949 err = tag_desc[TagType(t1)].coerce_to[TRAT](v1, &result.val); 950 if (err != PSUCCEED) return err; 951 Return_Numeric(v, t, result) 952} 953 954static int 955p_rationalize(value v1, type t1, value v, type t) 956{ 957 int err; 958 pword result; 959 if (IsRef(t1)) { Bip_Error(PDELAY_1) } 960 err = tag_desc[TagType(t1)].arith_op[ARITH_NICERAT](v1, &result); 961 if (err != PSUCCEED) return err; 962 Return_Numeric(v, t, result) 963} 964 965static int 966p_bignum2(value v1, type t1, value v, type t) 967{ 968 int err; 969 pword result; 970 if (IsRef(t1)) { Bip_Error(PDELAY_1) } 971 result.tag.kernel = TBIG; 972 err = tag_desc[TagType(t1)].coerce_to[TBIG](v1, &result.val); 973 /* We fail here instead of making a type error. Thus bignum/2 can be 974 * used more easily to test whether bignums are available or not */ 975 if (err != PSUCCEED) { Fail_; } 976 Return_Numeric(v, t, result) 977} 978 979static int 980p_breal2(value v1, type t1, value v, type t) 981{ 982 int err; 983 pword result; 984 if (IsRef(t1)) { Bip_Error(PDELAY_1) } 985 result.tag.kernel = TIVL; 986 err = tag_desc[TagType(t1)].coerce_to[TIVL](v1, &result.val); 987 if (err != PSUCCEED) return err; 988 Return_Numeric(v, t, result) 989} 990 991static int 992p_float2(value v1, type t1, value v, type t) 993{ 994 pword result; 995 int err; 996 if (IsRef(t1)) { Bip_Error(PDELAY_1) } 997 result.tag.kernel = TDBL; 998 err = tag_desc[TagType(t1)].coerce_to[TagType(result.tag)](v1, &result.val); 999 if (err != PSUCCEED) return(err); 1000 Return_Numeric(v, t, result) 1001} 1002 1003/* auxiliary function for power operation */ 1004 1005#if (SIZEOF_WORD == 4) 1006#define MAX_SQUAREABLE_WORD 46340 1007#else 1008#if (SIZEOF_WORD == 8) 1009#define MAX_SQUAREABLE_WORD 3037000499UL 1010#else 1011PROBLEM: Cannot deal with word size SIZEOF_WORD. 1012#endif 1013#endif 1014 1015static int 1016_int_pow(word x, 1017 uword y, /* y >= 0 */ 1018 word *r) 1019{ 1020 word result; 1021 int neg = 0; 1022 1023 if (y <= 1) { 1024 /* 0^0 is 1 to be consistent with float pow() */ 1025 *r = (y == 1) ? x : 1; /* x^0 x^1 */ 1026 return PSUCCEED; 1027 } 1028 if (x <= 1) { 1029 if (x >= -1) { 1030 *r = (x >= 0 || (y&1)) ? x : 1; /* -1^y 0^y 1^y */ 1031 return PSUCCEED; 1032 } 1033 if (x == MIN_S_WORD) return INTEGER_OVERFLOW; 1034 x = -x; 1035 if (y & 1) neg = 1; 1036 } 1037 1038 /* Now x>=2 and y>=2 */ 1039 while(!(y&1)) { 1040 if (x > MAX_SQUAREABLE_WORD) 1041 return INTEGER_OVERFLOW; 1042 x *= x; 1043 y >>= 1; 1044 } 1045 result = x; 1046 while(y>>=1) { 1047 if (x > MAX_SQUAREABLE_WORD) 1048 return INTEGER_OVERFLOW; 1049 x *= x; 1050 if (y&1) { 1051 word tmp = result*x; 1052 if (tmp/x != result) 1053 return INTEGER_OVERFLOW; 1054 result = tmp; 1055 } 1056 } 1057 *r = neg ? -result : result; 1058 return PSUCCEED; 1059} 1060 1061 1062static int 1063p_power(value v1, type t1, value v2, type t2, value v, type t) 1064{ 1065 pword result; 1066 int err; 1067 1068 if(IsInteger(t2)) 1069 { 1070 if (IsInteger(t1)) 1071 { 1072 if (v2.nint < 0) 1073 { 1074 if (GlobalFlags & PREFER_RATIONALS) 1075 { 1076 /* this will force bignum ^ int -> rational */ 1077 Bip_Error(INTEGER_OVERFLOW) 1078 } 1079 else 1080 { 1081 Make_Checked_Float(&result, 1082 SafePow((double)(v1.nint),(double)v2.nint)); 1083 } 1084 } 1085 else 1086 { 1087 result.tag.kernel = TINT; 1088 err = _int_pow(v1.nint, (uword) v2.nint, &result.val.nint); 1089 if (err) { Bip_Error(err); } 1090 } 1091 } 1092 else if (IsBignum(t1)) 1093 { 1094 err = tag_desc[TBIG].arith_op[ARITH_POW](v1, v2, &result); 1095 if (err) { Bip_Error(err); } 1096 } 1097 else if (IsRational(t1)) 1098 { 1099 err = tag_desc[TRAT].arith_op[ARITH_POW](v1, v2, &result); 1100 if (err) { Bip_Error(err); } 1101 } 1102 else 1103 { 1104 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_POW); 1105 } 1106 } 1107 else if(IsBignum(t2)) 1108 { 1109 Bip_Error(RANGE_ERROR) 1110 } 1111 else if (IsRational(t2) && !IsDouble(t1)) 1112 { 1113 Bip_Error(TYPE_ERROR) 1114 } 1115 else /* default also handles delay */ 1116 { 1117 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_POW); 1118 } 1119 Kill_DE; 1120 Return_Numeric(v, t, result) 1121} 1122 1123 1124static int 1125p_powm(value v1, type t1, value v2, type t2, value v3, type t3, value v, type t) 1126{ 1127 pword result; 1128 int err; 1129 if (IsRef(t1)) { Bip_Error(PDELAY_1) } 1130 if (IsRef(t2)) { Bip_Error(PDELAY_2) } 1131 if (IsRef(t3)) { Bip_Error(PDELAY_3) } 1132 err = tag_desc[TagType(t1)].coerce_to[TBIG](v1, &v1); 1133 if (err != PSUCCEED) return(err); 1134 err = tag_desc[TagType(t2)].coerce_to[TBIG](v2, &v2); 1135 if (err != PSUCCEED) return(err); 1136 err = tag_desc[TagType(t3)].coerce_to[TBIG](v3, &v3); 1137 if (err != PSUCCEED) return(err); 1138 err = tag_desc[TBIG].arith_op[ARITH_POWM](v1, v2, v3, &result); 1139 if (err != PSUCCEED) return(err); 1140 Return_Unify_Pw(v, t, result.val, result.tag); 1141} 1142 1143 1144/* 1145 * Shifts counts that exceed the word length are undefined in C. 1146 * If we shift left, it's an overflow, otherwise undefined (change this). 1147 * In the other cases, overflow is checked by shifting back and checking 1148 * if we obtain the original value. 1149 */ 1150 1151static int 1152p_lshift(value v1, type t1, value v2, type t2, value v, type t) 1153{ 1154 pword result; 1155 int err; 1156 1157 if (IsInteger(t2)) 1158 { 1159 if (IsInteger(t1)) 1160 { 1161 if (v2.nint > 0) /* shift left */ 1162 { 1163 if (v1.nint && v2.nint >= BITS_PER_WORD) 1164 { Bip_Error(INTEGER_OVERFLOW); } 1165 result.val.nint = v1.nint << v2.nint; 1166 if (result.val.nint >> v2.nint != v1.nint) 1167 { Bip_Error(INTEGER_OVERFLOW); } 1168 } 1169 else /* shift right */ 1170 { 1171 if (-v2.nint >= BITS_PER_WORD) 1172 result.val.nint = v1.nint >> BITS_PER_WORD-1; 1173 else 1174 result.val.nint = v1.nint >> -v2.nint; 1175 } 1176 result.tag.kernel = TINT; 1177 } 1178 else if (IsBignum(t1)) 1179 { 1180 err = tag_desc[TBIG].arith_op[ARITH_SHL](v1, v2, &result); 1181 if (err != PSUCCEED) { Bip_Error(err); } 1182 } 1183 else /* error */ 1184 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_SHL); 1185 } 1186 else if (IsBignum(t2)) 1187 { 1188 if (!BigNegative(v2.ptr)) 1189 { Bip_Error(RANGE_ERROR); } 1190 if (IsInteger(t1)) 1191 { 1192 Make_Integer(&result, v1.nint < 0 ? -1 : 0); 1193 } 1194 else if (IsBignum(t1)) 1195 { 1196 Make_Integer(&result, BigNegative(v1.ptr) ? -1 : 0); 1197 } 1198 else 1199 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_SHL); 1200 } 1201 else 1202 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_SHL); 1203 1204 Kill_DE; 1205 Return_Numeric(v, t, result) 1206} 1207 1208static int 1209p_rshift(value v1, type t1, value v2, type t2, value v, type t) 1210{ 1211 pword result; 1212 int err; 1213 1214 if (IsInteger(t2)) 1215 { 1216 if (IsInteger(t1)) 1217 { 1218 if (v2.nint >= 0) /* shift right */ 1219 { 1220 if (v2.nint >= BITS_PER_WORD) 1221 result.val.nint = v1.nint >> BITS_PER_WORD-1; 1222 else 1223 result.val.nint = v1.nint >> v2.nint; 1224 } 1225 else /* shift left */ 1226 { 1227 if (v1.nint && -v2.nint >= BITS_PER_WORD) 1228 { Bip_Error(INTEGER_OVERFLOW); } 1229 result.val.nint = v1.nint << -v2.nint; 1230 if (result.val.nint >> -v2.nint != v1.nint) 1231 { Bip_Error(INTEGER_OVERFLOW); } 1232 } 1233 result.tag.kernel = TINT; 1234 } 1235 else if (IsBignum(t1)) 1236 { 1237 err = tag_desc[TBIG].arith_op[ARITH_SHR](v1, v2, &result); 1238 if (err != PSUCCEED) { Bip_Error(err); } 1239 } 1240 else 1241 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_SHR); 1242 } 1243 else if (IsBignum(t2)) 1244 { 1245 if (BigNegative(v2.ptr)) 1246 { Bip_Error(RANGE_ERROR); } 1247 if (IsInteger(t1)) 1248 { 1249 Make_Integer(&result, v1.nint < 0 ? -1 : 0); 1250 } 1251 else if (IsBignum(t1)) 1252 { 1253 Make_Integer(&result, BigNegative(v1.ptr) ? -1 : 0); 1254 } 1255 else 1256 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_SHR); 1257 } 1258 else 1259 return binary_arith_op(v1, t1, v2, t2, v, t, ARITH_SHR); 1260 1261 Kill_DE; 1262 Return_Numeric(v, t, result) 1263} 1264 1265static int 1266p_minint(value v, type t) 1267{ 1268 pword result; 1269 Make_Integer(&result, MIN_S_WORD); 1270 Return_Numeric(v, t, result) 1271} 1272 1273static int 1274p_maxint(value v, type t) 1275{ 1276 pword result; 1277 Make_Integer(&result, MAX_S_WORD); 1278 Return_Numeric(v, t, result) 1279} 1280 1281 1282static int 1283p_setbit(value vi, type ti, value vn, type tn, value v, type t) /* argument order because of overflow handler */ 1284{ 1285 int err; 1286 pword result; 1287 1288 Check_Integer(tn); 1289 if (vn.nint < 0) 1290 { 1291 Bip_Error(RANGE_ERROR); 1292 } 1293 if (IsInteger(ti)) 1294 { 1295 if (vn.nint < BITS_PER_WORD-1) 1296 { 1297 Make_Integer(&result, vi.nint | ((word)1 << vn.nint)); 1298 } 1299 else if (vi.nint < 0) 1300 { 1301 Make_Integer(&result, vi.nint); 1302 } 1303 else 1304 { 1305 Bip_Error(INTEGER_OVERFLOW); 1306 } 1307 } 1308 else if (IsBignum(ti) || IsString(ti)) 1309 { 1310 err = tag_desc[TagType(ti)].arith_op[ARITH_SETBIT](vi, vn, &result); 1311 if (err != PSUCCEED) return err; 1312 } 1313 else 1314 return binary_arith_op(vi, ti, vn, tn, v, t, ARITH_SETBIT); 1315 1316 Kill_DE; 1317 Return_Numeric(v, t, result) 1318} 1319 1320static int 1321p_clrbit(value vi, type ti, value vn, type tn, value v, type t) 1322{ 1323 int err; 1324 pword result; 1325 1326 Check_Integer(tn); 1327 if (vn.nint < 0) 1328 { 1329 Bip_Error(RANGE_ERROR); 1330 } 1331 if (IsInteger(ti)) 1332 { 1333 if (vn.nint < BITS_PER_WORD-1) 1334 { 1335 Make_Integer(&result, vi.nint & ~((word)1 << vn.nint)); 1336 } 1337 else if (vi.nint >= 0) 1338 { 1339 Make_Integer(&result, vi.nint); 1340 } 1341 else 1342 { 1343 Bip_Error(INTEGER_OVERFLOW); 1344 } 1345 } 1346 else if (IsBignum(ti) || IsString(ti)) 1347 { 1348 err = tag_desc[TagType(ti)].arith_op[ARITH_CLRBIT](vi, vn, &result); 1349 if (err != PSUCCEED) return err; 1350 } 1351 else 1352 return binary_arith_op(vi, ti, vn, tn, v, t, ARITH_CLRBIT); 1353 1354 Kill_DE; 1355 Return_Numeric(v, t, result) 1356} 1357 1358static int 1359p_getbit(value vi, type ti, value vn, type tn, value v, type t) 1360{ 1361 int err; 1362 pword result; 1363 1364 Check_Integer(tn); 1365 if (vn.nint < 0) 1366 { 1367 Bip_Error(RANGE_ERROR); 1368 } 1369 if (IsInteger(ti)) 1370 { 1371 Make_Integer(&result, 1372 vn.nint < BITS_PER_WORD ? 1373 ((uword) vi.nint >> vn.nint) & 1 : 1374 vi.nint < 0 ? 1 : 0); 1375 } 1376 else if (IsBignum(ti) || IsString(ti)) 1377 { 1378 err = tag_desc[TagType(ti)].arith_op[ARITH_GETBIT](vi, vn, &result); 1379 if (err != PSUCCEED) return err; 1380 } 1381 else 1382 return binary_arith_op(vi, ti, vn, tn, v, t, ARITH_GETBIT); 1383 1384 Kill_DE; 1385 Return_Numeric(v, t, result) 1386} 1387 1388#if 0 1389static int 1390_strg_setbit(value v1, value v2, pword *pres) /* string x int -> string */ 1391{ 1392 int words_old = 2*(BufferPwords(v1.ptr)-1); 1393 int words_offset = v2.nint / BITS_PER_WORD; 1394 int words_new = (words_old > words_offset+1) ? words_old : words_offset+1; 1395 word *from, *to; 1396 int i; 1397 1398 pres->val.ptr = TG; 1399 Push_Buffer(words_new * SIZEOF_WORD); 1400 to = (word *) BufferStart(pres->val.ptr); 1401 from = (word *) BufferStart(v1.ptr); 1402 for (i = 0; i < words_old; i++) 1403 to[i] = from[i]; 1404 for (; i < words_offset; i++) 1405 to[i] = 0; 1406 to[words_offset] |= 1 << (v2.nint % BITS_PER_WORD); 1407 Succeed_; 1408} 1409#endif 1410 1411 1412/* 1413 * integer_list(+Integer, +ChunkSizeInBits, -ListOfChunks) 1414 * 1415 * Takes a (big) integer and splits it up into a list of small integers of 1416 * ChunkSizeInBits each. The first list element contains the least 1417 * significant bits. E.g. 1418 * ?- X is 8'1234567, sepia_kernel:integer_list(X,3,L). 1419 * X = 342391 1420 * L = [7, 6, 5, 4, 3, 2, 1] 1421 * ChunkSizeInBits must not exceed the wordsize. 1422 */ 1423 1424int 1425p_integer_list(value vi, type ti, value vsz, type tsz, value v, type t) 1426{ 1427 int err; 1428 pword result; 1429 1430 Check_Integer_Or_Bignum(ti) 1431 Check_Integer(tsz) 1432 err = tag_desc[TagType(ti)].coerce_to[TBIG](vi, &vi); 1433 if (err != PSUCCEED) return(err); 1434 err = ec_big_to_chunks(vi.ptr, vsz.nint, &result); 1435 Return_If_Error(err); 1436 Return_Unify_Pw(v, t, result.val, result.tag); 1437} 1438 1439 1440/*------------------------------------------------------------------------ 1441 * Generic auxiliary functions 1442 *-----------------------------------------------------------------------*/ 1443 1444int 1445un_arith_op( 1446 value v1, type t1, /* input */ 1447 pword *result, /* output */ 1448 int op, /* operation */ 1449 int top) /* the 'minimal' type for the result */ 1450{ 1451 int err; 1452 1453 if (tag_desc[TagType(t1)].numeric < tag_desc[top].numeric) 1454 { 1455 if (!IsNumber(t1)) 1456 { Bip_Error(ARITH_TYPE_ERROR); } 1457 err = tag_desc[TagType(t1)].coerce_to[top](v1, &v1); 1458 if (err != PSUCCEED) return err; 1459 } 1460 else 1461 /* CAUTION: must strip extra tag bits, e.g. PERSISTENT */ 1462 top = TagType(t1); 1463 return tag_desc[top].arith_op[op](v1, result); 1464} 1465 1466int 1467unary_arith_op( 1468 value v1, type t1, /* input */ 1469 value v, type t, /* output */ 1470 int op, /* operation */ 1471 int top) /* the 'minimal' type for the result */ 1472{ 1473 pword result; 1474 int err; 1475 if (IsRef(t1)) { Bip_Error(PDELAY_1) } 1476 err = un_arith_op(v1, t1, &result, op, top); 1477 if (err != PSUCCEED) return err; 1478 Return_Numeric(v, t, result) 1479} 1480 1481int 1482bin_arith_op(value v1, type t1, value v2, type t2, pword *pres, int op) 1483{ 1484 int err; 1485 1486 if (!SameType(t1,t2)) 1487 { 1488 if (tag_desc[TagType(t1)].numeric > tag_desc[TagType(t2)].numeric) 1489 { 1490 if (!IsNumber(t2)) 1491 { Bip_Error(ARITH_TYPE_ERROR); } 1492 err = tag_desc[TagType(t2)].coerce_to[TagType(t1)](v2, &v2); 1493 } 1494 else 1495 { 1496 if (!IsNumber(t1)) 1497 { Bip_Error(ARITH_TYPE_ERROR); } 1498 err = tag_desc[TagType(t1)].coerce_to[TagType(t2)](v1, &v1); 1499 t1.kernel = t2.kernel; 1500 } 1501 if (err != PSUCCEED) return err; 1502 } 1503 err = tag_desc[TagType(t1)].arith_op[op](v1, v2, pres); 1504 if (err != PSUCCEED) return err; /* return with untouched *pres! */ 1505 return PSUCCEED; 1506} 1507 1508int 1509binary_arith_op(value v1, type t1, value v2, type t2, value v, type t, int op) 1510{ 1511 pword result; 1512 int err; 1513 if (IsRef(t1)) { Bip_Error(PDELAY_1) } 1514 if (IsRef(t2)) { Bip_Error(PDELAY_2) } 1515 err = bin_arith_op(v1, t1, v2, t2, &result, op); 1516 if (err != PSUCCEED) return err; 1517 Kill_DE; /* in case it's a demon */ 1518 Return_Numeric(v, t, result) 1519} 1520 1521 1522/* 1523 * Arithmetically compare two numbers 1524 * 1525 * PRE: IsNumber(t1) && IsNumber(t2) 1526 * PRE: *res must be one of BIEq,BINe,BILt,BIGt,BILe,BIGe 1527 * (needed to adjust the comparison in the breal case) 1528 * 1529 * Returns: 1530 * PSUCCEED and result [-1,0,1] in *res 1531 * PDELAY not decidable 1532 * UNIMPLEMENTED from coercion 1533 */ 1534int 1535arith_compare(value v1, type t1, value v2, type t2, int *res) 1536{ 1537 int err; 1538 pword *old_tg = TG; /* coercion may create temporaries */ 1539 if (!SameType(t1,t2)) 1540 { 1541 if (tag_desc[TagType(t1)].numeric > tag_desc[TagType(t2)].numeric) 1542 { 1543 err = tag_desc[TagType(t2)].coerce_to[TagType(t1)](v2, &v2); 1544 t2.kernel = Tag(t1.kernel); 1545 } 1546 else 1547 { 1548 err = tag_desc[TagType(t1)].coerce_to[TagType(t2)](v1, &v1); 1549 t1.kernel = Tag(t2.kernel); 1550 } 1551 if (err != PSUCCEED) return err; 1552 } 1553 err = tag_desc[TagType(t1)].arith_compare(v1, v2, res); 1554 TG = old_tg; /* coercion may have created temporaries */ 1555 return err; 1556} 1557 1558static int 1559_int_nop(value v1, pword *pres) 1560{ 1561 Make_Integer(pres, v1.nint); 1562 Succeed_; 1563} 1564 1565static int 1566_dbl_nop(value v1, pword *pres) 1567{ 1568 pres->tag.kernel = TDBL; 1569 pres->val.all = v1.all; /* double or pointer */ 1570 Succeed_; 1571} 1572 1573static int 1574_noc(value in, value *out) 1575{ 1576 *out = in; 1577 Succeed_; 1578} 1579 1580/*------------------------------------------------------------------------ 1581 * Integer operations (most of them are expanded) 1582 *-----------------------------------------------------------------------*/ 1583 1584static int 1585_compare_int(value v1, value v2) 1586{ 1587 return (v1.nint > v2.nint ? 1 : v1.nint < v2.nint ? -1 : 0); 1588} 1589 1590static int 1591_arith_compare_int(value v1, value v2, int *res) 1592{ 1593 *res = (v1.nint > v2.nint ? 1 : v1.nint < v2.nint ? -1 : 0); 1594 Succeed_; 1595} 1596 1597/* CAUTION: The code for int_add and int_mul is duplicated in the emulator */ 1598 1599static int 1600_int_add(value v1, value v2, pword *pres) /* int x int -> int/big */ 1601{ 1602 word res = v1.nint + v2.nint; 1603 if (((v1.nint >= 0) == (v2.nint >= 0)) && (v1.nint >= 0) != (res >= 0)) 1604 { 1605 Bip_Error(INTEGER_OVERFLOW) 1606 } 1607 Make_Integer(pres, res); 1608 Succeed_; 1609} 1610 1611static int 1612_int_sub(value v1, value v2, pword *pres) /* int x int -> int/big */ 1613{ 1614 word res = v1.nint - v2.nint; 1615 if (((v1.nint >= 0) != (v2.nint >= 0)) && (v1.nint >= 0) != (res >= 0)) 1616 { 1617 Bip_Error(INTEGER_OVERFLOW) 1618 } 1619 Make_Integer(pres, res); 1620 Succeed_; 1621} 1622 1623static int 1624_int_mul(value v1, value v2, pword *pres) /* int x int -> int/big */ 1625{ 1626 word n1 = v1.nint; 1627 word n2 = v2.nint; 1628 word res; 1629 1630 if (n1 == 0) { 1631 Make_Integer(pres, 0); 1632 Succeed_ 1633 } 1634 if (n2 == MIN_S_WORD) { 1635 /* Not true if n1 == 1, but who cares... */ 1636 Bip_Error(INTEGER_OVERFLOW) 1637 } 1638 res = n1 * n2; 1639 if (res / n1 != n2) { 1640 Bip_Error(INTEGER_OVERFLOW) 1641 } 1642 Make_Integer(pres, res); 1643 Succeed_; 1644} 1645 1646static int 1647_int_neg(value v1, pword *pres) /* needed in the parser to evaluate signs */ 1648{ 1649 if (v1.nint == MIN_S_WORD) 1650 { Bip_Error(INTEGER_OVERFLOW); } 1651 Make_Integer(pres, -v1.nint); 1652 Succeed_; 1653} 1654 1655static int 1656_int_sgn(value v1, pword *pres) 1657{ 1658 Make_Integer(pres, v1.nint > 0 ? 1 : v1.nint < 0 ? -1: 0); 1659 Succeed_; 1660} 1661 1662static int 1663_int_min(value v1, value v2, pword *pres) 1664{ 1665 Make_Integer(pres, v1.nint > v2.nint ? v2.nint : v1.nint); 1666 Succeed_; 1667} 1668 1669static int 1670_int_max(value v1, value v2, pword *pres) 1671{ 1672 Make_Integer(pres, v1.nint < v2.nint ? v2.nint : v1.nint); 1673 Succeed_; 1674} 1675 1676static int 1677_int_abs(value v1, pword *pres) 1678{ 1679 if (v1.nint < 0) 1680 return _int_neg(v1, pres); 1681 Make_Integer(pres, v1.nint); 1682 Succeed_; 1683} 1684 1685static int 1686_int_gcd(value v1, value v2, pword *pres) 1687{ 1688 /* No special TINT implementation, we default to the TBIG one */ 1689 int err; 1690 err = tag_desc[TINT].coerce_to[TBIG](v1, &v1); 1691 if (err != PSUCCEED) return err; 1692 err = tag_desc[TINT].coerce_to[TBIG](v2, &v2); 1693 if (err != PSUCCEED) return err; 1694 return tag_desc[TBIG].arith_op[ARITH_GCD](v1, v2, pres); 1695} 1696 1697static int 1698_int_lcm(value v1, value v2, pword *pres) 1699{ 1700 /* No special TINT implementation, we default to the TBIG one */ 1701 int err; 1702 err = tag_desc[TINT].coerce_to[TBIG](v1, &v1); 1703 if (err != PSUCCEED) return err; 1704 err = tag_desc[TINT].coerce_to[TBIG](v2, &v2); 1705 if (err != PSUCCEED) return err; 1706 return tag_desc[TBIG].arith_op[ARITH_LCM](v1, v2, pres); 1707} 1708 1709static int 1710_int_atan2(value v1, value v2, pword *pres) 1711{ 1712 Make_Double(pres, Atan2((double)v1.nint, (double)v2.nint)); 1713 Succeed_; 1714} 1715 1716 1717/*---------------------------------------------------------------------------- 1718 * Doubles 1719 *----------------------------------------------------------------------------*/ 1720 1721static int 1722_compare_dbl(value v1, value v2) 1723{ 1724 return Dbl(v1) > Dbl(v2) ? 1 1725 : Dbl(v1) < Dbl(v2) ? -1 1726 : Dbl(v1) != 0.0 ? 0 1727 : PedanticZeroLess(Dbl(v1),Dbl(v2)) ? -1 1728 : PedanticZeroLess(Dbl(v2),Dbl(v1)) ? 1 1729 : 0; 1730} 1731 1732static int 1733_arith_compare_dbl(value v1, value v2, int *res) 1734{ 1735 *res = (Dbl(v1) > Dbl(v2) ? 1 1736 : Dbl(v1) < Dbl(v2) ? -1 : 0); 1737 Succeed_; 1738} 1739 1740#ifndef UNBOXED_DOUBLES 1741static int 1742_equal_dbl(pword *pw1, pword *pw2) 1743{ 1744 /* compare the doubles bitwise (as integers) in order to be 1745 * able to distinguish negative and positive zeros */ 1746#if SIZEOF_DOUBLE == SIZEOF_WORD 1747 return BufferStart(pw1)->val.all == BufferStart(pw2)->val.all; 1748#elif SIZEOF_DOUBLE == 2*SIZEOF_WORD 1749 return BufferStart(pw1)->val.all == BufferStart(pw2)->val.all 1750 && BufferStart(pw1)->tag.all == BufferStart(pw2)->tag.all; 1751#else 1752 PROBLEM: Cannot deal with word size SIZEOF_WORD. 1753#endif 1754} 1755#endif 1756 1757 1758/* 1759 * arithmetic operations on doubles 1760 */ 1761 1762static int 1763_dbl_neg(value v1, pword *pres) /* needed in the parser to evaluate signs */ 1764{ 1765 Make_Double(pres, -Dbl(v1)) 1766 Succeed_; 1767} 1768 1769static int 1770_dbl_sgn(value v1, pword *pres) 1771{ 1772 Make_Integer(pres, Dbl(v1) == 0.0 ? 0 : Dbl(v1) > 0.0 ? 1: -1); 1773 Succeed_; 1774} 1775 1776static int 1777_dbl_min(value v1, value v2, pword *pres) 1778{ 1779 double d1 = Dbl(v1); 1780 double d2 = Dbl(v2); 1781 pres->tag.kernel = TDBL; 1782 pres->val.all = PedanticLess(d1,d2) ? v1.all : v2.all; 1783 Succeed_; 1784} 1785 1786static int 1787_dbl_max(value v1, value v2, pword *pres) 1788{ 1789 double d1 = Dbl(v1); 1790 double d2 = Dbl(v2); 1791 pres->tag.kernel = TDBL; 1792 pres->val.all = PedanticGreater(d1,d2) ? v1.all : v2.all; 1793 Succeed_; 1794} 1795 1796static int 1797_dbl_add(value v1, value v2, pword *pres) 1798{ 1799 Make_Checked_Double(pres, Dbl(v1) + Dbl(v2)) 1800 Succeed_; 1801} 1802 1803static int 1804_dbl_sub(value v1, value v2, pword *pres) 1805{ 1806 Make_Checked_Double(pres, Dbl(v1) - Dbl(v2)) 1807 Succeed_; 1808} 1809 1810static int 1811_dbl_mul(value v1, value v2, pword *pres) 1812{ 1813 Make_Checked_Double(pres, Dbl(v1) * Dbl(v2)) 1814 Succeed_; 1815} 1816 1817static int 1818_dbl_div(value v1, value v2, pword *pres) 1819{ 1820 Make_Checked_Double(pres, Dbl(v1) / Dbl(v2)) 1821 Succeed_; 1822} 1823 1824static int 1825_dbl_abs(value v1, pword *pres) 1826{ 1827 if (Dbl(v1) < 0.0) 1828 { 1829 Make_Double(pres, -Dbl(v1)) 1830 } 1831 else if (Dbl(v1) == 0.0) /* for negative zero */ 1832 { 1833 Make_Double(pres, 0.0) 1834 } 1835 else 1836 { 1837 pres->tag.kernel = TDBL; 1838 pres->val.all = v1.all; 1839 } 1840 Succeed_; 1841} 1842 1843static int 1844_dbl_sin(value v1, pword *pres) 1845{ 1846 Make_Checked_Double(pres, sin(Dbl(v1))) 1847 Succeed_; 1848} 1849 1850static int 1851_dbl_cos(value v1, pword *pres) 1852{ 1853 Make_Checked_Double(pres, cos(Dbl(v1))) 1854 Succeed_; 1855} 1856 1857static int 1858_dbl_tan(value v1, pword *pres) 1859{ 1860 Make_Checked_Double(pres, tan(Dbl(v1))) 1861 Succeed_; 1862} 1863 1864static int 1865_dbl_asin(value v1, pword *pres) 1866{ 1867 double y = Dbl(v1); 1868 if (!OneMOne(y)) 1869 { Bip_Error(ARITH_EXCEPTION) } 1870 Make_Checked_Double(pres, asin(y)) 1871 Succeed_; 1872} 1873 1874static int 1875_dbl_acos(value v1, pword *pres) 1876{ 1877 double y = Dbl(v1); 1878 if (!OneMOne(y)) 1879 { Bip_Error(ARITH_EXCEPTION) } 1880 Make_Checked_Double(pres, acos(y)) 1881 Succeed_; 1882} 1883 1884static int 1885_dbl_atan(value v1, pword *pres) 1886{ 1887 Make_Checked_Double(pres, atan(Dbl(v1))) 1888 Succeed_; 1889} 1890 1891static int 1892_dbl_atan2(value v1, value v2, pword *pres) 1893{ 1894 Make_Checked_Double(pres, Atan2(Dbl(v1), Dbl(v2))) 1895 Succeed_; 1896} 1897 1898static int 1899_dbl_exp(value v1, pword *pres) 1900{ 1901 double d = Dbl(v1); 1902 /* Catch the special cases of raising 'e' to +Inf and -Inf as some 1903 platforms give incorrect results w.r.t. IEEE754 specs */ 1904 if ( d == -HUGE_VAL ) { 1905 d = 0.0; 1906 } else if ( d == HUGE_VAL ) { 1907 /* Do nothing as e^(+Inf) = +Inf */ 1908 } else { 1909 d = exp(d); 1910 } 1911 Make_Double(pres, d) 1912 Succeed_; 1913} 1914 1915static int 1916_dbl_ln(value v1, pword *pres) 1917{ 1918 double y = Dbl(v1); 1919 if (!NonNegative(y)) 1920 { Bip_Error(ARITH_EXCEPTION); } 1921 Make_Double(pres, log(y)) 1922 Succeed_; 1923} 1924 1925static int 1926_dbl_sqrt(value v1, pword *pres) 1927{ 1928 double y = Dbl(v1); 1929 if (!NonNegative(y)) 1930 { Bip_Error(ARITH_EXCEPTION); } 1931 Make_Checked_Double(pres, sqrt(y)) 1932 Succeed_; 1933} 1934 1935static int 1936_dbl_pow(value v1, value v2, pword *pres) 1937{ 1938 Make_Checked_Double(pres, SafePow(Dbl(v1), Dbl(v2))) 1939 Succeed_; 1940} 1941 1942static int 1943_dbl_round(value v1, pword *pres) 1944{ 1945 double x; 1946#if defined(HAVE_RINT) && !defined(HP_RINT) && !defined(HAVE_RINT_BUG) 1947 x = rint(Dbl(v1)); 1948#else 1949 /* 1950 * Round to even number we are exactly in the middle. 1951 * Make sure we round to -0.0 if between -0.5 and -0.0 1952 */ 1953 x = Ceil(Dbl(v1)); 1954 if (x - Dbl(v1) > 0.5 || 1955 (x - Dbl(v1) == 0.5 && ((word)x & 1))) 1956 x -= 1.0; 1957#endif /* rint */ 1958 Make_Checked_Double(pres, x) 1959 Succeed_; 1960} 1961 1962static int 1963_dbl_floor(value v1, pword *pres) 1964{ 1965 Make_Checked_Double(pres, floor(Dbl(v1))) 1966 Succeed_; 1967} 1968 1969static int 1970_dbl_ceil(value v1, pword *pres) 1971{ 1972 Make_Checked_Double(pres, Ceil(Dbl(v1))) 1973 Succeed_; 1974} 1975 1976static int 1977_dbl_truncate(value v1, pword *pres) 1978{ 1979#ifdef HAVE_TRUNC 1980 Make_Checked_Double(pres, trunc(Dbl(v1))) 1981#else 1982 double f = Dbl(v1); 1983 Make_Checked_Double(pres, f < 0 ? Ceil(f) : floor(f)); 1984#endif 1985 Succeed_; 1986} 1987 1988 1989/* This _dbl_fix() is a simplified version of the one in bigrat.c */ 1990 1991static int 1992_dbl_fix(value v1, pword *pres) 1993{ 1994 double f = Dbl(v1); 1995 if (MIN_S_WORD_DBL <= f && f < MAX_S_WORD_1_DBL) /* fits in word? */ 1996 { 1997 Make_Integer(pres, (word) f); 1998 } 1999 else if (finite(f)) 2000 { 2001 Bip_Error(INTEGER_OVERFLOW); 2002 } 2003 else 2004 { 2005 Bip_Error(ARITH_EXCEPTION); 2006 } 2007 Succeed_; 2008} 2009 2010static int 2011_dbl_int2(value v1, pword *pres) 2012{ 2013 double f = Dbl(v1); 2014 if (MIN_S_WORD_DBL <= f && f < MAX_S_WORD_1_DBL) /* fits in word? */ 2015 { 2016 double ipart; 2017 if (modf(f, &ipart) == 0.0) 2018 { 2019 Make_Integer(pres, (word) ipart); 2020 } 2021 else 2022 { 2023 Bip_Error(ARITH_EXCEPTION); 2024 } 2025 } 2026 else if (finite(f)) 2027 { 2028 Bip_Error(INTEGER_OVERFLOW); 2029 } 2030 else 2031 { 2032 Bip_Error(ARITH_EXCEPTION); 2033 } 2034 Succeed_; 2035} 2036 2037/*-------------------------------------------------------------------------- 2038 * type coercion functions 2039 *--------------------------------------------------------------------------*/ 2040 2041/*ARGSUSED*/ 2042static int 2043_arith_type_err(value in, value *out) /* CAUTION: we allow out == &in */ 2044{ 2045 return ARITH_TYPE_ERROR; 2046} 2047 2048/*ARGSUSED*/ 2049static int 2050_type_err(value in, value *out) /* CAUTION: we allow out == &in */ 2051{ 2052 return TYPE_ERROR; 2053} 2054 2055static int 2056_unimp_err(value in, value *out) /* CAUTION: we allow out == &in */ 2057{ 2058 return UNIMPLEMENTED; 2059} 2060 2061static int 2062_int_dbl(value in, value *out) /* CAUTION: we allow out == &in */ 2063{ 2064 Make_Double_Val(*out, (double) in.nint) 2065 Succeed_; 2066} 2067 2068 2069/*-------------------------------------------------------------------------- 2070 * Initialize the arithmetic system (tables and bips) 2071 *--------------------------------------------------------------------------*/ 2072 2073void 2074bip_arith_init(int flags) 2075{ 2076 int i, j; 2077 extern void bigrat_init(void); 2078 2079 /* Initialise rounding mode information for interval arithmetic. */ 2080 init_rounding_modes(); 2081 2082 for(i=0; i <= NTYPES; i++) 2083 tag_desc[i].numeric = 0; 2084 2085 tag_desc[TINT].numeric = 1; /* mark and order the numeric types */ 2086 tag_desc[TBIG].numeric = 2; 2087 tag_desc[TRAT].numeric = 3; 2088 tag_desc[TDBL].numeric = 4; 2089 tag_desc[TIVL].numeric = 5; 2090 2091 for(i=0; i <= NTYPES; i++) 2092 { 2093 for(j=0; j <= NTYPES; j++) 2094 tag_desc[i].coerce_to[j] = 2095 tag_desc[i].numeric ? _type_err : _arith_type_err; 2096 for(j=0; j < ARITH_OPERATIONS; j++) 2097 tag_desc[i].arith_op[j] = 2098 tag_desc[i].numeric ? _type_err : _arith_type_err; 2099 if (tag_desc[i].numeric) 2100 { 2101 tag_desc[i].coerce_to[TBIG] = _unimp_err; 2102 tag_desc[i].coerce_to[TRAT] = _unimp_err; 2103 } 2104 tag_desc[i].coerce_to[i] = _noc; 2105 } 2106 2107 tag_desc[TINT].compare = _compare_int; 2108 tag_desc[TINT].arith_compare = _arith_compare_int; 2109 tag_desc[TINT].coerce_to[TDBL] = _int_dbl; 2110 tag_desc[TINT].arith_op[ARITH_PLUS] = _int_nop; 2111 tag_desc[TINT].arith_op[ARITH_CHGSIGN] = 2112 tag_desc[TINT].arith_op[ARITH_NEG] = _int_neg; 2113 tag_desc[TINT].arith_op[ARITH_ADD] = _int_add; 2114 tag_desc[TINT].arith_op[ARITH_SUB] = _int_sub; 2115 tag_desc[TINT].arith_op[ARITH_MUL] = _int_mul; 2116 tag_desc[TINT].arith_op[ARITH_MIN] = _int_min; 2117 tag_desc[TINT].arith_op[ARITH_MAX] = _int_max; 2118 tag_desc[TINT].arith_op[ARITH_GCD] = _int_gcd; 2119 tag_desc[TINT].arith_op[ARITH_LCM] = _int_lcm; 2120 tag_desc[TINT].arith_op[ARITH_ABS] = _int_abs; 2121 tag_desc[TINT].arith_op[ARITH_ROUND] = _int_nop; 2122 tag_desc[TINT].arith_op[ARITH_FLOOR] = _int_nop; 2123 tag_desc[TINT].arith_op[ARITH_CEIL] = _int_nop; 2124 tag_desc[TINT].arith_op[ARITH_TRUNCATE] = _int_nop; 2125 tag_desc[TINT].arith_op[ARITH_FIX] = _int_nop; 2126 tag_desc[TINT].arith_op[ARITH_INT] = _int_nop; 2127 tag_desc[TINT].arith_op[ARITH_SGN] = _int_sgn; 2128 tag_desc[TINT].arith_op[ARITH_ATAN2] = _int_atan2; 2129 2130 tag_desc[TDBL].compare = _compare_dbl; 2131 tag_desc[TDBL].arith_compare = _arith_compare_dbl; 2132#ifndef UNBOXED_DOUBLES 2133 tag_desc[TDBL].equal = _equal_dbl; 2134#endif 2135 tag_desc[TDBL].arith_op[ARITH_PLUS] = _dbl_nop; 2136 tag_desc[TDBL].arith_op[ARITH_CHGSIGN] = 2137 tag_desc[TDBL].arith_op[ARITH_NEG] = _dbl_neg; 2138 tag_desc[TDBL].arith_op[ARITH_ADD] = _dbl_add; 2139 tag_desc[TDBL].arith_op[ARITH_SUB] = _dbl_sub; 2140 tag_desc[TDBL].arith_op[ARITH_MUL] = _dbl_mul; 2141 tag_desc[TDBL].arith_op[ARITH_DIV] = _dbl_div; 2142 tag_desc[TDBL].arith_op[ARITH_MIN] = _dbl_min; 2143 tag_desc[TDBL].arith_op[ARITH_MAX] = _dbl_max; 2144 tag_desc[TDBL].arith_op[ARITH_ABS] = _dbl_abs; 2145 tag_desc[TDBL].arith_op[ARITH_SIN] = _dbl_sin; 2146 tag_desc[TDBL].arith_op[ARITH_COS] = _dbl_cos; 2147 tag_desc[TDBL].arith_op[ARITH_TAN] = _dbl_tan; 2148 tag_desc[TDBL].arith_op[ARITH_ASIN] = _dbl_asin; 2149 tag_desc[TDBL].arith_op[ARITH_ACOS] = _dbl_acos; 2150 tag_desc[TDBL].arith_op[ARITH_ATAN] = _dbl_atan; 2151 tag_desc[TDBL].arith_op[ARITH_ATAN2] = _dbl_atan2; 2152 tag_desc[TDBL].arith_op[ARITH_EXP] = _dbl_exp; 2153 tag_desc[TDBL].arith_op[ARITH_LN] = _dbl_ln; 2154 tag_desc[TDBL].arith_op[ARITH_SQRT] = _dbl_sqrt; 2155 tag_desc[TDBL].arith_op[ARITH_POW] = _dbl_pow; 2156 tag_desc[TDBL].arith_op[ARITH_ROUND] = _dbl_round; 2157 tag_desc[TDBL].arith_op[ARITH_FLOOR] = _dbl_floor; 2158 tag_desc[TDBL].arith_op[ARITH_FIX] = _dbl_fix; 2159 tag_desc[TDBL].arith_op[ARITH_INT] = _dbl_int2; 2160 tag_desc[TDBL].arith_op[ARITH_CEIL] = _dbl_ceil; 2161 tag_desc[TDBL].arith_op[ARITH_TRUNCATE] = _dbl_truncate; 2162 tag_desc[TDBL].arith_op[ARITH_SGN] = _dbl_sgn; 2163 2164 bigrat_init(); /* implementation of bignums and rationals */ 2165 ec_intervals_init(); /* implementation of float intervals */ 2166 2167 if (!(flags & INIT_SHARED)) 2168 return; 2169 2170 /* plus/3 and times/3 have NONVAR because the bound argument is not known */ 2171 built_in(in_dict("succ", 2), p_succ, B_UNSAFE|U_SIMPLE) 2172 -> mode = BoundArg(1, NONVAR) | BoundArg(2, NONVAR); 2173 built_in(in_dict("plus", 3), p_plus, B_UNSAFE|U_SIMPLE|PROC_DEMON) 2174 -> mode = BoundArg(1, NONVAR) | BoundArg(2, NONVAR) | 2175 BoundArg(3, NONVAR); 2176 built_in(in_dict("times", 3), p_times, B_UNSAFE|U_SIMPLE|PROC_DEMON) 2177 -> mode = BoundArg(1, NONVAR) | BoundArg(2, NONVAR) | 2178 BoundArg(3, NONVAR); 2179 2180 (void) exported_built_in(in_dict("collect", 3), p_collect, B_UNSAFE); 2181 (void) exported_built_in(in_dict("collapse_linear", 2), p_collapse_linear, B_UNSAFE); 2182 (void) b_built_in(in_dict("between", 4), p_between, d_.kernel_sepia); 2183 2184 (void) built_in(in_dict("+", 2), p_uplus, B_UNSAFE|U_SIMPLE); 2185 (void) built_in(in_dict("abs", 2), p_abs, B_UNSAFE|U_SIMPLE); 2186 (void) built_in(in_dict("^", 3), p_power, B_UNSAFE|U_SIMPLE|PROC_DEMON); 2187 (void) built_in(in_dict("<<", 3), p_lshift, B_UNSAFE|U_SIMPLE|PROC_DEMON); 2188 (void) built_in(in_dict(">>", 3), p_rshift, B_UNSAFE|U_SIMPLE|PROC_DEMON); 2189 (void) built_in(in_dict("min", 3), p_min, B_UNSAFE|U_SIMPLE|PROC_DEMON); 2190 (void) built_in(in_dict("max", 3), p_max, B_UNSAFE|U_SIMPLE|PROC_DEMON); 2191 (void) built_in(in_dict("gcd", 3), p_gcd, B_UNSAFE|U_SIMPLE|PROC_DEMON); 2192 (void) built_in(in_dict("gcd", 5), p_gcd_ext, B_UNSAFE|U_GROUND|PROC_DEMON); 2193 (void) built_in(in_dict("lcm", 3), p_lcm, B_UNSAFE|U_SIMPLE|PROC_DEMON); 2194 (void) built_in(in_dict("setbit", 3), p_setbit, B_UNSAFE|U_SIMPLE); 2195 (void) built_in(in_dict("getbit", 3), p_getbit, B_UNSAFE|U_SIMPLE); 2196 (void) built_in(in_dict("clrbit", 3), p_clrbit, B_UNSAFE|U_SIMPLE); 2197 (void) built_in(in_dict("sin", 2), p_sin, B_UNSAFE|U_SIMPLE); 2198 (void) built_in(in_dict("cos", 2), p_cos, B_UNSAFE|U_SIMPLE); 2199 (void) built_in(in_dict("tan", 2), p_tan, B_UNSAFE|U_SIMPLE); 2200 (void) built_in(in_dict("asin", 2), p_asin, B_UNSAFE|U_SIMPLE); 2201 (void) built_in(in_dict("acos", 2), p_acos, B_UNSAFE|U_SIMPLE); 2202 (void) built_in(in_dict("atan", 2), p_atan, B_UNSAFE|U_SIMPLE); 2203 (void) built_in(in_dict("atan", 3), p_atan2, B_UNSAFE|U_SIMPLE|PROC_DEMON); 2204 (void) built_in(in_dict("exp", 2), p_exp, B_UNSAFE|U_SIMPLE); 2205 (void) built_in(in_dict("ln", 2), p_ln, B_UNSAFE|U_SIMPLE); 2206 (void) built_in(in_dict("sqrt", 2), p_sqrt, B_UNSAFE|U_SIMPLE); 2207 (void) built_in(in_dict("round", 2), p_round, B_UNSAFE|U_SIMPLE); 2208 (void) built_in(in_dict("floor", 2), p_floor, B_UNSAFE|U_SIMPLE); 2209 (void) built_in(in_dict("ceiling", 2), p_ceiling, B_UNSAFE|U_SIMPLE); 2210 (void) built_in(in_dict("truncate", 2), p_truncate, B_UNSAFE|U_SIMPLE); 2211 (void) built_in(in_dict("numerator", 2), p_numerator, B_UNSAFE|U_SIMPLE); 2212 (void) built_in(in_dict("denominator", 2), p_denominator,B_UNSAFE|U_SIMPLE); 2213 (void) built_in(in_dict("sgn", 2), p_sgn, B_UNSAFE|U_SIMPLE); 2214 (void) local_built_in(in_dict("pi", 1), p_pi, B_UNSAFE|U_SIMPLE); 2215 (void) local_built_in(in_dict("e", 1), p_e, B_UNSAFE|U_SIMPLE); 2216 2217 (void) built_in(in_dict("fix", 2), p_fix, B_UNSAFE|U_SIMPLE); 2218 (void) built_in(in_dict("integer", 2), p_integer2, B_UNSAFE|U_SIMPLE); 2219 (void) built_in(in_dict("float", 2), p_float2, B_UNSAFE|U_SIMPLE); 2220 (void) built_in(in_dict("rational", 2), p_rational2, B_UNSAFE|U_SIMPLE); 2221 (void) built_in(in_dict("rationalize", 2), p_rationalize, B_UNSAFE|U_SIMPLE); 2222 (void) exported_built_in(in_dict("minint", 1), p_minint, B_UNSAFE|U_SIMPLE); 2223 (void) exported_built_in(in_dict("maxint", 1), p_maxint, B_UNSAFE|U_SIMPLE); 2224 (void) exported_built_in(in_dict("bignum", 2), p_bignum2,B_UNSAFE|U_SIMPLE); 2225 (void) exported_built_in(in_dict("breal", 2), p_breal2,B_UNSAFE|U_SIMPLE); 2226 (void) exported_built_in(in_dict("is_zero", 1), p_is_zero,B_SAFE); 2227 (void) exported_built_in(in_dict("integer_list", 3), p_integer_list,B_UNSAFE|U_SIMPLE); 2228 2229 (void) exported_built_in(in_dict("powm", 4), p_powm, B_UNSAFE|U_SIMPLE|PROC_DEMON); 2230} 2231 2232 2233/* At least on SUNs, defining matherr() returning non-zero will 2234 * suppress error messages being printed by math library routines. 2235 */ 2236#ifdef HAVE_MATHERR 2237/*ARGSUSED*/ 2238int 2239matherr(struct exception *ex) 2240{ 2241 return 1; 2242} 2243#endif 2244