СНАУ методом Зэйделя

[Алекс2]

СНАУ методом Зэйделя

Исполнить «ручное» решение заданой системы уравнений методом Зэйделя с точностью ε = 0.01 (до 4-х итэраций);
Составить и отладить программу на Паскаль, которая реализует введение входных данных, решение системы методом Зэйделя, выведение результатов в удобной форме на экран и в файл; фрагменты программы оформить как процедуры и функции.
Решить СНАУ с помощью этой программы, сравнить с результатами «ручных» расчетов.
Построить графики сходимости итэррацыонного процесса
Вычислить корни системы уровнений при 5 разных значениях точности . Результаты оформить в виде таблицы


Я нашол программу для 3-х нелинейных уровнений с тремя неизвестными


но есть только фрагмент программы
Var E,x1,x2,x3,F1,F2,F3,
x11,x21,x31:real;
k: integer;
Label M,M1;
…
Begin
<введение данных >
K:=0;
M1: x11:=-x1*x1+2*x2*x3+0.1;
X21:=x2*x2-3*x11*x3-0.2;
X31:=-x3*x3-2*x11*x21+0.3;
F1:=x11+x11*x11-2*x21*x31-0.1;
F2:=x21-x21*x21+3*x11*x31+0.2;
F3:=x31+x31*x31+2*x11*x21-0.3;
if (abs(F1)<=E) and (abs(F2)<=E) and (abs(F3)<=E)
then goto M;
x1:=x11; x2:=x21; x3:=x31;
k:=k+1;
goto M1;
M: <выведенние результатов>
…
Я кое- что понял как составлять для 3-х уровнений.
Но как это применить к моей системе из 2-х уровнений я не понимаю.

Пожалуйста подскажите, а если можете то и код программы

Зараннее спасибо!

 

[Алекс2]

Я извиняюсь. забыл сказать добрый всем вечер.

 

[Vladimir_S]

Но как это применить к моей системе из 2-х уровнений я не понимаю.

А я Вам в сотый раз повторяю: в Вашей системе ОШИБКА, в таком виде она НЕСОВМЕСТНА, за исключением сингулярной точки x1=x2=0, каковую невозможно найти численными методами.
Просто поразительно: стараешься человеку помочь, объясняешь, доказываешь - впустую. Лоб и стенка. Решите систему, да и всё тут.
Еще раз: сверьте тот источник, откуда Вы взяли Вашу систему: наверняка при переписывании или еще как-то Вы допустили ошибку. В таком виде задача не-раз-ре-ши-ма!!! Потому что приводит к уравнению
112*(х2)² = (х2)²
И как прикажете такое решать?

 

[Алекс2]

Добрый вечер, уважаемый Vladimir S, вы на меня пожалуйста не злитесь, я только начинаю изучать программирование, и пытаюсь что-то понять. Я Вам очень признателен за помощь. А условие моей системы мне дал преподаватель и я не мог знать решается она или нет. ( и стоит это условие сразу к двум темам: метод Зэйделя и метод Ньютона-Раффсона). Пойду к преподавателю спрошу. Может быть у них какая-то ошибка при наборе произошла. Еще раз спасибо.

 

[Алекс2]

