Вы делаете это неправильно

По-моему, почти в каждый курс программирования входит задачка вроде «напишите программу, решающую квадратные уравнения». Обычно это второе или третье задание после «Hello, world!» — считается, что это хороший способ продемонстрировать нетривиальные инструкции ветвления. «Хорошее» решение сводится к вычислению дискриминанта и в случае, если он неотрицательный — вычислению корней квадратного трехчлена ax2+bx+c по формуле, известной из курса алгебры за седьмой, что ли, класс:

square1

Вот такой вариант решения обычно считается более-менее приемлемым (хотя его можно/нужно обвешать еще несколькими проверками — например, не равен ли коэффициент a нулю?):

int solve(double a, double b, double c, double *x1, double *x2){
	double d = b*b - 4.*a*c;
	if( d >= 0 ){
		d = sqrt(d);
		*x1 = (-b + d)/(2.*a);
		*x2 = (-b - d)/(2.*a);
		return 0;
	}
	return -1;
}

В чем проблема? На первый взгляд все более-менее хорошо, но… Давайте для тестирования будем подсовывать уравнения с известными корнями — используя для этого теорему Виета. А именно, зафиксируем коэффициент a=1, тогда уравнение с корнями x1 и x2 будет иметь коээффициенты b=-(x1+x2) и c=x1x2. Сравнивая корни, полученные при решении уравнения, с известными нам, оценим «качество» решения.

Если корни «нормальные» — те, с которыми справится шестиклассник — то все хорошо. Но что будет, если взять два «нехороших» корня — к примеру, x1=1, x2=10-14 (это для double; если вы пользуетесь float — то возьмите второй корень, равный 10-5)? Проверьте — не забыв включить вывод максимально возможного количества значащих цифр (в printf лучше всего использовать форматный спецификатор %g, при использовании вывода в стиле C++, через iostream, такой вывод включен по умолчанию). Ошибка при вычислении второго корня возникнет уже в четвертой значащей цифре, это, на самом деле, уже довольно неприлично.

В чем дело? Если обратить внимание на «обычную» формулу для вычисления корней квадратного уравнения, то мы увидим, что при вычислении меньшего корня корень из дискриминанта вычитается из b — а так как они в этом случае очень близки, то возникающая при этом вычислительная ошибка оказывается слишком большой.

Метод, разумеется, можно улучшить. Для начала — можно вспомнить о существовании еще одной формулы для корней квадратного уравнения:

square2

Выводится она абсолютно аналогичным образом, от «классической» отличается тем, что «не работает» при c=0.

Если переписать программу, чтобы она использовала эту формулу — то меньший корень «нехорошего» уравнения будет вычисляться точно, а проблемы возникнут с большим корнем. Причина та же самая — вычитание двух близких по величине чисел. Но ведь если вычислять больший корень по первой формуле, а меньший — по второй, то эта проблема исчезнет! Поэтому более правильный метод решения квадратного уравнения должен выглядеть так:

— вычисляем дискриминант D=b2-4ac
— если дискриминант неотрицателен, то вычисляем

q

— корни уравнения равны q/a и c/q.

Как реализовать это в программе — довольно очевидно, это особо ее не усложнит.

Какая здесь мораль? Численные методы и программирование — это две совершенно разных области человеческого знания. Математические задачки — вроде решения квадратного уравнения — могут показаться интересными с точки зрения обучения программированию, но это «чужая территория», и можно столкнуться с совершенно непредсказуемыми вещами. Признайтесь, многие ли слышали о сложностях, возникающих при решении на компьютере квадратных уравнений — хотя казалось бы, что может быть проще?

Вы делаете это неправильно: 4 комментария

  1. Шура, вы задрот. Ну какие численные методы? Какой код на сях? В тренде завно ООП и всякие интерпретируемые языки (каждый год новые), на которых пишут фронтэнды СУБД и встраиваемые интерфейсы всякой ерунды в гипертекстовых документах. Вы бы еще в машинных кодах пример разобрали. 95% аудитоиивсё это так далеко…

    1. 95% аудитории даже не знает, как это работает на уровне машинного кода. А на уровне прерываний тем более. Но это еще не говорит о том, что знать этого не нужно. Впрочем, пост статьи был совершенно не об этом, если кто-то не заметил )

    2. То же самое можно написать и в стиле ООП, и на любом интерпретируемом языке (Си я выбрал, как «нестареющую классику») — дело в другом. Задачка вроде бы элементарная, и якобы демонстрирует умение программиста написать несколько if … else … для обработки некоторых «нехороших» ситуаций (нулевой коэффициент при x^2, возникающие в вычислениях NaNы или бесконечности, …) — только к «правильному» решению все эти вещи никак не приближают, и более того — наверное, 95% «преподавателей» даже не задумываются о том, что выбранный ими алгоритм решения изначально не вполне годится. Ну а дальнейшее его «вылизывание», обвешивание тестами и так далее в целом бессмысленно.

Добавить комментарий

Ваш адрес email не будет опубликован. Обязательные поля помечены *