re PR fortran/32980 (Vendor extension: Intrinsic functions (D)GAMMA, LGAMMA (ALGAMA/DLGAMA))

2007-08-26  Tobias Burnus  <burnus@net-b.de>

	PR fortran/32980
	* intrinsic.h (gfc_simplify_gamma,gfc_simplify_lgamma,
	gfc_resolve_gamma,gfc_resolve_lgamma): New function declations.
	* mathbuiltins.def: Define GAMMA and LGAMMA.
	* intrinsic.c (add_functions): Add GAMMA, DGAMMA, LGAMMA, ALGAMA
	and DLGAMA.
	* simplify.c (gfc_simplify_gamma,gfc_simplify_lgamma): New functions.
	* iresolve.c (gfc_resolve_gamma,gfc_resolve_lgamma): New functions.
	* intrinsic.texi: Add documentation for GAMMA and LGAMMA.

2007-08-26  Tobias Burnus  <burnus@net-b.de>

	PR fortran/32980
	* gfortran.dg/gamma_1.f90: New.
	* gfortran.dg/gamma_2.f90: New.
	* gfortran.dg/gamma_3.f90: New.

From-SVN: r127809

7 files changed, 215 insertions, 0 deletions

diff --git a/gcc/fortran/ChangeLog b/gcc/fortran/ChangeLog
index 9adaf07..2597164 100644
--- a/gcc/fortran/ChangeLog
+++ b/gcc/fortran/ChangeLog
@@ -1,5 +1,17 @@
 2007-08-26  Tobias Burnus  <burnus@net-b.de>
 
+	PR fortran/32980
+	* intrinsic.h (gfc_simplify_gamma,gfc_simplify_lgamma,
+	gfc_resolve_gamma,gfc_resolve_lgamma): New function declations.
+	* mathbuiltins.def: Define GAMMA and LGAMMA.
+	* intrinsic.c (add_functions): Add GAMMA, DGAMMA, LGAMMA, ALGAMA
+	and DLGAMA.
+	* simplify.c (gfc_simplify_gamma,gfc_simplify_lgamma): New functions.
+	* iresolve.c (gfc_resolve_gamma,gfc_resolve_lgamma): New functions.
+	* intrinsic.texi: Add documentation for GAMMA and LGAMMA.
+
+2007-08-26  Tobias Burnus  <burnus@net-b.de>
+
 	PR fortran/33188
 	* parse.c (parse_derived): Support empty derived type
 	definitions for Fortran 2003.
diff --git a/gcc/fortran/intrinsic.c b/gcc/fortran/intrinsic.c
index 4128e44..d273f80 100644
--- a/gcc/fortran/intrinsic.c
+++ b/gcc/fortran/intrinsic.c
@@ -1453,6 +1453,16 @@ add_functions (void)
 
   make_generic ("fput", GFC_ISYM_FPUT, GFC_STD_GNU);
 
