Driving a GSL Root Finder

The GNU Scientific Library one-dimensional root finding routines are handy, but I never remember how to drive them properly. Here is what I essentially always end up re-distilling from the examples:

int fsolver_solve(
        gsl_root_fsolver * const s,
        gsl_function     * const f,
        const size_t             maxiter,
        const double             epsabs,
        const double             epsrel,
        double           * const lower,
        double           * const upper,
        double           * const root)
{
    // Initialize fsolver on [lower, upper]
    gsl_error_handler_t * h = gsl_set_error_handler_off();    // Push handler
    int status = gsl_root_fsolver_set(s, f, *lower, *upper);
    status = (GSL_SUCCESS == status) ? GSL_CONTINUE : status;

    // Proceed until success, failure, or maximum iterations reached
    // On success, overwrite *location as described in API documentation.
    *root = GSL_NAN;
    for (size_t iter = 1; iter <= maxiter && status == GSL_CONTINUE; ++iter) {
        if (GSL_SUCCESS != (status = gsl_root_fsolver_iterate(s))) break;
        *lower = gsl_root_fsolver_x_lower(s);
        *upper = gsl_root_fsolver_x_upper(s);
        status = gsl_root_test_interval(*lower, *upper, epsabs, epsrel);
    }
    if (status == GSL_SUCCESS) {
        *root = gsl_root_fsolver_root(s);
    }
    gsl_set_error_handler(h);                                 // Pop handler
    return status;
}