# Usage of Piecewise/RSolve to solve an equation based on condition

I am trying to compute price of call option on a stock using the finite difference method. I have a VBA code for it, and I am trying to find out how I can use Mathematica memoization functions and recursive function to code it efficiently.

The VBA code is as below.

optval(vol, intrate, expn, payoff, strike, Etype, NAS)
ReDim S(0 To NAS) ReDim vold(0 To NAS) ReDim Vnew(0 To NAS)
ReDim Dummy(0    To NAS, 1 To 3)
ds = 2 * strike / NAS dt = 0.9 / vol / vol / NAS / NAS
NTS = Int(expn / dt) + 1 dt = expn / NTS
q = 1 If payoff = "p" Then q = -1 For i = 0 To NAS
S(i) = i * ds Dummy(i, 1) = S(i) vold(i) = Application.Max(q * (S(i) -  Strike), 0)
Dummy(i, 2) = vold(i) Next i
For k = 1 To NTS For i = 1 To NAS - 1
Delta = (vold(i + 1) - vold(i - 1)) / 2 / ds
Gamma = (vold(i + 1) - 2 * vold(i) + vold(i - 1)) / ds / ds Theta = -0.5 *  vol * vol * S(i) * S(i) * Gamma - intrate * S(i) * Delta + intrate * vold(i)

Vnew(i) = vold(i) - Theta * dt
Next i
Vnew(0) = (1 - intrate * dt) * vold(0)
Vnew(NAS) = 2 * vold(NAS - 1) - vold(NAS - 2)
For i = 0 To NAS
vold(i) = Vnew(i)
Next i
If Etype = "Y" Then
For i = 0 To NAS
vold(i) = Application.Max(vold(i), Dummy(i, 2))
Next i
End If
Next k

For i = 0 To NAS
Dummy(i, 3) = vold(i)
Next i
optval = Dummy End Function


The computation method involves creation of array/table which uses three seperate equation. The above VBA code uses three variable dummy[[i+1, 1]], dummy[[I+1,2]] and dummy[[i+3]] to store final output, and in addition uses vold and vnew to compute one node from previous node backwards.

Can we use memoization functions/recursive functions to do this in Mathematica?

I uploaded the excel pictures of 4X4 matrix with expected value of individual element in it to explain the logic of computation. I am using Piecewise to create a array (table d using function fd[i,j]) using four conditions as in the below code to improve implementation the logic. The computation start from step j = ts+1 using step 1 below and recursively computes all the elements of array back to j = 1 using condions in step 2-4 below.

1. 1st equation is when fd[i,j] = Max[s[[i]] - strike1, 0] when j=ts+1 This condition is working in below code.

2. For ( i == 1 && 0 < j < ts + 1), the code is expected to use (1 - int)*fd[i, j + 1] to compute fd[I,j]. For example, with respect to uploaded matrix, element of Column 2 and 1st row = (1-int*dt)*element(1,3);

3. For (1 < i < as + 2 && 0 < j, ts + 1), the code is expected to use fd[i, j - 1] = del*dt.

For example with respect to uploaded matrix, I am expecting elements of row 2 and 3 and column 1 and 2 to be computed recursively using function del for example element of row3 and column2 is derived from elements of row 2,column 3 and row4 column 3 using function del. Hence element of 2nd row and column 1= ((element(3,2)-element(1,2))/ds)*dt, and element of 3rd row and column 2= ((element(4,3)-element(2,3))/ds)*dt

4. For (i = as + 2 && 0 < j < ts + 1), the code is expected to use 2*fd[i - 1, j] - fd[i - 2, j].

For example, with respect to uploaded matrix, element of row 4 and column 2 = 2*(elements of row 3 and column 2) -(elements of row2 and column 2). Similiarly element of row 4 and column 1 = 2*(elements of row 3 and column 1) -(elements of row2 and column 1)

The code can be run using parameter as recursive[4, 3, 3, 0.01, 0.02, 1]. The below code computes the 1st condition but gets into problems for rest of conditions.

I also tried to implement logic using a series of Do loops. Pl refers to code in my previous question - Matrix iteration failing after second step. There the code was incorrectly computing last row (i.e as+1) from last but one column.

I had a memoization issue in that part of loop which I didn't understand how to solve. Rest of the rows and columns were correctly computed for second iteration (i.e. last but one column). As successive columns always use data from previous column. The computation for successive columns gets incorrect from bottom (i.e. for row as, as-1 etc).

Can anybody help in implementing above logic using memoization or recursive functions like RSolve/Recurrence Table instead?

