octave-maintainers
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: gammaln(0) (octave-3.3.51+, MinGW)


From: Tatsuro MATSUOKA
Subject: Re: gammaln(0) (octave-3.3.51+, MinGW)
Date: Mon, 2 Aug 2010 09:45:40 +0900 (JST)

Hello 

In liboctave/lo-specfun,

I found as:

double
xlgamma (double x)
{
#if defined (HAVE_LGAMMA)
  return lgamma (x);
#else
  double result;
  double sgngam;

  if (xisnan (x))
    result = x;
  else if (xisinf (x))
    result = octave_Inf;
  else
    F77_XFCN (dlgams, DLGAMS, (x, result, sgngam));

  return result;
#endif
}

As I posted previously,

I found the below definition,

#define HAVE_LGAMMA 1
#define HAVE_LGAMMAF 1

However it seemed to be ignored, and therefore stalec was used.

Hmmmm ??

Apart from the failure in detecting the definition, the above description 
should be modified, because
stlatec gamma functions do not allow x=0.0, and beta.m assume that gammaln(0) 
is infinity as shown in
their own test as

%!test
%! a = 0.25 + (0:5) * 0.5;
%! tol = 10 * max (a) * eps;
%! assert (zeros (size (a)), beta (a, -a), tol)
%! assert (zeros (size (a)), beta (-a, a), tol)

Therefore in the case of slatec gamma function, perhaps x==0 should be treated

  if (xisnan (x))
    result = x;
  else if (xisinf (x))
    result = octave_Inf;
+  else if (x==0.0)
+    result = octave_Inf;
  else
    F77_XFCN (dlgams, DLGAMS, (x, result, sgngam));

  return result;

Is the above correct ?  It should be like 'abs(x)<tol' ?


Regards

Tatsuro

  

--- Tatsuro MATSUOKA wrote:

