# The challenge

The plastic number is a number related to the golden ratio, with many interesting mathematical properties. As such, there are many approaches that can be used to calculate the number.

In order to precisely specify the number for the purposes of this challenge, we'll use the following definition (although there are plenty of equivalent definitions, and you can use any definition you wish as long as it comes to the same number):

The plastic number is a real number ρ such that ρ³=ρ+1.

Your challenge is to write a program or function which takes an integer x as input (with x > 1), and produces an approximation to ρ as output, such that the larger the value of x gets, the closer the output gets to ρ (with at most finitely many exceptions; staying at the same value counts as "closer" for this purpose), and for any positive number δ, there's some input x to your program that produces an output that's within δ of ρ.

# Clarifications

• If you're outputting via a method that inherently outputs strings (e.g. the standard output stream), you can format output either in decimal (e.g. 1.3247179572), or as a ratio of two integers with a / character between them.
• If you're outputting as a value within your programming language (e.g. returning from a function), it must be of a fixed-point, floating-point, or rational type. (In particular, you can't use data types that store numbers symbolically, unless they're used only to hold the ratio of two integers. So if you're using Mathematica or a similar language, you'll need to include the extra code to actually generate the digits of the output.)
• Your answer must work in a hypothetical variant of your language in which integers can be arbitrarily large, and memory (including stack) is unlimited. You may not assume that floating-point arithmetic in your language is arbitrarily accurate, but must instead use its actual accuracy (meaning that outputting a floating-point number is only going to be possible in languages where the accuracy of floating-point numbers can be controlled at runtime).
• x can have any meaning you want (so long as increasing it gives more accurate outputs). I imagine that most submissions will have it control the number of digits of output to produce, or the number of iterations of the algorithm used by your program to converge on the plastic number, but other meanings are acceptable.

# Testcase

Here are the first few digits of the plastic number:

1.32471795724474602596090885


More digits are available on OEIS.

# Victory condition

As usual for , shorter is better, measured in bytes. However, feel free to post answers even if they don't win, so long as they add something (e.g. a different language, or a different algorithm) to the existing answers.

• hmm, (cbrt(108+12*sqrt(69))+cbrt(108-12*sqrt(69)))/6 this seems like a good time to use the Drake approximation': sqrt(69)=8.something bit.ly/2rCqedX ^_^ – DrQuarius Jun 16 '17 at 3:34
• Can we also assume the recursion/stack depth is unlimited? – xnor Jun 16 '17 at 3:40
• To clarify the second point, can we use arbitrary precision libraries (e.g. mpmath in Python)? They use an auxiliary data type, but do you count that as storing things "symbolically"? – Batman Jun 16 '17 at 4:08
• Well, at the very least I'd expect answers to converge to ρ. Also, an "honest" solution could easily fail the test x > y --> |ρx - ρ| > |ρy - ρ| for a finite number of (x, y) pairs. If that isn't acceptable, I think this should be made more explicit in the spec. – Dennis Jun 16 '17 at 7:58
• Many answerers have fallen into the trap(?) of computing an x digit approximation to ρ, the problem being that there are probably infinitely many x such that an (x + 1)-digit approximation is no better than an x digit approximation. You should probably clarify whether you intended this to be allowed. If you don’t, replace “closer” with “strictly closer”; if you do, “at least as close”, or something. You could also consider the looser requirement that the sequence converges to ρ, which would additionally allow xnor’s answer. – Anders Kaseorg Jun 16 '17 at 9:25

# Python 2, 49 bytes

n=x=input()
while n**3/x/x<n+x:n+=1
print n,'/',x


Try it online!

The idea is to express the ρ with ρ³=ρ+1 as a fraction n/x whose denominator x is the input accuracy parameter. We take (n/x)³=n/x+1 and clear denominators to get n³=x²(x+n).

Since the LHS increases in n faster than the RHS, we can approximate the equality point n as the smallest with n³≥x²(x+n). The code counts up n until this is the case, starting at x which is smaller.

