**rag@entropy.ms.washington.edu (David Ragozin)** (05/02/87)

Several notes about bc (dc) bugs. A) A simplified version of someones earlier sqrt bug program also shows the division problem, with integer arithmetic: a=10^18-10^12+1 #(12 can be replaced by 9 through 15) b=a*a b/a a sqrt(b) All three displayed numbers should be equal, but are not. First two should be equal even with scale=0 as the division should be exact. Failure of third to be correct could be attributed to scale problems, but see B) below. B) Since sqrt is computed using newton's method, it involves divisions, so can be expected to fail if division fails. Now for the "strange" part: If you prefix the previous program with scale=2n+1 (i.e. any odd integer) a-sqrt(b) computes as 0, while if scale=2n (even scale) then sqrt(b) computes as different from a! (b/a computes as unequal to a no matter what scale is set.) C) Above results hold on Vax11/750,785's under BSD4.3 and on 3b2/400 under SysVRel3 and ATT UNIX PC under 3.5. D) A TOTALLY different (trivial) bug occurs on all machines using /usr/lib/lib.b distributed by BSD labelled 4.1, dated 83/04/02. The code for a(x) (arctan of x) is missing braces, causing it to return pi/4 (up to 51 digits) independent of actual value of x. (This bug does not occur in the lib.b distributed with Sys V). A correct version of the offending lines: if(x==1) { if(scale<52) { return(.7853981633974483096156608458198757210492923498437764/1) }} simply adds the two sets of braces. Unfortunately, until the division bug is fixed, all bc library computations and sqrt's are suspect. Does anyone know of a domain of correct operation? David Ragozin, Department of Mathematics, GN-50, University of Washington, Seattle, WA. 98195. (206-543-1148) rag@uw.ms.entropy.edu or ...!uw-beaver!uw-entropy!rag

**woetzel@gmdzi.UUCP (Gerd Woetzel)** (05/05/87)

The division errors of bc(1) are caused by the dc(1) command, which is used indirectly by bc. Here are two fixes to dc.c made by Heinrich Schueth (schue%gmdzi@unido.UUCP). The first one addresses an error in the scale computation (line 60ff), the second is a correction to the divide algorithm (line 665ff). Please send only credits to Heinrich, and no flames to me. After the following diffs you'll find a little shell script to test your implementation of bc/dc. *** dc.c.OLD Tue May 5 13:27:28 1987 --- dc.c Tue May 5 14:13:14 1987 *************** *** 60,67 binop('*'); p = pop(); sunputc(p); ! savk = sk1+sk2; ! if(savk>k && savk>sk1 && savk>sk2){ sk = sk1; if(sk<sk2)sk = sk2; if(sk<k)sk = k; --- 60,67 ----- binop('*'); p = pop(); sunputc(p); ! savk = n = sk1+sk2; ! if(n>k && n>sk1 && n>sk2){ sk = sk1; if(sk<sk2)sk = sk2; if(sk<k)sk = k; *************** *** 65,71 sk = sk1; if(sk<sk2)sk = sk2; if(sk<k)sk = k; ! p = removc(p,savk-sk); savk = sk; } sputc(p,savk); --- 65,71 ----- sk = sk1; if(sk<sk2)sk = sk2; if(sk<k)sk = k; ! p = removc(p,n-sk); savk = sk; } sputc(p,savk); *************** *** 655,661 sbackc(p); salterc(p,dig); sbackc(p); ! if(--offset >= 0)divd->wt--; } if(divcarry != 0){ salterc(p,dig-1); --- 655,668 ----- sbackc(p); salterc(p,dig); sbackc(p); ! if(--offset >= 0){ ! if(d > 0){ ! sbackc(divd); ! dd=sbackc(divd); ! salterc(divd,dd+100); ! } ! divd->wt--; ! } } if(divcarry != 0){ salterc(p,dig-1); ----------------- Test script --- Cut here ------------------- # Error in division. Found by schue%gmdzi@unido.UUCP (Heinrich Schueth) # Fixed in dc.c line 665ff (see diffs above)" echo "Test 1 (division): Result of expression is" bc <<ENDSCRIPT1 scale=0 a=2799912455811601 b=99990001 (a/b)*b + (a%b) - a ENDSCRIPT1 echo "Result should be" echo "0" echo "" # Errors found by ronald@kulcs.UUCP (Ronald Cools) # Same fix echo "Test 2 (division): Results of divisions are" bc <<ENDSCRIPT2 scale=45 a=999.999000000000000000000167001 1/a 10/a 100/a 1000/a 10000/a ENDSCRIPT2 echo "Results should be" echo ".001000001000001000001000000832999665998498997" echo ".010000010000010000010000008329996659984989973" echo "0.100000100000100000100000083299966599849899733" echo "1.000001000001000001000000832999665998498997331" echo "10.000010000010000010000008329996659984989973319" echo "" # Error at computation of scale. Found by Heinrich Schueth" # Fix in dc.c line 60ff (see diffs above)" echo "Test 3 (scale): Result of multiplication is" bc <<ENDSCRIPT3 scale=64 a=1/9 a*a ENDSCRIPT3 echo "Result should be" echo ".0123456790123456790123456790123456790123456790123456790123456790"