use std::any::Any;
use crate::errors::QlResult;
use crate::exercise::{EuropeanExercise, Exercise, ExerciseType};
use crate::instruments::{OneAssetOptionResults, OptionArguments, StrikedTypePayoff};
use crate::math::randomnumbers::rngtraits::McRngTraits;
use crate::methods::montecarlo::{
ControlVariate, EarlyExercisePathPricer, LongstaffSchwartzPathPricer, LsmBasisSystem, Path,
PolynomialType,
};
use crate::patterns::observable::{AsObservable, Observable};
use crate::pricingengine::{Arguments, PricingEngine, Results};
use crate::pricingengines::McLongstaffSchwartzEngineBase;
use crate::pricingengines::vanilla::{AnalyticEuropeanEngine, EuropeanPathPricer};
use crate::processes::GeneralizedBlackScholesProcess;
use crate::shared::{Shared, shared};
use crate::stochasticprocess::StochasticProcess1D;
use crate::types::{Real, Size};
use crate::{fail, require};
pub struct AmericanPathPricer {
payoff: Shared<dyn StrikedTypePayoff>,
scaling_value: Real,
polynomial_order: Size,
polynomial_type: PolynomialType,
exercise_indices: Option<Vec<bool>>,
}
impl AmericanPathPricer {
pub fn new(
payoff: Shared<dyn StrikedTypePayoff>,
polynomial_order: Size,
polynomial_type: PolynomialType,
) -> AmericanPathPricer {
let scaling_value = 1.0 / payoff.strike();
AmericanPathPricer {
payoff,
scaling_value,
polynomial_order,
polynomial_type,
exercise_indices: None,
}
}
}
impl EarlyExercisePathPricer<Path> for AmericanPathPricer {
type State = Real;
fn value(&self, path: &Path, t: Size) -> Real {
if self
.exercise_indices
.as_ref()
.is_some_and(|indices| !indices[t])
{
return 0.0;
}
self.payoff.value(self.state(path, t) / self.scaling_value)
}
fn state(&self, path: &Path, t: Size) -> Real {
path[t] * self.scaling_value
}
fn basis_system(&self) -> Vec<Box<dyn Fn(Real) -> Real>> {
let mut v = LsmBasisSystem::path_basis_system(self.polynomial_order, self.polynomial_type);
let payoff = Shared::clone(&self.payoff);
let scaling_value = self.scaling_value;
v.push(Box::new(move |state: Real| {
payoff.value(state / scaling_value)
}));
v
}
}
pub struct MCAmericanEngine<RNG> {
base: McLongstaffSchwartzEngineBase<RNG>,
process: Shared<GeneralizedBlackScholesProcess>,
polynomial_order: Size,
polynomial_type: PolynomialType,
}
impl<RNG: McRngTraits> MCAmericanEngine<RNG> {
#[allow(clippy::too_many_arguments)]
pub fn new(
process: Shared<GeneralizedBlackScholesProcess>,
time_steps: Option<Size>,
time_steps_per_year: Option<Size>,
antithetic_variate: bool,
control_variate: bool,
required_samples: Option<Size>,
required_tolerance: Option<Real>,
max_samples: Option<Size>,
seed: u32,
polynomial_order: Size,
polynomial_type: PolynomialType,
n_calibration_samples: Option<Size>,
antithetic_variate_calibration: Option<bool>,
seed_calibration: Option<u32>,
) -> QlResult<MCAmericanEngine<RNG>> {
require!(
!matches!(
polynomial_type,
PolynomialType::Legendre | PolynomialType::Chebyshev
),
"insufficient polynomial type"
);
let base = McLongstaffSchwartzEngineBase::new(
Shared::clone(&process) as Shared<dyn StochasticProcess1D>,
time_steps,
time_steps_per_year,
false,
antithetic_variate,
control_variate,
required_samples,
required_tolerance,
max_samples,
seed,
n_calibration_samples,
antithetic_variate_calibration,
seed_calibration,
)?;
Ok(MCAmericanEngine {
base,
process,
polynomial_order,
polynomial_type,
})
}
pub fn lsm_base(&self) -> &McLongstaffSchwartzEngineBase<RNG> {
&self.base
}
pub fn polynomial_order(&self) -> Size {
self.polynomial_order
}
pub fn lsm_path_pricer(&self) -> QlResult<Shared<LongstaffSchwartzPathPricer>> {
let arguments = self.base.arguments();
let Some(exercise) = &arguments.exercise else {
fail!("no exercise given");
};
require!(
matches!(
exercise.exercise_type(),
ExerciseType::American | ExerciseType::Bermudan
),
"wrong exercise given"
);
require!(!exercise.payoff_at_expiry(), "payoff at expiry not handled");
let Some(payoff) = &arguments.payoff else {
fail!("no payoff given");
};
require!(
self.polynomial_type != PolynomialType::Chebyshev2nd
|| payoff.option_type() != crate::option::OptionType::Call,
"Chebyshev2nd basis is undefined for in-the-money call states"
);
let grid = self.base.time_grid()?;
let mut early = AmericanPathPricer::new(
Shared::clone(payoff),
self.polynomial_order,
self.polynomial_type,
);
if exercise.exercise_type() == ExerciseType::Bermudan {
let mut indices = vec![false; grid.size()];
for &time in grid.mandatory_times() {
indices[grid.index(time)?] = true;
}
early.exercise_indices = Some(indices);
}
Ok(shared(LongstaffSchwartzPathPricer::new(
&grid,
shared(early) as Shared<dyn EarlyExercisePathPricer<Path, State = Real>>,
&self.process.risk_free_rate(),
)?))
}
}
impl<RNG: McRngTraits> MCAmericanEngine<RNG> {
pub fn control_path_pricer(&self) -> QlResult<EuropeanPathPricer> {
let Some(payoff) = &self.base.arguments().payoff else {
fail!("no payoff given");
};
let grid = self.base.time_grid()?;
let Some(last_time) = grid.back() else {
fail!("empty time grid");
};
let discount = self
.process
.risk_free_rate()
.current_link()?
.discount(last_time, false)?;
EuropeanPathPricer::new(payoff.option_type(), payoff.strike(), discount)
}
pub fn control_variate_value(&self) -> QlResult<Real> {
let arguments = self.base.arguments();
let Some(exercise) = &arguments.exercise else {
fail!("no exercise given");
};
let mut engine = AnalyticEuropeanEngine::new(Shared::clone(&self.process));
let Some(control_arguments) =
(engine.arguments_mut() as &mut dyn Any).downcast_mut::<OptionArguments>()
else {
fail!("engine is using inconsistent arguments");
};
control_arguments.payoff = arguments.payoff.clone();
control_arguments.exercise =
Some(shared(EuropeanExercise::new(exercise.last_date())) as Shared<dyn Exercise>);
engine.calculate()?;
let Some(results) = (engine.results() as &dyn Any).downcast_ref::<OneAssetOptionResults>()
else {
fail!("engine is using inconsistent results");
};
let Some(value) = results.instrument.value else {
fail!("control variation pricing engine returned no value");
};
Ok(value)
}
}
impl<RNG: McRngTraits> AsObservable for MCAmericanEngine<RNG> {
fn observable(&self) -> &Observable {
self.base.observable()
}
}
impl<RNG: McRngTraits> PricingEngine for MCAmericanEngine<RNG> {
fn arguments_mut(&mut self) -> &mut dyn Arguments {
self.base.arguments_mut()
}
fn results(&self) -> &dyn Results {
self.base.results()
}
fn reset(&mut self) {
self.base.reset();
}
fn calculate(&mut self) -> QlResult<()> {
let pricer = self.lsm_path_pricer()?;
let control_variate = if self.base.control_variate() {
Some(ControlVariate {
value: self.control_variate_value()?,
path_pricer: Box::new(self.control_path_pricer()?),
})
} else {
None
};
self.base.calculate_with(pricer, control_variate)?;
if self.base.control_variate() {
let instrument = &mut self.base.results_mut().instrument;
instrument.value = instrument.value.map(|value| value.max(0.0));
}
Ok(())
}
}
pub struct MakeMcAmericanEngine<RNG> {
process: Shared<GeneralizedBlackScholesProcess>,
steps: Option<Size>,
steps_per_year: Option<Size>,
samples: Option<Size>,
max_samples: Option<Size>,
tolerance: Option<Real>,
antithetic: bool,
control_variate: bool,
seed: u32,
polynomial_order: Size,
calibration_samples: Option<Size>,
polynomial_type: PolynomialType,
_rng: std::marker::PhantomData<RNG>,
}
impl<RNG: McRngTraits> MakeMcAmericanEngine<RNG> {
pub fn new(process: Shared<GeneralizedBlackScholesProcess>) -> MakeMcAmericanEngine<RNG> {
MakeMcAmericanEngine {
process,
steps: None,
steps_per_year: None,
samples: None,
max_samples: None,
tolerance: None,
antithetic: false,
control_variate: false,
seed: 0,
polynomial_order: 2,
calibration_samples: None,
polynomial_type: PolynomialType::Monomial,
_rng: std::marker::PhantomData,
}
}
#[must_use]
pub fn with_steps(mut self, steps: Size) -> Self {
self.steps = Some(steps);
self
}
#[must_use]
pub fn with_steps_per_year(mut self, steps: Size) -> Self {
self.steps_per_year = Some(steps);
self
}
#[must_use]
pub fn with_samples(mut self, samples: Size) -> Self {
self.samples = Some(samples);
self
}
#[must_use]
pub fn with_absolute_tolerance(mut self, tolerance: Real) -> Self {
self.tolerance = Some(tolerance);
self
}
#[must_use]
pub fn with_max_samples(mut self, samples: Size) -> Self {
self.max_samples = Some(samples);
self
}
#[must_use]
pub fn with_seed(mut self, seed: u32) -> Self {
self.seed = seed;
self
}
#[must_use]
pub fn with_antithetic_variate(mut self, antithetic: bool) -> Self {
self.antithetic = antithetic;
self
}
#[must_use]
pub fn with_control_variate(mut self, control_variate: bool) -> Self {
self.control_variate = control_variate;
self
}
#[must_use]
pub fn with_polynomial_order(mut self, order: Size) -> Self {
self.polynomial_order = order;
self
}
#[must_use]
pub fn with_basis_system(mut self, family: PolynomialType) -> Self {
self.polynomial_type = family;
self
}
#[must_use]
pub fn with_calibration_samples(mut self, samples: Size) -> Self {
self.calibration_samples = Some(samples);
self
}
pub fn build(self) -> QlResult<MCAmericanEngine<RNG>> {
require!(
self.steps.is_some() || self.steps_per_year.is_some(),
"number of steps not given"
);
require!(
self.steps.is_none() || self.steps_per_year.is_none(),
"number of steps overspecified"
);
require!(
!(self.samples.is_some() && self.tolerance.is_some()),
"number of samples already set"
);
if self.tolerance.is_some() {
require!(
RNG::ALLOWS_ERROR_ESTIMATE,
"chosen random generator policy does not allow an error estimate"
);
}
MCAmericanEngine::new(
self.process,
self.steps,
self.steps_per_year,
self.antithetic,
self.control_variate,
self.samples,
self.tolerance,
self.max_samples,
self.seed,
self.polynomial_order,
self.polynomial_type,
self.calibration_samples,
None,
None,
)
}
}
#[cfg(test)]
mod tests {
use std::any::Any;
use super::*;
use crate::exercise::{AmericanExercise, EuropeanExercise, Exercise};
use crate::instruments::{OptionArguments, PlainVanillaPayoff};
use crate::math::array::Array;
use crate::math::randomnumbers::rngtraits::PseudoRandom;
use crate::math::timegrid::TimeGrid;
use crate::option::OptionType;
use crate::pricingengines::vanilla::test_market::{Market, market, time_to_days, today};
use crate::time::date::Date;
const STRIKE: Real = 40.0;
const TOL: Real = 1e-14;
fn american(expiry: Date) -> Shared<dyn Exercise> {
shared(AmericanExercise::over(today(), expiry).unwrap()) as Shared<dyn Exercise>
}
fn put_pricer(order: Size) -> AmericanPathPricer {
AmericanPathPricer::new(
shared(PlainVanillaPayoff::new(OptionType::Put, STRIKE))
as Shared<dyn StrikedTypePayoff>,
order,
PolynomialType::Monomial,
)
}
fn path(values: [Real; 3]) -> Path {
Path::new(TimeGrid::new(1.0, 2).unwrap(), Array::from(values)).unwrap()
}
fn flat_market() -> Market {
let market = market();
market.set(36.0, 0.0, 0.06, 0.20);
market
}
fn engine(market: &Market) -> MCAmericanEngine<PseudoRandom> {
MakeMcAmericanEngine::<PseudoRandom>::new(Shared::clone(&market.process))
.with_steps(4)
.with_samples(256)
.with_seed(42)
.build()
.unwrap()
}
fn set_option(engine: &mut MCAmericanEngine<PseudoRandom>, exercise: Shared<dyn Exercise>) {
let args = (engine.arguments_mut() as &mut dyn Any)
.downcast_mut::<OptionArguments>()
.unwrap();
args.payoff = Some(shared(PlainVanillaPayoff::new(OptionType::Put, STRIKE))
as Shared<dyn StrikedTypePayoff>);
args.exercise = Some(exercise);
}
#[test]
fn the_state_is_scaled_and_the_value_is_not() {
let pricer = put_pricer(3);
let p = path([36.0, 44.0, 38.0]);
assert!((pricer.state(&p, 0) - 0.9).abs() < TOL);
assert!((pricer.state(&p, 1) - 1.1).abs() < TOL);
assert!((pricer.value(&p, 0) - 4.0).abs() < TOL);
assert_eq!(pricer.value(&p, 1), 0.0, "an OTM put is worth 0");
assert!((pricer.value(&p, 2) - 2.0).abs() < TOL);
}
#[test]
fn the_exercise_value_round_trips_through_the_scaling() {
let pricer = AmericanPathPricer::new(
shared(PlainVanillaPayoff::new(OptionType::Put, 3.0)) as Shared<dyn StrikedTypePayoff>,
2,
PolynomialType::Monomial,
);
let p = path([0.89, 0.89, 0.89]);
assert_eq!(pricer.value(&p, 0), 2.1100000000000003);
assert_ne!(pricer.value(&p, 0), 3.0 - 0.89);
}
#[test]
fn the_payoff_is_the_extra_basis_function() {
for order in [2, 3] {
let basis = put_pricer(order).basis_system();
assert_eq!(
basis.len(),
order + 2,
"order {order}: monomials plus payoff"
);
assert!((basis[0](0.9) - 1.0).abs() < TOL);
assert!((basis[1](0.9) - 0.9).abs() < TOL);
assert!(
(basis[order + 1](0.9) - 4.0).abs() < TOL,
"the appended function must unscale before applying the payoff"
);
assert_eq!(
basis[order + 1](1.1),
0.0,
"and must stay a payoff, not a monomial"
);
}
}
#[test]
fn only_a_plain_american_exercise_builds_a_path_pricer() {
let market = flat_market();
let expiry = today() + time_to_days(1.0);
let mut accepted = engine(&market);
set_option(&mut accepted, american(expiry));
assert!(accepted.lsm_path_pricer().is_ok());
let mut rejected = engine(&market);
set_option(
&mut rejected,
shared(EuropeanExercise::new(expiry)) as Shared<dyn Exercise>,
);
assert_eq!(
rejected.lsm_path_pricer().err().unwrap().message(),
"wrong exercise given"
);
let mut at_expiry = engine(&market);
set_option(
&mut at_expiry,
shared(AmericanExercise::new(today(), expiry, true).unwrap()) as Shared<dyn Exercise>,
);
assert_eq!(
at_expiry.lsm_path_pricer().err().unwrap().message(),
"payoff at expiry not handled"
);
}
#[test]
fn a_missing_exercise_is_rejected() {
let market = flat_market();
let mut bare = engine(&market);
assert_eq!(
bare.lsm_path_pricer().err().unwrap().message(),
"no exercise given"
);
let args = (bare.arguments_mut() as &mut dyn Any)
.downcast_mut::<OptionArguments>()
.unwrap();
args.exercise = Some(american(today() + time_to_days(1.0)));
assert_eq!(
bare.lsm_path_pricer().err().unwrap().message(),
"no payoff given"
);
}
#[test]
fn the_builder_defaults_match_the_cpp_factory() {
let market = flat_market();
let built = engine(&market);
assert_eq!(built.polynomial_order(), 2);
assert_eq!(built.lsm_base().n_calibration_samples(), 2048);
assert!(!built.lsm_base().antithetic_variate_calibration());
let overridden = MakeMcAmericanEngine::<PseudoRandom>::new(Shared::clone(&market.process))
.with_steps(4)
.with_samples(256)
.with_polynomial_order(4)
.with_calibration_samples(64)
.build()
.unwrap();
assert_eq!(overridden.polynomial_order(), 4);
assert_eq!(overridden.lsm_base().n_calibration_samples(), 64);
}
#[test]
fn the_builder_validates_its_named_parameters() {
let market = flat_market();
let maker = || MakeMcAmericanEngine::<PseudoRandom>::new(Shared::clone(&market.process));
assert!(maker().with_samples(256).build().is_err(), "no steps given");
assert!(
maker()
.with_steps(4)
.with_steps_per_year(50)
.with_samples(256)
.build()
.is_err(),
"steps overspecified"
);
assert!(
maker()
.with_steps(4)
.with_samples(256)
.with_absolute_tolerance(0.02)
.build()
.is_err(),
"samples and tolerance are exclusive"
);
let antithetic = maker()
.with_steps(4)
.with_samples(256)
.with_antithetic_variate(true)
.build()
.unwrap();
assert!(
antithetic.lsm_base().antithetic_variate_calibration(),
"the pricing flag must reach the engine, not a hardcoded false"
);
assert!(
maker()
.with_steps(4)
.with_absolute_tolerance(0.02)
.with_max_samples(4_096)
.build()
.is_ok()
);
}
#[test]
fn bermudan_grid_retains_only_positive_mandatory_dates() {
use crate::exercise::BermudanExercise;
let market = flat_market();
let mut engine = engine(&market);
set_option(
&mut engine,
shared(
BermudanExercise::new(
vec![
today() - 3,
today(),
today() + 91,
today() + 91,
today() + 203,
today() + 365,
],
false,
)
.unwrap(),
),
);
let grid = engine.lsm_base().time_grid().unwrap();
assert_eq!(
grid.mandatory_times(),
&[91.0 / 360.0, 203.0 / 360.0, 365.0 / 360.0]
);
assert_eq!(grid.size(), 5);
assert_eq!(grid[1], 91.0 / 360.0);
assert_eq!(grid[2], 203.0 / 360.0);
assert!(engine.lsm_path_pricer().is_ok());
set_option(
&mut engine,
shared(BermudanExercise::new(vec![today()], false).unwrap()),
);
assert!(
engine
.lsm_path_pricer()
.err()
.unwrap()
.message()
.contains("no positive exercise time")
);
set_option(
&mut engine,
shared(BermudanExercise::new(vec![today() + 365], true).unwrap()),
);
assert_eq!(
engine.lsm_path_pricer().err().unwrap().message(),
"payoff at expiry not handled"
);
}
#[test]
fn bermudan_mask_forbids_exercise_between_contract_dates_in_both_passes() {
use crate::methods::montecarlo::PathPricer;
let market = flat_market();
let grid = TimeGrid::new(1.0, 2).unwrap();
for (allowed, expected_positive) in [(false, false), (true, true)] {
let mut early = put_pricer(2);
early.exercise_indices = Some(vec![false, allowed, true]);
let lsm = LongstaffSchwartzPathPricer::new(
&grid,
shared(early),
&market.process.risk_free_rate(),
)
.unwrap();
let sample = path([40.0, 1.0, 40.0]);
assert_eq!(lsm.price(&sample), 0.0);
lsm.calibrate().unwrap();
assert_eq!(lsm.price(&sample) > 0.0, expected_positive);
assert_eq!(
lsm.exercise_probability().unwrap(),
if allowed { 1.0 } else { 0.0 }
);
}
}
#[test]
fn unsupported_engine_bases_are_rejected_before_sampling() {
let market = flat_market();
for family in [PolynomialType::Legendre, PolynomialType::Chebyshev] {
let error = MakeMcAmericanEngine::<PseudoRandom>::new(Shared::clone(&market.process))
.with_steps(4)
.with_samples(16)
.with_basis_system(family)
.build()
.err()
.unwrap();
assert_eq!(error.message(), "insufficient polynomial type");
}
}
}
#[cfg(test)]
mod oracle {
use std::any::Any;
use super::{MCAmericanEngine, MakeMcAmericanEngine};
use crate::exercise::{AmericanExercise, EuropeanExercise, Exercise};
use crate::handle::Handle;
use crate::instrument::Instrument;
use crate::instruments::{OneAssetOptionResults, OptionArguments};
use crate::instruments::{PlainVanillaPayoff, StrikedTypePayoff, VanillaOption};
use crate::interestrate::Compounding;
use crate::math::randomnumbers::rngtraits::McRngTraits;
use crate::math::randomnumbers::rngtraits::PseudoRandom;
use crate::math::statistics::{GeneralStatistics, MeanStdDev, Statistics};
use crate::methods::montecarlo::{MonteCarloModel, PathGenerator};
use crate::option::OptionType;
use crate::pricingengine::PricingEngine;
use crate::pricingengines::vanilla::AnalyticEuropeanEngine;
use crate::processes::GeneralizedBlackScholesProcess;
use crate::quotes::{Quote, SimpleQuote};
use crate::settings::Settings;
use crate::shared::{Shared, SharedMut, shared, shared_mut};
use crate::stochasticprocess::StochasticProcess1D;
use crate::termstructures::volatility::{BlackConstantVol, BlackVolTermStructure};
use crate::termstructures::yields::FlatForward;
use crate::termstructures::yieldtermstructure::YieldTermStructure;
use crate::time::date::{Date, Month};
use crate::time::daycounters::actual365fixed::Actual365Fixed;
use crate::time::frequency::Frequency;
use crate::types::Size;
use crate::types::{Rate, Real, Volatility};
const UNDERLYING: Real = 36.0;
const DIVIDEND_YIELD: Rate = 0.0;
const RISK_FREE_RATE: Rate = 0.06;
const VOLATILITY: Volatility = 0.20;
const QUANTLIB_MC_VALUE: Real = 2.054422273006143;
const EXPECTED_EXERCISE_PROBABILITY: Real = 0.48013;
fn todays_date() -> Date {
Date::new(15, Month::May, 1998)
}
fn settlement_date() -> Date {
Date::new(17, Month::May, 1998)
}
fn maturity() -> Date {
Date::new(17, Month::May, 1999)
}
fn flat_curve(rate: Rate) -> Handle<dyn YieldTermStructure> {
Handle::new(shared(FlatForward::with_rate(
settlement_date(),
rate,
Actual365Fixed::new(),
Compounding::Continuous,
Frequency::Annual,
)) as Shared<dyn YieldTermStructure>)
}
fn process() -> Shared<GeneralizedBlackScholesProcess> {
process_with_vol(VOLATILITY)
}
fn process_with_vol(volatility: Real) -> Shared<GeneralizedBlackScholesProcess> {
let spot = Handle::new(shared(SimpleQuote::new(UNDERLYING)) as Shared<dyn Quote>);
let vol = Handle::new(shared(BlackConstantVol::new(
settlement_date(),
None,
volatility,
Actual365Fixed::new(),
)) as Shared<dyn BlackVolTermStructure>);
shared(GeneralizedBlackScholesProcess::new(
spot,
flat_curve(DIVIDEND_YIELD),
flat_curve(RISK_FREE_RATE),
vol,
))
}
fn payoff() -> Shared<dyn StrikedTypePayoff> {
shared(PlainVanillaPayoff::new(OptionType::Put, UNDERLYING))
as Shared<dyn StrikedTypePayoff>
}
#[test]
fn american_put_reproduces_the_quantlib_price_and_exercise_probability() {
let settings = shared(Settings::new());
settings.set_evaluation_date(todays_date());
let process = process();
let exercise = shared(AmericanExercise::over(settlement_date(), maturity()).unwrap())
as Shared<dyn Exercise>;
let mut american = VanillaOption::new(payoff(), exercise, Shared::clone(&settings));
american.base_mut().set_pricing_engine(shared_mut(
MakeMcAmericanEngine::<PseudoRandom>::new(Shared::clone(&process))
.with_steps(75)
.with_antithetic_variate(true)
.with_absolute_tolerance(0.02)
.with_seed(42)
.with_polynomial_order(3)
.build()
.unwrap(),
) as SharedMut<dyn PricingEngine>);
let calculated = american.npv().unwrap();
let error_estimate = american.error_estimate().unwrap();
let exercise_probability = american.result::<Real>("exerciseProbability").unwrap();
let mut european = VanillaOption::new(
payoff(),
shared(EuropeanExercise::new(maturity())) as Shared<dyn Exercise>,
settings,
);
european.base_mut().set_pricing_engine(
shared_mut(AnalyticEuropeanEngine::new(process)) as SharedMut<dyn PricingEngine>
);
let european_value = european.npv().unwrap();
assert!(
(calculated - QUANTLIB_MC_VALUE).abs() < 2.34 * error_estimate,
"price {calculated} +/- {error_estimate} misses QuantLib's \
{QUANTLIB_MC_VALUE} by {}, band {}",
(calculated - QUANTLIB_MC_VALUE).abs(),
2.34 * error_estimate
);
assert!(
(exercise_probability - EXPECTED_EXERCISE_PROBABILITY).abs() < 0.015,
"exercise probability {exercise_probability} vs \
{EXPECTED_EXERCISE_PROBABILITY}"
);
assert!(
calculated > european_value,
"an American put {calculated} cannot be worth less than its \
European twin {european_value}"
);
}
#[test]
fn five_remaining_american_cases_preserve_upstream_acceptance() {
let references = [
(36.0, 0.2 + 0.1, 3.430219464343261, 0.51678),
(36.0, 0.4, 4.785433883245219, 0.54598),
(40.0, 0.2, 4.485350467454564, 0.75549),
(40.0, 0.2 + 0.1, 5.736716960463202, 0.67569),
(40.0, 0.4, 7.107644070425648, 0.65562),
];
for (strike, volatility, fd, probability) in references {
let settings = shared(Settings::new());
settings.set_evaluation_date(todays_date());
let process = process_with_vol(volatility);
let exercise = shared(AmericanExercise::until(maturity(), false).unwrap());
let mut option = VanillaOption::new(
shared(PlainVanillaPayoff::new(OptionType::Put, strike)),
exercise,
settings,
);
option.base_mut().set_pricing_engine(shared_mut(
MakeMcAmericanEngine::<PseudoRandom>::new(process)
.with_steps(75)
.with_antithetic_variate(true)
.with_absolute_tolerance(0.02)
.with_seed(42)
.with_polynomial_order(3)
.build()
.unwrap(),
));
let value = option.npv().unwrap();
let error = option.error_estimate().unwrap();
assert!(
error.is_finite() && error > 0.0,
"invalid standard error {error}"
);
let exercise = option.result::<Real>("exerciseProbability").unwrap();
assert!(
(value - fd).abs() <= 2.34 * error,
"strike={strike}, vol={volatility}: {value} +/- {error} vs FD {fd}"
);
assert!(
(exercise - probability).abs() <= 0.015,
"strike={strike}, vol={volatility}: exercise {exercise} vs {probability}"
);
}
}
#[test]
fn bermudan_exercise_only_and_refined_grids_match_independent_oracles() {
use crate::exercise::BermudanExercise;
for (steps, expected) in [(3, 4.331587683099873), (75, 4.334758948830429)] {
let settings = shared(Settings::new());
settings.set_evaluation_date(todays_date());
let exercise = shared(
BermudanExercise::new(
vec![settlement_date() + 91, settlement_date() + 203, maturity()],
false,
)
.unwrap(),
);
let mut option = VanillaOption::new(
shared(PlainVanillaPayoff::new(OptionType::Put, 40.0)),
exercise,
settings,
);
option.base_mut().set_pricing_engine(shared_mut(
MakeMcAmericanEngine::<PseudoRandom>::new(process())
.with_steps(steps)
.with_samples(32768)
.with_calibration_samples(8192)
.with_seed(42)
.with_antithetic_variate(true)
.build()
.unwrap(),
));
let value = option.npv().unwrap();
let error = option.error_estimate().unwrap();
assert!(
error.is_finite() && error > 0.0,
"invalid standard error {error}"
);
assert!(
(value - expected).abs() <= 2.34 * error,
"steps={steps}: {value} +/- {error} vs {expected}"
);
}
}
const CV_SAMPLES: Size = 4096;
fn cv_engine(
strike: Real,
control_variate: bool,
samples: Size,
seed: u32,
) -> MCAmericanEngine<PseudoRandom> {
cv_engine_with(strike, true, control_variate, samples, seed)
}
fn cv_engine_with(
strike: Real,
antithetic_variate: bool,
control_variate: bool,
samples: Size,
seed: u32,
) -> MCAmericanEngine<PseudoRandom> {
let mut engine = MakeMcAmericanEngine::<PseudoRandom>::new(process())
.with_steps(75)
.with_antithetic_variate(antithetic_variate)
.with_control_variate(control_variate)
.with_samples(samples)
.with_seed(seed)
.with_polynomial_order(3)
.build()
.unwrap();
let args = (engine.arguments_mut() as &mut dyn Any)
.downcast_mut::<OptionArguments>()
.unwrap();
args.payoff = Some(shared(PlainVanillaPayoff::new(OptionType::Put, strike))
as Shared<dyn StrikedTypePayoff>);
args.exercise = Some(
shared(AmericanExercise::over(settlement_date(), maturity()).unwrap())
as Shared<dyn Exercise>,
);
engine.calculate().unwrap();
engine
}
fn outcome(engine: &MCAmericanEngine<PseudoRandom>) -> (Real, Real, Real) {
let results = (engine.results() as &dyn Any)
.downcast_ref::<OneAssetOptionResults>()
.unwrap();
let probability = *results.instrument.additional_results["exerciseProbability"]
.downcast_ref::<Real>()
.unwrap();
(
results.instrument.value.unwrap(),
results.instrument.error_estimate.unwrap(),
probability,
)
}
fn european_over_the_pricing_paths(
engine: &MCAmericanEngine<PseudoRandom>,
samples: Size,
seed: u32,
) -> (Real, Real) {
let grid = engine.lsm_base().time_grid().unwrap();
let generator = PseudoRandom::make_sequence_generator(grid.size() - 1, seed).unwrap();
let paths = PathGenerator::from_time_grid(
process() as Shared<dyn StochasticProcess1D>,
grid,
generator,
false,
)
.unwrap();
let mut model = MonteCarloModel::new(
paths,
engine.control_path_pricer().unwrap(),
GeneralStatistics::new(),
true,
)
.unwrap();
model.add_samples(samples).unwrap();
let accumulator = model.sample_accumulator();
(
accumulator.mean().unwrap(),
accumulator.error_estimate().unwrap(),
)
}
#[test]
fn control_variate_shifts_by_the_european_error_on_the_same_paths() {
const SEED: u32 = 42;
let plain = cv_engine(UNDERLYING, false, CV_SAMPLES, SEED);
let controlled = cv_engine(UNDERLYING, true, CV_SAMPLES, SEED);
let (v0, e0, p0) = outcome(&plain);
let (v1, e1, p1) = outcome(&controlled);
let analytic = controlled.control_variate_value().unwrap();
let (mean, mean_error) = european_over_the_pricing_paths(&controlled, CV_SAMPLES, SEED);
assert!(
((v1 - v0) - (analytic - mean)).abs() < 1e-12,
"cv {v1} - plain {v0} != analytic {analytic} - european mean {mean}"
);
assert!((analytic - mean).abs() > 1e-4, "the shift must be material");
assert_eq!(p1, p0, "the control must not touch the exercise decision");
assert!(e0.is_finite() && e1.is_finite() && e1 > 0.0);
let settings = shared(Settings::new());
settings.set_evaluation_date(todays_date());
let mut european = VanillaOption::new(
payoff(),
shared(EuropeanExercise::new(maturity())) as Shared<dyn Exercise>,
settings,
);
european.base_mut().set_pricing_engine(
shared_mut(AnalyticEuropeanEngine::new(process())) as SharedMut<dyn PricingEngine>
);
assert_eq!(analytic, european.npv().unwrap());
assert!(
(mean - analytic).abs() < 4.0 * mean_error,
"control pricer mean {mean} +/- {mean_error} vs analytic {analytic}"
);
}
#[test]
fn control_variate_reduces_the_error_where_the_twin_tracks_the_value() {
for (strike, antithetic) in [(32.0, true), (UNDERLYING, false)] {
let (_, plain, _) = outcome(&cv_engine_with(strike, antithetic, false, CV_SAMPLES, 42));
let (_, controlled, _) =
outcome(&cv_engine_with(strike, antithetic, true, CV_SAMPLES, 42));
assert!(
controlled < plain,
"K={strike}, antithetic={antithetic}: cv error {controlled} \
must fall below plain error {plain}"
);
}
}
#[test]
fn a_negative_control_variate_estimate_is_floored_at_zero() {
const STRIKE: Real = 22.0;
const SAMPLES: Size = 64;
const SEED: u32 = 16;
let plain = cv_engine(STRIKE, false, SAMPLES, SEED);
let controlled = cv_engine(STRIKE, true, SAMPLES, SEED);
let analytic = controlled.control_variate_value().unwrap();
let (mean, _) = european_over_the_pricing_paths(&controlled, SAMPLES, SEED);
let raw = outcome(&plain).0 + analytic - mean;
assert!(
raw < -0.01,
"fixture must produce a negative raw estimate: {raw}"
);
assert_eq!(outcome(&controlled).0, 0.0);
}
}