34 const ext::shared_ptr<FdmMesher>& mesher)
35 :
m_(mesher->layout()->size(), mesher->layout()->size()) {
37 const Integer hPoints = nPoints/2;
38 const bool isEven = (nPoints == 2*hPoints);
40 Array xValues = mesher->locations(direction);
41 std::set<Real> tmp(xValues.
begin(), xValues.
end());
42 xValues =
Array(tmp.begin(), tmp.end());
44 const Integer nx(mesher->layout()->dim()[direction]);
47 "inconsistent set of grid values in direction " << direction);
50 "inconsistent number of points");
52 Array xOffsets(nPoints);
53 const ext::function<
Real(
Real)> emptyFct;
55 for (
const auto& iter : *mesher->layout()) {
56 const auto ix =
Integer(iter.coordinates()[direction]);
57 const Integer offset = std::max(0, hPoints - ix)
58 - std::max(0, hPoints - (nx-((isEven)? 0 : 1) - ix));
60 const Integer ilx = ix - hPoints + offset;
62 for (
Integer j=0; j < nPoints; ++j) {
64 xOffsets[j] = xValues[idx] - xValues[ix];
70 const Size i = iter.index();
71 for (
Integer j=0; j < nPoints; ++j) {
72 const Size k = mesher->layout()->neighbourhood(iter, direction, ilx - ix + j);
74 m_(i, k) = weights[j];
1-D array used in linear algebra.
const_iterator end() const
Size size() const
dimension of the array
const_iterator begin() const
NthOrderDerivativeOp(Size direction, Size order, Integer nPoints, const ext::shared_ptr< FdmMesher > &mesher)
array_type apply(const array_type &r) const override
SparseMatrix toMatrix() const override
Numerical Differentiation on arbitrarily spaced grids.
const Array & weights() const
#define QL_REQUIRE(condition, message)
throw an error if the given pre-condition is not verified
memory layout of a fdm linear operator
QL_INTEGER Integer
integer number
std::size_t Size
size of a container
Array prod(const SparseMatrix &A, const Array &x)
boost::numeric::ublas::compressed_matrix< Real > SparseMatrix
n-th order derivative linear operator
numerical differentiation of arbitrary order and on irregular grids
ext::shared_ptr< YieldTermStructure > r