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_pde_1d_blackscholes_fd (d03nc)

Purpose

nag_pde_1d_blackscholes_fd (d03nc) solves the Black–Scholes equation for financial option pricing using a finite difference scheme.

Syntax

[s, t, f, theta, delta, gamma, lambda, rho, ifail] = d03nc(kopt, x, mesh, s, t, tdpar, r, q, sigma, ntkeep, 'ns', ns, 'nt', nt, 'alpha', alpha)
[s, t, f, theta, delta, gamma, lambda, rho, ifail] = nag_pde_1d_blackscholes_fd(kopt, x, mesh, s, t, tdpar, r, q, sigma, ntkeep, 'ns', ns, 'nt', nt, 'alpha', alpha)

Description

nag_pde_1d_blackscholes_fd (d03nc) solves the Black–Scholes equation (see Hull (1989) and Wilmott et al. (1995))
(f)/(t) + (rq)S(f)/(S) + (σ2S2)/2 (2f)/(S2) = rf
f t +(r-q)S f S +σ2S22 2f S2 =rf
(1)
Smin < S < Smax,  tmin < t < tmax,
Smin<S<Smax,  tmin<t<tmax,
(2)
for the value ff of a European or American, put or call stock option, with exercise price XX. In equation (1) tt is time, SS is the stock price, rr is the risk free interest rate, qq is the continuous dividend, and σσ is the stock volatility. According to the values in the array tdpar, the parameters rr, qq and σσ may each be either constant or functions of time. The function also returns values of various Greeks.
nag_pde_1d_blackscholes_fd (d03nc) uses a finite difference method with a choice of time-stepping schemes. The method is explicit for alpha = 0.0alpha=0.0 and implicit for nonzero values of alpha. Second order time accuracy can be obtained by setting alpha = 0.5alpha=0.5. According to the value of the parameter mesh the finite difference mesh may be either uniform, or user-defined in both SS and tt directions.

References

Hull J (1989) Options, Futures and Other Derivative Securities Prentice–Hall
Wilmott P, Howison S and Dewynne J (1995) The Mathematics of Financial Derivatives Cambridge University Press

Parameters

Compulsory Input Parameters

1:     kopt – int64int32nag_int scalar
Specifies the kind of option to be valued.
kopt = 1kopt=1
A European call option.
kopt = 2kopt=2
An American call option.
kopt = 3kopt=3
A European put option.
kopt = 4kopt=4
An American put option.
Constraint: kopt = 1kopt=1, 22, 33 or 44.
2:     x – double scalar
The exercise price XX.
3:     mesh – string (length ≥ 1)
Indicates the type of finite difference mesh to be used:
mesh = 'U'mesh='U'
Uniform mesh.
mesh = 'C'mesh='C'
Custom mesh supplied by you.
Constraint: mesh = 'U'mesh='U' or 'C''C'.
4:     s(ns) – double array
ns, the dimension of the array, must satisfy the constraint ns2ns2.
If mesh = 'C'mesh='C', s(i)si must contain the iith stock price in the mesh, for i = 1,2,,nsi=1,2,,ns. These values should be in increasing order, with s(1) = Smins1=Smin and s(ns) = Smaxsns=Smax.
If mesh = 'U'mesh='U', s(1)s1 must be set to SminSmin and s(ns)sns to SmaxSmax, but s(2),s(3),,s(ns1)s2,s3,,sns-1 need not be initialized, as they will be set internally by the function in order to define a uniform mesh.
Constraints:
  • if mesh = 'C'mesh='C', s(1)0.0s10.0 and s(i) < s(i + 1)si<si+1, for i = 1,2,,ns1i=1,2,,ns-1;
  • if mesh = 'U'mesh='U', 0.0s(1) < s(ns)0.0s1<sns.