recursive[as_, ts_, ds_, dt_, int_, strike1_] :=
Module[{s, d, del}, s = Table[i*ds, {i, 0, as + 1}];
fd[i_, j_] := fd[i, j] = Piecewise[{
{Max[s[[i]] - strike1, 0], j == ts + 1},
{(1 - int)*fd[i, j + 1], i == 1 && 1 < j < ts + 1},
{fd[i, j - 1] = del*dt, 1 < i < as + 2 && 0 < j< ts + 1},
{2*fd[i - 1, j] - fd[i - 2, j], i == as + 2 && 1 < j < ts + 1}
}];
del := ( fd[i + 1, j] - fd[i - 1, j])/ds;
d = Table[fd[i, j], {i, 1, as + 2}, {j, 1, ts + 1}];
Grid[d]
]

• Much better question. Thank you. I will study this and see what I can come come up with. – Bill May 7 '15 at 16:19

EDIT: Corrections VBA to Mathematica code translation

optval[vol_, intrate_, expn_, payoff_, strike_, Etype_, NAS_] :=
Module[{S, vold, Vnew, Dummy, ds, dt, NTS, q, gamma, Delta, Theta},
S = Table[0, {NAS + 1}];
vold = Table[0, {NAS + 1}];
Vnew = Table[0, {NAS + 1}];
Dummy = Table[0, {NAS + 1}, {3}];
ds = 2*strike/NAS;
dt = 0.9/vol^2/NAS^2;
NTS = Floor[expn/dt] + 1;
dt = expn/NTS;
q = 1;
If[payoff == "p", q = -1];(*End If payoff???*)
For[i = 0, i <= NAS, i++,
S[[i + 1]] = i*ds;
Dummy[[i + 1, 1]] = S[[i + 1]];
vold[[i + 1]] = Max[q*(S[[i + 1]] - strike), 0];
Dummy[[i + 1, 2]] = vold[[i + 1]]
];(*End Next i*)
For[k = 1, k <= NTS, k++,
For[i = 1, i <= NAS - 1, i++,
Delta = (vold[[i + 1 + 1]] - vold[[i - 1 + 1]])/2/ds;
gamma = (vold[[i + 1 + 1]] - 2*vold[[i + 1]] + vold[[i - 1+1]])/ds^2;
Theta = -0.5*vol^2*S[[i + 1]]^2*gamma - intrate*S[[i + 1]]*Delta +
intrate*vold[[i + 1]];
Vnew[[i + 1]] = vold[[i + 1]] - Theta*dt
];(*End Next i*)
Vnew[[0 + 1]] = (1 - intrate*dt)*vold[[0 + 1]];
Vnew[[NAS + 1]] = 2*vold[[NAS - 1 + 1]] - vold[[NAS - 2 + 1]];
For[i = 0, i <= NAS, i++,
vold[[i + 1]] = Vnew[[i + 1]]
];(*End Next i*)
If[Etype == "Y",
For[i = 0, i <= NAS, i++,
vold[[i + 1]] = Max[vold[[i + 1]], Dummy[[i + 1, 2]]]
](*End Next i*)
](*End If Etype*)
];(*End Next k*)
For[i = 0, i <= NAS, i++,
Dummy[[i + 1, 3]] = vold[[i + 1]]
];(*End Next i*)
Dummy
];
optval[0.2, 0.05, 1, "c", 100, "N", 10] // MatrixForm


EDIT: Simplified the first two columns, working on the third now.

optoptval[vol_, intrate_, expn_, payoff_, strike_, Etype_, NAS_] :=
Module[{q=If[payoff=="p",1,-1]},
Table[{2*i*strike/NAS, Max[0, q*(NAS-2*i)*strike/NAS]}, {i,0,NAS}]];

Print[optval[0.2, 0.05, 1, "c", 100, "N", 6] // MatrixForm,
optoptval[0.2, 0.05, 1, "c", 100, "N", 6] // MatrixForm]

Print[optval[0.2, 0.05, 1, "p", 100, "N", 10] // MatrixForm,
optoptval[0.2, 0.05, 1, "p", 100, "N", 10] // MatrixForm]

• I have simplified the first two columns. Please compare optoptval with optval for your test cases. If you find any errors then please show me the inputs. I am still working on the last column, but your code is much more complicated for that and I can't tell yet if I can make it much simpler. – Bill May 9 '15 at 16:29
• Yes 1st two column for output of optval is working fine. Agree last column is not easy to do. If you can please help in implementing recursive[4, 3, 3, 0.01, 0.02, 1]. I have simplified step3 for building recursive. Actual implementation in optval uses few more formulas than recursive. If you can implement recursive, I can use the concept to implement optval. Thanks for working on it. – Kausik May 10 '15 at 3:22
• Let me know where are you on implementation, particularly implementation of recursive functin. Do we need to flag it to moderator for getting more help. Sorry for following up. Thanks for your time and attention. – Kausik May 11 '15 at 17:54
• I'm sorry. I've tried. I suspect there was simplicity in the original mathematics, but trying to see that through your VB translation is just more than I am able to do. I hope it works out for you. I'm done. – Bill May 12 '15 at 17:27