|
| 1 | +/* |
| 2 | + Copyright (C) 2024 Quaternion Risk Management Ltd |
| 3 | + All rights reserved. |
| 4 | +
|
| 5 | + This file is part of ORE, a free-software/open-source library |
| 6 | + for transparent pricing and risk analysis - http://opensourcerisk.org |
| 7 | +
|
| 8 | + ORE is free software: you can redistribute it and/or modify it |
| 9 | + under the terms of the Modified BSD License. You should have received a |
| 10 | + copy of the license along with this program. |
| 11 | + The license is also available online at <http://opensourcerisk.org> |
| 12 | +
|
| 13 | + This program is distributed on the basis that it will form a useful |
| 14 | + contribution to risk analytics and model standardisation, but WITHOUT |
| 15 | + ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or |
| 16 | + FITNESS FOR A PARTICULAR PURPOSE. See the license for more details. |
| 17 | +*/ |
| 18 | + |
| 19 | +#include <qle/methods/fdmlgmop.hpp> |
| 20 | + |
| 21 | +#include <ql/methods/finitedifferences/meshers/fdmmesher.hpp> |
| 22 | +#include <ql/methods/finitedifferences/operators/fdmlinearoplayout.hpp> |
| 23 | +#include <ql/methods/finitedifferences/operators/secondderivativeop.hpp> |
| 24 | + |
| 25 | +namespace QuantExt { |
| 26 | + |
| 27 | +FdmLgmOp::FdmLgmOp(const ext::shared_ptr<FdmMesher>& mesher, const ext::shared_ptr<StochasticProcess1D>& process) |
| 28 | + : mesher_(mesher), process_(process), dxMap_(FirstDerivativeOp(0, mesher)), dxxMap_(SecondDerivativeOp(0, mesher)), |
| 29 | + mapT_(0, mesher) {} |
| 30 | + |
| 31 | +void FdmLgmOp::setTime(Time t1, Time t2) { |
| 32 | + Real v = process_->variance(t1, 0.0, t2 - t1) / (t2 - t1); |
| 33 | + mapT_.axpyb(Array(), dxMap_, dxxMap_.mult(0.5 * Array(mesher_->layout()->size(), v)), Array(1, 0)); |
| 34 | +} |
| 35 | + |
| 36 | +Size FdmLgmOp::size() const { return 1u; } |
| 37 | + |
| 38 | +Array FdmLgmOp::apply(const Array& u) const { return mapT_.apply(u); } |
| 39 | + |
| 40 | +Array FdmLgmOp::apply_direction(Size direction, const Array& r) const { |
| 41 | + if (direction == 0) |
| 42 | + return mapT_.apply(r); |
| 43 | + else { |
| 44 | + Array retVal(r.size(), 0.0); |
| 45 | + return retVal; |
| 46 | + } |
| 47 | +} |
| 48 | + |
| 49 | +Array FdmLgmOp::apply_mixed(const Array& r) const { |
| 50 | + Array retVal(r.size(), 0.0); |
| 51 | + return retVal; |
| 52 | +} |
| 53 | + |
| 54 | +Array FdmLgmOp::solve_splitting(Size direction, const Array& r, Real dt) const { |
| 55 | + if (direction == 0) |
| 56 | + return mapT_.solve_splitting(r, dt, 1.0); |
| 57 | + else { |
| 58 | + Array retVal(r); |
| 59 | + return retVal; |
| 60 | + } |
| 61 | +} |
| 62 | + |
| 63 | +Array FdmLgmOp::preconditioner(const Array& r, Real dt) const { return solve_splitting(0, r, dt); } |
| 64 | + |
| 65 | +#if !defined(QL_NO_UBLAS_SUPPORT) |
| 66 | +std::vector<QuantLib::SparseMatrix> FdmLgmOp::toMatrixDecomp() const { |
| 67 | + std::vector<QuantLib::SparseMatrix> retVal(1, mapT_.toMatrix()); |
| 68 | + return retVal; |
| 69 | +} |
| 70 | +#endif |
| 71 | +} // namespace QuantExt |
0 commit comments