5:     t(nt) – double array
nt, the dimension of the array, must satisfy the constraint nt2nt2.
If mesh = 'C'mesh='C' then t(j)tj must contain the jjth time in the mesh, for j = 1,2,,ntj=1,2,,nt. These values should be in increasing order, with t(1) = tmint1=tmin and t(nt) = tmaxtnt=tmax.
If mesh = 'U'mesh='U' then t(1)t1 must be set to tmintmin and t(nt)tnt to tmaxtmax, but t(2),t(3),,t(nt1)t2,t3,,tnt-1 need not be initialized, as they will be set internally by the function in order to define a uniform mesh.
Constraints:
  • if mesh = 'C'mesh='C', t(1)0.0t10.0 and t(j) < t(j + 1)tj<tj+1, for j = 1,2,,nt1j=1,2,,nt-1;
  • if mesh = 'U'mesh='U', 0.0t(1) < t(nt)0.0t1<tnt.
6:     tdpar(33) – logical array
Specifies whether or not various parameters are time-dependent. More precisely, rr is time-dependent if tdpar(1) = truetdpar1=true and constant otherwise. Similarly, tdpar(2)tdpar2 specifies whether qq is time-dependent and tdpar(3)tdpar3 specifies whether σσ is time-dependent.
7:     r( : :) – double array
Note: the dimension of the array r must be at least ntnt if tdpar(1) = truetdpar1=true, and at least 11 otherwise.
If tdpar(1) = truetdpar1=true then r(j)rj must contain the value of the risk-free interest rate r(t)r(t) at the jjth time in the mesh, for j = 1,2,,ntj=1,2,,nt.
If tdpar(1) = falsetdpar1=false then r(1)r1 must contain the constant value of the risk-free interest rate rr. The remaining elements need not be set.
8:     q( : :) – double array
Note: the dimension of the array q must be at least ntnt if tdpar(2) = truetdpar2=true, and at least 11 otherwise.
If tdpar(2) = truetdpar2=true then q(j)qj must contain the value of the continuous dividend q(t)q(t) at the jjth time in the mesh, for j = 1,2,,ntj=1,2,,nt.
If tdpar(2) = falsetdpar2=false then q(1)q1 must contain the constant value of the continuous dividend qq. The remaining elements need not be set.
9:     sigma( : :) – double array
Note: the dimension of the array sigma must be at least ntnt if tdpar(3) = truetdpar3=true, and at least 11 otherwise.
If tdpar(3) = truetdpar3=true then sigma(j)sigmaj must contain the value of the volatility σ(t)σ(t) at the jjth time in the mesh, for j = 1,2,,ntj=1,2,,nt.
If tdpar(3) = falsetdpar3=false then sigma(1)sigma1 must contain the constant value of the volatility σσ. The remaining elements need not be set.
10:   ntkeep – int64int32nag_int scalar
The number of solutions to be stored in the time direction. The function calculates the solution backwards from t(nt)tnt to t(1)t1 at all times in the mesh. These time solutions and the corresponding Greeks will be stored at times t(i)ti, for i = 1,2,,ntkeepi=1,2,,ntkeep, in the arrays f, theta, delta, gamma, lambda and rho. Other time solutions will be discarded. To store all time solutions set ntkeep = ntntkeep=nt.
Constraint: 1ntkeepnt1ntkeepnt.

Optional Input Parameters

1:     ns – int64int32nag_int scalar
Default: The dimension of the array s.
The number of stock prices to be used in the finite difference mesh.
Constraint: ns2ns2.
2:     nt – int64int32nag_int scalar
Default: The dimension of the array t.
The number of time-steps to be used in the finite difference method.
Constraint: nt2nt2.
3:     alpha – double scalar
The value of λλ to be used in the time-stepping scheme. Typical values include:
alpha = 0.0alpha=0.0
Explicit forward Euler scheme.
alpha = 0.5alpha=0.5
Implicit Crank–Nicolson scheme.
alpha = 1.0alpha=1.0
Implicit backward Euler scheme.
The value 0.50.5 gives second-order accuracy in time. Values greater than 0.50.5 give unconditional stability. Since 0.50.5 is at the limit of unconditional stability this value does not damp oscillations.
Default: 0.550.55
Constraint: 0.0alpha1.00.0alpha1.0.

