Ошибка в методе Ньютона

[felarl]

Ошибка в методе Ньютона

Пишет Вещественное деление на 0, сама программа выглядит так:

uses Crt;
var x,a,b,e: real;
function f1(d: real): real;
begin
f1:=(x-2)*(x-2)*2*x-1;
end;
function f2(d:real): real;
begin
f2:=6*x*x-16*x+8;
end;
begin
Clrscr;
a:=0.4;
b:=2;
e:=0.001;
if f1(a)*f2(a)>0 then x:=a
else x:=b;
while abs(f1(x))>e do
begin
x:=x-f1(x)/f2(x);
end;
Writeln ('x=',x,' f(x)=',f1(x));
Readln
end.

 

[felarl]

Буду благодарен за помощь!

 

[Vladimir_S]

Буду благодарен за помощь!

Легко. Всё правильно, кроме того, что в ТЕЛЕ обеих функций аргументом должно быть d, а не x. Раз уж в заголовке стоит d. С учетом этого и еще пары мелких поправок:

uses Crt;

const
 a=0.4;
 b=2;
 e=0.001;

var
 x: real;

function f1(d: real): real;
begin
 f1:=(d-2)*(d-2)*2*d-1;
end;

function f2(d:real): real;
begin
 f2:=6*d*d-16*d+8;
end;

begin
 Clrscr;
 if f1(a)*f2(a)>0 then x:=a else x:=b;
 while abs(f1(x))>e do
  begin
   x:=x-f1(x)/f2(x);
  end;
 Writeln ('x= ',x:0:5,'   f(x)= ',f1(x):0:12);
 Readln
end.

Результат работы программы:

x= 0.14536           f(x)= -0.000002649617

 

[felarl]

Легко. Всё правильно, кроме того, что в ТЕЛЕ обеих функций аргументом должно быть d, а не x. Раз уж в заголовке стоит d. С учетом этого и еще пары мелких поправок:

uses Crt;

const
 a=0.4;
 b=2;
 e=0.001;

var
 x: real;

function f1(d: real): real;
begin
 f1:=(d-2)*(d-2)*2*d-1;
end;

function f2(d:real): real;
begin
 f2:=6*d*d-16*d+8;
end;

begin
 Clrscr;
 if f1(a)*f2(a)>0 then x:=a else x:=b;
 while abs(f1(x))>e do
  begin
   x:=x-f1(x)/f2(x);
  end;
 Writeln ('x= ',x:0:5,'   f(x)= ',f1(x):0:12);
 Readln
end.

Результат работы программы:

x= 0.14536           f(x)= -0.000002649617

Тааак, извините, могли бы объяснить почему d? что-то я не понимаю...

 

[felarl]

И еще вопросик, почему х не попадает в заданный интервал?

 

[Vladimir_S]

Тааак, извините, могли бы объяснить почему d? что-то я не понимаю...

Потому что так Вы назвали (в шапке функции) формальный параметр. При обращении к функции он заменяется фактическим "х". Но в теле функции должен стоять именно формальный.

И еще вопросик, почему х не попадает в заданный интервал?

Знаете - я вообще не понимаю смысла интервала, равно как и метода выбора стартовой точки. В принципе, стартовой может быть любая точка не очень далеко расположенная от ожидаемого значения решения, лишь бы там не было корней производной.

 

[felarl]

Потому что так Вы назвали (в шапке функции) формальный параметр. При обращении к функции он заменяется фактическим "х". Но в теле функции должен стоять именно формальный.Знаете - я вообще не понимаю смысла интервала, равно как и метода выбора стартовой точки. В принципе, стартовой может быть любая точка не очень далеко расположенная от ожидаемого значения решения, лишь бы там не было корней производной.

С d разобрался... А почему он пишет ответ 0.14, а не 1.45 например, на графике видно что там есть как этот корень, так и другой, да и метод дихотомии выдает ответ 1.45.
P.S. Извините если вопросы кажутся абсурдными, просто только учусь, и хотелось бы понимать что написанно.

 

[Vladimir_S]

А почему он пишет ответ 0.14, а не 1.45 например, на графике видно что там есть как этот корень, так и другой, да и метод дихотомии выдает ответ 1.45.

Точнее

x= 1.40303      f(x)= 0.000003728748

А есть и третий корень - чай, уравнение-то кубическое!

x= 2.45161      f(x)= 0.000003331817

Все три корня можно получить, задавая исходную точку где-то в окрестности каждого из них. Вот поэтому прежде, чем численно решать уравнение, нужно его исследовать, в простейшем случае - построить график да и посмотреть, где корень, единствен ли он, а если их несколько - какой нужен и т.п. В нашем случае график выглядит так:

 

[felarl]

Точнее

x= 1.40303      f(x)= 0.000003728748

А есть и третий корень - чай, уравнение-то кубическое!

x= 2.45161      f(x)= 0.000003331817

Все три корня можно получить, задавая исходную точку где-то в окрестности каждого из них. Вот поэтому прежде, чем численно решать уравнение, нужно его исследовать, в простейшем случае - построить график да и посмотреть, где корень, единствен ли он, а если их несколько - какой нужен и т.п. В нашем случае график выглядит так:

А что считать исходной точкой в программе? я так понимаю точку а=0,4?

 

[Vladimir_S]

А что считать исходной точкой в программе? я так понимаю точку а=0,4?

Да почему же? ЛЮБУЮ точку где-то поблизости от предполагаемого корня (находится, например, по графику). Вот Вам вариант программы, выдающий все три корня - и причем вообще безо всяких a и b:

uses Crt;
const
 e=0.001;
var
 x: real;

function f1(d: real): real;
begin
 f1:=(d-2)*(d-2)*2*d-1;
end;

