White R.E. Computational Mathematics: Models, Methods, and Analysis with MATLAB and MPI
Подождите немного. Документ загружается.
230 CHAPTER 5. EPIDEMICS, IMAGES AND MONEY
5.6.4 Method
The explicit time discretization with projection to the payo obstacle will be
used. Again replace the time variable by W w= Now abuse notation a
little and write S (V
1
> V
2
> ) in place of S(V
1
> V
2
> w) = S (V
1
> V
2
> W ) so that
S
replaces S
w
. With the initial and boundary conditions given one may
apply the explicit finite di
erence method as used for heat diusion in two
directions to obtain
S
n+1
l>m
approximations for S (lV
1
> mV
2
> (n + 1)). For
the inequality condition simply choose the S
n+1
l>m
= max(H lV
1
mV
2
> 0)
if
S
n+1
l>m
? max(H lV
1
mV
2
> 0)=
Explicit Method with Projection for (5.6.1)-(5.6.8).
Let
1
(@(V
1
)
2
)(
2
1
@2)>
2
(@(V
2
)
2
)(
2
2
@2)>
12
(@(2V
1
2V
2
))
12
(
1
2
@2)
V
1l
lV
1
and V
2m
mV
2
S
n+1
l>m
= S
n
l>m
+
1
V
2
1
l
(S
n
l
1>m
2S
n
l>m
+ S
n
l
+1>m
)
+2V
1l
V
2m
12
((S
l+1>m+1
S
l1>m+1
) (S
l+1>m1
S
l1>m1
))
+
2
V
2
2
m
(S
n
l>m
1
2S
n
l>m
+ S
n
l>m
+1
)
+ (
@V
1
) uV
1l
(S
n
l
+1>m
S
n
l>m
)
+ (
@V
2
) uV
2m
(S
n
l>m
+1
S
n
l>m
)
uS
n
l>m
S
n+1
lm
= max(S
n+1
lm
> max(H lV
1
mV
2
> 0))=
5.6.5 Implementation
The two underlying assets model contains four one dimensional models along
the boundary of the two underlying assets domain. In the interior of the domain
the two dimensional Black-Scholes equation must be solved. In the MATLAB
cod e bs2d.m the input is given in lines 1-35, and two possible payo functions
are defined in lines 26-35. The time loop is done in lines 36-124. The one
dimensional Black-Scholes equation along the boundaries are solved in lines 44-
57 for | = 0, 58-72 for | = O, 73-87 for { = 0 and 88-102 for { = O. For
| = 0 the projection to the payo obstacle is done in lines 50-52, and the times
when the put value separates from the payo
obstacle is found in lines 53-56.
The two dimensional Black-Scholes equation is s olved for the interior nodes in
lines 103-123. Lines 124-137 generate four graphs depicting the value of the put
MATLAB Code bs2d.m
1. % Program bs2d
2. % Black-Scholes Equation
3. % Two underlying assets
© 2004 by Chapman & Hall/CRC
contract and the optimal exercise times, see Figures 5.6.2 and 5.6.3.
5.6. BLACK-SCHOLES MODEL FOR TWO ASSETS 231
4. % Explicit time with projection
5. clear;
6. n = 40;
7. maxk = 1000;
8. f = .00;
9. T=.5;
10. dt = T/maxk;
11. L=4;
12. dx = L/n;
13. dy = L/n;
14. sig1 = .4;
15. sig2 = .4;
16. rho12 = .3;
17. E1 = 1;
18. E2 = 1.5;
19. total = E1 + E2;
20. alpha1 = .5*sig1^2*dt/(dx*dx);
21. alpha2 = .5*sig2^2*dt/(dy*dy);
22. alpha12 = .5*sig1*sig2*rho12*dt/(2*dx*2*dy);
23. r = .12;
24. sur = zeros(n+1);
25. % Insert Initial Condition
26. for j = 1:n+1
27. y(j) = dy*(j-1);
28. for i = 1:n+1
29. x(i) = dx*(i-1);
30. %Define the payo
function
31. u(i,j,1) = max(E1-x(i),0.0) + max(E2-y(j),0.0);
32. % u(i,j,1) = max(total -(x(i) + y(j)),0.0);
33. suro(i,j) = u(i,j,1);
34. end;
35. end;
36. % Begin Time Steps
37. for k = 1:maxk
38. % Insert Boundary Conditions
39. u(n+1,1,k+1) = E2;
40. u(n+1,n+1,k+1) = 0.0;
41. u(1,n+1,k+1) = E1 ;
42. u(1,1,k+1) = total;
43. % Do y = 0.
44. j=1;
45. for i = 2:n
46. u(i,j,k+1) = dt*f+x(i)*x(i)*alpha1*...
47. (u(i-1,j,k) + u(i+1,j,k)-2.*u(i,j,k)) ...
48. +u(i,j,k)*(1 -r*dt)- ...
© 2004 by Chapman & Hall/CRC
232 CHAPTER 5. EPIDEMICS, IMAGES AND MONEY
49. r*x(i)*dt/dx*(u(i,j,k)-u(i+1,j,k));
50. if (u(i,j,k+1)
?suro(i,j))
51. u(i,j,k+1) = suro(i,j);
52. end
53. if ((u(i,j,k+1)
Asuro(i,j))&...
54. (u(i,j,k)==suro(i,j)))
55. sur(i,j)= k+.5;
56. end
57. end
58. % Do y = L.
59. j=n+1;
60. for i = 2:n
61. u(i,j,k+1) = dt*f+x(i)*x(i)*alpha1*...
62. (u(i-1,j,k) + u(i+1,j,k)-2.*u(i,j,k)) ...
63. +u(i,j,k)*(1 -r*dt)- ...
64. r*x(i)*dt/dx*(u(i,j,k)-u(i+1,j,k));
65. if (u(i,j,k+1)
?suro(i,j))
66. u(i,j,k+1) = suro(i,j);
67. end
68. if ((u(i,j,k+1)
Asuro(i,j))&...
69. (u(i,j,k)==suro(i,j)))
70. sur(i,j)= k+.5;
71. end
72. end
73. % Do x = 0.
74. i=1;
75. for j = 2:n
76. u(i,j,k+1) = dt*f+y(j)*y(j)*alpha2*...
77. (u(i,j-1,k) + u(i,j+1,k)-2.*u(i,j,k))...
78. +u(i,j,k)*(1 -r*dt)-...
79. r*y(j)*dt/dy*(u(i,j,k)-u(i,j+1,k));
80. if (u(i,j,k+1)
?suro(i,j))
81. u(i,j,k+1) = suro(i,j);
82. end
83. if ((u(i,j,k+1)
Asuro(i,j)) &...
84. (u(i,j,k)==suro(i,j)))
85. sur(i,j)= k+.5;
86. end
87. end
88. % Do x = L.
89. i=n+1;
90. for j = 2:n
91. u(i,j,k+1) = dt*f+y(j)*y(j)*alpha2*...
92. (u(i,j-1,k) + u(i,j+1,k)-2.*u(i,j,k))...
93. +u(i,j,k)*(1 -r*dt)-...
© 2004 by Chapman & Hall/CRC
5.6. BLACK-SCHOLES MODEL FOR TWO ASSETS 233
94. r*y(j)*dt/dy*(u(i,j,k)-u(i,j+1,k));
95. if (u(i,j,k+1)
?suro(i,j))
96. u(i,j,k+1) = suro(i,j);
97. end
98. if ((u(i,j,k+1)
Asuro(i,j))&...
99. (u(i,j,k)==suro(i,j)))
100. sur(i,j)= k+.5;
101. end
102. end
103. % Solve for Interior Nodes
104. for j= 2:n
105. for i = 2:n
106. u(i,j,k+1) = dt*f+x(i)*x(i)*alpha1*...
107. (u(i-1,j,k) + u(i+1,j,k)-2.*u(i,j,k))...
108. +u(i,j,k)*(1 -r*dt)...
109. -r*x(i)*dt/dx*(u(i,j,k)-u(i+1,j,k))...
110. +y(j)*y(j)*alpha2*...
111. (u(i,j-1,k) + u(i,j+1,k)-2.*u(i,j,k)) ...
112. -r*y(j)*dt/dy*(u(i,j,k)-u(i,j+1,k)) ...
113. +2.0*x(i)*y(j)*alpha12*...
114. (u(i+1,j+1,k)-u(i-1,j+1,k)
-u(i+1,j-1,k)+u(i-1,j-1,k));
115. if (u(i,j,k+1)
?suro(i,j))
116. u(i,j,k+1) = suro(i,j);
117. end
118. if ((u(i,j,k+1)
Asuro(i,j)) &...
119. (u(i,j,k)==suro(i,j)))
120. sur(i,j)= k+.5;
121. end
122. end
123. end
124. end
125. figure(1)
126. subplot(2,2,1)
127. mesh(x,y,suro’)
128. title(’Payo
Value’)
129. subplot(2,2,2)
130. mesh(x,y,u(:,:,maxk+1)’)
131. title(’Put Value’)
132. subplot(2,2,3)
133. mesh(x,y,u(:,:,maxk+1)’-suro’)
134. title(’Put Minus Payo
’)
135. subplot(2,2,4)
136. mesh(x,y,sur’*dt)
137. title(’Optimal Exercise Times’)
© 2004 by Chapman & Hall/CRC
234 CHAPTER 5. EPIDEMICS, IMAGES AND MONEY
Figure 5.6.2: max(H
1
+ H
2
V
1
V
2
> 0)
Figure 5.6.2 is for the lumped payo
max(H
1
+ H
2
V
1
V
2
> 0) with two
assets with di
erent volatilities
1
= 0=4 and
2
= 0=1= All four graphs have
the underlying asset values on the { and | axes. The upper left graph is for
the payo
function, and the upper right graph is for the value of the American
put option contract at time equal to 0
=5 before the expiration date, W . Note
the di
erence in the graphs along the axes, which can be attributed to the
larger volatility for the first underlying asset. The lower left graph depicts the
dierence in the put value at time equal to 0=5 and the value of the payo. The
lower right graph depicts time of optimal exercise versus the two underlying
assets. Here the vertical axis has
= W w> and the interface for the optimal
exercise of the put moves towards the vertical axis as
increases.
max(H
1
V
1
> 0) + max(H
2
V
2
> 0)
with two assets and with equal volatilities
1
= 0=4 and
2
= 0=4= and with
dierent exercise values H
1
= 1=0 and H
2
= 1=5= All four graphs have the
underlying asset values on the
{ and | axes. The upper left graph is for the
payo
function, and the upper right graph is for the value of the American
put option contract at time equal to 0=5 before the expiration date. Note the
di
erence in the graphs along the axes, which can be attributed to the dierent
exercise values. The lower left graph depicts the di erence in the put value at
time equal to 0
=5 and the value of the payo. The lower right graph has time
of optimal exercise versus the two underlying assets. Here the vertical axis has
= W w where W is the expiration date. There are three interface curves for
possible optimal exercise times. One is inside the underlying asset region [0
H
1
]× [0 H
2
] and moves towards the } axis as increases. The other two are in
© 2004 by Chapman & Hall/CRC
Figure 5.6.3 is for the distributedpayo
5.6. BLACK-SCHOLES MODEL FOR TWO ASSETS 235
Figure 5.6.3: max(H
1
V
1
> 0) + max(H
2
V
2
> 0)
the regions [
H
1
4] × [0 H
2
] and [0 H
1
] × [H
2
4] and move towards the vertical
lines containing the points (4> 0> 0) and (0> 4> 0), resp ectively.
5.6.6 Assessment
The parameters in the Black-Scholes model may depend on time and may be
di!cult to estimate. The precise assumptions under which the Black-Scholes
equation models option contracts should be carefully studied. The question of
convergence of the discrete model to the continuous model needs to be exam-
ined.
The jump from one to multiple asset models presents some interesting dif-
ferences from heat and mass transfer problems. The boundary conditions are
implied from the solution of the one asset Black-Scholes equation. When going
from one to two assets, there may be a variety of payo
or obstacle functions
that will result in multiple interfaces for the possible optimal exercise opportu-
nities. There may be many assets as contrasted to only three space directions
for heat and mass transfer.
5.6.7 Exercises
1. Experiment with variable time and asset space steps.
2.
3. Experiment with other payo
functions=
4. Experiment with variable interest rates u.
© 2004 by Chapman & Hall/CRC
Duplicate the calculations in Figures 5.6.2 and 5.6.3.
236 CHAPTER 5. EPIDEMICS, IMAGES AND MONEY
5. Experiment with variable exercise values
H=
6. Experiment with variable expiration time W=
7. Experiment with variable volatility =
© 2004 by Chapman & Hall/CRC
Chapter 6
High Performance
Computing
Because many applications are very complicated such as weather prediction,
there will often be a large number of unknowns. For heat di
usion in one direc-
tion the long thin wire was broken into
q smaller cells and the temperature was
approximated for all the time steps in each segment. Typical values for
q could
range from 50 to 100. If one needs to model temperature in a tropical storm,
then there will be di
usion in three directions. So, if there are q = 100 segments
in each direction, then there will be
q
3
= 100
3
= 10
6
cells each with unknown
temperatures, velocity components and pressures. Three dimensional problems
present strong challenges in memory, computational units and s torage of data.
The first three sections are a very brief description of serial, vector and mul-
tiprocessing computer architectures. The last three sections illustrate the use
of the IBM/SP and MPI for the parallel computation of matrix products and
the two and three space variable models of heat di
usion and pollutant trans-
fer.
subroutines for communication among the pro cessors. Additional introductory
materials on parallel computations can be found in P. S. Pacheco [21] and more
advanced topics in Dongarra, Du
, Sorensen and van der Vorst [6].
6.1 Vector Computers and Matrix Products
6.1.1 Intro duction
In this section we consider the components of a computer and the various ways
they are connected. In particular, the idea behind a vector pipeline is in-
troduced, and a model for speedup is presented. Applications to matrix-vector
products and to heat and mass transfer in two directions will be presented. The
sizes of matrix models substantially increase when the heat and mass transfer
237
© 2004 by Chapman & Hall/CRC
Chapter 7 contains amore detailed description of MPI and the essential
238 CHAPTER 6. HIGH PERFORMANCE COMPUTING
Figure 6.1.1: von Neumann Computer
in two or three directions are modeled. This is the reason for considering vector
and multiprocessor computers.
The von Neumann definition of a computer contains three parts: main
memory, input-output device and central proce ssing unit (CPU). The CPU
has three components: the arithmetic logic unit, the control unit and the local
memory. The arithmetic logic unit does the floating point calculations while
the control unit governs the instructions and data. Figure 6.1.1 illustrates a
von Neumann computer with the three basic component s . The lo cal memory
is small compared to the main memory, but moving data within the CPU is
usually very fast. Hence, it is important to move data from the main memory to
the local memory and do as much computation with this data as possible before
moving it back to the main memory. Algorithms that have been optimized for
a particular computer will take these facts into careful consideration.
Another way of describing a computer is the hierarchical classification of its
components. There are three levels: the processor level with wide band commu-
nication paths, the register level with several bytes (8 bits per byte) pathways
and the gate or logic level with several bits in its pathways.
processor level depiction of a multiprocessor computer with four CPUs. The
CPUs communicate with each other via the shared memory. The switch con-
trols access to the shared memory, and here there is a potential for a bottleneck.
The purpose of multiprocessors is to do more computation in less time. This is
critical in many applications such as weather prediction.
Within the CPU is the arithmetic logic u nit with many floating point adders,
which are register level devices. A floating point add can be described in four
distinct steps each requiring a distinct hardware segment. For example, use
four digits to do a floating point add 100.1 + (-3.6):
CE: compare expressions .1001 · 10
3
and =36 · 10
1
AE: mantissa alignment =1001 · 10
3
and =0036 · 10
3
AD: mantissa add 1001 0036 = 0965
NR: normalization
=9650 · 10
2
=
ways with several bytes of data. The data moves from left to right in time
intervals equal to the clock cycle time of the particular computer. If each step
© 2004 by Chapman & Hall/CRC
Figure 6.1.2 is a
This is depicted by Figure 6.1.3 where the lines indicate communication path-
6.1. VECTOR COMPUTERS AND MATRIX PRODUCTS 239
Figure 6.1.2: Shared Memory Multiprocessor
Figure 6.1.3: Floating Point Add
takes one clock cycle and the clock cycle time is 6 nanoseconds, then a floating
point add will take 24 nanoseconds (10
9
sec.).
Within a floating point adder there are many devices that add integers.
These devices typically deal with just a few bits of information and are examples
of gate or logic level devices. Here the integer adds can be done by b ase two
numbers. These devices are combinations of a small number of transistors
designed to simulate truth tables that reflect basic binary operations. Table
6.1.1 indicates how one digit of a base two number with 64 digits (or bits) can
the input is
{, | and f
1
(the carry from the previous
digit) and the output is
} and f
0
(the carry to the next digit).
Table 6.1.1: Truth Table for Bit Adder
x y c
1
z c
0
0 0 0 0 0
0 0 1 1 0
0 1 0 1 0
0 1 1 0 1
1 0 0 1 0
1 0 1 0 1
1 1 0 0 1
1 1 1 1 1
© 2004 by Chapman & Hall/CRC
be added. In Figure 6.1.4