Feigenbaum constant calculation: Difference between revisions

From Rosetta Code
Content added Content deleted
m (→‎{{header|REXX}}: added the display of the true value (to the number of specified digits).)
m (→‎{{header|REXX}}: corrected misspellings, changed precision defaults and number of iterations for I.)
Line 465: Line 465:
=={{header|REXX}}==
=={{header|REXX}}==
{{trans|Sidef}}
{{trans|Sidef}}
<lang rexx>/*REXX pgm calculates the Feigenbaurn bifurcation velocity using any # of decimal digits*/
<lang rexx>/*REXX pgm calculates the Feigenbaum bifurcation velocity using any # of decimal digits.*/
parse arg digs maxi maxj . /*obtain optional argument from the CL.*/
parse arg digs maxi maxj . /*obtain optional argument from the CL.*/
if digs=='' | digs=="," then digs=40 /* " " " " " " */
if digs=='' | digs=="," then digs=30 /* " " " " " " */
if maxi=='' | maxi=="," then maxi=15 /*Not specified? Then use the default.*/
if maxi=='' | maxi=="," then maxi=20 /*Not specified? Then use the default.*/
if maxj=='' | maxj=="," then maxj=10 /* " " " " " " */
if maxj=='' | maxj=="," then maxj=10 /* " " " " " " */
numeric digits digs /*use the specified # of decimal digits*/
numeric digits digs /*use the specified # of decimal digits*/
Line 484: Line 484:
x= 0; y= 0
x= 0; y= 0
do k=1 for 2**i
do k=1 for 2**i
y= 1 - 2 * y * x
y= 1 - 2 * x * y
x= a - x**2
x= a - x**2
end /*k*/
end /*k*/
Line 498: Line 498:
#= 4.669201609102990671853203820466201617258185577475768632745651343004134330211314737138,
#= 4.669201609102990671853203820466201617258185577475768632745651343004134330211314737138,
|| 68974402394801381716
|| 68974402394801381716
say ' true value= ' # / 1 /*true value of Feigenbaurn's constant.*/</lang>
say ' true value= ' # / 1 /*true value of Feigenbaum's constant. */</lang>
{{out|output|text=&nbsp; when using the default inputs:}}
{{out|output|text=&nbsp; when using the default inputs:}}
<pre>
<pre>
Using 10 iterations for J, with 40 decimal digits:
Using 10 iterations for J, with 30 decimal digits:


────i──── ───────────────d───────────────
────i──── ────────────────────d────────────────────
2 3.218511422038087912270504530742813256018
2 3.21851142203808791227050453077
3 4.385677598568339085744948568775522346173
3 4.3856775985683390857449485682
4 4.600949276538075357811694698623834984934
4 4.60094927653807535781169469969
5 4.655130495391980136486254995856898818963
5 4.65513049539198013648625498649
6 4.666111947828571388331213696711776471107
6 4.66611194782857138833121364654
7 4.668548581446840948044543680148146102083
7 4.66854858144684094804454708811
8 4.669060660648268239132599822630273970875
8 4.66906066064826823913257549468
9 4.669171555379511388886004609897560033836
9 4.6691715553795113888859465442
10 4.669195156030017174021108801191558304938
10 4.66919515603001717402161720542
11 4.669200229086856497938353781003810044639
11 4.66920022908685649793393149233
12 4.66920131329420417116475494118414885682
12 4.66920131329420417113719511412
13 4.669201545780906707506058109960038118631
13 4.66920154578090670783369507315
14 4.669201595537493910292470639266101619701
14 4.66920159553749390966169074155
15 4.669201606198152157723831098067070167449
15 4.66920160619815215840788706632
16 4.66920160848080435144581223484
17 4.66920160896974538458267849027
18 4.66920160907444981238909862845
19 4.66920160909687888294310165196
20 4.66920160910169069039564432665


true value= 4.669201609102990671853203820466201617258
true value= 4.66920160910299067185320382047
</pre>
</pre>


Revision as of 05:29, 19 September 2018

Feigenbaum constant calculation is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.


Calculate the Feigenbaum constant. See details: Feigenbaum constant



C

Translation of: Ring

<lang c>#include <stdio.h>