Добрый вечер, уважаемый Vladimir_S, я узнал у преподавателя за это задание. оказывается мне нужно с помощью програмы доказать что х1 =х2 =0. и построить график на котором (как я понял)будут две параболы которые имеют одну точку пересечения. но как это сделать ума не приложу. ( нам рассказывали только о системе из 3-х уровнений которие решаются.
Посмотрите пожалуйста, что можно сделать. Заранее спасибо

 

[Vladimir_S]

Добрый вечер, уважаемый Vladimir_S,

Добрый!

я узнал у преподавателя за это задание. оказывается мне нужно с помощью програмы доказать что х1 =х2 =0.

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

и построить график на котором (как я понял)будут две параболы которые имеют одну точку пересечения. но как это сделать ума не приложу.

Ну это-то как раз просто: переписываем оба уравнения так, чтобы получились зависимости х1(х2) и строим соответствующие кривые:

 

[Vladimir_S]

Ну, значит, так, сударь!
Поисследовал я малость эту задачку.
Начнем с метода Ньютона (при чем тут Рафсон, вот хоть убейте, не понимаю).
По идее, тут надо иметь дело с такими двумя функциями:
f1(х1,х2) = -√(x1) + 4(x2) = 0
f2(х1,х2) = -7(x1) + (x2)² = 0,
и при этом обе функции должны быть дифференцируемы по обеим переменным на всей числовой оси, включая нулевые и отрицательные значения. К сожалению, наличие корня в первом уравнении нарушает это требование, и если формально запрограммировать алгоритм Ньютона, то программа попросту вылетит с ошибкой (проверено). Поэтому попробуем обойти этот момент, введя такие переобозначения:
t = √(x1)
x1 = t
то есть вместо х1 будем работать с корнем из него, а в конце полученное значение возведем в квадрат, т.е. исходный набор функций преобразуем к виду:
f1(х1,х2) = -(x1) + 4(x2) = 0
f2(х1,х2) = -7(x1)² + (x2)² = 0,
Тогда задача с грехом пополам решается.

Const
 x1_0=0.1;
 x2_0=0.1;
 Eps=0.00001;

Type
 Matr=Array[1..2,1..2] of Real;

Var
 x1,x2,dx1,dx2,Det_0,Det_X1,Det_X2:Real;
 M:Matr;

Function Determ(D:Matr):Real;
begin
 Determ:=D[1,1]*D[2,2]-D[1,2]*D[2,1];
end;

Function f1(z1,z2:real):Real;
begin
 f1:=-z1+z2*4;
end;

Function f2(z1,z2:real):Real;
begin
 f2:=-Sqr(z1)*7+Sqr(z2);
end;

Function Der_x1_f1(z1,z2:real):Real;
begin
  Der_x1_f1:=-1.0;
end;

Function Der_x2_f1(z1,z2:real):Real;
begin
  Der_x2_f1:=4.0;
end;

Function Der_x1_f2(z1,z2:real):Real;
begin
 Der_x1_f2:=-z1*14;
end;

Function Der_x2_f2(z1,z2:real):Real;
begin
 Der_x2_f2:=2.0*z2;
end;

Begin
 x1:=x1_0;
 x2:=x2_0;
 Repeat
  M[1,1]:=Der_x1_f1(x1,x2);
  M[1,2]:=Der_x2_f1(x1,x2);
  M[2,1]:=Der_x1_f2(x1,x2);
  M[2,2]:=Der_x2_f2(x1,x2);
  Det_0:=Determ(M);
  M[1,1]:=-f1(x1,x2);
  M[2,1]:=-f2(x1,x2);
  Det_X1:=Determ(M);
  M[1,1]:=Der_x1_f1(x1,x2);
  M[1,2]:=-f1(x1,x2);
  M[2,1]:=Der_x1_f2(x1,x2);
  M[2,2]:=-f2(x1,x2);
  Det_X2:=Determ(M);
  dx1:=Det_X1/Det_0;
  dx2:=Det_X2/Det_0;
  x1:=x1+dx1;
  x2:=x2+dx2; {
 Writeln('x1 = ',x1:0:6);
 Writeln('x2 = ',x2:0:6);
 Readln; Halt; }
 Until (Abs(dx1)<Sqrt(Eps)) and (Abs(dx2)<Eps);
 Writeln('x1 = ',Sqr(x1):0:6);
 Writeln('x2 = ',x2:0:6);
 Readln
End.

Окончательный результат (уже после возведения в квадрат х1, то есть для ИСХОДНОЙ системы):
х1 = 0.000000
х2 = 0.000005
Ладно, с этим кое-как справились.

Теперь насчет Гаусса-Зейделя. Тут всё совсем тухло.
Прежде всего, это метод решения не СНАУ, а СЛАУ, то есть системы ЛИНЕЙНЫХ алгебраических уравнений, поэтому надо как-то привести задачу к СЛАУ. Для этого проще всего возвести обе части первого уравнения в квадрат:
х1 = 16(х2)²,
затем переобозначив
(х2)² → (х2),
пытаться решить уже линейную систему
(х1) - 16(х2) = 0
7(х1) - (х2) = 0.
Ну так вот, ничего из этой затеи не выйдет: даже простейший анализ "на бумаге" показывает, что для такой системы итерационная процедура пойдет "вразнос", т.е. каждая следующая пара значений (х1) и (х2) будет убегать от нулевой точки, а не приближаться к ней.

Пожалуй, это всё, чем могу помочь.

P.S. А препод ваш - идиот и извращенец. Таково моё мнение.

 

[Алекс2]

Добрый вечер, или скорей всего ночи, Vladimir_S.
Хочу сказать Вам огромнейшее спасибо. Вы меня очень выручили. К такому решению я бы сам не додумался вовек.
А на щет преподавателя я с Вами совершенно согласен.

 

[Алекс2]

Добрый вечер, уважаемый Vladimir_S. Не поможете ли еще немного.

Вопрос: нет ли у Вас программы для графиков функций?
Мне нужно составить отдельно программу которая нарисует мне эти 2 графика.
Помогите, пожалуйста, если можно, а то у меня не получается.

 

[Vladimir_S]

Добрый вечер, уважаемый Vladimir_S. Не поможете ли еще немного.

Вопрос: нет ли у Вас программы для графиков функций?
Мне нужно составить отдельно программу которая нарисует мне эти 2 графика.
Помогите, пожалуйста, если можно, а то у меня не получается.

Тут так. Я пользуюсь DOS Turbo Pascal с его Borland-графикой. Вряд ли это Вас устроит. Например, если у Вас АВС, то там графпакет другой, в котором я не пишу.

 

[Алекс2]

Добрый вечер,Vladimir_S. Спасибо за ответ. У меня Турбо Паскаль. Буду искать пример.

 

[Vladimir_S]

Добрый вечер,Vladimir_S. Спасибо за ответ. У меня Турбо Паскаль. Буду искать пример.

Не очень понял. Но если Borland-пакет устраивает, то:

Uses Graph, CRT;

Var
 x,y,step,Kx,Ky,Dx,Dy:Real;
 i:Integer;
 S:String;

Begin
 InitGraph(.........); [COLOR=Red][B]{Здесь - самостоятельно}[/B][/COLOR]
 SetLineStyle(SolidLn,0,ThickWidth);
 Line(20,GetMaxY-60,GetMaxX-20,GetMaxY-60);
 Line((GetMaxX div 2),GetMaxY-60,(GetMaxX div 2),20);
 SetTextJustify(CenterText,TopText);
 OutTextXY(GetMaxX-30,GetMaxY-35,'x');
 OutTextXY(GetMaxX-20,GetMaxY-30,'2');
 SetTextJustify(LeftText,CenterText);
 OutTextXY((GetMaxX div 2)+10,25,'x');
 OutTextXY((GetMaxX div 2)+20,30,'1');
 Kx:=(GetMaxX-40)/0.1;
 Ky:=(GetMaxY-80)/0.001;
 SetLineStyle(SolidLn,0,NormWidth);
 step:=0.01;
 SetTextJustify(CenterText,TopText);
 for i:=-5 to 5 do
  begin
   Line(Round(GetMaxX/2+step*i*Kx),GetMaxY-64,
        Round(GetMaxX/2+step*i*Kx),GetMaxY-56);
   Str(step*i:5:2,S);
   OutTextXY(Round(GetMaxX/2+step*i*Kx),GetMaxY-52,S);
  end;
 step:=0.0001;
 SetTextJustify(RightText,CenterText);
 for i:=1 to 10 do
  begin
   Line((GetMaxX div 2)-4,Round(GetMaxY-60-step*i*Ky),
        (GetMaxX div 2)+4,Round(GetMaxY-60-step*i*Ky));
   Str(step*i:6:4,S);
   OutTextXY((GetMaxX div 2)-8,Round(GetMaxY-60-step*i*Ky),S);
  end;
 SetColor(LightRed);
 SetLineStyle(SolidLn,0,ThickWidth);
 MoveTo(Round(GetMaxX/2-0.008*Kx),Round(GetMaxY-60-Sqr(0.008)*16*Ky));
 for i:=1 to 16 do
  LineTo(Round(GetMaxX/2+(-0.008+0.001*i)*Kx),
         Round(GetMaxY-60-Sqr(-0.008+0.001*i)*16*Ky));
 SetColor(LightBlue);
 MoveTo(Round(GetMaxX/2-0.05*Kx),Round(GetMaxY-60-Sqr(0.05)/7*Ky));
 for i:=1 to 100 do
  LineTo(Round(GetMaxX/2+(-0.05+0.001*i)*Kx),
         Round(GetMaxY-60-Sqr(-0.05+0.001*i)/7*Ky));
 ReadKey
End.

О конфигурировании графического режима и параметров InitGraph см.
http://www.tehnari.ru/f43/t95841/index2.html#post1026837 и дальше, а также
http://www.tehnari.ru/f43/t95841/#post1026107

 

[Алекс2]

Есть программа которая рисует 2 графика. но она на FreePaskal, помогите, пожалуйста переделать на Турбо Паскаль 7.

программа построения дух графиков на одной оси координат}
program P;