A small byte save is to divide both sides by x² to write n³/x²≥x+n (negated in the while condition). This is floor division in the code, but the fractional part lost is negligible.

A same-length alternative instead puts x as the numerator:

Python 2, 49 bytes

n=x=input()
while x**3/n/n<n+x:n-=1
print x,'/',n


Try it online!

• Although this output converges to ρ (∀ε > 0 ∃x₀ ∀x ≥ x₀ |f(x) − ρ| < ε), it doesn’t satisfy “the larger the value of x gets, the closer the output gets to ρ (with at most finitely many exceptions)” (∃x₀ ∀x ≥ x₀ |f(x + 1) − ρ| < |f(x) − ρ|). – Anders Kaseorg Jun 16 '17 at 9:00
• This problem can be fixed by using 2**input() rather than just input(); then, each approximation will be as least as accurate as the one before. – user62131 Jun 16 '17 at 17:46

# Mathematica, 20 bytes

#^3-#-1&~Root~1~N~#&


Mathematica's builtin Root function gives the solutions to a polynomial equation f[x] == 0.

### Explanation

#^3-#-1&~Root~1~N~#&
&  (* Function *)
#^3-#-1&              (* A pure-function polynomial, x^3-x-1 *)
~Root~1       (* Find the first root *)
~N~#   (* approximate to (input) digits *)


### Sample I/O

In:= f=#^3-#-1&~Root~1~N~#&;
f

Out= 1.

In:= f

Out= 1.32471796

In:= f

Out= 1.324717957244746025960908854478097340734404056901733364534015050302827851245547594054699347981787280

• PS: Root[x^3-x-1,1]~N~#& works fine (despite not telling it that x is a variable) for the same byte count. – Greg Martin Jun 16 '17 at 7:43
• @AndersKaseorg: I changed that rule because it was clearly broken. No valid answers were invalidated, but some answers (like this one) became valid. – user62131 Jun 16 '17 at 17:43

# Mathematica, 27 bytes

x/.Solve[x^3==x+1>2,x]~N~#&


-1 byte from Martin
-2 bytes from ovs

input



output

{1.32471795724474602596090885}

• Solve[x^3==x+1>2,x]~N~#& for 24 bytes – ovs Jun 16 '17 at 6:02
• The result of this is {{x -> 1.32...}} though. You might want to check with ais whether that's a valid output format. – Martin Ender Jun 16 '17 at 6:10
• ok.. all fixed i guess – J42161217 Jun 16 '17 at 6:21
• It's still {1.32...} actually, but that format is probably less contentious. – Martin Ender Jun 16 '17 at 6:23
• I made the challenge a bit more general so that this would be valid, it wasn't meant to disallow "first x digits" solutions. So this is valid now, even though it wasn't before. – user62131 Jun 16 '17 at 17:45

# sed, 67 60 (59+1) bytes

