[yorick-ygsl] 01/02: New functions gsl_poly_solve_quadratic and gsl_poly_solve_cubic.

[Thibaut Jean-Claude Paumard]

This is an automated email from the git hooks/post-receive script.

thibaut pushed a commit to annotated tag upstream/1.2.0
in repository yorick-ygsl.

commit c073d6dc66a9986b7b882de61d314fa6d69ef498
Author: Éric Thiébaut <eric.thiebaut at univ-lyon1.fr>
Date:   Fri May 16 13:03:05 2014 +0200

    New functions gsl_poly_solve_quadratic and gsl_poly_solve_cubic.
---
 NEWS   |  4 ++++
 gsl.i  | 27 ++++++++++++++++++++++
 ygsl.c | 83 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 3 files changed, 114 insertions(+)

diff --git a/NEWS b/NEWS
index 466c0e6..747359f 100644
--- a/NEWS
+++ b/NEWS
@@ -1,3 +1,7 @@
+ * 2014-05-16: Release 1.1.1
+   New functions gsl_poly_solve_quadratic and gsl_poly_solve_cubic to find
+   the roots of quadratic and cubic polynomials.
+
  * 2013-07-09: Release 1.1.1
    README.md used in place of README.
 
diff --git a/gsl.i b/gsl.i
index a4c4be2..7f48062 100644
--- a/gsl.i
+++ b/gsl.i
@@ -857,6 +857,33 @@ extern gsl_sf_eta;
  * SEE ALSO: gsl_sf.
  */
 
+extern gsl_poly_solve_quadratic;
+extern gsl_poly_solve_cubic;
+/* DOCUMENT x = gsl_poly_solve_quadratic(a, b, c);
+         or x = gsl_poly_solve_quadratic(v);
+         or x = gsl_poly_solve_cubic(a, b, c);
+         or x = gsl_poly_solve_cubic(v);
+
+      These functions return the real roots of a quadratic or cubic polynomials
+      with real coefficients A, B and C.  When called with a single argument,
+      it must be a vector of coefficients: V = [A,B,C].
+
+      If there are no roots, an empty result is returned otherwise a vector of
+      1, 2, or 3 roots is returned.  The roots are sorted in ascending order.
+      The case of coincident roots is not considered special.  Therefore,
+      either 0 or 1 or 2 roots are returned for a quadratic polynomial (a
+      single root only occurs if A=0) and either 0 or 1 or 3 roots are
+      returned for a cubic polynomial.
+
+      The roots X are such that:
+
+        A*X^2 + B*X + C = 0
+        X^3 + A*X^2 + B*X + C = 0
+
+      for a quadratic and a cubic polynomial respectively.
+
+*/
+
 /*
  * Local Variables:
  * mode: Yorick
diff --git a/ygsl.c b/ygsl.c
index 17f6c1c..004b073 100644
--- a/ygsl.c
+++ b/ygsl.c
@@ -96,6 +96,7 @@
 #include <yapi.h>
 
 #include <gsl/gsl_sf.h>
+#include <gsl/gsl_poly.h>
 
 /* Define some macros to get rid of some GNU extensions when not compiling
    with GCC. */
@@ -492,6 +493,87 @@ FN(gsl_sf_bessel_lnKnu)
 FN(gsl_sf_bessel_Knu_scaled)
 #undef FN
 
+/*---------------------------------------------------------------------------*/
+/* POLYNOMIAL ROOTS */
+
+static void push_vector_d(long n, const double inp[])
+{
+  long dims[2];
+  double* out;
+  long j;
+
+  if (n > 0) {
+    dims[0] = 1;
+    dims[1] = n;
+    out = ypush_d(dims);
+    for (j = 0; j < n; ++j) {
+      out[j] = inp[j];
+    }
+  } else {
+    ypush_nil();
+  }
+}
+
+void Y_gsl_poly_solve_quadratic(int argc)
+{
+  double a, b, c;
+  long dims[Y_DIMSIZE];
+  double x[2];
+  const double* coef;
+  long n;
+
+  if (argc == 1) {
+    coef = ygeta_d(0, &n, dims);
+    if (dims[0] == 1 && dims[1] == 3) {
+      a = coef[0];
+      b = coef[1];
+      c = coef[2];
+    } else {
+      goto bad_args;
+    }
+  } else if (argc == 3) {
+    a = ygets_d(2);
+    b = ygets_d(1);
+    c = ygets_d(0);
+  } else {
+    bad_args:
+    y_error("expecting a 3-element vector or 3 arguments");
+    return;
+  }
+  n = gsl_poly_solve_quadratic(a, b, c, &x[0], &x[1]);
+  push_vector_d(n, x);
+}
+
+void Y_gsl_poly_solve_cubic(int argc)
+{
+  double a, b, c;
+  long dims[Y_DIMSIZE];
+  double x[3];
+  const double* coef;
+  long n;
+
+  if (argc == 1) {
+    coef = ygeta_d(0, &n, dims);
+    if (dims[0] == 1 && dims[1] == 3) {
+      a = coef[0];
+      b = coef[1];
+      c = coef[2];
+    } else {
+      goto bad_args;
+    }
+  } else if (argc == 3) {
+    a = ygets_d(2);
+    b = ygets_d(1);
+    c = ygets_d(0);
+  } else {
+    bad_args:
+    y_error("expecting a 3-element vector or 3 arguments");
+    return;
+  }
+  n = gsl_poly_solve_cubic(a, b, c, &x[0], &x[1], &x[2]);
+  push_vector_d(n, x);
+}
+
 /*
  * Local Variables:
  * mode: C
@@ -500,5 +582,6 @@ FN(gsl_sf_bessel_Knu_scaled)
  * indent-tabs-mode: nil
  * fill-column: 78
  * coding: utf-8
+ * ispell-local-dictionary: "american"
  * End:
  */

-- 
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/yorick-ygsl.git