uses
Graph,
CRT;

function f1(X: real): real;
begin
f1 := 16 * sqr(X);
end;

function f2(X: real): real;
begin
f2 := 1 / 7 * sqr(X);
end;

var
GraphDevice, GraphMode: integer;
ScreenXstart, ScreenXfinish, ScreenYstart, ScreenYfinish: integer;
ScaleX, ScaleY: real;
Xstart, Xfinish, Ystart, Yfinish: real;
X, Y, dX: real;
Xscreen, Yscreen: integer;
IsNextVisiblePixel: boolean;
begin
{параметры графиков}
Xstart := -0.05;
Xfinish := 0.05;
Ystart := 0.0;
Yfinish := 0.0010;
{инициализация графики}
GraphDevice := Detect;
GraphMode := Detect;
InitGraph(GraphDevice, GraphMode, '');
if GraphResult <> grOk then
begin
writeln('Error');
Exit;
end;
{от края экрана отступим по 10% экрана}
ScreenXstart := GetMaxX div 10;
ScreenXfinish := GetMaxX - ScreenXstart;
ScreenYstart := GetMaxY div 10;
ScreenYfinish := GetMaxY - ScreenYstart;
{расчёт масштабных коэффициентов}
ScaleX := (ScreenXfinish - ScreenXstart) / (Xfinish - Xstart);
ScaleY := (ScreenYfinish - ScreenYstart) / (Yfinish - Ystart);
{построение осей координат}
SetColor(7);
MoveTo(ScreenXstart, ScreenYstart);
LineTo(ScreenXstart, ScreenYfinish);
LineTo(ScreenXfinish, ScreenYfinish);
LineTo(ScreenXfinish, ScreenYstart);
LineTo(ScreenXstart, ScreenYstart);

