Find The Real Solutions of a Cubic

Question

Description

All cubic equations can be solved, and every cubic has at least one solution. The goal of this challenge is to find the real solutions to a given cubic using inputs, and (obviously) the smallest program size.

Rules

The solution must find the real solutions of a cubic correctly, using any method. The cubic must be in the standard form $$ax^3 + bx^2 + cx + d = 0$$

The first coefficient cannot be zero, as then it would be a quadratic. There must a different output if it is zero, and outputting nothing at all is acceptable.

Scoring

The score is based on the size of the program alone.

Inputs

The four coefficients of the cubic, a, b, c, and d.

The input may be formatted however you like.

Outputs

Every real solution to the cubic. There will either be three, or one for every cubic.

Example

With 3 solutions

Input:

1
3
0
-1

Output:

-2.879
-0.653
0.532

With 1 solution

Input:

1
1
1
1

Output:

-1

Additional Info

While you can use any method, you can find an image with the cubic equation (similar to the quadratic equation, can solve any cubic every time) equation can be found here.

You can learn more about cubics here.

Answer 1

R, 42 41 bytes

function(c,a=polyroot(c))a[!Im(a)]/!!c[4]

Try it online!

Input is vector of d,c,b,a. Outputs Inf if the cubic coefficient (a) is zero.

Using a[!Im(a)] to select only real solutions is very susceptible to floating-point rounding errors; the TIO header rounds values less than 0.0000000001 to zero to prevent this. Including rounding in the code costs 5 3 more bytes.

---

R, 111 bytes

function(e,z=1i^(1:3)){for(i in 1:99)z=z-e%*%rbind(z^3,z^2,z,1)/(z-(x=c(z,z))[2:4])/(z-x[3:5]);if(e)z[!Im(z)]}	

Try it online!

Non-builtin Durant-Kerner strategy copied from 97.100.97.109's answer (upvote that!). Using this approach to find all the roots in parallel works very nicely with R's vectorized operations.

Answer 2

Python, 143 134 133 bytes

f=lambda a,b,c,d,R=[-1,1j,1]:[R:=[(z:=R[2])-(((a*z+b)*z+c)*z+d)/(z-R[0])/(z-R[1])]+R for i in" "*99]and[x for x in R[:3]if x.imag==0]

Attempt This Online!

This is the first solution that doesn't use a builtin, and it suffers greatly for it. Instead, I'm using the Durant-Kerner method to find all approximations, then filtering for whether they're complex.

---

-9 bytes from @Steffan for various small changes.

Answer 3

MATL, 12 bytes

1)?G&ZQ&Zj~)

Try it online!

How it works

       % Implicit input
1)     % Get the first element
?      % If non-zero
  G    %   Push input again
  &ZQ  %   Implicit input. Roots of polynomial. Gives a numeric vector
  &Zj  %   Push real and imaginary parts
  ~    %   Negate. Gives "true" if imaginary part is zero, or "false" otherwise
  )    %   Keep only real parts corresponding to "true".
       % Implicit end
       % Implicit display

Answer 4

Vyxal P, 4 bytes

hß∆P

Try It Online!

Good old built-ins. This may not work on some interpreters until the SymPy ACE vulnerabilities are fixed. Input as a list [a,b,c,d].

Edit: Changed output to the input if leading coefficient is 0.

Answer 5

Both functions below take a polynomial in x.

PARI/GP, 27 bytes

