contrib/bc/gen/lib2.bc

   1 /*
   2  * *****************************************************************************
   3  *
   4  * SPDX-License-Identifier: BSD-2-Clause
   5  *
   6  * Copyright (c) 2018-2021 Gavin D. Howard and contributors.
   7  *
   8  * Redistribution and use in source and binary forms, with or without
   9  * modification, are permitted provided that the following conditions are met:
  10  *
  11  * * Redistributions of source code must retain the above copyright notice, this
  12  *   list of conditions and the following disclaimer.
  13  *
  14  * * Redistributions in binary form must reproduce the above copyright notice,
  15  *   this list of conditions and the following disclaimer in the documentation
  16  *   and/or other materials provided with the distribution.
  17  *
  18  * THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
  19  * AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
  20  * IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
  21  * ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT HOLDER OR CONTRIBUTORS BE
  22  * LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR
  23  * CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF
  24  * SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS
  25  * INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN
  26  * CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE)
  27  * ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
  28  * POSSIBILITY OF SUCH DAMAGE.
  29  *
  30  * *****************************************************************************
  31  *
  32  * The second bc math library.
  33  *
  34  */
  35
  36 define p(x,y){
  37         auto a
  38         a=y$
  39         if(y==a)return (x^a)@scale
  40         return e(y*l(x))
  41 }
  42 define r(x,p){
  43         auto t,n
  44         if(x==0)return x
  45         p=abs(p)$
  46         n=(x<0)
  47         x=abs(x)
  48         t=x@p
  49         if(p<scale(x)&&x-t>=5>>p+1)t+=1>>p
  50         if(n)t=-t
  51         return t
  52 }
  53 define ceil(x,p){
  54         auto t,n
  55         if(x==0)return x
  56         p=abs(p)$
  57         n=(x<0)
  58         x=abs(x)
  59         t=(x+((x@p<x)>>p))@p
  60         if(n)t=-t
  61         return t
  62 }
  63 define f(n){
  64         auto r
  65         n=abs(n)$
  66         for(r=1;n>1;--n)r*=n
  67         return r
  68 }
  69 define perm(n,k){
  70         auto f,g,s
  71         if(k>n)return 0
  72         n=abs(n)$
  73         k=abs(k)$
  74         f=f(n)
  75         g=f(n-k)
  76         s=scale
  77         scale=0
  78         f/=g
  79         scale=s
  80         return f
  81 }
  82 define comb(n,r){
  83         auto s,f,g,h
  84         if(r>n)return 0
  85         n=abs(n)$
  86         r=abs(r)$
  87         s=scale
  88         scale=0
  89         f=f(n)
  90         h=f(r)
  91         g=f(n-r)
  92         f/=h*g
  93         scale=s
  94         return f
  95 }
  96 define log(x,b){
  97         auto p,s
  98         s=scale
  99         if(scale<K)scale=K
 100         if(scale(x)>scale)scale=scale(x)
 101         scale*=2
 102         p=l(x)/l(b)
 103         scale=s
 104         return p@s
 105 }
 106 define l2(x){return log(x,2)}
 107 define l10(x){return log(x,A)}
 108 define root(x,n){
 109         auto s,m,r,q,p
 110         if(n<0)sqrt(n)
 111         n=n$
 112         if(n==0)x/n
 113         if(x==0||n==1)return x
 114         if(n==2)return sqrt(x)
 115         s=scale
 116         scale=0
 117         if(x<0&&n%2==0)sqrt(x)
 118         scale=s+2
 119         m=(x<0)
 120         x=abs(x)
 121         p=n-1
 122         q=10^ceil((length(x$)/n)$,0)
 123         while(r!=q){
 124                 r=q
 125                 q=(p*r+x/r^p)/n
 126         }
 127         if(m)r=-r
 128         scale=s
 129         return r@s
 130 }
 131 define cbrt(x){return root(x,3)}
 132 define pi(s){
 133         auto t,v
 134         if(s==0)return 3
 135         s=abs(s)$
 136         t=scale
 137         scale=s+1
 138         v=4*a(1)
 139         scale=t
 140         return v@s
 141 }
 142 define t(x){
 143         auto s,c,l
 144         l=scale
 145         scale+=2
 146         s=s(x)
 147         c=c(x)
 148         scale=l
 149         return s/c
 150 }
 151 define a2(y,x){
 152         auto a,p
 153         if(!x&&!y)y/x
 154         if(x<=0){
 155                 p=pi(scale+2)
 156                 if(y<0)p=-p
 157         }
 158         if(x==0)a=p/2
 159         else{
 160                 scale+=2
 161                 a=a(y/x)+p
 162                 scale-=2
 163         }
 164         return a@scale
 165 }
 166 define sin(x){return s(x)}
 167 define cos(x){return c(x)}
 168 define atan(x){return a(x)}
 169 define tan(x){return t(x)}
 170 define atan2(y,x){return a2(y,x)}
 171 define r2d(x){
 172         auto r,i,s
 173         s=scale
 174         scale+=5
 175         i=ibase
 176         ibase=A
 177         r=x*180/pi(scale)
 178         ibase=i
 179         scale=s
 180         return r@s
 181 }
 182 define d2r(x){
 183         auto r,i,s
 184         s=scale
 185         scale+=5
 186         i=ibase
 187         ibase=A
 188         r=x*pi(scale)/180
 189         ibase=i
 190         scale=s
 191         return r@s
 192 }
 193 define frand(p){
 194         p=abs(p)$
 195         return irand(10^p)>>p
 196 }
 197 define ifrand(i,p){return irand(abs(i)$)+frand(p)}
 198 define srand(x){
 199         if(irand(2))return -x
 200         return x
 201 }
 202 define brand(){return irand(2)}
 203 define void output(x,b){
 204         auto c
 205         c=obase
 206         obase=b
 207         x
 208         obase=c
 209 }
 210 define void hex(x){output(x,G)}
 211 define void binary(x){output(x,2)}
 212 define ubytes(x){
 213         auto p,b,i
 214         b=ibase
 215         ibase=A
 216         x=abs(x)$
 217         i=2^8
 218         for(p=1;i-1<x;p*=2){i*=i}
 219         ibase=b
 220         return p
 221 }
 222 define sbytes(x){
 223         auto p,b,n,z
 224         z=(x<0)
 225         x=abs(x)
 226         x=x$
 227         n=ubytes(x)
 228         b=ibase
 229         ibase=A
 230         p=2^(n*8-1)
 231         if(x>p||(!z&&x==p))n*=2
 232         ibase=b
 233         return n
 234 }
 235 define void output_byte(x,i){
 236         auto j,p,y,b
 237         j=ibase
 238         ibase=A
 239         s=scale
 240         scale=0
 241         x=abs(x)$
 242         b=x/(2^(i*8))
 243         b%=2^8
 244         y=log(256,obase)
 245         if(b>1)p=log(b,obase)+1
 246         else p=b
 247         for(i=y-p;i>0;--i)print 0
 248         if(b)print b
 249         scale=s
 250         ibase=j
 251 }
 252 define void output_uint(x,n){
 253         auto i,b
 254         b=ibase
 255         ibase=A
 256         for(i=n-1;i>=0;--i){
 257                 output_byte(x,i)
 258                 if(i)print" "
 259                 else print"\n"
 260         }
 261         ibase=b
 262 }
 263 define void hex_uint(x,n){
 264         auto o
 265         o=obase
 266         obase=G
 267         output_uint(x,n)
 268         obase=o
 269 }
 270 define void binary_uint(x,n){
 271         auto o
 272         o=obase
 273         obase=2
 274         output_uint(x,n)
 275         obase=o
 276 }
 277 define void uintn(x,n){
 278         if(scale(x)){
 279                 print"Error: ",x," is not an integer.\n"
 280                 return
 281         }
 282         if(x<0){
 283                 print"Error: ",x," is negative.\n"
 284                 return
 285         }
 286         if(x>=2^(n*8)){
 287                 print"Error: ",x," cannot fit into ",n," unsigned byte(s).\n"
 288                 return
 289         }
 290         binary_uint(x,n)
 291         hex_uint(x,n)
 292 }
 293 define void intn(x,n){
 294         auto t
 295         if(scale(x)){
 296                 print"Error: ",x," is not an integer.\n"
 297                 return
 298         }
 299         t=2^(n*8-1)
 300         if(abs(x)>=t&&(x>0||x!=-t)){
 301                 print "Error: ",x," cannot fit into ",n," signed byte(s).\n"
 302                 return
 303         }
 304         if(x<0)x=2^(n*8)-(-x)
 305         binary_uint(x,n)
 306         hex_uint(x,n)
 307 }
 308 define void uint8(x){uintn(x,1)}
 309 define void int8(x){intn(x,1)}
 310 define void uint16(x){uintn(x,2)}
 311 define void int16(x){intn(x,2)}
 312 define void uint32(x){uintn(x,4)}
 313 define void int32(x){intn(x,4)}
 314 define void uint64(x){uintn(x,8)}
 315 define void int64(x){intn(x,8)}
 316 define void uint(x){uintn(x,ubytes(x))}
 317 define void int(x){intn(x,sbytes(x))}