MoveTo((ScreenXfinish + ScreenXstart) div 2, ScreenYstart);
LineTo((ScreenXfinish + ScreenXstart) div 2, ScreenYfinish);
{построение 1-го графика}
SetColor(5);
IsNextVisiblePixel := False;
X := Xstart;
dX := (Xfinish - Xstart) / 1000;
while X <= Xfinish do
begin
Y := f1(X);
Xscreen := round((X - Xstart) * ScaleX) + ScreenXstart;
Yscreen := -round((Y - Ystart) * ScaleY) + ScreenYfinish;
if (Xscreen >= ScreenXstart) and (Xscreen <= ScreenXfinish) and
(Yscreen > ScreenYstart) and (Yscreen < ScreenYfinish) then
begin
if IsNextVisiblePixel then
LineTo(Xscreen, Yscreen)
else
MoveTo(Xscreen, Yscreen);
IsNextVisiblePixel := True;
end
else
IsNextVisiblePixel := False;
X := X + dX;
end;
{построение 2-го графика}
SetColor(4);
IsNextVisiblePixel := False;
X := Xstart;
dX := (Xfinish - Xstart) / 1000;
while X <= Xfinish do
begin
Y := f2(X);
Xscreen := round((X - Xstart) * ScaleX) + ScreenXstart;
Yscreen := -round((Y - Ystart) * ScaleY) + ScreenYfinish;
if (Xscreen >= ScreenXstart) and (Xscreen <= ScreenXfinish) and
(Yscreen > ScreenYstart) and (Yscreen < ScreenYfinish) then
begin
if IsNextVisiblePixel then
LineTo(Xscreen, Yscreen)
else
MoveTo(Xscreen, Yscreen);
IsNextVisiblePixel := True;
end
else
IsNextVisiblePixel := False;
X := X + dX;
end;
ReadKey;
end.

 

[Vladimir_S]

Есть программа которая рисует 2 графика. но она на FreePaskal, помогите, пожалуйста переделать на Турбо Паскаль 7.

Без разницы. Ничего переделывать не нужно.
И если Вам больше по сердцу такой слепой график

чем тот, что предложил Вам я - с разметкой осей - то что же, как говорится, "насильно мил не будешь".
Успехов!

 

[Алекс2]

Добрый день,Vladimir_S. Хочу еще раз сказать спасибо именно за вашу программу. Проблема была в том, что я не мог увидеть результат. Теперь все нормально. график есть. Огромное спасибо еще раз!!!!!!!!!