void feigenbaum() { 

   int i, j, k, max_it = 13, max_it_j = 10;
   double a, x, y, d, a1 = 1.0, a2 = 0.0, d1 = 3.2;
   printf(" i       d\n");
   for (i = 2; i <= max_it; ++i) {
       a = a1 + (a1 - a2) / d1;
       for (j = 1; j <= max_it_j; ++j) {
           x = 0.0;
           y = 0.0;
           for (k = 1; k <= 1 << i; ++k) {
                y = 1.0 - 2.0 * y * x;
                x = a - x * x;
           }
           a -= x / y;
       }
       d = (a1 - a2) / (a - a1);
       printf("%2d    %.8f\n", i, d);
       d1 = d;
       a2 = a1;
       a1 = a;
   }

}

int main() { 

   feigenbaum();
   return 0;

}</lang>

Output:
 i       d
 2    3.21851142
 3    4.38567760
 4    4.60094928
 5    4.65513050
 6    4.66611195
 7    4.66854858
 8    4.66906066
 9    4.66917155
10    4.66919515
11    4.66920026
12    4.66920098
13    4.66920537

C#

Translation of: Kotlin

<lang csharp>using System;

namespace FeigenbaumConstant { 

   class Program {
       static void Main(string[] args) {
           var maxIt = 13;
           var maxItJ = 10;
           var a1 = 1.0;
           var a2 = 0.0;
           var d1 = 3.2;
           Console.WriteLine(" i       d");
           for (int i = 2; i <= maxIt; i++) {
               var a = a1 + (a1 - a2) / d1;
               for (int j = 1; j <= maxItJ; j++) {
                   var x = 0.0;
                   var y = 0.0;
                   for (int k = 1; k <= 1<<i; k++) {
                       y = 1.0 - 2.0 * y * x;
                       x = a - x * x;
                   }
                   a -= x / y;
               }
               var d = (a1 - a2) / (a - a1);
               Console.WriteLine("{0,2:d}    {1:f8}", i, d);
               d1 = d;
               a2 = a1;
               a1 = a;
           }
       }
   }

}</lang>

Output:
 i       d
 2    3.21851142
 3    4.38567760
 4    4.60094928
 5    4.65513050
 6    4.66611195
 7    4.66854858
 8    4.66906066
 9    4.66917155
10    4.66919515
11    4.66920026
12    4.66920098
13    4.66920537

D

<lang d>import std.stdio;

void main() { 

   int max_it = 13;
   int max_it_j = 10;
   double a1 = 1.0;
   double a2 = 0.0;
   double d1 = 3.2;
   double a;

   writeln(" i       d");
   for (int i=2; i<=max_it; i++) {
       a = a1 + (a1 - a2) / d1;
       for (int j=1; j<=max_it_j; j++) {
           double x = 0.0;
           double y = 0.0;
           for (int k=1; k <= 1<<i; k++) {
               y = 1.0 - 2.0 * y * x;
               x = a - x * x;
           }
           a -= x / y;
       }
       double d = (a1 - a2) / (a - a1);
       writefln("%2d    %.8f", i, d);
       d1 = d;
       a2 = a1;
       a1 = a;
   }

}</lang>

Output:
 i       d
 2    3.21851142
 3    4.38567760
 4    4.60094928
 5    4.65513050
 6    4.66611195
 7    4.66854858
 8    4.66906066
 9    4.66917155
10    4.66919515
11    4.66920028
12    4.66920099
13    4.66920555

Go

Translation of: Ring

<lang go>package main

import "fmt"

func feigenbaum() { 

   maxIt, maxItJ := 13, 10
   a1, a2, d1 := 1.0, 0.0, 3.2
   fmt.Println(" i       d")
   for i := 2; i <= maxIt; i++ {
       a := a1 + (a1-a2)/d1
       for j := 1; j <= maxItJ; j++ {
           x, y := 0.0, 0.0
           for k := 1; k <= 1<<uint(i); k++ {
               y = 1.0 - 2.0*y*x
               x = a - x*x
           }
           a -= x / y
       }
       d := (a1 - a2) / (a - a1)
       fmt.Printf("%2d    %.8f\n", i, d)
       d1, a2, a1 = d, a1, a
   }

}

func main() { 

   feigenbaum()

}</lang>

Output:
 i       d
 2    3.21851142
 3    4.38567760
 4    4.60094928
 5    4.65513050
 6    4.66611195
 7    4.66854858
 8    4.66906066
 9    4.66917155
