← Back to branch summary

~ubuntu-branches/ubuntu/karmic/grace/karmic

« back to all changes in this revision

Viewing changes to cephes/yn.c

  • Committer: Bazaar Package Importer
  • Author(s): Torsten Werner
  • Date: 2002-03-19 14:19:58 UTC
  • Revision ID: james.westby@ubuntu.com-20020319141958-5gxna6vo1ek3zjml
Tags: upstream-5.1.7
Import upstream version 5.1.7

Show diffs side-by-side

added added

removed removed

Lines of Context:
cephes/yn.c
 
1
/*                                                      yn.c
 
2
 *
 
3
 *      Bessel function of second kind of integer order
 
4
 *
 
5
 *
 
6
 *
 
7
 * SYNOPSIS:
 
8
 *
 
9
 * double x, y, yn();
 
10
 * int n;
 
11
 *
 
12
 * y = yn( n, x );
 
13
 *
 
14
 *
 
15
 *
 
16
 * DESCRIPTION:
 
17
 *
 
18
 * Returns Bessel function of order n, where n is a
 
19
 * (possibly negative) integer.
 
20
 *
 
21
 * The function is evaluated by forward recurrence on
 
22
 * n, starting with values computed by the routines
 
23
 * y0() and y1().
 
24
 *
 
25
 * If n = 0 or 1 the routine for y0 or y1 is called
 
26
 * directly.
 
27
 *
 
28
 *
 
29
 *
 
30
 * ACCURACY:
 
31
 *
 
32
 *
 
33
 *                      Absolute error, except relative
 
34
 *                      when y > 1:
 
35
 * arithmetic   domain     # trials      peak         rms
 
36
 *    DEC       0, 30        2200       2.9e-16     5.3e-17
 
37
 *    IEEE      0, 30       30000       3.4e-15     4.3e-16
 
38
 *
 
39
 *
 
40
 * ERROR MESSAGES:
 
41
 *
 
42
 *   message         condition      value returned
 
43
 * yn singularity   x = 0              MAXNUM
 
44
 * yn overflow                         MAXNUM
 
45
 *
 
46
 * Spot checked against tables for x, n between 0 and 100.
 
47
 *
 
48
 */
 
49
 
 
50
/*
 
51
Cephes Math Library Release 2.1:  December, 1988
 
52
Copyright 1984, 1987 by Stephen L. Moshier
 
53
Direct inquiries to 30 Frost Street, Cambridge, MA 02140
 
54
*/
 
55
 
 
56
#include "mconf.h"
 
57
#include "cephes.h"
 
58
 
 
59
#ifndef HAVE_YN
 
60
 
 
61
extern double MAXNUM, MAXLOG;
 
62
 
 
63
double yn( n, x )
 
64
int n;
 
65
double x;
 
66
{
 
67
double an, anm1, anm2, r;
 
68
int k, sign;
 
69
 
 
70
if( n < 0 )
 
71
        {
 
72
        n = -n;
 
73
        if( (n & 1) == 0 )      /* -1**n */
 
74
                sign = 1;
 
75
        else
 
76
                sign = -1;
 
77
        }
 
78
else
 
79
        sign = 1;
 
80
 
 
81
 
 
82
if( n == 0 )
 
83
        return( sign * y0(x) );
 
84
if( n == 1 )
 
85
        return( sign * y1(x) );
 
86
 
 
87
/* test for overflow */
 
88
if( x <= 0.0 )
 
89
        {
 
90
        mtherr( "yn", SING );
 
91
        return( -MAXNUM );
 
92
        }
 
93
 
 
94
/* forward recurrence on n */
 
95
 
 
96
anm2 = y0(x);
 
97
anm1 = y1(x);
 
98
k = 1;
 
99
r = 2 * k;
 
100
do
 
101
        {
 
102
        an = r * anm1 / x  -  anm2;
 
103
        anm2 = anm1;
 
104
        anm1 = an;
 
105
        r += 2.0;
 
106
        ++k;
 
107
        }
 
108
while( k < n );
 
109
 
 
110
 
 
111
return( sign * an );
 
112
}
 
113
 
 
114
#endif /* HAVE_YN */
 
115