Алгоритм получение корней квадратного уравнения

Получение корней квадратного уравнения относится к элементарным понятиям алгебры; формулы для корней известны очень давно. В то же время, при численной реализации этих формул могут встретиться неприятности, обход одной из которых и будет здесь разобран. Программа написана для случая действительных коэффициентов (корни могут быть комплексными). 

Квадратное уравнение записывается в виде:
a*x²+b*x+c=0,   a!=0.
Классическая формула для нахождения его корней (как в действительном, так и в комплексном случае):

x₁=(-b-sqrt(b²-4ac))/(2a),    x₂=(-b+sqrt(b²-4ac))/(2a).

Существует и альтернативный способ получения корней: x₁=2c/(-b-sqrt(b²-4ac)),    x₂=2c/(-b+sqrt(b²-4ac)).

Трудность, о которой будет идти здесь речь (и она не единственная), связана с тем, что в том случае, когда b² значительно превышает |4ac|, то при пользовании классическими формулами один из корней будет формироваться с использованием вычитания двух близких чисел и это может привести к потере точности. Чтобы этого избежать, по классическим формулам вычисляется только один корень (при вычислении которого этого не происходит), а второй корень находится по формуле Виета.

Случай, когда корни комплексные, для действительных корней уравнения указанной трудности не представляет.

 

Программа для решения квадратного уравнения с действительными коэффициентами

Ниже расположена программа для решения квадратного уравнения с действительными коэффициентами.

/* Решение квадратного уравнения. Случай с действительными коэффициентами.

   int Quadratic(double *x,double a,double b,double c);
   Параметры:
   x - массив решений (размера 2). 
   На выходе:
       2 действительных корня -> тогда x заполняется ими;
       2 сложно-соединенных корня -> x[0] - действительная часть, 
         x[1] - неотрицательная мнимая часть. 
       другие случаи -> x[0] - уникальный корень, если возвращено (-1),
                     иначе нет верных корней.
   a, b, c - коэффициенты: ax^2 + bx + c = 0. 
   Результаты: 2 - 2 действительных корня;
            1 - 1 действительный корень (x[0]=x[1]);
            0 - 2 комплексных корня;
            -1 - один действительный корень в случае a==0;
            -2 - нет корней в случае a=0, b=0;
            -3 - бесконечное число корней (a=b=c=0). 
*/

#include <math.h>   /* for sqrt() */

int Quadratic(double *x,double a,double b,double c) {
  double d;
  /* вырожденные случаи */
  if(a==0.) {
    if(b==0.) {
      if(c==0.) return(-3);
      return(-2);
    }
    x[0]=-c/b; return(-1);
  }
  /* главный случай */
  d=b*b-4.*a*c;       /* дискриминант */
  /* единственный корень */
  if(d==0.) {
    x[0]=x[1]=-b/(2.*a); return(1);
  }
  /* сложно-соединенные корни */
  if(d<0.) {
    double t=0.5/a;
    x[0]=-b*t; x[1]=sqrt(-d)*t;
    return(0);
  }
  /* 2 действительных корня */
  if(b>=0.) d=(-0.5)*(b+sqrt(d));
  else d=(-0.5)*(b-sqrt(d));
  x[0]=d/a; x[1]=c/d;
  return(2);
}