+  add_sym_1 ("gamma", GFC_ISYM_GAMMA, CLASS_ELEMENTAL, ACTUAL_NO, BT_REAL, dr, GFC_STD_GNU,
+	     gfc_check_fn_r, gfc_simplify_gamma, gfc_resolve_gamma,
+	     x, BT_REAL, dr, REQUIRED);
+
+  add_sym_1 ("dgamma", GFC_ISYM_GAMMA, CLASS_ELEMENTAL, ACTUAL_NO, BT_REAL, dr, GFC_STD_GNU,
+	     gfc_check_fn_r, gfc_simplify_gamma, gfc_resolve_gamma,
+	     x, BT_REAL, dr, REQUIRED);
+
+  make_generic ("gamma", GFC_ISYM_GAMMA, GFC_STD_GNU);
+
   /* Unix IDs (g77 compatibility)  */
   add_sym_1 ("getcwd", GFC_ISYM_GETCWD, NO_CLASS, ACTUAL_NO, BT_INTEGER, di,  GFC_STD_GNU,
 	     NULL, NULL, gfc_resolve_getcwd,
@@ -1690,6 +1700,21 @@ add_functions (void)
 
   make_generic ("len_trim", GFC_ISYM_LEN_TRIM, GFC_STD_F95);
 
+  add_sym_1 ("lgamma", GFC_ISYM_LGAMMA, CLASS_ELEMENTAL, ACTUAL_NO, BT_REAL, dr, GFC_STD_GNU,
+	     gfc_check_fn_r, gfc_simplify_lgamma, gfc_resolve_lgamma,
+	     x, BT_REAL, dr, REQUIRED);
+
+  add_sym_1 ("algama", GFC_ISYM_LGAMMA, CLASS_ELEMENTAL, ACTUAL_NO, BT_REAL, dr, GFC_STD_GNU,
+	     gfc_check_fn_r, gfc_simplify_lgamma, gfc_resolve_lgamma,
+	     x, BT_REAL, dr, REQUIRED);
+
+  add_sym_1 ("dlgama", GFC_ISYM_LGAMMA, CLASS_ELEMENTAL, ACTUAL_NO, BT_REAL, dr, GFC_STD_GNU,
+	     gfc_check_fn_r, gfc_simplify_lgamma, gfc_resolve_lgamma,
+	     x, BT_REAL, dr, REQUIRED);
+
+  make_generic ("lgamma", GFC_ISYM_LGAMMA, GFC_STD_GNU);
+
+
   add_sym_2 ("lge", GFC_ISYM_LGE, CLASS_ELEMENTAL, ACTUAL_NO, BT_LOGICAL, dl, GFC_STD_F77,
 	     NULL, gfc_simplify_lge, NULL,
 	     sta, BT_CHARACTER, dc, REQUIRED, stb, BT_CHARACTER, dc, REQUIRED);
diff --git a/gcc/fortran/intrinsic.h b/gcc/fortran/intrinsic.h
index cf242b8..e284a6c 100644
--- a/gcc/fortran/intrinsic.h
+++ b/gcc/fortran/intrinsic.h
@@ -221,6 +221,7 @@ gfc_expr *gfc_simplify_exponent (gfc_expr *);
 gfc_expr *gfc_simplify_float (gfc_expr *);
 gfc_expr *gfc_simplify_floor (gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_fraction (gfc_expr *);
+gfc_expr *gfc_simplify_gamma (gfc_expr *);
 gfc_expr *gfc_simplify_huge (gfc_expr *);
 gfc_expr *gfc_simplify_iachar (gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_iand (gfc_expr *, gfc_expr *);
@@ -243,6 +244,7 @@ gfc_expr *gfc_simplify_kind (gfc_expr *);
 gfc_expr *gfc_simplify_lbound (gfc_expr *, gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_len (gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_len_trim (gfc_expr *, gfc_expr *);
+gfc_expr *gfc_simplify_lgamma (gfc_expr *);
 gfc_expr *gfc_simplify_lge (gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_lgt (gfc_expr *, gfc_expr *);
 gfc_expr *gfc_simplify_lle (gfc_expr *, gfc_expr *);
@@ -354,6 +356,7 @@ void gfc_resolve_fget (gfc_expr *, gfc_expr *);
 void gfc_resolve_fputc (gfc_expr *, gfc_expr *, gfc_expr *);
 void gfc_resolve_fput (gfc_expr *, gfc_expr *);
 void gfc_resolve_g77_math1 (gfc_expr *, gfc_expr *);
+void gfc_resolve_gamma (gfc_expr *, gfc_expr *);
 void gfc_resolve_getcwd (gfc_expr *, gfc_expr *);
 void gfc_resolve_getgid (gfc_expr *);
 void gfc_resolve_getpid (gfc_expr *);
@@ -384,6 +387,7 @@ void gfc_resolve_kill (gfc_expr *, gfc_expr *, gfc_expr *);
 void gfc_resolve_lbound (gfc_expr *, gfc_expr *, gfc_expr *, gfc_expr *);
 void gfc_resolve_len (gfc_expr *, gfc_expr *, gfc_expr *);
 void gfc_resolve_len_trim (gfc_expr *, gfc_expr *, gfc_expr *);
+void gfc_resolve_lgamma (gfc_expr *, gfc_expr *);
 void gfc_resolve_link (gfc_expr *, gfc_expr *, gfc_expr *);
 void gfc_resolve_loc (gfc_expr *, gfc_expr *);
 void gfc_resolve_log (gfc_expr *, gfc_expr *);
diff --git a/gcc/fortran/intrinsic.texi b/gcc/fortran/intrinsic.texi
index f09246c..e94a7e3 100644
--- a/gcc/fortran/intrinsic.texi
+++ b/gcc/fortran/intrinsic.texi
@@ -121,6 +121,7 @@ Some basic guidelines for editing this document:
 * @code{FSEEK}:         FSEEK,     Low level file positioning subroutine
 * @code{FSTAT}:         FSTAT,     Get file status
 * @code{FTELL}:         FTELL,     Current stream position
+* @code{GAMMA}:         GAMMA,     Gamma function
 * @code{GERROR}:        GERROR,    Get last system error message
 * @code{GETARG}:        GETARG,    Get command line arguments
 * @code{GET_COMMAND}:   GET_COMMAND, Get the entire command line
@@ -161,6 +162,7 @@ Some basic guidelines for editing this document:
 * @code{LBOUND}:        LBOUND,    Lower dimension bounds of an array
 * @code{LEN}:           LEN,       Length of a character entity
 * @code{LEN_TRIM}:      LEN_TRIM,  Length of a character entity without trailing blank characters
+* @code{LGAMMA}:        LGAMMA,    Logarithm of the Gamma function
 * @code{LGE}:           LGE,       Lexical greater than or equal
 * @code{LGT}:           LGT,       Lexical greater than
 * @code{LINK}:          LINK,      Create a hard link
@@ -4484,6 +4486,65 @@ END PROGRAM
 
 
 
+@node GAMMA
+@section @code{GAMMA} --- Gamma function
+@fnindex GAMMA
+@fnindex DGAMMA
+@cindex Gamma function
+@cindex Factorial function
+
+@table @asis
+@item @emph{Description}:
+@code{GAMMA(X)} computes Gamma (@math{\Gamma}) of @var{X}. For positive,
+integer values of @var{X} the Gamma function simplifies to the factorial
+function @math{\Gamma(x)=(x-1)!}.
+
+@tex
+$$
+\Gamma(x) = \int_0^\infty t^{x-1}{\rm e}^{-t}\,{\rm d}t
+$$
+@end tex
+
+@item @emph{Standard}:
+GNU Extension
+
+@item @emph{Class}:
+Elemental function
+
+@item @emph{Syntax}:
+@code{X = GAMMA(X)}
+
+@item @emph{Arguments}:
+@multitable @columnfractions .15 .70
+@item @var{X} @tab Shall be of type @code{REAL} and neither zero
+nor a negative integer.
+@end multitable
+
+@item @emph{Return value}:
+The return value is of type @code{REAL} of the same kind as @var{X}.
+
+@item @emph{Example}:
+@smallexample
+program test_gamma
+  real :: x = 1.0
+  x = gamma(x) ! returns 1.0
+end program test_gamma
+@end smallexample
+
+@item @emph{Specific names}:
+@multitable @columnfractions .20 .20 .20 .25
+@item Name             @tab Argument         @tab Return type       @tab Standard
+@item @code{GAMMA(X)}  @tab @code{REAL(4) X} @tab @code{REAL(4)}    @tab GNU Extension
+@item @code{DGAMMA(X)} @tab @code{REAL(8) X} @tab @code{REAL(8)}    @tab GNU Extension
+@end multitable
+
+@item @emph{See also}:
+Logarithm of the Gamma function: @ref{LGAMMA}
+
+@end table
+
+
+
 @node GERROR
 @section @code{GERROR} --- Get last system error message
 @fnindex GERROR
@@ -6230,6 +6291,60 @@ The return value is of type @code{INTEGER} and of kind @var{KIND}. If
 
 
 
+@node LGAMMA
+@section @code{LGAMMA} --- Logarithm of the Gamma function
+@fnindex GAMMA
+@fnindex ALGAMA
+@fnindex DLGAMA
+@cindex Gamma function, logarithm of
+@cindex 
+
+@table @asis
+@item @emph{Description}:
+@code{GAMMA(X)} computes the natural logrithm of the absolute value of the
+Gamma (@math{\Gamma}) function.
+
+@item @emph{Standard}:
+GNU Extension
+
+@item @emph{Class}:
+Elemental function
+
+@item @emph{Syntax}:
+@code{X = LGAMMA(X)}
+
+@item @emph{Arguments}:
+@multitable @columnfractions .15 .70
+@item @var{X} @tab Shall be of type @code{REAL} and neither zero
+nor a negative integer.
+@end multitable
+
+@item @emph{Return value}:
+The return value is of type @code{REAL} of the same kind as @var{X}.
+
+@item @emph{Example}:
+@smallexample
+program test_log_gamma
+  real :: x = 1.0
+  x = lgamma(x) ! returns 0.0
+end program test_log_gamma
+@end smallexample
+
+@item @emph{Specific names}:
+@multitable @columnfractions .20 .20 .20 .25
+@item Name             @tab Argument         @tab Return type       @tab Standard
+@item @code{LGAMMA(X)} @tab @code{REAL(4) X} @tab @code{REAL(4)}    @tab GNU Extension
+@item @code{ALGAMA(X)} @tab @code{REAL(4) X} @tab @code{REAL(4)}    @tab GNU Extension
+@item @code{DLGAMA(X)} @tab @code{REAL(8) X} @tab @code{REAL(8)}    @tab GNU Extension
+@end multitable
+
+@item @emph{See also}:
+Gamma function: @ref{GAMMA}
+
+@end table
+
+
+
 @node LGE
 @section @code{LGE} --- Lexical greater than or equal
 @fnindex LGE
diff --git a/gcc/fortran/iresolve.c b/gcc/fortran/iresolve.c
index fc837e1..7948b14 100644
--- a/gcc/fortran/iresolve.c
+++ b/gcc/fortran/iresolve.c
@@ -757,6 +757,15 @@ gfc_resolve_g77_math1 (gfc_expr *f, gfc_expr *x)
 
 
 void
+gfc_resolve_gamma (gfc_expr *f, gfc_expr *x)
+{
+  f->ts = x->ts;
+  f->value.function.name
+    = gfc_get_string ("__gamma_%d", x->ts.kind);
+}
+
+
+void
 gfc_resolve_getcwd (gfc_expr *f, gfc_expr *n ATTRIBUTE_UNUSED)
 {
   f->ts.type = BT_INTEGER;
@@ -1114,6 +1123,15 @@ gfc_resolve_len_trim (gfc_expr *f, gfc_expr *string, gfc_expr *kind)
 
 
 void
+gfc_resolve_lgamma (gfc_expr *f, gfc_expr *x)
+{
+  f->ts = x->ts;
+  f->value.function.name
+    = gfc_get_string ("__lgamma_%d", x->ts.kind);
+}
+
+
+void
 gfc_resolve_link (gfc_expr *f, gfc_expr *p1 ATTRIBUTE_UNUSED,
 		  gfc_expr *p2 ATTRIBUTE_UNUSED)
 {
diff --git a/gcc/fortran/mathbuiltins.def b/gcc/fortran/mathbuiltins.def
index 0fc7368..33e87d1 100644
--- a/gcc/fortran/mathbuiltins.def
+++ b/gcc/fortran/mathbuiltins.def
@@ -30,3 +30,5 @@ DEFINE_MATH_BUILTIN   (Y1,    "y1",     0)
 DEFINE_MATH_BUILTIN   (YN,    "yn",     2)
 DEFINE_MATH_BUILTIN   (ERF,   "erf",    0)
 DEFINE_MATH_BUILTIN   (ERFC,  "erfc",   0)
+DEFINE_MATH_BUILTIN   (GAMMA, "gamma",  0)
+DEFINE_MATH_BUILTIN   (LGAMMA,"lgamma", 0)
diff --git a/gcc/fortran/simplify.c b/gcc/fortran/simplify.c
index a395b04..5e129ae 100644
--- a/gcc/fortran/simplify.c
+++ b/gcc/fortran/simplify.c
@@ -1183,6 +1183,24 @@ gfc_simplify_fraction (gfc_expr *x)
 
 
 gfc_expr *
+gfc_simplify_gamma (gfc_expr *x)
+{
+  gfc_expr *result;
+
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+
+  gfc_set_model_kind (x->ts.kind);
+
+  mpfr_gamma (result->value.real, x->value.real, GFC_RND_MODE);
+
+  return range_check (result, "GAMMA");
+}
+
+
+gfc_expr *
 gfc_simplify_huge (gfc_expr *e)
 {
   gfc_expr *result;
@@ -2212,6 +2230,27 @@ gfc_simplify_len_trim (gfc_expr *e, gfc_expr *kind)
   return range_check (result, "LEN_TRIM");
 }
 
+gfc_expr *
+gfc_simplify_lgamma (gfc_expr *x __attribute__((unused)))
+{
+#if MPFR_VERSION >= MPFR_VERSION_NUM(2,3,0)
+  gfc_expr *result;
+
+  if (x->expr_type != EXPR_CONSTANT)
+    return NULL;
+
+  result = gfc_constant_result (x->ts.type, x->ts.kind, &x->where);
+
+  gfc_set_model_kind (x->ts.kind);
+
+  mpfr_lgamma (result->value.real, x->value.real, GFC_RND_MODE);
+
+  return range_check (result, "LGAMMA");
+#else
+  return NULL;
+#endif
+}
+
 
 gfc_expr *
 gfc_simplify_lge (gfc_expr *a, gfc_expr *b)