ESS_GUIESS_GUI_DocsESS_GUI_batch_fittingESS_GUI_bumps_abstractionESS_GUI_iss1116ESS_GUI_iss879ESS_GUI_iss959ESS_GUI_openclESS_GUI_orderingESS_GUI_sync_sascalccostrafo411magnetic_scattrelease-4.1.1release-4.1.2release-4.2.2release_4.0.1ticket-1009ticket-1094-headlessticket-1242-2d-resolutionticket-1243ticket-1249ticket885unittest-saveload
Last change
on this file since 74f419d was
79492222,
checked in by krzywon, 10 years ago
|
Changed the file and folder names to remove all SANS references.
|
-
Property mode set to
100644
|
File size:
1.8 KB
|
Line | |
---|
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.8: June, 2000 |
---|
52 | Copyright 1984, 1987, 2000 by Stephen L. Moshier |
---|
53 | */ |
---|
54 | |
---|
55 | #include "mconf.h" |
---|
56 | #ifdef ANSIPROT |
---|
57 | extern double y0 ( double ); |
---|
58 | extern double y1 ( double ); |
---|
59 | extern double log ( double ); |
---|
60 | #else |
---|
61 | double y0(), y1(), log(); |
---|
62 | #endif |
---|
63 | extern double MAXNUM, MAXLOG; |
---|
64 | |
---|
65 | double yn( n, x ) |
---|
66 | int n; |
---|
67 | double x; |
---|
68 | { |
---|
69 | double an, anm1, anm2, r; |
---|
70 | int k, sign; |
---|
71 | |
---|
72 | if( n < 0 ) |
---|
73 | { |
---|
74 | n = -n; |
---|
75 | if( (n & 1) == 0 ) /* -1**n */ |
---|
76 | sign = 1; |
---|
77 | else |
---|
78 | sign = -1; |
---|
79 | } |
---|
80 | else |
---|
81 | sign = 1; |
---|
82 | |
---|
83 | |
---|
84 | if( n == 0 ) |
---|
85 | return( sign * y0(x) ); |
---|
86 | if( n == 1 ) |
---|
87 | return( sign * y1(x) ); |
---|
88 | |
---|
89 | /* test for overflow */ |
---|
90 | if( x <= 0.0 ) |
---|
91 | { |
---|
92 | mtherr( "yn", SING ); |
---|
93 | return( -MAXNUM ); |
---|
94 | } |
---|
95 | |
---|
96 | /* forward recurrence on n */ |
---|
97 | |
---|
98 | anm2 = y0(x); |
---|
99 | anm1 = y1(x); |
---|
100 | k = 1; |
---|
101 | r = 2 * k; |
---|
102 | do |
---|
103 | { |
---|
104 | an = r * anm1 / x - anm2; |
---|
105 | anm2 = anm1; |
---|
106 | anm1 = an; |
---|
107 | r += 2.0; |
---|
108 | ++k; |
---|
109 | } |
---|
110 | while( k < n ); |
---|
111 | |
---|
112 | |
---|
113 | return( sign * an ); |
---|
114 | } |
---|
Note: See
TracBrowser
for help on using the repository browser.