Input Parameters Omitted from the MATLAB Interface

ldf work iwork

Output Parameters

1:     s(ns) – double array
If mesh = 'U'mesh='U', the elements of s define a uniform mesh over [Smin,Smax][Smin,Smax].
If mesh = 'C'mesh='C', the elements of s are unchanged.
2:     t(nt) – double array
If mesh = 'U'mesh='U', the elements of t define a uniform mesh over [tmin,tmax][tmin,tmax].
If mesh = 'C'mesh='C', the elements of t are unchanged.
3:     f(ldf,ntkeep) – double array
ldf = ns ldf=ns .
f(i,j)fij, for i = 1,2,,nsi=1,2,,ns and j = 1,2,,ntkeepj=1,2,,ntkeep, contains the value ff of the option at the iith mesh point s(i)si at time t(j)tj.
4:     theta(ldf,ntkeep) – double array
5:     delta(ldf,ntkeep) – double array
6:     gamma(ldf,ntkeep) – double array
7:     lambda(ldf,ntkeep) – double array
8:     rho(ldf,ntkeep) – double array
ldf = ns ldf=ns .
The values of various Greeks at the iith mesh point s(i)si at time t(j)tj, as follows:
theta(i,j) = (f)/(t), delta(i,j) = (f)/(S), gamma(i,j) = (2f)/(S2),
lambda(i,j) = (f)/(σ), rho(i,j) = (f)/(r).
thetaij= f t , deltaij= f S , gammaij= 2f S2 , lambdaij= f σ , rhoij= f r .
9:     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:
  ifail = 1ifail=1
On entry,kopt < 1kopt<1,
orkopt > 4kopt>4,
ormesh'U'mesh'U' or 'C''C',
orns < 2ns<2,
ornt < 2nt<2,
ors(1) < 0.0s1<0.0,
ort(1) < 0.0t1<0.0,
oralpha < 0.0alpha<0.0,
oralpha > 1.0alpha>1.0,
orntkeep < 1ntkeep<1,
orntkeep > ntntkeep>nt,
orldf < nsldf<ns.
  ifail = 2ifail=2
mesh = 'U'mesh='U' and the constraints: are violated. Thus the end points of the uniform mesh are not in order.
  ifail = 3ifail=3
mesh = 'C'mesh='C' and the constraints:
  • s(i) < s(i + 1)si<si+1, for i = 1,2,,ns1i=1,2,,ns-1,
  • t(i) < t(i + 1)ti<ti+1, for i = 1,2,,nt1i=1,2,,nt-1 
are violated. Thus the mesh points are not in order.

Accuracy

The accuracy of the solution ff and the various derivatives returned by the function is dependent on the values of ns and nt supplied, the distribution of the mesh points, and the value of alpha chosen. For most choices of alpha the solution has a truncation error which is second-order accurate in SS and first order accurate in tt. For alpha = 0.5alpha=0.5 the truncation error is also second-order accurate in tt.
The simplest approach to improving the accuracy is to increase the values of both ns and nt.

Further Comments

Timing

Each time-step requires the construction and solution of a tridiagonal system of linear equations. To calculate each of the derivatives lambda and rho requires a repetition of the entire solution process. The time taken for a call to the function is therefore proportional to ns × ntns×nt.

Algorithmic Details

nag_pde_1d_blackscholes_fd (d03nc) solves equation (1) using a finite difference method. The solution is computed backwards in time from tmaxtmax to tmintmin using a λλ scheme, which is implicit for all nonzero values of λλ, and is unconditionally stable for values of λ > 0.5λ>0.5. For each time-step a tridiagonal system is constructed and solved to obtain the solution at the earlier time. For the explicit scheme (λ = 0λ=0) this tridiagonal system degenerates to a diagonal matrix and is solved trivially. For American options the solution at each time-step is inspected to check whether early exercise is beneficial, and amended accordingly.
To compute the arrays lambda and rho, which are derivatives of the stock value ff with respect to the problem parameters σσ and rr respectively, the entire solution process is repeated with perturbed values of these parameters.