s,^,1/1/1 ,
:;s,(1*/(1*)/(1*).*)1$,\2\3/\1, t s,(/1*).*,\1,  Try it online! +1 for the -E flag (ERE instead of BRE). Input and output are both unary: input 11111 for x=5 e.g. Output is a fraction of two unary numbers: the aforementioned 11111 input yields output 11111/1111 (5/4 in decimal). Approximates the plastic number as a fraction between to consecutive elements of the Padovan sequence. • FWIW you don't need a space after the b command, but you can make it shorter still by using the empty label (: and b with no argument). tio.run/#%23K05N@f@/… – Jordan Jun 17 '17 at 3:13 • Oh excellent. And I can save yet another 4 bytes by using t instead of b, so that's a pretty nice save. Thanks :) – FireFly Jun 17 '17 at 10:38 ## Mathematica, 27 bytes Nest[(1+#)^(1/3)&,1,#]~N~#&  Uses a truncated approximation of the nested cubic radical form ³√(1+³√(1+³√(1+...))). While the output will always have x-1 decimal places, the result is actually less accurate than that, because the expression converges more slowly than one digit per iteration (x is also used as the number of nested radicals that are computed). For example x = 100 gives _________________________________________________________________________ 1.324717957244746025960908854478097340734404056901733364534015050302827850993693624204577670741656151  where the overlined part is correct. • I was planning to write this algorithm in dc, but got stymied because it turns out that it doesn't have a cube root operation, and raising a number to the power ⅓ doesn't work either :-( At least you can always count on Mathematica to have appropriate builtins… – user62131 Jun 16 '17 at 6:10 • @ais523 There's actually CubeRoot but ain't nobody got bytes for that. – Martin Ender Jun 16 '17 at 6:11 # Octave, 50 bytes @(n)char(digits(n)*0+vpasolve(sym('r^3-r-1'))(1));  Try it online! Defines an anonymous function, withn the desired number of digits of output. This answer abuses that digits returns the current setting for the number of digits in variable precision arithmetic. This means we can just use it in an anonymous function without errors about 'Too many output arguments'. Other than that, it's really straightforward: vpasolve is short for Variable-Precision Arithmetic Solve, with the precision set by the last call of digits. Since vpa is a Symbolic data type in Octave, which is banned per the spec, we just wrap the whole function in char(...) to get string output. Note that in solve and vpasolve, the f==0 is implied, so r^3==r+1 has been replaced by r^3-r-1 (==0) • I went and changed the question so that it doesn't disallow answers like this (it wasn't meant to). – user62131 Jun 16 '17 at 17:44 • @ais523 Thanks for the notification! – Sanchises Jun 16 '17 at 20:47 # MATL (27 28 bytes) 7BG:"t@)y@Q)+h]tG3+)yG2+)/  My first solution (27 bytes) Try it online! It's certainly not optimal, I'm still getting used to MATL. ## Explanation: I create a Padovan sequence up to input+3 then find the ratio of the last two numbers. 7B % Turn 7 into binary to give 1 1 1 G:" % For k=1:input do... t@) % Existing sequence member k y@1+) % Existing sequence member k+1 +h % Add them together and concatenate to the sequence array ] % End loop tG3+) % Final sequence member yG2+) % Second last sequence member / % Divide to approximate ρ  ## Proper fraction output (35 bytes) (28 bytes, @Sanchises): However, the first solution doesn't fulfill the need for arbitrary precision being the floating point limit of default MATL settings. So rather than adding several bytes to extend this precision, it's simpler to take the proper fraction route and write a fraction of the final two integers in the (N-1)th and Nth elements of the truncated Padovan sequence. e.g "114/86" 7BG:"t@)y@1+)+h]tG3+)V'/'YcyG2+)VYc 7BG:"t@tQh)sh]tJ)V47hyJq)Vh&  Courtesy of user @Sanchises. :) Try it online! ## Non-iterative evaluation: Notably, my shortest code for the 'exact' version is (23 bytes): 1-1h69X^12**108+1I/^6/s  Try it online! ...but doesn't give arbitrary precision. I wonder if anyone can adjust this to fulfill the rules (use the input etc) and still add less than 5 bytes? :P • 1+ can be shortened to Q.With that in mind, you can replace @)y@1+)+ with just @tQh)s. Furthermore, you can use J to indicate the end of an array; and finally, MATL doesn't distinguish between normal arrays and character arrays, so you can replace Yc by h (you don't need the extra functionality of Yc). This gives just 28 bytes: 7BG:"t@tQh)sh]tJ)V47hyJq)Vh& (notice the & to prevent superfluous output, and replacing '/' by 47). – Sanchises Jun 16 '17 at 8:30 • Kudos for the 7B though, much better than naively pushing lllv – Sanchises Jun 16 '17 at 8:33 • @DrQuarius The latest version can always be found in this GitHub link – Luis Mendo Jun 17 '17 at 14:17 • @DrQuarius No, this behaviour is present in the rather old MATL spec I generally use. You should really check out Table 3. Not only does clipboard J by default contain 1j, but the clipboard L also contains many useful indexing functions (note that 1j equals end in MATL). – Sanchises Jun 17 '17 at 15:07 • Also, don't worry, I'm an mechanical engineer. I think MATL(AB) has little use outside of a scientific environment, so I would guess that the majority of MATL(AB)/Octave golfers are from outside CS. – Sanchises Jun 17 '17 at 15:17 # M, 15 14 bytes ²×3’ *3Ḥ‘÷Ç Ç¡  Try it online! ### Algorithm This uses rationals and Newton's method. Specifically, for input x, the first x iterations with starting value x are applied. We're trying to find a specific root of the polynomial p(t) = t³ - t - 1. Newton's method achieves this by taking a starting value t0 – sufficiently close to ρ – and recursively defining a sequence by tn+1 = tn - p(tn) / p'(tn). Since p'(t) = 3t² -1, we get tn+1 = tn - (tn³ - tn - 1)/(3tn² - 1) = (3tn³ - tn - tn³ + tn + 1) / (3tn² - 1) = (2tn³ + 1) / (3tn² - 1). Note that the initial approximation x gets progressively worse as x increases. While the output for x = 3 is slightly less precise than the output for x = 2, since Newton's method converges quadratically to ρ, this shouldn't be an issue for large values of x. ### How it works Ç¡ Main link. Argument: x Ç¡ Call the second helper link x times, which initial argument x. *3Ḥ‘÷Ç Second helper link. Argument: t *3 Compute t³. Ḥ Unhalve; yield 2t³. ‘ Increment; yield 2t³+1. Ç Call the first helper link with argument t. ÷ Divide the left result by the right one. ²×3’ First helper link. Argument: t ² Compute t². ×3 Compute 3t². ’ Decrement; yield 3t²-1.  • Too bad you can't use ...µ¡... – Erik the Outgolfer Jun 26 '17 at 7:24 # Julia 0.5, 44 40 bytes |(x,r=1)=x<1?r:~-x|big(2r^3+1)//(3r^2-1)  Uses rationals and Newton's method. Try it online! # 05AB1E, 23 bytes U[X3m¹/¹/¹X+›#X>U]X'/¹J  Try it online! Direct port of https://codegolf.stackexchange.com/a/126822/59376 by xnor. # Charcoal, 28 bytes ＡＩθθＡθνＷ‹∕∕Ｘν³θθ⁺νθＡ⁺ν¹νＩ∕νθ  Try it online! Link to verbose mode. Also I apparently messed up Divide and IntDivide :| Uses the same method as the Python and JavaScript answers. # NewStack, 14 bytes ¹Fᵢ{E2x³⁺÷3x²⁻  ### Break down: ¹ Add arbitrary number 1 to the stack. Fᵢ{ Define for loop with a user's input amount of itterations. E Define new edit for element 0 (element 0 being the 1 added. earlier). 2x³⁺÷3x²⁻ update x to equal (2x^3+1)/(3x^2-1). (x = element 0).  ### How it works: The formula (2x3+1)/(3x2-1) comes from the simplification of Newton's method for the equasion x3=x+1. You can find it here. Repeating this process an infinite amoune of times converges to the plastic number. It's rate of convergence is rather quick at around 2.6 decimals per iteration. INPUT ITERATION >> VALUE 0 >> 1 1 >> 1.5 2 >> 1.3478260869565217 3 >> 1.325200398950907 4 >> 1.3247181739990537 5 >> 1.3247179572447898 6 >> 1.324717957244746 <- 16 decimal precision in 6 iterations! ... 100 >> 1.324717957244746  # Padovan sequence alternative, 2725 17 bytes ¹Fᵢ{[ƨ2+ƨ3]ℲƤƨ/ƨ2  ### Break down: ¹ Append first element of Padovan sequence. Fᵢ{ Ⅎ Define for loop of user's input amount of iterations. [ƨ2+ƨ3] Append sum second and third to last elements. Ƥƨ/ƨ2 Print ratio of last two elements.  -2 bytes by choosing better print strategy -8 bytes by choosing better way to index stack ### How it works: As the Padovan sequence continues, the ratio of the last two elements converge to the plastic number. INPUT ITERATION >> VALUE 0 >> 1 1 >> 2 ... 10 >> 1.3157894736842106 ... 89 >> 1.324717957244746 <- 16 decimal precision in 89 iterations ... 100> > 1.324717957244746  # Ruby, 35 bytes ->x{n=1;n+=0.1**x while n**3-n<1;n}  Try it online! ## Clojure, 46 bytes #(nth(iterate(fn[i](Math/pow(inc i)(/ 3)))1)%)  Uses the iterated cube-root formula. This is a bit more interesting but longer: (def f #(apply comp(repeat %(fn[i](Math/pow(inc i)(/ 3)))))) ((f 10)1) 1.3247179361449652  • “You may not assume that floating-point arithmetic in your language is arbitrarily accurate, but must instead use its actual accuracy (meaning that outputting a floating-point number is only going to be possible in languages where the accuracy of floating-point numbers can be controlled at runtime).” – Anders Kaseorg Jun 16 '17 at 9:30 • Ooh I did not notice that, what a bummer. And implementing cubic root with BigDecimal seems quite tricky. – NikoNyrh Jun 16 '17 at 10:48 ## Javascript, 36 bytes f=(x,n=x)=>n**3/x/x<n+x?f(x,++n):n/x  Works the same as the top python answer. No console.log was included because if you run f(x) in console it will be logged automatically. f=(x,n=x)=>n**3/x/x<n+x?f(x,++n):n/x console.log(f(300)) # ><>, 38+3 = 41 bytes 11\n; ?!\}2,:01{::::}**-+0(?$-{+{1-:}