p->if(#p>3,polrootsreal(p))

Attempt This Online!

Returns 0 for quadratics and below.

PARI/GP, 101 bytes

p->k=[Pi,I,2];for(q=1,99,k-=[subst(p,x,k[i])/vecprod([k[i]-z|z<-k[^i]])|i<-[1..3]]);[r|r<-k,!imag(r)]

Attempt This Online!

Yet another implementation of the Durand–Kerner method.

Answer 6

Jelly, 8? 11 bytes

ṛ/ȧÆr,N+AƲƇ

A monadic Link that accepts a list of the coefficients ([d,c,b,a]) and yields a list of the (one or three) real roots, unless a is zero in which case it yields an empty list.

Try it online! Or see the test-suite.

...\$8\$ bytes if we can handle a being zero (in which case it works for all polynomials):

Ær,NċAƲƇ

Try it online!

How?

The heavy lifting is done by a built-in...

ṛ/ȧÆr,N+AƲƇ - Link: coefficients, P = [d,b,c,a]
 /          - reduce (P) by:
ṛ           -   right
               -> a
   Ær       - polynomial roots (of P)
  ȧ         - logical AND
               -> [d,c,b,a] or 0
          Ƈ - keep those for which (with 0, keep those of range(0)=[]):
         Ʋ  -   last four links as a monad - f(x):
      N     -     negate (x)
     ,      -     (P) paired with (that)
        A   -     absolute value (of x)
       ċ    -     count occurrences (of that in the pair)
                    -> Truthy (1 or 2*) if x is real, falsey (0) if not
                                  * when x=0

Answer 7

Python 3, 99 bytes

import numpy
f=lambda a,b,c,d:[float(i)for i in numpy.roots([a,b,c,d])if round(float(i))!=0if a!=0]

Try it online!

Answer 8

J, 32 bytes

[:{:^:([:-.*./@(=+)){:@:>@p.^:{:

I am pretty bad at conditionals in J, but after some finagling, this was the best I could come up with. Expects coefficients as d,c,b,a

Attempt This Online!

[:{:^:([:-.*./@(=+)){:@:>@p.^:{:
                    {:@:>@p.^:{:  : u^:v executes u if v is true, else return x
                              {:  : is last item 0?
                    {:@:>@p.      : if non-zero
                          p.      : converts coefficient form to multiplier-roots form
                        >         : unbox the results
                    {:            : take the roots, not the multiplier
  {:^:([:-.*./@(=+))              : The check for only one real
      ([:-.*./@(=+))              : u(v(x)) where u is -. and v is *./@(=+)
               (=+)               : each x equal to its complex conjugate?
           *./                    : AND reduce result
         -.                       : negate result
  {:                              : if roots are not all reals, take the last value

extra:
u@v -> apply u monadically to v(x) for every application of v
u@:v -> apply u to the entire result of v(x)
[:u v -> strictly u(v(x))

Answer 9

J, 20 bytes

1(#~]=+)@{::[:p.]*{:

Attempt This Online!

Takes d c b a, and returns the array of real roots, or an empty array if a is 0.

1(#~]=+)@{::[:p.]*{:
                ]*{:    multiply the entire polynomial with `a`
                        which doesn't change the roots if `a` is nonzero
                        but gives 0 0 0 0 otherwise, which has zero roots according to `p.`
            [:p.    convert to multipler-roots form
1        {::        extract the roots
 (     )@      filter real roots:
  #~]=+        filter the list by "complex conjugate equals self"

Answer 10

Wolfram Language (Mathematica), 50 bytes

Try it online!

g[k_]/;k[[4]]!=0:=Solve[x^Range[0,3].k==0,x,Reals]

Input is the number of coefficients in reverse order {d, c, b, a}. This defines a function g only if the 4th argument is unequal 0. Otherwise this function remains undefined and returns the function call as output. Uses the built-in Solve function over the Reals domain.

Answer 11

Wolfram Language (Mathematica), 35 bytes

#^0NSolve[x#+#2&~Fold~!##,x,Reals]&

Try it online!

Input [a, b, c, d]. If \$a=0\$, returns a list of as many {Indeterminate}s as there would be solutions.

When given non-equation expressions, NSolve finds their roots.

          x#+#2&~Fold~!##           convert to a polynomial in x
   NSolve[               ,x,Reals]  get real solutions
#^0                                 invalidate solutions if a=0

Answer 12

MATLAB, 29 bytes

syms a b c d;roots([a b c d])