Example

function nag_pde_1d_blackscholes_fd_example
kopt = int64(4);
x = 50;
mesh = 'U';
s = [0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     100];
t = [0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0.4166667];
tdpar = [false; false; false];
r = [0.1];
q = [0];
sigma = [0.4];
ntkeep = int64(4);
[sOut, tOut, f, theta, delta, gamma, lambda, rho, ifail] = ...
    nag_pde_1d_blackscholes_fd(kopt, x, mesh, s, t, tdpar, r, q, sigma, ntkeep,'alpha',1)
 

sOut =

     0
     5
    10
    15
    20
    25
    30
    35
    40
    45
    50
    55
    60
    65
    70
    75
    80
    85
    90
    95
   100


tOut =

         0
    0.0417
    0.0833
    0.1250
    0.1667
    0.2083
    0.2500
    0.2917
    0.3333
    0.3750
    0.4167


f =

   50.0000   50.0000   50.0000   50.0000
   45.0000   45.0000   45.0000   45.0000
   40.0000   40.0000   40.0000   40.0000
   35.0000   35.0000   35.0000   35.0000
   30.0000   30.0000   30.0000   30.0000
   25.0000   25.0000   25.0000   25.0000
   20.0000   20.0000   20.0000   20.0000
   15.0000   15.0000   15.0000   15.0000
   10.1543   10.0958   10.0464   10.0117
    6.5848    6.4424    6.2916    6.1306
    4.0672    3.8785    3.6729    3.4463
    2.4264    2.2423    2.0454    1.8336
    1.4174    1.2662    1.1096    0.9481
    0.8195    0.7072    0.5953    0.4851
    0.4724    0.3941    0.3190    0.2484
    0.2726    0.2202    0.1717    0.1282
    0.1573    0.1233    0.0929    0.0667
    0.0897    0.0685    0.0501    0.0347
    0.0484    0.0363    0.0259    0.0175
    0.0211    0.0156    0.0110    0.0073
         0         0         0         0


theta =

         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
   -1.4043   -1.1857   -0.8328   -0.2806
   -3.4185   -3.6183   -3.8646   -4.1880
   -4.5285   -4.9339   -5.4387   -6.0796
   -4.4165   -4.7277   -5.0821   -5.4821
   -3.6294   -3.7585   -3.8748   -3.9632
   -2.6946   -2.6860   -2.6441   -2.5561
   -1.8790   -1.8018   -1.6941   -1.5505
   -1.2578   -1.1621   -1.0461   -0.9097
   -0.8154   -0.7282   -0.6301   -0.5231
   -0.5084   -0.4411   -0.3689   -0.2943
   -0.2928   -0.2484   -0.2024   -0.1566
   -0.1324   -0.1108   -0.0888   -0.0674
         0         0         0         0


delta =

   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -0.9846   -0.9904   -0.9954   -0.9988
   -0.8415   -0.8558   -0.8708   -0.8869
   -0.6087   -0.6217   -0.6373   -0.6565
   -0.4158   -0.4200   -0.4246   -0.4297
   -0.2650   -0.2612   -0.2563   -0.2498
   -0.1607   -0.1535   -0.1450   -0.1348
   -0.0945   -0.0872   -0.0791   -0.0700
   -0.0547   -0.0487   -0.0424   -0.0357
   -0.0315   -0.0271   -0.0226   -0.0182
   -0.0183   -0.0152   -0.0122   -0.0093
   -0.0109   -0.0087   -0.0067   -0.0049
   -0.0069   -0.0053   -0.0039   -0.0027
   -0.0048   -0.0036   -0.0026   -0.0017
   -0.0042   -0.0031   -0.0022   -0.0015