Expects the input to be present on the stack at program start, so +3 bytes for the -v flag.

Try it online!

Effectively performs a binary search to narrow in on the output value. Increasing x increases the number of iterations to perform.

Edit: refactored calculation slightly to save 1 byte, previous version:

11\n;
?!\}2,:0{::::}**$-1-0)?$-{+{1-:}


# k, 27 bytes

{"/"/$1_|x({y,z,x+y}.)/3#1}  Try it online! This assumes infinite integers (which, alas, is not true). It uses the Padovan sequence. # TI-BASIC, 21 bytes :Prompt X //Prompt for input, 3 bytes :While X //While X, 3 bytes :³√(1+Y→Y //Calculate cube root of 1+Y and store to Y, 7 bytes :DS<(X,0 //Decrement X and skip next command (which doesn't do anything), 5 bytes :End //End While loop, 2 bytes :Y //Display Y, 1 byte  Uses this recursive formula. Interestingly, hard-coding the number and rounding it gives the same byte-count: # TI-BASIC, 21 bytes :Prompt X //Prompt for input, 3 bytes :.5√(3 //Store √(3)/2 to Ans, 5 bytes :Ansֿ¹cosh(3ֿ¹coshֿ¹(3Ans //Store the plastic number to Ans, 9 bytes :round(Ans,X //Round the plastic number X decimal digits, 4 bytes  Uses this trigonometric formula. • I don't think you can use TI-BASIC's floats here: Your answer must work in a hypothetical variant of your language in which integers can be arbitrarily large, and memory (including stack) is unlimited. You may not assume that floating-point arithmetic in your language is arbitrarily accurate, but must instead use its actual accuracy (meaning that outputting a floating-point number is only going to be possible in languages where the accuracy of floating-point numbers can be controlled at runtime). – lirtosiast Jul 11 '17 at 1:28 ## C#, 317 bytes using m=System.Math;a=x=>{if(a==0)return "1/1";var d=a(x-1).Split('/');var b=int.Parse(d);var c=int.Parse(d);return string.Format("{0}/{1}",(2*m.Pow(b,3)+m.Pow(c,3)).ToString(new string('#',int.MaxValue.ToString().Length)),(3*m.Pow(b,2)*c-m.Pow(c,3)).ToString(new string('#',int.MaxValue.ToString().Length)));};  It returns the result as a fraction. ### Explanation It uses the Newton's method with x iterations for finding the root of the polynomial p^3-p-1=0. The formula is x_n=1-(f(x_(n-1)))/(f'(x_(n-1))), and x_0 is a starting point. The polynomials derivative is 3p^2-1, and let's say x_(n-1)=b/c. Then, by using the above formula we get, that x_n=(2 b^3+c^3)/(3 b^2 c-c^3). Let's also say, that we start from 1, this will happen, when x=2, because x>1, and is an integer. Idented, and commented code: using System; string PlasticNumber(int x) { if (x == 2) return "1/1"; //If x=2, we return our starting value, but we need to return it as a fraction var d = PlasticNumber(x - 1).Split('/'); var b = System.Convert.ToInt32(d); var c = int.Parse(d); //We parse the previous value of the fraction, and put it into two variables return string.Format("{0}/{1}", (2 * Math.Pow(b, 3) + Math.Pow(c, 3)) .ToString(new string('#', int.MaxValue.ToString().Length)), (3 * Math.Pow(b, 2) * c - Math.Pow(c, 3)) .ToString(new string('#', int.MaxValue.ToString().Length))); //We return the result as a fraction, but it's important not to return it in scientific notation, because that will cause issues in the parsing process }  # PHP, 86 bytes for($p=[1,1,1];$i++<$argn;)$p[]=bcadd($p[$i],$p[$i-1]);echo bcdiv($p[$i+1],$p[$i],$i);


PHP Sandbox Online

Creates the Padovan Spiral and print the ratio of the last two numbers.

# Axiom, 96 bytes

h(n:NNI):Float==(n>1.E5=>-1;n:=n+1;j:=digits(n::PI);r:=solve(x^3-x=1,10.^-n);digits(j);rhs(r.1))


results

(31) -> [h(i) for i in 0..10]
(31)
[1.0, 1.3, 1.33, 1.325, 1.3247, 1.32472, 1.324718, 1.324718, 1.32471796,
1.324717957, 1.3247179572]
Type: List Float


how you can see h(2) should be 1.32 and not 1.33 so there is some error in last digits

Then there would be this one of 110 bytes

g(n:NNI):Float==(n>1.E5=>-1;n:=n+1;j:=digits(n::PI);x:=sqrt(23./108);r:=(.5+x)^(1/3)+(.5-x)^(1/3);digits(j);r)


It use the formula for resolve equation of III grade of type x^3-3*p*x-2*q=0 in the case q^2-p^3>=0 that is m=sqrt(q^2-p^3) and x=(q+m)^(1/3)+(q-m)^(1/3)

In our case r^3-r-1=0 this can be written as r^3-3*(1/3)r-2*(1/2)=0 so p=1/3 q=1/2 m=1/4-1/27=23/108 x=(0.5+m)^(1/3)+(0.5-m)^(1/3)

this one that use Newton iteration with start point r=1

f(n:NNI):Float==(n>1.E5=>-1;n:=n+1;j:=digits(n::PI);e:=10^-n;r:=1.;repeat(v:=(r^3-r-1)/(3*r^2-1);abs(v)<e=>break;r:=r-v);digits(j);r)


it change in the function, digits value for obtain one obj of n+1 digits afther the float point. At end the digits() value is assigned back to preciding value.

# APL (Dyalog Extended), 19 bytes

{{(⍺+⍵*÷3)*÷3}/⍵/1}


Try it online!

Uses the cube root approximation method, which very, very slowly creeps towards the value of the plastic number.

## Explanation

{{(⍺+⍵*÷3)*÷3}/⍵/1} ⍵ → input
⍵/1} replicate 1 ⍵ times
{           }/     reduce by inner fn:
*÷3       cube root of
⍵*÷3           cube root of right arg
⍺+               plus left arg
`