gcc/testsuite/go.test/test/chan/powser1.go

   1 // run
   2
   3 // Copyright 2009 The Go Authors. All rights reserved.
   4 // Use of this source code is governed by a BSD-style
   5 // license that can be found in the LICENSE file.
   6
   7 // Test concurrency primitives: power series.
   8
   9 // Power series package
  10 // A power series is a channel, along which flow rational
  11 // coefficients.  A denominator of zero signifies the end.
  12 // Original code in Newsqueak by Doug McIlroy.
  13 // See Squinting at Power Series by Doug McIlroy,
  14 //   http://www.cs.bell-labs.com/who/rsc/thread/squint.pdf
  15
  16 package main
  17
  18 import "os"
  19
  20 type rat struct  {
  21         num, den  int64 // numerator, denominator
  22 }
  23
  24 func (u rat) pr() {
  25         if u.den==1 {
  26                 print(u.num)
  27         } else {
  28                 print(u.num, "/", u.den)
  29         }
  30         print(" ")
  31 }
  32
  33 func (u rat) eq(c rat) bool {
  34         return u.num == c.num && u.den == c.den
  35 }
  36
  37 type dch struct {
  38         req chan  int
  39         dat chan  rat
  40         nam int
  41 }
  42
  43 type dch2 [2] *dch
  44
  45 var chnames string
  46 var chnameserial int
  47 var seqno int
  48
  49 func mkdch() *dch {
  50         c := chnameserial % len(chnames)
  51         chnameserial++
  52         d := new(dch)
  53         d.req = make(chan int)
  54         d.dat = make(chan rat)
  55         d.nam = c
  56         return d
  57 }
  58
  59 func mkdch2() *dch2 {
  60         d2 := new(dch2)
  61         d2[0] = mkdch()
  62         d2[1] = mkdch()
  63         return d2
  64 }
  65
  66 // split reads a single demand channel and replicates its
  67 // output onto two, which may be read at different rates.
  68 // A process is created at first demand for a rat and dies
  69 // after the rat has been sent to both outputs.
  70
  71 // When multiple generations of split exist, the newest
  72 // will service requests on one channel, which is
  73 // always renamed to be out[0]; the oldest will service
  74 // requests on the other channel, out[1].  All generations but the
  75 // newest hold queued data that has already been sent to
  76 // out[0].  When data has finally been sent to out[1],
  77 // a signal on the release-wait channel tells the next newer
  78 // generation to begin servicing out[1].
  79
  80 func dosplit(in *dch, out *dch2, wait chan int ) {
  81         both := false   // do not service both channels
  82
  83         select {
  84         case <-out[0].req:
  85
  86         case <-wait:
  87                 both = true
  88                 select {
  89                 case <-out[0].req:
  90
  91                 case <-out[1].req:
  92                         out[0], out[1] = out[1], out[0]
  93                 }
  94         }
  95
  96         seqno++
  97         in.req <- seqno
  98         release := make(chan  int)
  99         go dosplit(in, out, release)
 100         dat := <-in.dat
 101         out[0].dat <- dat
 102         if !both {
 103                 <-wait
 104         }
 105         <-out[1].req
 106         out[1].dat <- dat
 107         release <- 0
 108 }
 109
 110 func split(in *dch, out *dch2) {
 111         release := make(chan int)
 112         go dosplit(in, out, release)
 113         release <- 0
 114 }
 115
 116 func put(dat rat, out *dch) {
 117         <-out.req
 118         out.dat <- dat
 119 }
 120
 121 func get(in *dch) rat {
 122         seqno++
 123         in.req <- seqno
 124         return <-in.dat
 125 }
 126
 127 // Get one rat from each of n demand channels
 128
 129 func getn(in []*dch) []rat {
 130         n := len(in)
 131         if n != 2 { panic("bad n in getn") }
 132         req := new([2] chan int)
 133         dat := new([2] chan rat)
 134         out := make([]rat, 2)
 135         var i int
 136         var it rat
 137         for i=0; i<n; i++ {
 138                 req[i] = in[i].req
 139                 dat[i] = nil
 140         }
 141         for n=2*n; n>0; n-- {
 142                 seqno++
 143
 144                 select {
 145                 case req[0] <- seqno:
 146                         dat[0] = in[0].dat
 147                         req[0] = nil
 148                 case req[1] <- seqno:
 149                         dat[1] = in[1].dat
 150                         req[1] = nil
 151                 case it = <-dat[0]:
 152                         out[0] = it
 153                         dat[0] = nil
 154                 case it = <-dat[1]:
 155                         out[1] = it
 156                         dat[1] = nil
 157                 }
 158         }
 159         return out
 160 }
 161
 162 // Get one rat from each of 2 demand channels
 163
 164 func get2(in0 *dch, in1 *dch) []rat {
 165         return getn([]*dch{in0, in1})
 166 }
 167
 168 func copy(in *dch, out *dch) {
 169         for {
 170                 <-out.req
 171                 out.dat <- get(in)
 172         }
 173 }
 174
 175 func repeat(dat rat, out *dch) {
 176         for {
 177                 put(dat, out)
 178         }
 179 }
 180
 181 type PS *dch    // power series
 182 type PS2 *[2] PS // pair of power series
 183
 184 var Ones PS
 185 var Twos PS
 186
 187 func mkPS() *dch {
 188         return mkdch()
 189 }
 190
 191 func mkPS2() *dch2 {
 192         return mkdch2()
 193 }
 194
 195 // Conventions
 196 // Upper-case for power series.
 197 // Lower-case for rationals.
 198 // Input variables: U,V,...
 199 // Output variables: ...,Y,Z
 200
 201 // Integer gcd; needed for rational arithmetic
 202
 203 func gcd (u, v int64) int64 {
 204         if u < 0 { return gcd(-u, v) }
 205         if u == 0 { return v }
 206         return gcd(v%u, u)
 207 }
 208
 209 // Make a rational from two ints and from one int
 210
 211 func i2tor(u, v int64) rat {
 212         g := gcd(u,v)
 213         var r rat
 214         if v > 0 {
 215                 r.num = u/g
 216                 r.den = v/g
 217         } else {
 218                 r.num = -u/g
 219                 r.den = -v/g
 220         }
 221         return r
 222 }
 223
 224 func itor(u int64) rat {
 225         return i2tor(u, 1)
 226 }
 227
 228 var zero rat
 229 var one rat
 230
 231
 232 // End mark and end test
 233
 234 var finis rat
 235
 236 func end(u rat) int64 {
 237         if u.den==0 { return 1 }
 238         return 0
 239 }
 240
 241 // Operations on rationals
 242
 243 func add(u, v rat) rat {
 244         g := gcd(u.den,v.den)
 245         return  i2tor(u.num*(v.den/g)+v.num*(u.den/g),u.den*(v.den/g))
 246 }
 247
 248 func mul(u, v rat) rat {
 249         g1 := gcd(u.num,v.den)
 250         g2 := gcd(u.den,v.num)
 251         var r rat
 252         r.num = (u.num/g1)*(v.num/g2)
 253         r.den = (u.den/g2)*(v.den/g1)
 254         return r
 255 }
 256
 257 func neg(u rat) rat {
 258         return i2tor(-u.num, u.den)
 259 }
 260
 261 func sub(u, v rat) rat {
 262         return add(u, neg(v))
 263 }
 264
 265 func inv(u rat) rat {   // invert a rat
 266         if u.num == 0 { panic("zero divide in inv") }
 267         return i2tor(u.den, u.num)
 268 }
 269
 270 // print eval in floating point of PS at x=c to n terms
 271 func evaln(c rat, U PS, n int) {
 272         xn := float64(1)
 273         x := float64(c.num)/float64(c.den)
 274         val := float64(0)
 275         for i:=0; i<n; i++ {
 276                 u := get(U)
 277                 if end(u) != 0 {
 278                         break
 279                 }
 280                 val = val + x * float64(u.num)/float64(u.den)
 281                 xn = xn*x
 282         }
 283         print(val, "\n")
 284 }
 285
 286 // Print n terms of a power series
 287 func printn(U PS, n int) {
 288         done := false
 289         for ; !done && n>0; n-- {
 290                 u := get(U)
 291                 if end(u) != 0 {
 292                         done = true
 293                 } else {
 294                         u.pr()
 295                 }
 296         }
 297         print(("\n"))
 298 }
 299
 300 // Evaluate n terms of power series U at x=c
 301 func eval(c rat, U PS, n int) rat {
 302         if n==0 { return zero }
 303         y := get(U)
 304         if end(y) != 0 { return zero }
 305         return add(y,mul(c,eval(c,U,n-1)))
 306 }
 307
 308 // Power-series constructors return channels on which power
 309 // series flow.  They start an encapsulated generator that
 310 // puts the terms of the series on the channel.
 311
 312 // Make a pair of power series identical to a given power series
 313
 314 func Split(U PS) *dch2 {
 315         UU := mkdch2()
 316         go split(U,UU)
 317         return UU
 318 }
 319
 320 // Add two power series
 321 func Add(U, V PS) PS {
 322         Z := mkPS()
 323         go func() {
 324                 var uv []rat
 325                 for {
 326                         <-Z.req
 327                         uv = get2(U,V)
 328                         switch end(uv[0])+2*end(uv[1]) {
 329                         case 0:
 330                                 Z.dat <- add(uv[0], uv[1])
 331                         case 1:
 332                                 Z.dat <- uv[1]
 333                                 copy(V,Z)
 334                         case 2:
 335                                 Z.dat <- uv[0]
 336                                 copy(U,Z)
 337                         case 3:
 338                                 Z.dat <- finis
 339                         }
 340                 }
 341         }()
 342         return Z
 343 }
 344
 345 // Multiply a power series by a constant
 346 func Cmul(c rat,U PS) PS {
 347         Z := mkPS()
 348         go func() {
 349                 done := false
 350                 for !done {
 351                         <-Z.req
 352                         u := get(U)
 353                         if end(u) != 0 {
 354                                 done = true
 355                         } else {
 356                                 Z.dat <- mul(c,u)
 357                         }
 358                 }
 359                 Z.dat <- finis
 360         }()
 361         return Z
 362 }
 363
 364 // Subtract
 365
 366 func Sub(U, V PS) PS {
 367         return Add(U, Cmul(neg(one), V))
 368 }
 369
 370 // Multiply a power series by the monomial x^n
 371
 372 func Monmul(U PS, n int) PS {
 373         Z := mkPS()
 374         go func() {
 375                 for ; n>0; n-- { put(zero,Z) }
 376                 copy(U,Z)
 377         }()
 378         return Z
 379 }
 380
 381 // Multiply by x
 382
 383 func Xmul(U PS) PS {
 384         return Monmul(U,1)
 385 }
 386
 387 func Rep(c rat) PS {
 388         Z := mkPS()
 389         go repeat(c,Z)
 390         return Z
 391 }
 392
 393 // Monomial c*x^n
 394
 395 func Mon(c rat, n int) PS {
 396         Z:=mkPS()
 397         go func() {
 398                 if(c.num!=0) {
 399                         for ; n>0; n=n-1 { put(zero,Z) }
 400                         put(c,Z)
 401                 }
 402                 put(finis,Z)
 403         }()
 404         return Z
 405 }
 406
 407 func Shift(c rat, U PS) PS {
 408         Z := mkPS()
 409         go func() {
 410                 put(c,Z)
 411                 copy(U,Z)
 412         }()
 413         return Z
 414 }
 415
 416 // simple pole at 1: 1/(1-x) = 1 1 1 1 1 ...
 417
 418 // Convert array of coefficients, constant term first
 419 // to a (finite) power series
 420
 421 /*
 422 func Poly(a []rat) PS {
 423         Z:=mkPS()
 424         begin func(a []rat, Z PS) {
 425                 j:=0
 426                 done:=0
 427                 for j=len(a); !done&&j>0; j=j-1)
 428                         if(a[j-1].num!=0) done=1
 429                 i:=0
 430                 for(; i<j; i=i+1) put(a[i],Z)
 431                 put(finis,Z)
 432         }()
 433         return Z
 434 }
 435 */
 436
 437 // Multiply. The algorithm is
 438 //      let U = u + x*UU
 439 //      let V = v + x*VV
 440 //      then UV = u*v + x*(u*VV+v*UU) + x*x*UU*VV
 441
 442 func Mul(U, V PS) PS {
 443         Z:=mkPS()
 444         go func() {
 445                 <-Z.req
 446                 uv := get2(U,V)
 447                 if end(uv[0])!=0 || end(uv[1]) != 0 {
 448                         Z.dat <- finis
 449                 } else {
 450                         Z.dat <- mul(uv[0],uv[1])
 451                         UU := Split(U)
 452                         VV := Split(V)
 453                         W := Add(Cmul(uv[0],VV[0]),Cmul(uv[1],UU[0]))
 454                         <-Z.req
 455                         Z.dat <- get(W)
 456                         copy(Add(W,Mul(UU[1],VV[1])),Z)
 457                 }
 458         }()
 459         return Z
 460 }
 461
 462 // Differentiate
 463
 464 func Diff(U PS) PS {
 465         Z:=mkPS()
 466         go func() {
 467                 <-Z.req
 468                 u := get(U)
 469                 if end(u) == 0 {
 470                         done:=false
 471                         for i:=1; !done; i++ {
 472                                 u = get(U)
 473                                 if end(u) != 0 {
 474                                         done = true
 475                                 } else {
 476                                         Z.dat <- mul(itor(int64(i)),u)
 477                                         <-Z.req
 478                                 }
 479                         }
 480                 }
 481                 Z.dat <- finis
 482         }()
 483         return Z
 484 }
 485
 486 // Integrate, with const of integration
 487 func Integ(c rat,U PS) PS {
 488         Z:=mkPS()
 489         go func() {
 490                 put(c,Z)
 491                 done:=false
 492                 for i:=1; !done; i++ {
 493                         <-Z.req
 494                         u := get(U)
 495                         if end(u) != 0 { done= true }
 496                         Z.dat <- mul(i2tor(1,int64(i)),u)
 497                 }
 498                 Z.dat <- finis
 499         }()
 500         return Z
 501 }
 502
 503 // Binomial theorem (1+x)^c
 504
 505 func Binom(c rat) PS {
 506         Z:=mkPS()
 507         go func() {
 508                 n := 1
 509                 t := itor(1)
 510                 for c.num!=0 {
 511                         put(t,Z)
 512                         t = mul(mul(t,c),i2tor(1,int64(n)))
 513                         c = sub(c,one)
 514                         n++
 515                 }
 516                 put(finis,Z)
 517         }()
 518         return Z
 519 }
 520
 521 // Reciprocal of a power series
 522 //      let U = u + x*UU
 523 //      let Z = z + x*ZZ
 524 //      (u+x*UU)*(z+x*ZZ) = 1
 525 //      z = 1/u
 526 //      u*ZZ + z*UU +x*UU*ZZ = 0
 527 //      ZZ = -UU*(z+x*ZZ)/u
 528
 529 func Recip(U PS) PS {
 530         Z:=mkPS()
 531         go func() {
 532                 ZZ:=mkPS2()
 533                 <-Z.req
 534                 z := inv(get(U))
 535                 Z.dat <- z
 536                 split(Mul(Cmul(neg(z),U),Shift(z,ZZ[0])),ZZ)
 537                 copy(ZZ[1],Z)
 538         }()
 539         return Z
 540 }
 541
 542 // Exponential of a power series with constant term 0
 543 // (nonzero constant term would make nonrational coefficients)
 544 // bug: the constant term is simply ignored
 545 //      Z = exp(U)
 546 //      DZ = Z*DU
 547 //      integrate to get Z
 548
 549 func Exp(U PS) PS {
 550         ZZ := mkPS2()
 551         split(Integ(one,Mul(ZZ[0],Diff(U))),ZZ)
 552         return ZZ[1]
 553 }
 554
 555 // Substitute V for x in U, where the leading term of V is zero
 556 //      let U = u + x*UU
 557 //      let V = v + x*VV
 558 //      then S(U,V) = u + VV*S(V,UU)
 559 // bug: a nonzero constant term is ignored
 560
 561 func Subst(U, V PS) PS {
 562         Z:= mkPS()
 563         go func() {
 564                 VV := Split(V)
 565                 <-Z.req
 566                 u := get(U)
 567                 Z.dat <- u
 568                 if end(u) == 0 {
 569                         if end(get(VV[0])) != 0 {
 570                                 put(finis,Z)
 571                         } else {
 572                                 copy(Mul(VV[0],Subst(U,VV[1])),Z)
 573                         }
 574                 }
 575         }()
 576         return Z
 577 }
 578
 579 // Monomial Substition: U(c x^n)
 580 // Each Ui is multiplied by c^i and followed by n-1 zeros
 581
 582 func MonSubst(U PS, c0 rat, n int) PS {
 583         Z:= mkPS()
 584         go func() {
 585                 c := one
 586                 for {
 587                         <-Z.req
 588                         u := get(U)
 589                         Z.dat <- mul(u, c)
 590                         c = mul(c, c0)
 591                         if end(u) != 0 {
 592                                 Z.dat <- finis
 593                                 break
 594                         }
 595                         for i := 1; i < n; i++ {
 596                                 <-Z.req
 597                                 Z.dat <- zero
 598                         }
 599                 }
 600         }()
 601         return Z
 602 }
 603
 604
 605 func Init() {
 606         chnameserial = -1
 607         seqno = 0
 608         chnames = "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz"
 609         zero = itor(0)
 610         one = itor(1)
 611         finis = i2tor(1,0)
 612         Ones = Rep(one)
 613         Twos = Rep(itor(2))
 614 }
 615
 616 func check(U PS, c rat, count int, str string) {
 617         for i := 0; i < count; i++ {
 618                 r := get(U)
 619                 if !r.eq(c) {
 620                         print("got: ")
 621                         r.pr()
 622                         print("should get ")
 623                         c.pr()
 624                         print("\n")
 625                         panic(str)
 626                 }
 627         }
 628 }
 629
 630 const N=10
 631 func checka(U PS, a []rat, str string) {
 632         for i := 0; i < N; i++ {
 633                 check(U, a[i], 1, str)
 634         }
 635 }
 636
 637 func main() {
 638         Init()
 639         if len(os.Args) > 1 {  // print
 640                 print("Ones: "); printn(Ones, 10)
 641                 print("Twos: "); printn(Twos, 10)
 642                 print("Add: "); printn(Add(Ones, Twos), 10)
 643                 print("Diff: "); printn(Diff(Ones), 10)
 644                 print("Integ: "); printn(Integ(zero, Ones), 10)
 645                 print("CMul: "); printn(Cmul(neg(one), Ones), 10)
 646                 print("Sub: "); printn(Sub(Ones, Twos), 10)
 647                 print("Mul: "); printn(Mul(Ones, Ones), 10)
 648                 print("Exp: "); printn(Exp(Ones), 15)
 649                 print("MonSubst: "); printn(MonSubst(Ones, neg(one), 2), 10)
 650                 print("ATan: "); printn(Integ(zero, MonSubst(Ones, neg(one), 2)), 10)
 651         } else {  // test
 652                 check(Ones, one, 5, "Ones")
 653                 check(Add(Ones, Ones), itor(2), 0, "Add Ones Ones")  // 1 1 1 1 1
 654                 check(Add(Ones, Twos), itor(3), 0, "Add Ones Twos") // 3 3 3 3 3
 655                 a := make([]rat, N)
 656                 d := Diff(Ones)
 657                 for i:=0; i < N; i++ {
 658                         a[i] = itor(int64(i+1))
 659                 }
 660                 checka(d, a, "Diff")  // 1 2 3 4 5
 661                 in := Integ(zero, Ones)
 662                 a[0] = zero  // integration constant
 663                 for i:=1; i < N; i++ {
 664                         a[i] = i2tor(1, int64(i))
 665                 }
 666                 checka(in, a, "Integ")  // 0 1 1/2 1/3 1/4 1/5
 667                 check(Cmul(neg(one), Twos), itor(-2), 10, "CMul")  // -1 -1 -1 -1 -1
 668                 check(Sub(Ones, Twos), itor(-1), 0, "Sub Ones Twos")  // -1 -1 -1 -1 -1
 669                 m := Mul(Ones, Ones)
 670                 for i:=0; i < N; i++ {
 671                         a[i] = itor(int64(i+1))
 672                 }
 673                 checka(m, a, "Mul")  // 1 2 3 4 5
 674                 e := Exp(Ones)
 675                 a[0] = itor(1)
 676                 a[1] = itor(1)
 677                 a[2] = i2tor(3,2)
 678                 a[3] = i2tor(13,6)
 679                 a[4] = i2tor(73,24)
 680                 a[5] = i2tor(167,40)
 681                 a[6] = i2tor(4051,720)
 682                 a[7] = i2tor(37633,5040)
 683                 a[8] = i2tor(43817,4480)
 684                 a[9] = i2tor(4596553,362880)
 685                 checka(e, a, "Exp")  // 1 1 3/2 13/6 73/24
 686                 at := Integ(zero, MonSubst(Ones, neg(one), 2))
 687                 for c, i := 1, 0; i < N; i++ {
 688                         if i%2 == 0 {
 689                                 a[i] = zero
 690                         } else {
 691                                 a[i] = i2tor(int64(c), int64(i))
 692                                 c *= -1
 693                         }
 694                 }
 695                 checka(at, a, "ATan")  // 0 -1 0 -1/3 0 -1/5
 696 /*
 697                 t := Revert(Integ(zero, MonSubst(Ones, neg(one), 2)))
 698                 a[0] = zero
 699                 a[1] = itor(1)
 700                 a[2] = zero
 701                 a[3] = i2tor(1,3)
 702                 a[4] = zero
 703                 a[5] = i2tor(2,15)
 704                 a[6] = zero
 705                 a[7] = i2tor(17,315)
 706                 a[8] = zero
 707                 a[9] = i2tor(62,2835)
 708                 checka(t, a, "Tan")  // 0 1 0 1/3 0 2/15
 709 */
 710         }
 711 }