@@ -35,7 +35,8 @@ namespace QuantLib {
3535 Root (Handle<DefaultProbabilityTermStructure> dts, Real pd)
3636 : dts_(std::move(dts)), pd_(pd) {}
3737 Real operator ()(Real t) const {
38- QL_REQUIRE (t >= 0.0 , " t < 0" );
38+ QL_REQUIRE (t >= 0.0 , " GaussianRandomDefaultModel: internal error, t < 0 ("
39+ << t << " ) during root searching." );
3940 return dts_->defaultProbability (t, true ) - pd_;
4041 }
4142 private:
@@ -73,17 +74,20 @@ namespace QuantLib {
7374 Real p = CumulativeNormalDistribution ()(y);
7475
7576 if (dts->defaultProbability (tmax) < p)
76- pool_->setTime (name, tmax+ 1 );
77+ pool_->setTime (name, tmax + 1 );
7778 else {
78- try {
79- Brent brent;
80- brent.setLowerBound (0.0 );
81- brent.setUpperBound (tmax);
82- pool_->setTime (name, brent.solve (Root (dts,p),accuracy_, tmax/2.0 , 1 ));
83- }
84- catch (...){
85- pool_->setTime (name, Bisection ().solve (Root (dts,p), accuracy_, tmax/2.0 , 0.0 , tmax));
86- }
79+ // we know there is a zero of f(t) = dts->defaultProbability(t) - p in [0, tmax]
80+ try {
81+ // try bracketing the root and find it with Brent
82+ Brent brent;
83+ brent.setLowerBound (0.0 );
84+ brent.setUpperBound (tmax);
85+ pool_->setTime (name, brent.solve (Root (dts, p), accuracy_, tmax / 2.0 , 1.0 ));
86+ } catch (...) {
87+ // if Brent fails, use Bisection, this is guaranteed to find the root
88+ pool_->setTime (
89+ name, Bisection ().solve (Root (dts, p), accuracy_, tmax / 2.0 , 0.0 , tmax));
90+ }
8791 }
8892 }
8993 }
0 commit comments