← Back to branch summary

~ubuntu-fr-scripts/ufrs-math/math-dev

« back to all changes in this revision

Viewing changes to bc-lib/cf.bc

  • Committer: Jean-Michel Juzan
  • Date: 2009-04-01 12:03:01 UTC
  • Revision ID: jean-michel@juzan.org-20090401120301-j6g6xj1z0561zdg3
copie de la branche principale

Show diffs side-by-side

added added

removed removed

Lines of Context:
bc-lib/cf.bc
 
1
### CF.BC - Continued fraction experimentation for GNU bc
 
2
 
 
3
# Create a continued fraction representation of x in the global array cf[]
 
4
define cf(x) {
 
5
 auto sign, i, s3
 
6
 s3=int(scale/2);
 
7
 sign=sgn(x);x*=sign
 
8
 for(i=0;i<=s3;i++) {
 
9
  cf[i]=int(x)
 
10
  if(cf[i]>10^s3)break;
 
11
  x-=cf[i]
 
12
  if(x==0){i+=1;break}
 
13
  x=1/x
 
14
 }
 
15
 cf[i]=0
 
16
}
 
17
 
 
18
# Create a continued fraction representation of x in the global array cf[]
 
19
# using alternating signs
 
20
define cf2(x) {
 
21
 auto sign, i, s3
 
22
 s3=int(scale/2);
 
23
 sign=sgn(x);x*=sign
 
24
 for(i=0;i<=s3;i++) {
 
25
  cf[i]=1+int(x)
 
26
  if(cf[i]>10^s3)break;
 
27
  x-=cf[i]
 
28
  if(i!=2*int(i/2))cf[i]*=-1
 
29
  if(x==0){i+=1;break}
 
30
  x=1/-x
 
31
 }
 
32
 cf[i]=0
 
33
}
 
34
 
 
35
# Output the global array cf[] formatted as a continued fraction
 
36
define printcf() {
 
37
 auto i;
 
38
 print "[", cf[0], "; ";
 
39
 i=1;if(cf[i]!=0)for(i=1;cf[i+1]!=0;i++)print cf[i], ", ";
 
40
 print cf[i], "]\n";
 
41
}
 
42
 
 
43
# Convert global array cf[] back into a number
 
44
define cf2num() {
 
45
 auto n, i;
 
46
 for(i=0;cf[i]!=0;i++){}
 
47
 n=cf[--i]
 
48
 for(;i>=0;i--)n=1/n+cf[i]
 
49
 return(n);
 
50
}
 
51
 
 
52
# Turn the alternating signed continued fraction representation of x back
 
53
# into a number by first taking the absolute value of the chain
 
54
define cfflip(x) {
 
55
 cf2(x)
 
56
 for(i=0;cf[i]!=0;i++)cf[i]=abs(cf[i])
 
57
 return(cf2num());
 
58
}
 
59
 
 
60
# Turn the alternating signed continued fraction representation of x back
 
61
# into a number by first taking the absolute value less one of the chain
 
62
define cfflip1(x) {
 
63
 cf2(x)
 
64
 for(i=0;cf[i]!=0;i++)cf[i]=abs(cf[i])-1
 
65
 return(cf2num());
 
66
}
 
67
 
 
68
# Turn the binary representation of x into a continued fraction-like 
 
69
# structure using global array cf[].
 
70
# e.g. 0.1001000011111101110111 -> [0;1,2,1,4,6,1,3,1,3]
 
71
define bincf(x) {
 
72
 auto i,b,bb,n,j;
 
73
 x=abs(x);if(x>1)x=frac(x)
 
74
 x*=2;b=int(x);x-=b;n=1;j=1;cf[0]=0
 
75
 for(i=0;i<scale;i++){
 
76
  x*=2;bb=int(x)
 
77
  if(bb==b){n+=1}else{cf[j]=n;j+=1;n=1}
 
78
  b=bb;x-=b
 
79
 }
 
80
 if(n!=0){cf[j]=n;j+=1}
 
81
 cf[j-1]*=b
 
82
}
 
83
 
 
84
# Read a continued fraction as a representation of alternating groups of 
 
85
# bits in a binary number.
 
86
# e.g. [1;4,1,4,1,4,1,4,1...] -> 0.11110111101111011110...
 
87
define cfbin() {
 
88
 auto n,i,b,x;
 
89
 #cf[0]=0;
 
90
 for(i=1;cf[i]!=0;i++){}
 
91
 n=cf[--i];b=1;x=0
 
92
 for(;i>=0;i--){
 
93
  for(j=0;j<n;j++)x=x/2+b
 
94
  b=1-b
 
95
 }
 
96
 return(x)
 
97
}