ba237a9d91
Files correlati : Ricompilazione Demo : [ ] Commento : Aggiunti i sorgenti per Greenleaf Math Library (gfm.dll) git-svn-id: svn://10.65.10.50/trunk@10079 c028cbd2-c16b-5b4b-a496-9718f37d4682
101 lines
2.5 KiB
C
Executable File
101 lines
2.5 KiB
C
Executable File
/* int _LnDec80Bit(pDst,pSrc)
|
|
* ARGUMENT
|
|
* pDst is a DEC pointer to the destination.
|
|
* pSrc is a DEC pointer to the source.
|
|
*
|
|
* DESCRIPTION
|
|
* Takes the natural logarithm of pSrc, storing result in pDst.
|
|
*
|
|
* SIDE EFFECTS
|
|
* Possible errors:
|
|
* GM_NULLPOINTER
|
|
* GM_IMAG if pSrc <= 0
|
|
* GM_OVERFLOW if pSrc >= 0 and overflow occurs
|
|
*
|
|
* RETURNS
|
|
* GM_SUCCESS if successful, error code otherwise.
|
|
*
|
|
* ALGORTIHM
|
|
* pSrc = 10^m * 2^(-d) * x .7 <= x < 1.4
|
|
* z = (x - 1) / (x + 1)
|
|
* ln pSrc = m * ln 10 - d * ln 2 + ln x
|
|
* ln x = 2 * SUM [z^(2 * n + 1) / (2 * n + 1)]
|
|
* n = 0, 1, 2, 3,...
|
|
*
|
|
* AUTHOR
|
|
* Jared Levy 4 / 8 / 1987
|
|
* Copyright (C) 1987-1990 Greenleaf Software Inc. All rights reserved.
|
|
*
|
|
* MODIFICATIONS
|
|
*
|
|
*/
|
|
|
|
#include <stdio.h>
|
|
#include "gm.h"
|
|
#include "gmsystem.h"
|
|
|
|
int _LnDec80Bit(pDst,pSrc)
|
|
DEC *pDst, *pSrc;
|
|
{
|
|
|
|
int m;
|
|
DEC *pm, dm, *nsrc, dnsrc, *one, *num, dnum;
|
|
DEC *den, dden, *z, dz, *zsq, dzsq, *pdt, dterm, *fact, dfact;
|
|
DEC *powz, dpowz;
|
|
|
|
if (!_MacIsDecP(pSrc))
|
|
return(GM_IMAG);
|
|
|
|
/* find order of magnitude and normalize to range .1 < nsrc < 1 */
|
|
m = MagnitudeOfDecimal(pSrc) + 1;
|
|
nsrc = &dnsrc;
|
|
_MacDCopy(nsrc, pSrc);
|
|
nsrc->dc.id += m;
|
|
if (nsrc->dc.id > GM_IMAXID) {
|
|
_DivUnsArrByPwrOf10(nsrc->dc.sl, 5, nsrc->dc.id-GM_IMAXID);
|
|
nsrc->dc.id = GM_IMAXID;
|
|
}
|
|
/* calculate m * ln(10) */
|
|
pm = ConvLongToDecimal(&dm, (long) m);
|
|
(void) _MulDec80Bit(pDst, &decLn10, pm);
|
|
|
|
/* normalize to .7 < nsrc <= 1.4, adjust by factors of log 2 */
|
|
while (CompareDecimal(nsrc, &decPoint7) <= 0) {
|
|
/* double nsrc */
|
|
(void) _AddDec80Bit(nsrc, nsrc, nsrc);
|
|
/* subtract a ln 2 */
|
|
(void) _SubDec80Bit(pDst, pDst, &decLn2);
|
|
}
|
|
|
|
/* calculate z = (nsrc - 1) / (nsrc + 1) */
|
|
one = &decOne;
|
|
num = &dnum;
|
|
(void) _SubDec80Bit(&dnum, nsrc, one);
|
|
den = &dden;
|
|
(void) _AddDec80Bit(&dden, nsrc, one);
|
|
z = &dz;
|
|
(void) _DivRndDec80Bit(&dz, num, den, 23);
|
|
|
|
/* calculate log nsrc as an expansion of z */
|
|
zsq = &dzsq;
|
|
(void) _MulDec80Bit(&dzsq, z, z);
|
|
fact = ConvLongToDecimal(&dfact, 3L);
|
|
powz = &dpowz;
|
|
/* first term is 2 * z */
|
|
(void) _AddDec80Bit(powz, z, z);
|
|
(void) _AddDec80Bit(pDst, pDst, powz);
|
|
pdt = &dterm;
|
|
do {
|
|
/* powz = z^(2n+1) */
|
|
(void) _MulDec80Bit(powz, powz, zsq);
|
|
/* term = powz/(2n+1) */
|
|
(void) _DivRndDec80Bit(&dterm, powz, fact,23);
|
|
/* then add new term */
|
|
(void) _AddDec80Bit(pDst, pDst, pdt);
|
|
/* and increase factorial count */
|
|
fact->dc.sl[0] += 2;
|
|
} while (!(_MacIsDecZ(pdt)));
|
|
|
|
return(GM_SUCCESS);
|
|
}
|