final double theta = 0.5;
final double yL = Math.log(SPOT / 6);
final double yH = Math.log(6 * SPOT);
final ConvectionDiffusionPDESolver solver = new ThetaMethodFiniteDifference(theta, false);
final BoundaryCondition lower1 = new NeumannBoundaryCondition(1.0, yL, true);
final BoundaryCondition upper1 = new NeumannBoundaryCondition(1.0, yH, false);
final MeshingFunction timeMesh1 = new ExponentialMeshing(0, EXPIRY - DIVIDEND_DATE - 1e-6, 50, 0.0);
final MeshingFunction timeMesh2 = new ExponentialMeshing(EXPIRY - DIVIDEND_DATE + 1e-6, EXPIRY, 50, 0.0);
final MeshingFunction spaceMesh = new ExponentialMeshing(yL, yH, 101, 0.0);
final PDEGrid1D grid1 = new PDEGrid1D(timeMesh1, spaceMesh);
final double[] sNodes1 = grid1.getSpaceNodes();
//run the PDE solver backward to the dividend date
final PDE1DDataBundle<ConvectionDiffusionPDE1DCoefficients> db1 = new PDE1DDataBundle<>(pde, initialCon, lower1, upper1, grid1);
final PDETerminalResults1D res1 = (PDETerminalResults1D) solver.solve(db1);
//Map the spot nodes after (in calendar time) the dividend payment to nodes before
final int nSNodes = sNodes1.length;
final double[] sNodes2 = new double[nSNodes];
final double lnBeta = Math.log(1 - BETA);
for (int i = 0; i < nSNodes; i++) {
final double temp = sNodes1[i];
if (temp < 0) {
sNodes2[i] = Math.log(Math.exp(temp) + ALPHA) - lnBeta;
}
else {
sNodes2[i] = temp + Math.log(1 + ALPHA * Math.exp(-temp)) - lnBeta;
}
}
final PDEGrid1D grid2 = new PDEGrid1D(timeMesh2.getPoints(), sNodes2);
final BoundaryCondition lower2 = new NeumannBoundaryCondition(1.0, sNodes2[0], true);
final BoundaryCondition upper2 = new NeumannBoundaryCondition(1.0, sNodes2[nSNodes - 1], false);
//run the PDE solver backward from the dividend date to zero
final PDE1DDataBundle<ConvectionDiffusionPDE1DCoefficients> db2 = new PDE1DDataBundle<>(pde, res1.getTerminalResults(), lower2, upper2, grid2);
final PDETerminalResults1D res2 = (PDETerminalResults1D) solver.solve(db2);