Лекция 7.

Задачи численного анализа

Поиск корня функции методом деления отрезка пополам

Условие на функцию

f(a)*f(b) <= 0

Бинарный поиск, использованный для поиска в отсортированном массиве, — универсальный метод, применяемый в разных областях.

Так, его можно применять для нахождения корней функции на отрезке, если функция принимает на концах отрезка значения разного знака. Рассмотрим значение функции в середине данного отрезка. Если это значение имеет разный знак со значением функции на левой границе отрезка, то имеет смысл искать корень функции в левой половине отрезка. В противном случае имеет смысл искать корень функции в правой половине отрезка. Каждая такая итерация сокращает диапазон поиска корней вдвое, что ведет к быстрому уменьшению длины рассматриваемого отрезка. Как только длина рассматриваемого отрезка становится меньше некоего эпсилон, заданного пользователем, поиск можно считать законченным, так как любая точка этого отрезка приближает корень функции с заданной точностью эпсилон.

Программа поиска корня функции f(x) на отрезке [a, b] с точностью eps:

#include <stdio.h>

double f(double x);

int main(void) {

double a, b, c, eps;

printf("Input interval (a, b) and epsilon\n");

scanf("%lf%lf%lf", &a, &b, &eps);

if (f(a)*f(b) > 0) {

printf("Function values at the interval ends are the same sign\n");

return -1;

}

while ((b-a) >= eps) {

c = (a + b) / 2;

if (f(a)*f(c) <=0 ) {

b = c;

} else {

a = c;

}

printf("Function root: %lf\n", (a + b)/2);

return 0;

}

double f(double x) {

return x*x – 5;

}

Поиск корня функции методом Ньютона (касательных)

Метод Ньютона поиска корней функции основан на построении рекуррентной последовательности точек, все более точно приближающих искомый корень функции.

Он не является абсолютно надежным, как метод деления отрезка пополам, так как может оказаться, что последовательность расходится.

Однако в некоторых случаях он оказывается полезным.

Каждая точка в последовательности x₀→x₁→x₂→x₃→… (кроме точки x₀, задаваемой пользователем) находится из предыдущей следующим образом: в точке x_i проводится касательная к функции и место пересечения этой касательной и оси X выбирается в качестве точки x_i+1.

Программа поиска корня функции f(x) методом Ньютона:

#include <stdio.h>

#include <math.h>

double f(double x);

double df(double x);

int main(void) {

double a, eps;

printf("Input starting point a and epsilon\n");

scanf("%lf%lf", &a, &eps);

while (fabs(f(a)) > eps) {

/* if (N>1000) { … } */

if (df(a) == 0) {

printf("The method is inappropriate\n");

return -1;

}

a = a – f(a)/df(a);

}

printf("Function root: %lf\n", a);

return 0;

}

double f(double x) {

return x*x - 5;

}

double df(double x) {

return 2*x;

}

Поиск корня функции методом хорд (секущих)

Метод хорд поиска корней функции также основан на построении рекуррентной последовательности точек, все более точно приближающих искомый корень функции.

Он также не является абсолютно надежным, как метод деления отрезка пополам, так как может оказаться, что последовательность расходится.

Однако в некоторых случаях он также оказывается полезным.

Каждая точка в последовательности x₀→x₁→x₂→x₃→… (кроме точек x₀ и x₁, задаваемых пользователем) находится из двух предыдущей следующим образом: через точки (x_i ,f(x_i)) и (x_i+1,f(x_i+1)) проводится прямая (секущая) и место пересечения этой прямой и оси X выбирается в качестве точки x_i+₂.

Составим уравнение для поиска точки c пересечения хорды, проходящей через точки (a,f(a)) и (b,f(b)):

f(b)/f(a) = (b-c)/(a-c)

Отсюда находим точку с = (f(a)*b-f(b)*a)/(f(a)-f(b)).

Программа поиска корня функции f(x) методом хорд:

#include <stdio.h>

#include <math.h>

double f(double x);

int main(void) {

double a, b, c, eps;

printf("Input starting points a and b and epsilon\n");

scanf("%lf%lf%lf", &a, &b, &eps);

while (fabs(f(a)) > eps) {

if (f(a) == f(b)) {

printf("The method is inappropriate\n");

return -1;

}

c = (f(a)*b-f(b)*a)/(f(a)-f(b));

a = b;

b = c;

}

printf("Function root: %lf\n", b);

return 0;

}

double f(double x) {

return x*x – 4;

}

Поиск машинного эпсилон

Целое число представляется в памяти компьютера в двоичном виде:

Например, число 13:

31	30	29	28	27	26	25	...	7	6	5	4	3	2	1	0
0 (знак)	0	0	0	0	0	0		0	0	0	0	1	1	0	1

Например, число -13:

Прямое представление:

31	30	29	28	27	26	25	...	7	6	5	4	3	2	1	0
1 (знак)	0	0	0	0	0	0		0	0	0	0	1	1	0	1

Дополнительный код (используется чаще всего):

31	30	29	28	27	26	25	...	7	6	5	4	3	2	1	0
1 (знак)	1	1	1	1	1	1		1	1	1	1	0	0	1	1

0 (знак)

1+1 = 0 (1 в уме)

1+0 (+1) = 0 (1 в уме)

…

1+0 (+1) = 0 (1 в уме — идет в бит процессора (ARM бит Overflow))

Старший бит — пропадает из представления числа (но не из процессора)

Вещественные числа:

1000000000 = 0.1 * 10^10

0.00000000000001 = 0.1 * 10^(-13)

1000000000.00000000000001 = 0.100000000000000000000001 * 10^10

31	30	29	28	27	26	25	...	7	6	5	4	3	2	1	0
0 (знак)	0	0	0	0	0	0		0	0	0	0	1	1	0	1

Мантисса (выделено красным) — значащие цифры

Ордината (порядок) (выделено синим) — степень 10

Если мантисса не помещается в отведеннные под нее биты, то младшие значащие цифры пропадают (округляются)

Пример 1.
1+x = 1

x = 0.1 * 10*^(-22)

1 + 0.00000000000000000000001 = 1.00000000000000000000001 = 0.100000000000000000000001 * 10^(-1) =

= 0.10000000000000000000000 * 10^(-1)

Если вся мантисса суммы не поместится в отведенные под нее биты, то младшие цифра отбросятся =>

результат будет равен 1

Машинным эпсилон для данной вычислительной системы является наибольшее среди вещественнных чисел x, представимых в этой системе и удовлетворяющих в ней соотношению 1+x=1. Такое соотношение является возможным для некоторых очень маленьких чисел, так как в вычислительной системе значащие цифры числа (мантисса) хранятся с некоторой точностью, то есть округляются так, чтобы поместиться в ячейке памяти некоторого размера, определенного для данной вычислительной системы. При прибавлении такого очень маленького числа к единице количество значащих цифр увеличивается за счет появления большого количества нулей после старшей 1. Например, у суммы 1 + 0.000001 = 1.000001 количество значащих цифр 7, в то время как у исходных чисел было по одной значащей цифре.

При помещении результата суммирования в памяти вычислительной системы последние значащие цифры могут урезаться, чтобы результат поместился в ячейке памяти. Если таким образом урезаются все значащие цифры маленького числа, то результатом, помещенным в ячейку памяти, будет 1.

Программа поиска машинного эпсилон

#include <stdio.h>

int main(void) {

double x = 1;

while (1 + x > 1) {

x /= 2;

}

printf("Mach epsilon: %lf\n", x);

return 0;

}

Задача: Найти 1000 чисел, которые при сложении в программе в прямом и обратном порядке дают разный результат (подсказка — использовать машинное эпсилон)