d01ja Method
関数リスト一覧 　　NagLibrary Namespaceへ　　ライブラリイントロダクション　　本ヘルプドキュメントのchm形式版

d01ja attempts to evaluate an integral over an n-dimensional sphere (n=2, 3, or 4), to a user-specified absolute or relative accuracy, by means of a modified Sag–Szekeres method. The method can handle singularities on the surface or at the centre of the sphere, and returns an error estimate.

Syntax

C#
public static void d01ja(
	D01..::.D01JA_F f,
	int ndim,
	double radius,
	double epsa,
	double epsr,
	int method,
	int icoord,
	out double result,
	out double esterr,
	out int nevals,
	out int ifail
)
Visual Basic (Declaration)
Public Shared Sub d01ja ( _
	f As D01..::.D01JA_F, _
	ndim As Integer, _
	radius As Double, _
	epsa As Double, _
	epsr As Double, _
	method As Integer, _
	icoord As Integer, _
	<OutAttribute> ByRef result As Double, _
	<OutAttribute> ByRef esterr As Double, _
	<OutAttribute> ByRef nevals As Integer, _
	<OutAttribute> ByRef ifail As Integer _
)
Visual C++
public:
static void d01ja(
	D01..::.D01JA_F^ f, 
	int ndim, 
	double radius, 
	double epsa, 
	double epsr, 
	int method, 
	int icoord, 
	[OutAttribute] double% result, 
	[OutAttribute] double% esterr, 
	[OutAttribute] int% nevals, 
	[OutAttribute] int% ifail
)
F#
static member d01ja : 
        f:D01..::.D01JA_F * 
        ndim:int * 
        radius:float * 
        epsa:float * 
        epsr:float * 
        method:int * 
        icoord:int * 
        result:float byref * 
        esterr:float byref * 
        nevals:int byref * 
        ifail:int byref -> unit

Parameters

f
Type: NagLibrary..::.D01..::.D01JA_F
f must return the value of the integrand f at a given point.

A delegate of type D01JA_F.

See also [Further Comments].
ndim
Type: System..::.Int32
On entry: n, the dimension of the sphere.
Constraint: 2≤ndim≤4.
radius
Type: System..::.Double
On entry: α, the radius of the sphere.
Constraint: radius≥0.0.
epsa
Type: System..::.Double
On entry: the requested absolute tolerance. If epsa<0.0, its absolute value is used. See [Accuracy].
epsr
Type: System..::.Double
On entry: the requested relative tolerance.
epsr<0.0
Its absolute value is used.
epsr<10×machine precision
The latter value is used as epsr by the method. See [Accuracy].
method
Type: System..::.Int32
On entry: must specify the transformation to be used by the method. The choice depends on the behaviour of the integrand and on the required accuracy.
For well-behaved functions and functions with mild singularities on the surface of the sphere only:
method=1
Low accuracy required.
method=2
High accuracy required.
For functions with severe singularities on the surface of the sphere only:
method=3
Low accuracy required.
method=4
High accuracy required.
(in this case icoord must be set to 2, and the function defined in special spherical coordinates).
For functions with a singularity at the centre of the sphere (and possibly with singularities on the surface as well):
method=5
Low accuracy required.
method=6
High accuracy required.
method=0 can be used as a default value and is equivalent to method=1 if epsr>10-6, and to method=2 if epsr≤10-6.
The distinction between low and high required accuracies, as mentioned above, depends also on the behaviour of the function. Roughly one may assume the critical value of epsa and epsr to be 10-6, but the critical value will be smaller for a well-behaved integrand and larger for an integrand with severe singularities.
Suggested value: method=0.
Constraint: method=0, 1, 2, 3, 4, 5 or 6.
If icoord=2, method=3 or 4
icoord
Type: System..::.Int32
On entry: must specify which kind of coordinates are used in f.
icoord=0
Cartesian coordinates xi, for i=1,2,…,n.
icoord=1
Spherical coordinates (see [Spherical Polar Co-ordinates]): x[0]=ρ; x[i-1]=θi-1, for i=2,3,…,n.
icoord=2,
Special spherical polar coordinates (see [Machine Dependencies]), with the additional transformation ρ=α-λ: x[0]=λ=α-ρ; x[i-1]=θi-1, for i=2,3,…,n.
Constraint: icoord=0, 1 or 2.
If method=3 or 4, icoord=2
result
Type: System..::.Double %
On exit: the approximation to the integral I.
esterr
Type: System..::.Double %
On exit: an estimate of the modulus of the absolute error.
nevals
Type: System..::.Int32 %
On exit: the number of function evaluations used.
ifail
Type: System..::.Int32 %
On exit: ifail=0 unless the method detects an error (see [Error Indicators and Warnings]).

Description

d01ja calculates an approximation to the n-dimensional integral
I=∫⋯∫SFx1,…,xndx1⋯dxn,  2≤n≤4,
where S is the hypersphere
x12+⋯+xn2≤α<∞
(the integrand function may also be defined in spherical coordinates). The algorithm is based on the Sag–Szekeres method (see Sag and Szekeres (1964)), applying the product trapezoidal formula after a suitable radial transformation. An improved transformation technique is developed: depending on the behaviour of the function and on the required accuracy, different transformations can be used, some of which are ‘double exponential’, as defined by Takahasi and Mori (1974). The resulting technique allows the method to deal with integrand singularities on the surface or at the centre of the sphere. When the estimated error of the approximation with mesh size h is larger than the tolerated error, the trapezoidal formula with mesh size h/2 is calculated. A drawback of this method is the exponential growth of the number of function evaluations in the successive approximations (this number grows with a factor ≈2n). This introduces the restriction n≤4. Because the convergence rate of the successive approximations is normally better than linear, the error estimate is based on the linear extrapolation of the difference between the successive approximations (see Robinson and de Doncker (1981) and Roose and de Doncker (1981)). For further details of the algorithm, see Roose and de Doncker (1981).

References

Error Indicators and Warnings

Accuracy

You can specify an absolute and/or a relative tolerance, setting epsa and epsr. The method attempts to calculate an approximation result such that
I-result≤maxepsa,epsr×I.
If 0≤ifail≤3, esterr returns an estimate of, but not necessarily a bound for, I-result.

Further Comments

Timing

Spherical Polar Co-ordinates

Machine Dependencies

As a consequence of the transformation technique, the severity of the singularities which can be handled by d01ja depends on the precision and range of real numbers on the machine. method=3 or 4 must be used when the singularity on the surface is ‘severe’ in view of the requested accuracy and machine precision. In practice one has to set method=3 or 4 if d01ja terminates with ifail=3 when called with method=0, 1 or 2.
When integrating a function with a severe singular behaviour on the surface of the sphere, the additional transformation ρ=α-λ helps to avoid the loss of significant figures due to round-off error in the calculation of the integration nodes which are very close to the surface. For these points, the value of λ can be computed more accurately than the value of ρ. Naturally, care must be taken that the function subprogram does not contain expressions of the form α-λ, which could cause a large round-off error in the calculation of the integrand at the boundary of the sphere.
Care should be taken to avoid underflow and/or overflow problems in the function subprogram, because some of the integration nodes used by d01ja may be very close to the surface or to the centre of the sphere.
Example:
  • suppose the function
    fρ=1-ρ2-0.7
    is to be integrated over the unit sphere, with method=3 or 4. Then icoord should be set to 2; the transformation ρ=1-λ gives fρ=2λ-λ2 -0.7; and f could be coded thus:
Note that d01ja ensures that λ=x[0]>x02am, but underflow could occur in the computation of λ2.

Example

See Also