gamma =

         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
    0.0062    0.0038    0.0019    0.0005
    0.0510    0.0500    0.0480    0.0443
    0.0421    0.0436    0.0454    0.0479
    0.0351    0.0371    0.0396    0.0429
    0.0253    0.0264    0.0277    0.0291
    0.0164    0.0167    0.0169    0.0169
    0.0100    0.0098    0.0095    0.0091
    0.0059    0.0056    0.0052    0.0047
    0.0034    0.0031    0.0027    0.0024
    0.0019    0.0017    0.0014    0.0012
    0.0011    0.0009    0.0007    0.0006
    0.0006    0.0005    0.0004    0.0003
    0.0003    0.0002    0.0002    0.0001
         0         0         0         0


lambda =

         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
    6.3243    5.1893    3.8089    2.1118
   10.7215    9.9718    9.2140    8.4953
   12.3808   11.8073   11.2277   10.6365
   11.4834   10.8366   10.1417    9.3795
    9.3227    8.5840    7.7870    6.9211
    6.9621    6.2206    5.4412    4.6264
    4.9268    4.2651    3.5937    2.9227
    3.3602    2.8204    2.2920    1.7866
    2.2221    1.8126    1.4248    1.0683
    1.4122    1.1240    0.8586    0.6225
    0.8269    0.6459    0.4825    0.3408
    0.3789    0.2925    0.2155    0.1498
         0         0         0         0


rho =

         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
   -7.1918   -6.0114   -4.5204   -2.5855
   -8.4541   -7.6378   -6.8479   -6.1657
   -7.5988   -6.9323   -6.2879   -5.6707
   -5.8905   -5.2837   -4.6809   -4.0772
   -4.1854   -3.6547   -3.1306   -2.6135
   -2.8221   -2.3904   -1.9743   -1.5775
   -1.8437   -1.5137   -1.2055   -0.9228
   -1.1812   -0.9407   -0.7233   -0.5316
   -0.7451   -0.5768   -0.4292   -0.3038
   -0.4591   -0.3466   -0.2506   -0.1716
   -0.2655   -0.1966   -0.1389   -0.0927
   -0.1228   -0.0898   -0.0626   -0.0410
         0         0         0         0


ifail =

                    0


function d03nc_example
kopt = int64(4);
x = 50;
mesh = 'U';
s = [0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     100];
t = [0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0;
     0.4166667];
tdpar = [false; false; false];
r = [0.1];
q = [0];
sigma = [0.4];
ntkeep = int64(4);
[sOut, tOut, f, theta, delta, gamma, lambda, rho, ifail] = ...
    d03nc(kopt, x, mesh, s, t, tdpar, r, q, sigma, ntkeep,'alpha',1)
 

sOut =

     0
     5
    10
    15
    20
    25
    30
    35
    40
    45
    50
    55
    60
    65
    70
    75
    80
    85
    90
    95
   100


tOut =

         0
    0.0417
    0.0833
    0.1250
    0.1667
    0.2083
    0.2500
    0.2917
    0.3333
    0.3750
    0.4167


f =

   50.0000   50.0000   50.0000   50.0000
   45.0000   45.0000   45.0000   45.0000
   40.0000   40.0000   40.0000   40.0000
   35.0000   35.0000   35.0000   35.0000
   30.0000   30.0000   30.0000   30.0000
   25.0000   25.0000   25.0000   25.0000
   20.0000   20.0000   20.0000   20.0000
   15.0000   15.0000   15.0000   15.0000
   10.1543   10.0958   10.0464   10.0117
    6.5848    6.4424    6.2916    6.1306
    4.0672    3.8785    3.6729    3.4463
    2.4264    2.2423    2.0454    1.8336
    1.4174    1.2662    1.1096    0.9481
    0.8195    0.7072    0.5953    0.4851
    0.4724    0.3941    0.3190    0.2484
    0.2726    0.2202    0.1717    0.1282
    0.1573    0.1233    0.0929    0.0667
    0.0897    0.0685    0.0501    0.0347
    0.0484    0.0363    0.0259    0.0175
    0.0211    0.0156    0.0110    0.0073
         0         0         0         0


