hide long namesshow long names
hide short namesshow short names
Integer type:  int32  int64  nag_int  show int32  show int32  show int64  show int64  show nag_int  show nag_int

PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

NAG Toolbox: nag_specfun_1f1_real_scaled (s22bb)

Purpose

nag_specfun_1f1_real_scaled (s22bb) returns a value for the confluent hypergeometric function 1F1 (a ; b ; x) F 1 1 (a;b;x)  with real parameters aa, bb and real argument xx in the scaled form 1F1 (a ; b ; x) = mf × 2ms F 1 1 (a;b;x) = mf × 2ms . This function is sometimes also known as Kummer's function M(a,b,x)M(a,b,x).

Syntax

[frm, scm, ifail] = s22bb(ani, adr, bni, bdr, x)
[frm, scm, ifail] = nag_specfun_1f1_real_scaled(ani, adr, bni, bdr, x)

Description

nag_specfun_1f1_real_scaled (s22bb) returns a value for the confluent hypergeometric function 1F1 (a ; b ; x) F1 1 (a;b;x)  with real parameters aa, bb and xx in the scaled form 1F1 (a ; b ; x) = mf × 2ms F1 1 (a;b;x) = mf × 2 ms , where mfmf is the real scaled component and msms is the integer power of two scaling. This function is unbounded or not uniquely defined for bb equal to zero or a negative integer.
The confluent hypergeometric function is defined by the confluent series,
1F1(a ; b ; x) = M(a,b,x) = ( (a)s xs )/( (b)s s! ) = 1 + (a)/b x + ( a(a + 1) )/( b(b + 1) 2! )x2 +
s = 0
F1 1 (a;b;x) = M(a,b,x) = s=0 (a)s xs (b)s s! = 1 + a b x + a(a+1) b(b+1) 2! x2 +
where (a)s = 1 (a) (a + 1) (a + 2) (a + s1) (a)s = 1 (a) (a+1) (a+2) (a+s-1)  is the rising factorial of aa. M(a,b,x) M(a,b,x)  is a solution to the second order ODE (Kummer's Equation):
x (d2M)/(dx2) + (bx) (dM)/(dx) a M = 0 .
x d2M dx2 + (b-x) dM dx - a M = 0 .
(1)
Given the parameters and argument (a,b,x) (a,b,x) , this function determines a set of safe values {(αi,βi,ζi)i2} { (αi,βi,ζi) i2 }  and selects an appropriate algorithm to accurately evaluate the functions Mi (αi,βi,ζi) Mi (αi,βi,ζi) . The result is then used to construct the solution to the original problem M(a,b,x) M(a,b,x)  using, where necessary, recurrence relations and/or continuation.
For improved precision in the final result, this function accepts aa and bb split into an integral and a decimal fractional component. Specifically a = ai + ara=ai+ar, where |ar|0.5|ar|0.5 and ai = aarai=a-ar is integral. bb is similarly deconstructed.
Additionally, an artificial bound, arbndarbnd is placed on the magnitudes of aiai, bibi and xx to minimize the occurrence of overflow in internal calculations. arbnd = 0.0001 × Imax arbnd = 0.0001 × Imax , where Imax = x02bbImax=x02bb. It should, however, not be assumed that this function will produce an accurate result for all values of aiai, bibi and xx satisfying this criterion.
Please consult the NIST Digital Library of Mathematical Functions or the companion (2010) for a detailed discussion of the confluent hypergeometric function including special cases, transformations, relations and asymptotic approximations.

References

NIST Handbook of Mathematical Functions (2010) (eds F W J Olver, D W Lozier, R F Boisvert, C W Clark) Cambridge University Press
Pearson J (2009) Computation of hypergeometric functions MSc Dissertation, Mathematical Institute, University of Oxford

Parameters

Compulsory Input Parameters

1:     ani – double scalar
aiai, the nearest integer to aa, satisfying ai = aar ai = a-ar.
Constraints:
  • ani = aniani=ani;
  • |ani|arbnd|ani|arbnd.
2:     adr – double scalar
arar, the signed decimal remainder satisfying ar = aai ar = a-ai  and |ar| 0.5 |ar| 0.5.
Constraint: |adr|0.5|adr|0.5.
Note: if |adr| < 100.0ε|adr|<100.0ε, ar = 0.0ar=0.0 will be used, where εε is the machine precision as returned by nag_machine_precision (x02aj).
3:     bni – double scalar
bibi, the nearest integer to bb, satisfying bi = bbrbi=b-br.
Constraints:
  • bni = bnibni=bni;
  • |bni|arbnd|bni|arbnd;
  • if bdr = 0.0bdr=0.0, bni > 0bni>0.
4:     bdr – double scalar
brbr, the signed decimal remainder satisfying br = bbi br = b-bi and |br| 0.5 |br| 0.5.
Constraint: |bdr|0.5|bdr|0.5.
Note: if |bdradr| < 100.0ε|bdr-adr|<100.0ε, ar = brar=br will be used, where εε is the machine precision as returned by nag_machine_precision (x02aj).
5:     x – double scalar
The argument xx of the function.
Constraint: |x|arbnd|x|arbnd.

Optional Input Parameters

None.

Input Parameters Omitted from the MATLAB Interface

None.

Output Parameters

1:     frm – double scalar
mfmf, the scaled real component of the solution satisfying mf = M(a,b,x) × 2msmf=M(a,b,x)×2-ms.
Note: if overflow occurs upon completion, as indicated by ifail = 2ifail=2, the value of mfmf returned may still be correct. If overflow occurs in a subcalculation, as indicated by ifail = 5ifail=5, this should not be assumed.
2:     scm – int64int32nag_int scalar
msms, the scaling power of two, satisfying ms = log2((M(a,b,x))/(mf))ms= log2( M(a,b,x) mf ).
Note: if overflow occurs upon completion, as indicated by ifail = 2ifail=2, then msImaxmsImax, where ImaxImax is the largest representable integer (see nag_machine_integer_max (x02bb)). If overflow occurs during a subcalculation, as indicated by ifail = 5ifail=5, msms may or may not be greater than ImaxImax. In either case, scm = x02bbscm=x02bb will have been returned.
3:     ifail – int64int32nag_int scalar
ifail = 0ifail=0 unless the function detects an error (see [Error Indicators and Warnings]).

Error Indicators and Warnings

Errors or warnings detected by the function:

Cases prefixed with W are classified as warnings and do not generate an error of type NAG:error_n. See nag_issue_warnings.

W ifail = 1ifail=1
Underflow occurred during the evaluation of M(a,b,x)M(a,b,x).
The returned value may be inaccurate.
W ifail = 2ifail=2
On completion, overflow occurred in the evaluation of M(a,b,x)M(a,b,x).
W ifail = 3ifail=3
All approximations have completed, and the final residual estimate indicates some precision may have been lost.
  ifail = 4ifail=4
All approximations have completed, and the final residual estimate indicates no accuracy can be guaranteed.
  ifail = 5ifail=5
Overflow occurred in a subcalculation of M(a,b,x)M(a,b,x).
The answer may be completely incorrect.
  ifail = 11ifail=11
Constraint: .
  ifail = 13ifail=13
Constraint: ani = aniani=ani.
  ifail = 21ifail=21
Constraint: |adr|0.5|adr|0.5.
  ifail = 31ifail=31
Constraint: .
  ifail = 32ifail=32
On entry.
M(a,b,x)M(a,b,x) is undefined when bb is zero or a negative integer.
  ifail = 33ifail=33
Constraint: bni = bnibni=bni.
  ifail = 41ifail=41
Constraint: |bdr|0.5|bdr|0.5.
  ifail = 51ifail=51
Constraint: .

Accuracy

In general, if ifail = 0ifail=0, the value of MM may be assumed accurate, with the possible loss of one or two decimal places. Assuming the result does not under or overflow, an error estimate resres is made internally using equation (1). If the magnitude of resres is sufficiently large a nonzero ifail will be returned. Specifically,
ifail = 0ifail=0 res1000εres1000ε
ifail = 3ifail=3 1000 ε < res0.11000ε<res0.1
ifail = 4ifail=4 res > 0.1res>0.1
A further estimate of the residual can be constructed using equation (1), and the differential identity,
( d M(a,b,x) )/(dx) = a/b M (a + 1,b + 1,x) ,
( d2 M(a,b,x) )/(dx2) = (a(a + 1))/(b(b + 1)) M (a + 2,b + 2,x) .
d M(a,b,x) dx = ab M (a+1,b+1,x) , d2 M(a,b,x) dx2 = a(a+1) b(b+1) M (a+2,b+2,x) .
This estimate is however dependent upon the error involved in approximating M (a + 1,b + 1,x) M (a+1,b+1,x)  and M (a + 2,b + 2,x) M (a+2,b+2,x) .

Further Comments

The values returned in frm (mfmf) and scm (msms) may be used to explicitly evaluate M(a,b,x)M(a,b,x), and may also be used to evaluate products and ratios of multiple values of MM as follows,
M(a,b,x) = mf × 2ms
M (a1,b1,x1) × M (a2,b2,x2) = (mf1 × mf2) × 2(ms1 + ms2)
( M (a1,b1,x1) )/( M (a2,b2,x2) ) = (mf1)/(mf2) × 2(ms1ms2)
ln|M(a,b,x)| = ln|mf| + ms × ln(2)
.
M(a,b,x) = mf × 2ms M (a1,b1,x1) × M (a2,b2,x2) = ( mf1 × mf2 ) × 2 ( ms1 + ms2 ) M (a1,b1,x1) M (a2,b2,x2) = mf1 mf2 × 2 ( ms1 - ms2 ) ln| M (a,b,x) | = ln|mf| + ms × ln(2) .

Example

function nag_specfun_1f1_real_scaled_example
ai = -10;
bi = 30;
delta = 1e-4;
ar = delta;
br = -delta;
x = 25;
frmv = zeros(2, 1);
scmv = zeros(2, 1, 'int64');
fprintf('\n         a          b          x          frm    scm     M(a,b,x)\n');
[frmv(1), scmv(1), ifail] = nag_specfun_1f1_real_scaled(ai,  ar, bi,  br, x);
print_result(ai+ar, bi+br, x, frmv(1), scmv(1), ifail);
[frmv(2), scmv(2), ifail] = nag_specfun_1f1_real_scaled(ai, -ar, bi, -br, x);
print_result(ai-ar, bi-br, x, frmv(2), scmv(2), ifail);

% Calculate the product M1*M2
frm = prod(frmv);
scm = scmv(1)+scmv(2);

if scm < nag_machine_model_maxexp
  scale = frm*double(nag_machine_model_base)^double(scm);
  fprintf('\nSolution product %12.4e %6d %12.4e\n', frm, scm, scale);
else
  fprintf('\nSolution product %12.4e %6d Not representable\n', frm, scm);
end

% Calculate the ratio M1/M2
if frmv(2) ~= 0
  frm = frmv(1)/frmv(2);
  scm = scmv(1) - scmv(2);
  if scm < nag_machine_model_maxexp
    scale = frm*double(nag_machine_model_base)^double(scm);
    fprintf('\nSolution ratio   %12.4e %6d %12.4e\n', frm, scm, scale);
  else
    fprintf('\nSolution ratio   %12.4e %6d Not representable\n', frm, scm);
  end
end

function print_result(a, b, x, frm, scm, ifail)
  if ifail ==2 || ifail > 3
    % Either the result has overflowed, no accuracy may be assumed,
    % or an input error has been detected.
    fprintf('%10.4f %10.4f %10.4f                      FAILED\n', a, b, x);
  elseif scm < nag_machine_model_maxexp
    scale = frm*double(nag_machine_model_base)^double(scm);
    fprintf('%10.4f %10.4f %10.4f %12.4e %6d %12.4e\n', a, b, x, frm, scm, scale);
  else
    fprintf('%10.4f %10.4f %10.4f %12.4e %6d Not representable\n', a, b, x, frm, scm);
  end
 

         a          b          x          frm    scm     M(a,b,x)
   -9.9999    29.9999    25.0000  -7.7329e-01    -15  -2.3599e-05
  -10.0001    30.0001    25.0000  -7.7318e-01    -15  -2.3596e-05

Solution product   5.9789e-01    -30   5.5683e-10

Solution ratio     1.0001e+00      0   1.0001e+00

function s22bb_example
ai = -10;
bi = 30;
delta = 1e-4;
ar = delta;
br = -delta;
x = 25;
frmv = zeros(2, 1);
scmv = zeros(2, 1, 'int64');
fprintf('\n         a          b          x          frm    scm     M(a,b,x)\n');
[frmv(1), scmv(1), ifail] = s22bb(ai,  ar, bi,  br, x);
print_result(ai+ar, bi+br, x, frmv(1), scmv(1), ifail);
[frmv(2), scmv(2), ifail] = s22bb(ai, -ar, bi, -br, x);
print_result(ai-ar, bi-br, x, frmv(2), scmv(2), ifail);

% Calculate the product M1*M2
frm = prod(frmv);
scm = scmv(1)+scmv(2);

if scm < x02bl
  scale = frm*double(x02bh)^double(scm);
  fprintf('\nSolution product %12.4e %6d %12.4e\n', frm, scm, scale);
else
  fprintf('\nSolution product %12.4e %6d Not representable\n', frm, scm);
end

% Calculate the ratio M1/M2
if frmv(2) ~= 0
  frm = frmv(1)/frmv(2);
  scm = scmv(1) - scmv(2);
  if scm < x02bl
    scale = frm*double(x02bh)^double(scm);
    fprintf('\nSolution ratio   %12.4e %6d %12.4e\n', frm, scm, scale);
  else
    fprintf('\nSolution ratio   %12.4e %6d Not representable\n', frm, scm);
  end
end

function print_result(a, b, x, frm, scm, ifail)
  if ifail ==2 || ifail > 3
    % Either the result has overflowed, no accuracy may be assumed,
    % or an input error has been detected.
    fprintf('%10.4f %10.4f %10.4f                      FAILED\n', a, b, x);
  elseif scm < x02bl
    scale = frm*double(x02bh)^double(scm);
    fprintf('%10.4f %10.4f %10.4f %12.4e %6d %12.4e\n', a, b, x, frm, scm, scale);
  else
    fprintf('%10.4f %10.4f %10.4f %12.4e %6d Not representable\n', a, b, x, frm, scm);
  end
 

         a          b          x          frm    scm     M(a,b,x)
   -9.9999    29.9999    25.0000  -7.7329e-01    -15  -2.3599e-05
  -10.0001    30.0001    25.0000  -7.7318e-01    -15  -2.3596e-05

Solution product   5.9789e-01    -30   5.5683e-10

Solution ratio     1.0001e+00      0   1.0001e+00


PDF version (NAG web site, 64-bit version, 64-bit version)
Chapter Contents
Chapter Introduction
NAG Toolbox

© The Numerical Algorithms Group Ltd, Oxford, UK. 2009–2013