Skip to content

Commit

Permalink
Local vol grid fixes
Browse files Browse the repository at this point in the history
  • Loading branch information
gavbrennan committed Sep 17, 2024
1 parent 248b3d6 commit 9b150c5
Show file tree
Hide file tree
Showing 2 changed files with 125 additions and 26 deletions.
69 changes: 67 additions & 2 deletions src/Qwack.Options/VolSurfaces/LocalVol.cs
Original file line number Diff line number Diff line change
@@ -1,5 +1,6 @@
using System;
using System.Collections.Generic;
using System.Diagnostics;
using System.Linq;
using Qwack.Core.Basic;
using Qwack.Core.Basic.Correlation;
Expand Down Expand Up @@ -80,7 +81,7 @@ public static double[][] ComputeLocalVarianceOnGrid(this IVolSurface VanillaSurf
return lvGrid;
}

public static double[][] ComputeLocalVarianceOnGridFromCalls(this IVolSurface VanillaSurface, double[][] strikes, double[] timeSteps, Func<double, double> forwardFunc)
public static double[][] ComputeLocalVarianceOnGridFromCalls(this IVolSurface VanillaSurface, double[][] strikes, double[] timeSteps, Func<double, double> forwardFunc, int firstTimeStep = 0)
{
var numberOfTimesteps = timeSteps.Length;
var deltaK = 0.001 * forwardFunc(timeSteps[0]);
Expand All @@ -89,7 +90,7 @@ public static double[][] ComputeLocalVarianceOnGridFromCalls(this IVolSurface Va
var fwds = timeSteps.Select(t => forwardFunc(t)).ToArray();

//ParallelUtils.Instance.For(1, numberOfTimesteps, 1, it =>
for (var it = 1; it < numberOfTimesteps; it++)
for (var it = firstTimeStep + 1; it < numberOfTimesteps; it++)
{
var T = timeSteps[it];
var T1 = timeSteps[it - 1];
Expand Down Expand Up @@ -131,6 +132,70 @@ public static double[][] ComputeLocalVarianceOnGridFromCalls(this IVolSurface Va
return lvGrid;
}

public static double[][] ComputeLocalVarianceOnGridFromCalls2(this IVolSurface VanillaSurface, double[][] strikes, double[] timeSteps, Func<double, double> forwardFunc, int firstTimeStep = 0)
{
var numberOfTimesteps = timeSteps.Length;
var deltaK = 0.001 * forwardFunc(timeSteps[0]);
var lvGrid = new double[numberOfTimesteps - 1][];

var fwds = timeSteps.Select(t => forwardFunc(t)).ToArray();

//ParallelUtils.Instance.For(1, numberOfTimesteps, 1, it =>
for (var it = firstTimeStep + 1; it < numberOfTimesteps; it++)
{
var T = timeSteps[it];
var T1 = timeSteps[it - 1];
var fwd = fwds[it];
var bump = fwd / 10000;
var fwdtm1 = fwds[it - 1];
var rmq = Log(fwd / fwdtm1) / (T - T1);
var numberOfStrikes = strikes[it - 1].Length;

lvGrid[it - 1] = new double[numberOfStrikes];

if (numberOfStrikes > 1)
{
for (var ik = 0; ik < numberOfStrikes; ik++)
{

var K = strikes[it][ik];
var V = VanillaSurface.GetVolForAbsoluteStrike(K, T, fwd);
var C = BlackFunctions.BlackPV(fwd, K, 0.0, T, V, OptionType.C);
var Vtm1 = VanillaSurface.GetVolForAbsoluteStrike(K, T1, fwdtm1);

var Cp1 = BlackFunctions.BlackPV(fwd, K + bump, 0.0, T, V, OptionType.C);
var Cp2 = BlackFunctions.BlackPV(fwd, K + 2*bump, 0.0, T, V, OptionType.C);
var Cm1 = BlackFunctions.BlackPV(fwd, K - bump, 0.0, T, V, OptionType.C);
var Cm2 = BlackFunctions.BlackPV(fwd, K - 2 * bump, 0.0, T, V, OptionType.C);


//var dcdt = -BlackFunctions.BlackTheta(fwd, K, 0.0, T, V, OptionType.C);
var dcdt = -(BlackFunctions.BlackPV(fwdtm1, K, 0.0, T1, Vtm1, OptionType.C) - C) / (T - T1);
var dcdk = (Cp1 - Cm1) / (2 * bump);
var dcdkp = (Cp2 - Cp1) / bump;
var dcdkm = (Cm1 - Cm2) / bump;
var d2cdk2 = (dcdkp - dcdkm) / bump;

var localVariance = d2cdk2 == 0 ? V * V : (dcdt - rmq * (C - K * dcdk)) / (0.5 * K * K * d2cdk2);
//if(double.IsNaN(localVariance))
//{
// Debug.Assert(true);
//}
lvGrid[it - 1][ik] = localVariance;
}
}
else
{
var K = strikes[it][0];
var V = VanillaSurface.GetVolForAbsoluteStrike(K, T, fwd);
lvGrid[it - 1][0] = V * V;
}
}//, false).Wait();

return lvGrid;
}


public static List<double[][]> LocalCorrelationRaw(this IAssetFxModel model, double[] times)
{
var o = new List<double[][]>();
Expand Down
82 changes: 58 additions & 24 deletions src/Qwack.Paths/Processes/LVSingleAsset.cs
Original file line number Diff line number Diff line change
Expand Up @@ -32,6 +32,7 @@ public class LVSingleAsset : IPathProcess, IRequiresFinish
private bool _isComplete;
private Vector<double>[] _drifts;
private IInterpolator1D[] _lvInterps;
private Vector<double>[] _fixings;

private readonly Vector<double> _two = new(2.0);

Expand All @@ -42,7 +43,7 @@ public LVSingleAsset(IVolSurface volSurface, DateTime startDate, DateTime expiry
_surface = volSurface;
_startDate = startDate;
_expiryDate = expiryDate;
_numberOfSteps = nTimeSteps;
_numberOfSteps = nTimeSteps == 0 ? 100 : nTimeSteps;
_name = name;
_pastFixings = pastFixings ?? (new Dictionary<DateTime, double>());
_forwardCurve = forwardCurve;
Expand All @@ -66,18 +67,41 @@ public void Finish(IFeatureCollection collection)
return;
}

//fixings first
var dates = collection.GetFeature<ITimeStepsFeature>();
var fixings = new List<Vector<double>>();
for (var d = 0; d < dates.Dates.Length; d++)
{
var date = dates.Dates[d];
if (date >= _startDate) break;
try
{
var vect = new Vector<double>(_pastFixings[date]);
fixings.Add(vect);
}
catch (Exception e)
{
}

}
_fixings = [.. fixings];
var firstTime = _timesteps.Times[_fixings.Length];

//drifts and vols...
_drifts = new Vector<double>[_timesteps.TimeStepCount];
_lvInterps = new IInterpolator1D[_timesteps.TimeStepCount - 1];

var strikes = new double[_timesteps.TimeStepCount][];
var atmVols = new double[_timesteps.TimeStepCount];
for (var t = 0; t < strikes.Length; t++)
var adjTimes = new double[_timesteps.TimeStepCount];
for (var t = _fixings.Length; t < strikes.Length; t++)
{
var fwd = _forwardCurve(_timesteps.Times[t]);
atmVols[t] = _surface.GetVolForDeltaStrike(0.5, _timesteps.Times[t], fwd);
var time = _timesteps.Times[t] - firstTime;
adjTimes[t] = time;
var fwd = _forwardCurve(time);
atmVols[t] = _surface.GetVolForDeltaStrike(0.5, time, fwd);

if (_timesteps.Times[t] == 0)
if (time == 0)
{
strikes[t] = new double[] { fwd };
continue;
Expand All @@ -86,28 +110,30 @@ public void Finish(IFeatureCollection collection)
{
var nStrikes = 200;
var strikeStep = 1.0 / nStrikes;
strikes[t] = new double[nStrikes];
strikes[t] = new double[nStrikes-1];
for (var k = 0; k < strikes[t].Length; k++)
{
var deltaK = -(strikeStep + strikeStep * k);
strikes[t][k] = Options.BlackFunctions.AbsoluteStrikefromDeltaKAnalytic(fwd, deltaK, 0, _timesteps.Times[t], atmVols[t]);
strikes[t][k] = Options.BlackFunctions.AbsoluteStrikefromDeltaKAnalytic(fwd, deltaK, 0, time, atmVols[t]);
}
}
}

var lvSurface = Options.LocalVol.ComputeLocalVarianceOnGridFromCalls(_surface, strikes, _timesteps.Times, _forwardCurve);
var lvSurface = Options.LocalVol.ComputeLocalVarianceOnGridFromCalls2(_surface, strikes, adjTimes, _forwardCurve, _fixings.Length);

for (var t = 0; t < _lvInterps.Length; t++)
for (var t = _fixings.Length; t < _lvInterps.Length; t++)
{
_lvInterps[t] = InterpolatorFactory.GetInterpolator(strikes[t], lvSurface[t], t == 0 ? Interpolator1DType.DummyPoint : Interpolator1DType.LinearFlatExtrap);
_lvInterps[t] = InterpolatorFactory.GetInterpolator(strikes[t], lvSurface[t], t == _fixings.Length ? Interpolator1DType.DummyPoint : Interpolator1DType.LinearFlatExtrap);
}

var prevSpot = _forwardCurve(0);
for (var t = 1; t < _drifts.Length; t++)

for (var t = _fixings.Length + 1; t < _drifts.Length; t++)
{
var fxAtmVol = _adjSurface == null ? 0.0 : _adjSurface.GetForwardATMVol(0, _timesteps.Times[t]);
var driftAdj = _adjSurface == null ? 1.0 : Exp(atmVols[t] * fxAtmVol * _timesteps.Times[t] * _correlation);
var spot = _forwardCurve(_timesteps.Times[t]) * driftAdj;
var time = _timesteps.Times[t] - firstTime;
var fxAtmVol = _adjSurface == null ? 0.0 : _adjSurface.GetForwardATMVol(0, time);
var driftAdj = _adjSurface == null ? 1.0 : Exp(atmVols[t] * fxAtmVol * time * _correlation);
var spot = _forwardCurve(time) * driftAdj;

_drifts[t] = new Vector<double>(Log(spot / prevSpot) / _timesteps.TimeSteps[t]);

Expand All @@ -122,14 +148,17 @@ public void Process(IPathBlock block)
{
var previousStep = new Vector<double>(_forwardCurve(0));
var steps = block.GetStepsForFactor(path, _factorIndex);
var c = 0;
foreach (var kv in _pastFixings.Where(x => x.Key < _startDate).OrderBy(x => x.Key))
{
steps[c] = new Vector<double>(kv.Value);
c++;
}
var c = _fixings.Length;
_fixings.AsSpan().CopyTo(steps);
steps[c] = previousStep;

//foreach (var kv in _pastFixings.Where(x => x.Key < _startDate).OrderBy(x => x.Key))
//{
// steps[c] = new Vector<double>(kv.Value);
// c++;
//}
//steps[c] = previousStep;

if (_siegelInvert)
{
for (var step = c + 1; step < block.NumberOfSteps; step++)
Expand Down Expand Up @@ -175,11 +204,16 @@ public void SetupFeatures(IFeatureCollection pathProcessFeaturesCollection)
_timesteps = pathProcessFeaturesCollection.GetFeature<ITimeStepsFeature>();
_timesteps.AddDates(_pastFixings.Keys.Where(x => x < _startDate));

var stepSize = (_expiryDate - _startDate).TotalDays / _numberOfSteps;
var simDates = new List<DateTime> { _startDate };
for (var i = 0; i < _numberOfSteps; i++)
var periodSize = (_expiryDate - _startDate).TotalDays;
var stepSizeLinear = periodSize / _numberOfSteps;
var simDates = new List<DateTime>();

for (double i = 0; i < _numberOfSteps; i++)
{
simDates.Add(_startDate.AddDays(i * stepSize).Date);
var linStep = i * 1.0/_numberOfSteps;
var coshStep = (Cosh(linStep) - 1) / (Cosh(1) - 1);
var stepDays = coshStep * periodSize;
simDates.Add(_startDate.AddDays(stepDays));
}

_timesteps.AddDates(simDates.Distinct());
Expand Down

0 comments on commit 9b150c5

Please sign in to comment.