theta =

         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
   -1.4043   -1.1857   -0.8328   -0.2806
   -3.4185   -3.6183   -3.8646   -4.1880
   -4.5285   -4.9339   -5.4387   -6.0796
   -4.4165   -4.7277   -5.0821   -5.4821
   -3.6294   -3.7585   -3.8748   -3.9632
   -2.6946   -2.6860   -2.6441   -2.5561
   -1.8790   -1.8018   -1.6941   -1.5505
   -1.2578   -1.1621   -1.0461   -0.9097
   -0.8154   -0.7282   -0.6301   -0.5231
   -0.5084   -0.4411   -0.3689   -0.2943
   -0.2928   -0.2484   -0.2024   -0.1566
   -0.1324   -0.1108   -0.0888   -0.0674
         0         0         0         0


delta =

   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -1.0000   -1.0000   -1.0000   -1.0000
   -0.9846   -0.9904   -0.9954   -0.9988
   -0.8415   -0.8558   -0.8708   -0.8869
   -0.6087   -0.6217   -0.6373   -0.6565
   -0.4158   -0.4200   -0.4246   -0.4297
   -0.2650   -0.2612   -0.2563   -0.2498
   -0.1607   -0.1535   -0.1450   -0.1348
   -0.0945   -0.0872   -0.0791   -0.0700
   -0.0547   -0.0487   -0.0424   -0.0357
   -0.0315   -0.0271   -0.0226   -0.0182
   -0.0183   -0.0152   -0.0122   -0.0093
   -0.0109   -0.0087   -0.0067   -0.0049
   -0.0069   -0.0053   -0.0039   -0.0027
   -0.0048   -0.0036   -0.0026   -0.0017
   -0.0042   -0.0031   -0.0022   -0.0015


gamma =

         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
    0.0062    0.0038    0.0019    0.0005
    0.0510    0.0500    0.0480    0.0443
    0.0421    0.0436    0.0454    0.0479
    0.0351    0.0371    0.0396    0.0429
    0.0253    0.0264    0.0277    0.0291
    0.0164    0.0167    0.0169    0.0169
    0.0100    0.0098    0.0095    0.0091
    0.0059    0.0056    0.0052    0.0047
    0.0034    0.0031    0.0027    0.0024
    0.0019    0.0017    0.0014    0.0012
    0.0011    0.0009    0.0007    0.0006
    0.0006    0.0005    0.0004    0.0003
    0.0003    0.0002    0.0002    0.0001
         0         0         0         0


lambda =

         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
    6.3243    5.1893    3.8089    2.1118
   10.7215    9.9718    9.2140    8.4953
   12.3808   11.8073   11.2277   10.6365
   11.4834   10.8366   10.1417    9.3795
    9.3227    8.5840    7.7870    6.9211
    6.9621    6.2206    5.4412    4.6264
    4.9268    4.2651    3.5937    2.9227
    3.3602    2.8204    2.2920    1.7866
    2.2221    1.8126    1.4248    1.0683
    1.4122    1.1240    0.8586    0.6225
    0.8269    0.6459    0.4825    0.3408
    0.3789    0.2925    0.2155    0.1498
         0         0         0         0


rho =

         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
         0         0         0         0
   -7.1918   -6.0114   -4.5204   -2.5855
   -8.4541   -7.6378   -6.8479   -6.1657
   -7.5988   -6.9323   -6.2879   -5.6707
   -5.8905   -5.2837   -4.6809   -4.0772
   -4.1854   -3.6547   -3.1306   -2.6135
   -2.8221   -2.3904   -1.9743   -1.5775
   -1.8437   -1.5137   -1.2055   -0.9228
   -1.1812   -0.9407   -0.7233   -0.5316
   -0.7451   -0.5768   -0.4292   -0.3038
   -0.4591   -0.3466   -0.2506   -0.1716
   -0.2655   -0.1966   -0.1389   -0.0927
   -0.1228   -0.0898   -0.0626   -0.0410
         0         0         0         0


ifail =

                    0



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