PRO solve_realroots_quadratic, a, b, c, x1, x2, nroots, status
; Description: This module solves for the real root(s) of a quadratic equation, of which there
; may be zero, one or two. A quadratic equation has the form:
;
; a*x^2 + b*x + c = 0
;
; where a, b and c are the polynomial coefficients, which also happen to be the
; input parameters "a", "b" and "c" respectively. The module returns the distinct
; real roots "x1" and "x2" in ascending order, and the number "nroots" of distinct
; real roots of the equation. If there are no real roots, then "x1" and "x2" are
; set to "0.0", and if there is only one distinct real root, then "x2" is set to
; "0.0".
; Note that where at least one real root exists, the module makes use of the
; quadratic formula to calculate the real root with the greatest magnitude, and
; then employs a useful relationship between the real roots to obtain the other
; root. The quadratic formula is:
;
; x = ( -b +- sqrt( b^2 - (4*a*c) )) / (2*a)
;
; The useful relationship is:
;
; x2 = c / (a*x1)
;
; This approach avoids the risk of "catastrophic cancellation", where substantial
; loss of precision for one of the roots can occur when the magnitudes of the
; terms "b" and "sqrt( b^2 - (4*a*c) )" in the quadratic equation are similar.
; To allow many similar equations to be solved using a single call to this
; module, any one of the input parameters "a", "b" and "c" may be supplied as a
; vector rather than as a scalar. The remaining input parameters must be supplied
; as scalars. In this case, the output parameters are returned as vectors of the
; same length as the only vector input parameter.
;
; Input Parameters:
;
; a - FLOAT/DOUBLE SCALAR/VECTOR - The polynomial coefficient of the quadratic term in the
; quadratic equation. This parameter must be non-zero (in
; order for the equation to be a quadratic equation). If
; this parameter is a vector, then all of the values must
; be non-zero.
; b - FLOAT/DOUBLE SCALAR/VECTOR - The polynomial coefficient of the linear term in the
; quadratic equation.
; c - FLOAT/DOUBLE SCALAR/VECTOR - The constant term in the quadratic equation.
;
; Output Parameters:
;
; x1 - DOUBLE SCALAR/VECTOR - The smallest distinct real root of the quadratic equation. If
; the equation does not have any real roots, then "x1" is
; assigned a value of "0.0".
; x2 - DOUBLE SCALAR/VECTOR - The largest distinct real root of the quadratic equation. If
; the equation only has one distinct real root, or no real roots,
; then "x2" is assigned a value of "0.0".
; nroots - INTEGER SCALAR/VECTOR - The number of distinct real roots of the quadratic equation
; which may take one of the values "0", "1" or "2".
; status - INTEGER - If the module successfully determines the real root(s) of the quadratic
; equation(s), then "status" is returned with a value of "1", otherwise it
; is returned with a value of "0".
;
; N.B: If all of the input parameters are scalars, then all of the output parameters are
; scalars. If one of the input parameters is a vector, then all of the output parameters
; are vectors of the same length as the input vector. No more than ONE of the input
; parameters is allowed to be supplied as a vector.
;
; Author: Dan Bramich (dan.bramich@hotmail.co.uk)
;
; History:
;
; 27/10/2009 - Module created (dmb).
;Set the default output parameter values
;Determine the dimensions of the input parameters
;Check that at most only one of the input parameters is a vector
;If all of the input parameters are scalars
;Check that the input parameters are of the correct number type
;Check that the polynomial coefficient of the quadratic term in the quadratic equation is not
;zero
;Convert the quadratic equation to a monic quadratic equation by dividing all the polynomial
;coefficients by the coefficient of the quadratic term
;Calculate the discriminant of the quadratic equation
;If the quadratic equation has two distinct real roots, then calculate the values of the real
;roots, and sort them into ascending order
;If "b" is non-zero, then use the quadratic formula to calculate the real root with the
;greatest magnitude, and subsequently calculate the value of the remaining root using the
;useful relationship between the roots.
;Calculate the values of the two distinct real roots
;Sort the roots into ascending order
;If "b" is zero, then use the quadratic formula to calculate the values of both real roots
;Calculate the values of the two distinct real roots
;If the quadratic equation has one distinct real root, then calculate the value of the real
;root
;Calculate the value of the one distinct real root
;Set "status" to "1"
;Return with the roots of the quadratic equation
;If the polynomial coefficient of the quadratic term is a vector
;Check that "a" is of the correct number type, and determine its length
;Check that the remaining parameters are scalars of the correct number type
;Check that "a" has no zero values
;Convert the quadratic equations to monic quadratic equations by dividing all the polynomial
;coefficients by the coefficient of the quadratic term
;If the polynomial coefficient of the linear term is a vector
;Check that "b" is of the correct number type, and determine its length
;Check that the remaining parameters are scalars of the correct number type
;Check that "a" is not zero
;Convert the quadratic equations to monic quadratic equations by dividing all the polynomial
;coefficients by the coefficient of the quadratic term
;If the constant term is a vector
;Check that "c" is of the correct number type, and determine its length
;Check that the remaining parameters are scalars of the correct number type
;Check that "a" is not zero
;Convert the quadratic equations to monic quadratic equations by dividing all the polynomial
;coefficients by the coefficient of the quadratic term
;Set up the output parameters
;Calculate the discriminant of the quadratic equation
;Determine where the quadratic equation(s) have one distinct real root, and where they
;have two distinct real roots. In the cases where the quadratic equation(s) have two
;distinct real roots, distinguish between the cases where "b" is non-zero and where "b"
;is equal to zero.
;If there is at least one quadratic equation with one distinct real root, then for each
;of the relevant quadratic equations calculate the value of the real root
;Calculate the value of the real root
;If there is at least one quadratic equation with two distinct real roots and with "b"
;non-zero, then for each of the relevant quadratic equations calculate the values of the
;real roots, and sort them into ascending order.
;Use the quadratic formula to calculate the real root with the greatest magnitude, and
;subsequently calculate the value of the remaining root using the useful relationship
;between the roots.
;Sort the roots into ascending order
;If there is at least one quadratic equation with two distinct real roots and with "b"
;equal to zero, then for each of the relevant quadratic equations calculate the values
;of the real roots.
;Calculate the values of the two distinct real roots
;Set "status" to "1"