10    4.66919515
11    4.66920026
12    4.66920098
13    4.66920537

Haskell

<lang haskell>import Data.List (mapAccumL)

feigenbaumApprox :: Int -> [Double] feigenbaumApprox mx = snd $ mitch mx 10 

 where
   mitch :: Int -> Int -> ((Double, Double, Double), [Double])
   mitch mx mxj =
     mapAccumL
       (\(a1, a2, d1) i ->
           let a =
                 iterate
                   (\a ->
                       let (x, y) =
                             iterate
                               (\(x, y) -> (a - (x * x), 1.0 - ((2.0 * x) * y)))
                               (0.0, 0.0) !!
                             (2 ^ i)
                       in a - (x / y))
                   (a1 + (a1 - a2) / d1) !!
                 mxj
               d = (a1 - a2) / (a - a1)
           in ((a, a1, d), d))
       (1.0, 0.0, 3.2)
       [2 .. (1 + mx)]

-- TEST ------------------------------------------------------------------ main :: IO () main = 

 (putStrLn . unlines) $
 zipWith
   (\i s -> justifyRight 2 ' ' (show i) ++ '\t' : s)
   [1 ..]
   (show <$> feigenbaumApprox 13)
 where
   justifyRight n c s = drop (length s) (replicate n c ++ s)</lang>
Output:
 1    3.2185114220380866
 2    4.3856775985683365
 3    4.600949276538056
 4    4.6551304953919646
 5    4.666111947822846
 6    4.668548581451485
 7    4.66906066077106
 8    4.669171554514976
 9    4.669195154039278
10    4.669200256503637
11    4.669200975097843
12    4.669205372040318
13    4.669207514010413

Kotlin

Translation of: Ring

<lang scala>// Version 1.2.40

fun feigenbaum() { 

   val maxIt = 13
   val maxItJ = 10
   var a1 = 1.0
   var a2 = 0.0
   var d1 = 3.2
   println(" i       d")
   for (i in 2..maxIt) {
       var a = a1 + (a1 - a2) / d1
       for (j in 1..maxItJ) {
           var x = 0.0
           var y = 0.0
           for (k in 1..(1 shl i)) {
                y = 1.0 - 2.0 * y * x
                x = a - x * x
           }
           a -= x / y
       }
       val d = (a1 - a2) / (a - a1)
       println("%2d    %.8f".format(i,d))
       d1 = d
       a2 = a1
       a1 = a
   }

}

fun main(args: Array<String>) { 

   feigenbaum()

}</lang>

Output:
 i       d
 2    3.21851142
 3    4.38567760
 4    4.60094928
 5    4.65513050
 6    4.66611195
 7    4.66854858
 8    4.66906066
 9    4.66917155
10    4.66919515
11    4.66920026
12    4.66920098
13    4.66920537

Modula-2

<lang modula2>MODULE Feigenbaum; FROM FormatString IMPORT FormatString; FROM LongStr IMPORT RealToStr; FROM Terminal IMPORT WriteString,WriteLn,ReadChar;

VAR 

   buf : ARRAY[0..63] OF CHAR;
   i,j,k,max_it,max_it_j : INTEGER;
   a,x,y,d,a1,a2,d1 : LONGREAL;

BEGIN 

   max_it := 13;
   max_it_j := 10;

   a1 := 1.0;
   a2 := 0.0;
   d1 := 3.2;

   WriteString(" i       d");
   WriteLn;
   FOR i:=2 TO max_it DO
       a := a1 + (a1 - a2) / d1;
       FOR j:=1 TO max_it_j DO
           x := 0.0;
           y := 0.0;
           FOR k:=1 TO INT(1 SHL i) DO
               y := 1.0 - 2.0 * y * x;
               x := a - x * x
           END;
           a := a - x / y
       END;
       d := (a1 - a2) / (a - a1);
       FormatString("%2i    ", buf, i);
       WriteString(buf);
       RealToStr(d, buf);
       WriteString(buf);
       WriteLn;
       d1 := d;
       a2 := a1;
       a1 := a
   END;

   ReadChar

END Feigenbaum.</lang>

Perl 6

Works with: Rakudo version 2018.04.01
Translation of: Ring

<lang perl6>my $a1 = 1; my $a2 = 0; my $d = 3.2;

say ' i d';