> Hello
> 
> > In gfortran DLGAMA() is intrinsically defined.  
> > 
> > If it is used, perhaps the error will be disappeared.
> > http://gcc.gnu.org/onlinedocs/gfortran/LOG_005fGAMMA.html#LOG_005fGAMMA
> > 
> > In the case slatec gamma function, we cannot use it at x=0 so that special 
> > treatment is
> > required, for
> > gammaln in octave, I think.
> > 
> > Am I wrong?
> 
> In gcc, lgamma, i.e.log(abs(gamma)) is defined. (also lgammaf) 
> I have looked into "config.log" and found
> configure:47767: checking for lgamma
> <snip>
> configure:47767: result: yes
> configure:47767: checking for lgammaf
> <snip>
> configure:47767: result: yes
> 
> And in config.h, there exist
> 
> /* Define to 1 if you have the `lgamma' function. */
> #define HAVE_LGAMMA 1
> 
> /* Define to 1 if you have the `lgammaf' function. */
> #define HAVE_LGAMMAF 1
> 
> 
> Perhaps for gamma() in the octave-3.3.51+ built by me, lgamma is used.  
> However, for gammaln() in, subroutine dlgams in slatec is used.
> 
> Because I can see in the error message,
>  error: lgamma: exception encountered in Fortran subroutine dlgams_
> 
> Is there any reason that buid-in lgamma in gcc is not used but the slatec 
> dlgams is used,
> although
> 'configure' detects it ?
> 
> Regards
> 
> Tatsuro
> 
> --- Tatsuro MATSUOKA  wrote:
> 
> > Hello
> > 
> > > I cannot find source of xstopx in libcruft source files.
> > > 
> > > However when the code link against libcruft.dll.a, undefined symbol error 
> > > messages for
> > 'xstopx'
> > > were
> > > disappeared.
> > 
> > I found the reason. xstopx is defined in   misc/f77-fcn.h in libcruft.
> > 
> > > However
> > > 
> > >       program test
> > >       DOUBLE PRECISION X, DLGAM, SGNGAM, DLNGAM
> > >       X=0
> > >       write(*,*) DGAMMA(X)
> > >       end
> > > 
> > > The above code gives
> > > 
> > > $ ./test
> > >                  +Infinity
> > 
> > This is reasonable.  In this case, dgamma intrinsic in gfortran  is used.
> > If I add 'EXTERNAL DGAMMA', slatac dgamma is used and the same error as 
> > that by the dlgams
> > subroutime
> > is caused.
> > 
> > In gfortran DLGAMA() is intrinsically defined.  
> > 
> > If it is used, perhaps the error will be disappeared.
> > http://gcc.gnu.org/onlinedocs/gfortran/LOG_005fGAMMA.html#LOG_005fGAMMA
> > 
> > In the case slatec gamma function, we cannot use it at x=0 so that special 
> > treatment is
> > required, for
> > gammaln in octave, I think.
> > 
> > Am I wrong?
> > 
> > Regards
> > 
> > Tatsuro
> > 
> > 
> > 
> > --- Tatsuro MATSUOKA  wrote:
> > 
> > > Hello
> > > 
> > > I have reported error in test beta.m.
> > > Development source on July 2010 (platform MinGW, GCC-4.5.0)
> > > 
> > > 
> > > http://octave.1599824.n4.nabble.com/error-in-test-beta-m-tc2303186.html#a2303407
> > > 
> > > This problem reduced to 'gammaln(0)'
> > > 
> > > octave:5> gammaln(0)
> > >  ***MESSAGE FROM ROUTINE DGAMMA IN LIBRARY SLATEC.
> > >  ***FATAL ERROR, PROG ABORTED, TRACEBACK REQUESTED
> > >  *  X IS 0
> > >  *  ERROR NUMBER = 4
> > >  *
> > >  ***END OF MESSAGE
> > > 
> > >  ***JOB ABORT DUE TO FATAL ERROR.
> > > 0          ERROR MESSAGE SUMMARY
> > >  LIBRARY    SUBROUTINE MESSAGE START             NERR     LEVEL     COUNT
> > >  SLATEC     DGAMMA     X IS 0                       4         2         2
> > > 
> > > error: lgamma: exception encountered in Fortran subroutine dlgams_
> > > 
> > > I can fined symbol 'dlgams'
> > > 
> > > I have been trying to build simple test file using DLGAMS in slatec-fn 
> > > like
> > > 
> > >       program test
> > >       DOUBLE PRECISION X, DLGAM, SGNGAM, DLNGAM
> > >       X=0
> > >       call DLGAMS (X, DLGAM, SGNGAM)
> > >       write (*,*) DLGAM, SGNGAM
> > >       end
> > > 
> > > In the slatec-err routines (or other libcruft source files), I found 'call
> xstopx('messeage')'
> > > but
> > > I cannot find source of xstopx in libcruft source files.
> > > 
> > > However when the code link against libcruft.dll.a, undefined symbol error 
> > > messages for
> > 'xstopx'
> > > were
> > > disappeared.
> > > 
> > > And test result is 
> > > 
> > > $ gdb test
> > > GNU gdb (GDB) 7.0.1
> > > Copyright (C) 2009 Free Software Foundation, Inc.
> > > License GPLv3+: GNU GPL version 3 or later 
> > > <http://gnu.org/licenses/gpl.html>
> > > This is free software: you are free to change and redistribute it.
> > > There is NO WARRANTY, to the extent permitted by law.  Type "show copying"
> > > and "show warranty" for details.
> > > This GDB was configured as "mingw32".
> > > For bug reporting instructions, please see:
> > > <http://www.gnu.org/software/gdb/bugs/>...
> > > Reading symbols from 
> > > c:\usr\Tatsu\mingwhome\octaves\dgammatest/test.exe...done.
> > > (gdb) br test
> > > Breakpoint 1 at 0x4013c9: file test.f, line 3.
> > > (gdb) run
> > > Starting program: c:\usr\Tatsu\mingwhome\octaves\dgammatest/test.exe
> > > [New Thread 1048.0x8fc]
> > > 
> > > Breakpoint 1, test () at test.f:3
> > > 3             X=0
> > > Current language:  auto
> > > The current source language is "auto; currently fortran".
> > > (gdb) step
> > > 4             call DLGAMS (X, DLGAM, SGNGAM)
> > > (gdb) step
> > >  ***MESSAGE FROM ROUTINE DGAMMA IN LIBRARY SLATEC.
> > >  ***FATAL ERROR, PROG ABORTED, TRACEBACK REQUESTED
> > >  *  X IS 0
> > >  *  ERROR NUMBER = 4
> > >  *
> > >  ***END OF MESSAGE
> > > 
> > >  ***JOB ABORT DUE TO FATAL ERROR.
> > > 0          ERROR MESSAGE SUMMARY
> > >  LIBRARY    SUBROUTINE MESSAGE START             NERR     LEVEL     COUNT
> > >  SLATEC     DGAMMA     X IS 0                       4         2         1
> > > 
> > > gdb: unknown target exception 0xc0000027 at 0x77be5464
> > > 
> > > Program received signal ?, Unknown signal.
> > > 0x77be5464 in msvcrt!_global_unwind2 () from 
> > > C:\WINDOWS\system32\msvcrt.dll
> > > (gdb) bt
> > > #0  0x77be5464 in msvcrt!_global_unwind2 () from 
> > > C:\WINDOWS\system32\msvcrt.dll
> > > #1  0x77be6d8c in msvcrt!longjmp () from C:\WINDOWS\system32\msvcrt.dll
> > > #2  0x00000000 in ?? ()
> > > 
> > > I have obtained the same error gcc-4.4.0 / MinGW
> > > 
> > > Unfortunately I did not install gdb on cygwin but I could the following 
> > > message in cygwin. 
> > > 
> > >  ***MESSAGE FROM ROUTINE DGAMMA IN LIBRARY SLATEC.
> > >  ***FATAL ERROR, PROG ABORTED, TRACEBACK REQUESTED
> > >  *  X IS 0
> > >  *  ERROR NUMBER = 4
> > >  *
> > >  ***END OF MESSAGE
> > > 
> > >  ***JOB ABORT DUE TO FATAL ERROR.
> > > 0          ERROR MESSAGE SUMMARY
> > >  LIBRARY    SUBROUTINE MESSAGE START             NERR     LEVEL     COUNT
> > >  SLATEC     DGAMMA     X IS 0                       4         2         1
> > > 
> > > ***********************************
> > > 
> > > The error message claim that X is 0 in DGAMMA.
> > > 
> > > However
> > > 
> > >       program test
> > >       DOUBLE PRECISION X, DLGAM, SGNGAM, DLNGAM
> > >       X=0
> > >       write(*,*) DGAMMA(X)
> > >       end
> > > 
> > > The above code gives
> > > 
> > > $ ./test
> > >                  +Infinity
> > > 
> > > ???????
> > > 
> > > Any suggestions?
> > > 
> > > Regards
> > > 
> > > Tatsuro
> > > 
> > > 
> > > 
> > > 
> > > --------------------------------------
> > > Get the new Internet Explorer 8 optimized for Yahoo! JAPAN
> > > http://pr.mail.yahoo.co.jp/ie8/
> > > 
> > 
> > 
> > --------------------------------------
> > Get the new Internet Explorer 8 optimized for Yahoo! JAPAN
> > http://pr.mail.yahoo.co.jp/ie8/
> > 
> 
> 
> --------------------------------------
> Get the new Internet Explorer 8 optimized for Yahoo! JAPAN
> http://pr.mail.yahoo.co.jp/ie8/
> 


--------------------------------------
Get the new Internet Explorer 8 optimized for Yahoo! JAPAN
http://pr.mail.yahoo.co.jp/ie8/


reply via email to

[Prev in Thread] Current Thread [Next in Thread]