function f2(d:real): real;
begin
 f2:=6*d*d-16*d+8;
end;

begin
 Clrscr;

 x:=0.2;
 while abs(f1(x))>e do
  begin
   x:=x-f1(x)/f2(x);
  end;
 Writeln ('x= ',x:0:5,'   f(x)= ',f1(x):0:12);

 x:=1.6;
 while abs(f1(x))>e do
  begin
   x:=x-f1(x)/f2(x);
  end;
 Writeln ('x= ',x:0:5,'   f(x)= ',f1(x):0:12);

 x:=2.5;
 while abs(f1(x))>e do
  begin
   x:=x-f1(x)/f2(x);
  end;
 Writeln ('x= ',x:0:5,'   f(x)= ',f1(x):0:12);

 Readln
end.

 

[felarl]

Снова нужна помощь... Решаю уравнение методом Ньютона, вроде все получается, найденый мною корень хоть и существует, но не входит в заданный интервал, можете указать на ошибку? Заранее благодарен)
Само уравнение: ln(2-x)-cos2x=0, интервал [1;2]
Сама программа:

uses crt;
var x,a,b,e: real;
function f1(x: real): real;
begin
f1:=ln(2-x)-cos(2*x);
end;
function f2(x:real): real;
begin
f2:=(-1/2-x)+2*sin(2*x);
end;
begin
clrscr;
a:=1;
b:=2;
e:=0.001;
if f1(a)*f2(a)>0 then x:=a
else x:=b;
while abs(f1(x))>e do
begin
x:=x-f1(x)/f2(x);
end;
writeln ('x=',x,' f(x)=',f1(x));
end.

 

[Vladimir_S]

Прежде всего, правильно сосчитайте производную (первое слагаемое), а там дальше посмотрим. Пока что в качестве f2 написана полная чушь.

 

[felarl]

Прежде всего, правильно сосчитайте производную (первое слагаемое), а там дальше посмотрим. Пока что в качестве f2 написана полная чушь.

Не понимаю где ошибка... вроде правильно производную взял, может записал в строчку не правильно...

 

[Vladimir_S]

Не понимаю где ошибка... вроде правильно производную взял, может записал в строчку не правильно...

Охо-хо, грехи наши тяжкие... Объясняю.
У Вас:
f2:=(-1/2-x)+2*sin(2*x);
Надо:
f2:=-1/(2-x)+2*sin(2*x);

 

[felarl]

Охо-хо, грехи наши тяжкие... Объясняю.
У Вас:
f2:=(-1/2-x)+2*sin(2*x);
Надо:
f2:=-1/(2-x)+2*sin(2*x);

Да, понял свою ошибку, Спасибо) часто ошибки бывают именно в таких местах...

 

[felarl]

Если возможно, подскажите почему здесь происходит деление на ноль, и можно ли это исправить, не меняя интервал [0;2]

uses crt;
var x,a,b,e: real;
function f1(x: real): real;
begin
f1:=(x-2)*(x-2)*2*x-1;
end;
function f2(x:real): real;
begin
f2:=6*x*x-16*x+8;
end;
begin
clrscr;
a:=0;
b:=2;
e:=0.001;
if f1(a)*f2(a)>0 then x:=a
else x:=b;
while abs(f1(x))>e do
begin
x:=x-f1(x)/f2(x);
end;
writeln ('x=',x,' f(x)=',f1(x));
end.

 

[Vladimir_S]

Если возможно

Возможно.

подскажите почему здесь происходит деление на ноль

Потому что в указанном интервале производная проходит через 0 (горизонтальная касательная, см. график выше).

и можно ли это исправить, не меняя интервал [0;2]

Нет.

 

[felarl]

Возможно.Потому что в указанном интервале производная проходит через 0 (горизонтальная касательная, см. график выше).Нет.

Огромное спасибо за ответы)

 

[Vladimir_S]

И вообще, хотел бы я понять, откуда берется вот эта, вставляемая Вами во все программы, ахинея:
if f1(a)*f2(a)>0 then x:=a
else x:=b;
Ни малейшего смысла в ней я, извините, не вижу. Требуется проверка условия положительности произведения производных на концах интервала:
if f2(a)*f2(b)>0 then
Если это условие НЕ выполняется, то задача заведомо методом Ньютона не решается. Если выполняется, то тоже не факт, что решается: функция может оказаться сильно осциллирующей и производная внутри интервала может несколько раз пройти через ноль.
Поэтому, как я Вам уже писал и повторю снова: исходная точка должна выбираться "вручную", например, по графику функции. Алгоритмов выбора ее не существует. Мне, по крайней мере, они не известны.
Эх, разбирали-разбирали задачу, график рисовали, про выбор исходных точек жевали-переживали - да всё, как видно, впустую...

 

[felarl]

И вообще, хотел бы я понять, откуда берется вот эта, вставляемая Вами во все программы, ахинея:
if f1(a)*f2(a)>0 then x:=a
else x:=b;
Ни малейшего смысла в ней я, извините, не вижу. Требуется проверка условия положительности произведения производных на концах интервала:
if f2(a)*f2(b)>0 then
Если это условие НЕ выполняется, то задача заведомо методом Ньютона не решается. Если выполняется, то тоже не факт, что решается: функция может оказаться сильно осциллирующей и производная внутри интервала может несколько раз пройти через ноль.
Поэтому, как я Вам уже писал и повторю снова: исходная точка должна выбираться "вручную", например, по графику функции. Алгоритмов выбора ее не существует. Мне, по крайней мере, они не известны.
Эх, разбирали-разбирали задачу, график рисовали, про выбор исходных точек жевали-переживали - да всё, как видно, впустую...

Нас так учат, к сожалению, и поэтому я не могу взять в ручную начальную точку, т.к. это будет "не то что от нас хотят"