for 2 .. 13 -> $exp { 

   my $a = $a1 + ($a1 - $a2) / $d;
   do {
       my $x = 0;
       my $y = 0;
       for ^2 ** $exp {
           $y = 1 - 2 * $y * $x;
           $x = $a - $x²;
       }
       $a -= $x / $y;
   } xx 10;
    $d = ($a1 - $a2) / ($a - $a1);
    ($a2, $a1) = ($a1, $a);
    printf "%2d %.8f\n", $exp, $d;

}</lang>

Output:
 i d
 2 3.21851142
 3 4.38567760
 4 4.60094928
 5 4.65513050
 6 4.66611195
 7 4.66854858
 8 4.66906066
 9 4.66917155
10 4.66919515
11 4.66920026
12 4.66920098
13 4.66920537

Phix

Translation of: Ring

<lang Phix>constant maxIt = 13, 

       maxItJ = 10

atom a1 = 1.0, 

    a2 = 0.0,
    d1 = 3.2

puts(1," i d\n") for i=2 to maxIt do 

    atom a = a1 + (a1 - a2) / d1
    for j=1 to maxItJ do
         atom x = 0, y = 0
         for k=1 to power(2,i) do
              y = 1 - 2*y*x
              x = a - x*x 
         end for
         a = a - x/y
    end for
    atom d = (a1-a2)/(a-a1)
    printf(1,"%2d %.8f\n",{i,d})
    d1 = d
    a2 = a1
    a1 = a

end for</lang>

Output:
 i d
 2 3.21851142
 3 4.38567760
 4 4.60094928
 5 4.65513050
 6 4.66611195
 7 4.66854858
 8 4.66906066
 9 4.66917155
10 4.66919515
11 4.66920026
12 4.66920098
13 4.66920537

Racket

Translation of: C

