Численная градиентная расходимость в GSL

Я использую метод сопряженных градиентов

gsl_multimin_fdfminimizer_conjugate_fr

в пакете GSL, чтобы минимизировать функцию f. Этот метод требует производных, поэтому я передал производную двумя способами:

static void dfun(const gsl_vector *v, void *params, gsl_vector *df) {
return fun(x+epsilon)-f(x)/epsilon
}

и другой способ, который менее прямой, но более быстрый. Обе функции возвращают одинаковые значения для градиента во всех переменных с точностью до ~ 10 цифр. Однако, когда я запускаю минимизацию с использованием последней, я не достигаю минимума, и минимизатор возвращает 27: неспособен получить лучший результат. При использовании dfun в качестве производной я сходлюсь к правильному минимуму.

Что может быть причиной такого поведения, и какие шаги я могу предпринять, чтобы исправить это?