1 files changed, 130 insertions, 0 deletions
diff --git a/libm/double/gdtr.c b/libm/double/gdtr.c
new file mode 100644
index 000000000..6b27d9abb
--- /dev/null
+++ b/libm/double/gdtr.c
@@ -0,0 +1,130 @@
+/*							gdtr.c
+ *
+ *	Gamma distribution function
+ *
+ *
+ *
+ * SYNOPSIS:
+ *
+ * double a, b, x, y, gdtr();
+ *
+ * y = gdtr( a, b, x );
+ *
+ *
+ *
+ * DESCRIPTION:
+ *
+ * Returns the integral from zero to x of the gamma probability
+ * density function:
+ *
+ *
+ *                x
+ *        b       -
+ *       a       | |   b-1  -at
+ * y =  -----    |    t    e    dt
+ *       -     | |
+ *      | (b)   -
+ *               0
+ *
+ *  The incomplete gamma integral is used, according to the
+ * relation
+ *
+ * y = igam( b, ax ).
+ *
+ *
+ * ACCURACY:
+ *
+ * See igam().
+ *
+ * ERROR MESSAGES:
+ *
+ *   message         condition      value returned
+ * gdtr domain         x < 0            0.0
+ *
+ */
+/*							gdtrc.c
+ *
+ *	Complemented gamma distribution function
+ *
+ *
+ *
+ * SYNOPSIS:
+ *
+ * double a, b, x, y, gdtrc();
+ *
+ * y = gdtrc( a, b, x );
+ *
+ *
+ *
+ * DESCRIPTION:
+ *
+ * Returns the integral from x to infinity of the gamma
+ * probability density function:
+ *
+ *
+ *               inf.
+ *        b       -
+ *       a       | |   b-1  -at
+ * y =  -----    |    t    e    dt
+ *       -     | |
+ *      | (b)   -
+ *               x
+ *
+ *  The incomplete gamma integral is used, according to the
+ * relation
+ *
+ * y = igamc( b, ax ).
+ *
+ *
+ * ACCURACY:
+ *
+ * See igamc().
+ *
+ * ERROR MESSAGES:
+ *
+ *   message         condition      value returned
+ * gdtrc domain         x < 0            0.0
+ *
+ */
+
+/*							gdtr()  */
+
+
+/*
+Cephes Math Library Release 2.8:  June, 2000
+Copyright 1984, 1987, 1995, 2000 by Stephen L. Moshier
+*/
+
+#include <math.h>
+#ifdef ANSIPROT
+extern double igam ( double, double );
+extern double igamc ( double, double );
+#else
+double igam(), igamc();
+#endif
+
+double gdtr( a, b, x )
+double a, b, x;
+{
+
+if( x < 0.0 )
+	{
+	mtherr( "gdtr", DOMAIN );
+	return( 0.0 );
+	}
+return(  igam( b, a * x )  );
+}
+
+
+
+double gdtrc( a, b, x )
+double a, b, x;
+{
+
+if( x < 0.0 )
+	{
+	mtherr( "gdtrc", DOMAIN );
+	return( 0.0 );
+	}
+return(  igamc( b, a * x )  );
+}