<lang racket>#lang racket (define (feigenbaum #:max-it (max-it 13) #:max-it-j (max-it-j 10)) 

 (displayln " i       d" (current-error-port))
 (define-values (_a _a1 d)
   (for/fold ((a 1) (a1 0) (d 3.2))
             ((i (in-range 2 (add1 max-it))))      
     (let* ((a′ (for/fold ((a (+ a (/ (- a a1) d))))
                          ((j (in-range max-it-j)))
                  (let-values (([x y] (for/fold ((x 0) (y 0))
                                                ((k (expt 2 i)))
                                        (values (- a (* x x))
                                                (- 1 (* 2 y x))))))
                    (- a (/ x y)))))
            (d′ (/ (- a a1) (- a′ a))))
       (eprintf "~a   ~a\n" (~a i #:width 2) (real->decimal-string d′ 8))
       (values a′ a d′))))
 d)

(module+ main 

 (feigenbaum))</lang>
Output:
 i       d
2    3.21851142
3    4.38567760
4    4.60094928
5    4.65513050
6    4.66611195
7    4.66854858
8    4.66906066
9    4.66917155
10   4.66919515
11   4.66920026
12   4.66920098
13   4.66920537
4.669205372040318

REXX

Translation of: Sidef

<lang rexx>/*REXX pgm calculates the Feigenbaum bifurcation velocity using any # of decimal digits.*/ parse arg digs maxi maxj . /*obtain optional argument from the CL.*/ if digs== | digs=="," then digs=30 /* " " " " " " */ if maxi== | maxi=="," then maxi=20 /*Not specified? Then use the default.*/ if maxj== | maxj=="," then maxj=10 /* " " " " " " */ numeric digits digs /*use the specified # of decimal digits*/ 

   a1=  1
   a2=  0
   d1=  3.2
  pad= left(, 3)                              /*literal for spacing between values.  */

say 'Using ' maxj " iterations for " j', with ' digits() " decimal digits:" say say center('i',9,"─") pad center('d',digs+1,"─") /*show a centered title for I and D.*/ 

  do i=2  for maxi-1
  a= a1  +  (a1 - a2) / d1
                              do j=1  for maxj
                              x= 0;   y= 0
                                                   do k=1  for 2**i
                                                   y= 1  -  2 * x * y
                                                   x= a  -  x**2
                                                   end   /*k*/
                              a= a  -  x / y
                              end   /*j*/
  d= (a1 - a2)  /  (a - a1)
  say center(i, 9)  pad  d                      /*display values for  I & D ──►terminal*/
  d1= d
  a2= a1
  a1= a
  end   /*i*/                                   /*stick a fork in it,  we're all done. */

say

  1. = 4.669201609102990671853203820466201617258185577475768632745651343004134330211314737138,
  || 68974402394801381716

say ' true value= ' # / 1 /*true value of Feigenbaum's constant. */</lang>

output   when using the default inputs:
Using  10  iterations for  J,  with  30  decimal digits:

────i────     ───────────────d───────────────
    2         3.21851142203808791227050453077
    3         4.3856775985683390857449485682
    4         4.60094927653807535781169469969
    5         4.65513049539198013648625498649
    6         4.66611194782857138833121364654
    7         4.66854858144684094804454708811
    8         4.66906066064826823913257549468
    9         4.6691715553795113888859465442
   10         4.66919515603001717402161720542
   11         4.66920022908685649793393149233
   12         4.66920131329420417113719511412
   13         4.66920154578090670783369507315
   14         4.66920159553749390966169074155
   15         4.66920160619815215840788706632
   16         4.66920160848080435144581223484
   17         4.66920160896974538458267849027
   18         4.66920160907444981238909862845
   19         4.66920160909687888294310165196
   20         4.66920160910169069039564432665

 true value=  4.66920160910299067185320382047

Ring

<lang ring>

  1. Project : Feigenbaum constant calculation

decimals(8) see "Feigenbaum constant calculation:" + nl maxIt = 13 maxItJ = 10 a1 = 1.0 a2 = 0.0 d1 = 3.2 see "i " + "d" + nl for i = 2 to maxIt 

    a = a1 + (a1 - a2) / d1
    for j = 1 to maxItJ
         x = 0
         y = 0
         for k = 1 to pow(2,i)
              y = 1 - 2 * y * x
              x = a - x * x 
         next   
         a = a - x / y
    next
    d = (a1 - a2) / (a - a1)
    if i < 10
       see "" + i + "    " + d + nl
    else
       see "" + i + "  " + d + nl
    ok
    d1 = d
    a2 = a1
    a1 = a

next </lang> Output: 

Feigenbaum constant calculation:
i  d
2  3.21851142
3  4.38567760
4  4.60094928
5  4.65513050
6  4.66611195
7  4.66854858
8  4.66906066
9  4.66917155
10 4.66919515
11 4.66920026
12 4.66920098
13 4.66920537

Sidef

Translation of: Perl 6

<lang ruby>var a1 = 1 var a2 = 0 var δ = 3.2.float

say " i\tδ"

for i in (2..15) { 

   var a0 = ((a1 - a2)/δ + a1)
   10.times {
       var (x, y) = (0, 0)
       2**i -> times {
           y = (1 - 2*x*y)
           x = (a0 - x²)
       }
       a0 -= x/y
   }
   δ = ((a1 - a2) / (a0 - a1))
   (a2, a1) = (a1, a0)
   printf("%2d %.8f\n", i, δ)

}</lang>

Output:
 i	δ
 2 3.21851142
 3 4.38567760
 4 4.60094928
 5 4.65513050
 6 4.66611195
 7 4.66854858
 8 4.66906066
 9 4.66917156
10 4.66919516
11 4.66920023
12 4.66920131
13 4.66920155
14 4.66920160
15 4.66920161

zkl

Translation of: Kotlin

<lang zkl>fcn feigenbaum{ 

  maxIt,maxItJ,a1,a2,d1,a,d := 13, 10, 1.0, 0.0, 3.2, 0, 0;
  println(" i       d");
  foreach i in ([2..maxIt]){
     a=a1 + (a1 - a2)/d1;
     foreach j in ([1..maxItJ]){
        x,y := 0.0, 0.0;

foreach k in ([1..(1).shiftLeft(i)]){ y,x = 1.0 - 2.0*y*x, a - x*x; } a-=x/y 

     }
     d=(a1 - a2)/(a - a1);
     println("%2d    %.8f".fmt(i,d));
     d1,a2,a1 = d,a1,a;
  }

}();</lang>

Output:
 i       d
 2    3.21851142
 3    4.38567760
 4    4.60094928
 5    4.65513050
 6    4.66611195
 7    4.66854858
 8    4.66906066
 9    4.66917155
10    4.66919515
11    4.66920026
12    4.66920098
13    4.66920537
Retrieved from "https://rosettacode.org/wiki/Feigenbaum_constant_calculation?oldid=280862"