Математическое моделирование и численные методы в решении технических задач

1. Операции с файлами

2. Решение иррациональных уравнений

2.1 Метод хорд

2.2 Метод половинного деления

3. Вычисление определенного интеграла

3.1 Метод парабол (Симпсона)

3.2 Метод прямоугольников

3.3 Метод трапеций

4. Решение систем линейных алгебраических уравнений

4.1 Метод Гаусса

4.2 Метод Ньютона

4.3 Метод Зейделя

5. Решение дифференциальных уравнений

5.1 Метод Эйлера

5.2 Метод Рунге-Кутта

6. Ряды Фурье

1. Операции с файлами

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

В данном случае мы каждый элемент таблицы умножили на 3.

Текст задачи:

const c=3;

var t, p: text; a, b: array [1.10,1.10] of integer;

i, j, n, m, k:integer;

begin

assign (t,'H:Phoenix1.txt');

assign (p,'H: Phoenix2. txt');

reset (t); i:=0;

while not eof (t) do

begin i:=i+1; j:=0;

while not eoln (t) do

begin j:=j+1; read (t, a[i, j]);

end; readln (t); end;

n:=i; m:=j;

for i:=1 to n do

begin for j:=1 to m do

write (a[i, j]: 4); writeln; end; writeln;

for i:=1 to n do begin

for j:=1 to m do begin

b[i, j]: =a[i, j]*3;

write (b[i, j]:4); end;

writeln; end;

rewrite (p);

for i:=1 to n do begin

for j:=1 to m do

write (p, b[i, j],' ');

writeln (p); end; close (p);

end.

Исходные данные:

Полученные данные:

Вывод:

Данный способ работы с файлами достаточно прост и удобен. В данном случае выполняется всего одно действие (умножение на 3). Данные файла преобразуются не вручную, а автоматически.

2. Решение иррациональных уравнений

2.1 Метод хорд

Теоретические сведения:

Будем искать корень функции. Выберем две начальные точки (;) и (;) и проведем через них прямую. Она пересечет ось абсцисс в точке (;0). Теперь найдем значение функции с абсциссой. Временно будем считать корнем на отрезке [;]. Пусть точка имеет абсцисcу и лежит на графике. Теперь вместо точек и мы возьмём точку и точку. Теперь с этими двумя точками проделаем ту же операцию и так далее, то есть будем получать две точки и и повторять операцию с ними. Отрезок, соединяющий последние 2 точки, пересекает ось абсцисс в точке, значение абсциссы которой можно приближённо считать корнем. Эти действия нужно повторять до тех пор, пока не получим значение корня с нужным приближением.

Рис. 1 Иллюстрация метода хорд

Рисунок иллюстрирует работу метода хорд. В данном случае вторая производная функции положительна, поэтому в качестве начального приближения выбрана точка хо = b. Как видно из рисунка, метод имеет очень быструю сходимость среди всех методов решения нелинейных уравнений: обычно заданная точность достигается за 2−3 итерации.

Блок схема

Текст задачи

var a, b, e, c, cc, ee: real;

k:integer;

function f (x: real):real;

begin

f:=x*x-2*x;

end;

begin

a:=1;

b:=4;

ee:=0.001;

repeat

cc:=c;

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

k:=k+1;

if f (a)*f (c)<0 then b:=c else a:=c;

e:=abs ((c-cc)/c);

writeln (k,' ', c:0:5,' ', e:0:5,' ', f (c):0:5);

until e<=ee;

writeln;

end.

Результаты

а=1, b=4

Вывод

Достоинство метода хорд: очень быстрая сходимость по сравнению с методом половинного деления к заданной точности.

2.2 Метод половинного деления

Теоретические сведения:

Метод половинного деления при нахождении корня уравнения f (x)=0 состоит в делении пополам отрезка [a; b], где находится корень. Затем анализируется изменение знака функции на половинных отрезках, и одна из границ отрезка [a; b] переносится в его середину. Переносится та граница, со стороны которой функция на половине отрезка знака не меняет. Далее процесс повторяется. Итерации прекращаются, когда длина интервала [a; b] становится меньше заданной погрешности нахождения корня.

Разделим отрезок [a, b] пополам точкой. Если (что практически наиболее вероятно), то возможны два случая: либо f(x) меняет знак на отрезке [a, c] (Рис. 2), либо на отрезке [c, b] (Рис. 3)

Рис. 2

Рис. 3

Блок схема

Текст задачи

var a, b, e, c, cc, ee: real;

k:integer;

function f (x: real):real;

begin

F:=x*x-2*x;

end;

begin

a:=1;

b:=4;

ee:=0.001;

repeat

cc:=c;

c:=(a+b)/2;

k:=k+1;

if f (a)*f (c)<0 then b:=c else a:=c;

e:=abs ((c-cc)/c);

writeln (k,' ', c:0:5,' ', e:0:5,' ', f (c):0:5,' ');

until e<=ee;

writeln;

end.

Результаты:

a=1;b=4

Вывод:

Недостаток метода половинного деления: медленная сходимость по сравнению с методом хорд к заданной точности. Достоинство: небольшой алгоритм вычисления значения функции.

3. Вычисление определенного интеграла

3.1 Метод парабол

Теоретические сведения

Метод парабол — метод вычисления корней многочлена с комплексными коэффициентами, и позволяет найти все корни многочлена без предварительной информации о начальном приближении.

Пусть функция y = f (x) непрерывна на отрезке [a; b] и нам требуется вычислить определенный интеграл. Разобьем отрезок [a; b] на n элементарных отрезков длины 2h=(b-a)/n точками [x (2i-2);x (2i)]. Пусть точки x (2i-2) (i=1,2,3.n) являются серединами отрезков [x (2i-2);x (2i)] соответственно. В этом случае узлы определяются равенством x (i)=a+h*i. Суть метода Симпсона (парабол) заключается в том, что на каждом интервале [x (2i-2);x (2i)] подынтегральная функция приближается квадратичной параболой y=a (i)*x*x+b (i)*x+c (i), проходящей через точки .(x (2i-2);f (x (2i-2)), .(x (2i-1);f (x (2i-1)), .(x (2i);f (x (2i)). Отсюда и название метода. Геометрически это выглядит так:

Рис. 5 Рис. 6 Рис. 7

Графическая иллюстрация метода Симпсона. Красной линией изображен график функции y = f (x), синей линией показано приближение графика функции y = f (x) квадратичными параболами на каждом элементарном отрезке разбиения.

Рис. 8

Блок схема

begin

Текст задачи

const

a=1;

b=2;

ee=0.1;

var e, z, h, ss, s, x:real;

k3:integer;

function y (x:real):real;

begin

y:=exp (-x)-sqr (x-1)+1;

end;

begin

h:=b-a;

repeat

k3:=k3+1;

ss:=s;

h:=h/2;

s:=y (a)+y (b);

x:=a+h;

z:=-1;

while x

begin

z:=-z;

s:=s+y (x)*(3+z);

x:=x+h;

end;

s:=s*h/3;

e:=abs ((ss-s)/s);

writeln ('k3=', k3,' ','s=', s:0:7,' ','e=', e:0:7);

until

e<=ee;

end.

Результаты

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

3.2 Метод прямоугольников

Теоретические сведения Пусть функция y = f (x) непрерывна на отрезке [a; b]. Нам требуется вычислить определенный интеграл. Так же как в методе парабол разбиваем отрезки.

Суть метода прямоугольников заключается в том, что в качестве приближенного значения определенного интеграла берут интегральную сумму. Если отрезок интегрирования [a; b] разбить на равные части длины h точками, и в качестве точек выбрать середины элементарных отрезков, то приближенное равенство

можно записать в виде

.

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

называют шагом разбиения отрезка [a; b]. Приведем графическую иллюстрацию.

Рис. 9 Из чертежа видно, что подынтегральная функция y = f (x) приближается кусочной ступенчатой функцией

на отрезке интегрирования. С геометрической точки зрения для неотрицательной функции y = f (x) на отрезке [a; b] точное значение определенного интеграла представляет собой площадь криволинейной трапеции, а приближенное значение по методу прямоугольников — площадь ступенчатой фигуры.

Рис. 10 Рис. 11 Перейдем к оценке абсолютной погрешности метода прямоугольников. Сначала оценим погрешность на элементарном интервале. Погрешность метода прямоугольников в целом будет равна сумме абсолютных погрешностей на каждом элементарном интервале.

Блок схема

Текст задачи

const

aa=1;

bb=2;

nn=1;

ee=0.1;

var a, b, s, ss, e, deltax, x: real;

k, n: integer;

function f (x:real):real;

begin

f:=exp (-x)-sqr (x-1)+1;

end;

begin

a:=aa; b:=bb; n:=nn;

repeat begin ss:=s;

k:=k+1;

s:=0;

n:=n*2;

deltax:=(b-a)/n;

x:=a;

while x

s:=s+f (x)*deltax;

x:=x+deltax;

end;

e:=abs ((ss-s)/s);

writeln ('k= ', k,' ','s= ', s:0:7,' ','e= ', e:0:7);

end;

until e

end.

Результаты

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

3.3 Метод Трапеций

Теоретические сведения Пусть нам требуется вычислить определенный интеграл, где y = f (x) непрерывна на отрезке [a; b]. Разобьем отрезок [a; b] на n равных интервалов длины h точками. В этом случае шаг разбиения определяется так же как в методе парабол. Рассмотрим функцию на элементарных отрезках [x (i-1);x (i)]. Возможны четыре случая (на рисунке показаны простейшие из них, к которым все сводится при увеличении n):

На каждом отрезке [x (i-1);x (i)] заменим функцию y = f (x) отрезком прямой, проходящей через точки и. На рисунке показаны синими линиями:

В качестве приближенного значения интеграла возьмем выражение

Давайте выясним, что означает в геометрическом смысле записанное приближенное равенство. В первом случае площадь криволинейной трапеции приближенно равна площади трапеции с основаниями f (x (i-1)), f (x (i)) и высотой h, в последнем случае определенный интеграл приближенно равен площади трапеции с основаниямиf (x (i-1)), -f (x (i))и высотой h, взятой со знаком минус. Во втором и третьем случаях приближенное значение определенного интеграла равно разности площадей красной и синей областей, изображенных на рисунке ниже.

Рис. 11

Рис. 12

Блок схема

Текст задачи

Const aa=1; bb=2; ee=0.1; nn=1;

var a, b, s, ss, e, deltax, x: real;

k, n: integer;

function f (x:real):real;

begin

f:=exp (-x)-sqr (x-1)+1;

end;

begin

a:=aa; b:=bb; n:=nn;

repeat ss:=s;

k:=k+1;

s:=0;

n:=n*2;

DELTAx:=(b-a)/n;

x:=a;

while x

s:=s+deltax*(f (x)+f (x+deltax))/2; x:=x+deltax;

end;

e:=abs ((ss-s)/s);

writeln ('k=', k,' ','s=', s:0:7,' ','e=', e:0:7);

until e

end.

Результаты

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

4. Вывод о трех методах

Каждый способ приближённого решения определённого интеграла имеет свои преимущества и недостатки, в зависимости от поставленной задачи следует использовать конкретные методы. Если необходимо быстро получить решение, но нет необходимости в большой точности ответа, следует воспользоваться одним из методов прямоугольника. Если же необходимо получить наиболее точный результат, идеально подходит метод Симпсона. Метод трапеций даёт ответ более точный, чем метод прямоугольников, но методу Симпсона он сильно уступает, этот метод можно назвать «золотой серединой» между двумя другими.

Решение систем линейных алгебраических уравнений

4.1 Метод Гаусса

Теоретические сведения Метод Гаусса прекрасно подходит для решения систем линейных алгебраических уравнений (СЛАУ). Он обладает рядом преимуществ по сравнению с другими методами:

· во-первых, нет необходимости предварительно исследовать систему уравнений на совместность;

· во-вторых, методом Гаусса можно решать не только СЛАУ, в которых число уравнений совпадает с количеством неизвестных переменных и основная матрица системы невырожденная, но и системы уравнений, в которых число уравнений не совпадает с количеством неизвестных переменных или определитель основной матрицы равен нулю;

· в-третьих, метод Гаусса приводит к результату при сравнительно небольшом количестве вычислительных операций.

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

Записываем эту совокупность в матричном виде, и можно производить некоторые действия для получения неизвестных переменных:

поменять местами два уравнения,

умножить обе части какого-либо уравнения на произвольное и отличное от нуля действительное (или комплексное) число k,

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

Блок схема

Текст задачи

var

t:text;

a:array[1.1000,1.1000] of real;

c:real;

n, m, i, j, k:integer;

x:array[1.1000] of real;

begin

assign (t, 'H:Phoenix3.txt');

reset (t);

while not eof (t) do

begin

i:=i+1;

j:=0;

while not eoln (t) do

begin

j:=j+1;

read (t, a[i, j]);

end;

readln (t);

end;

n:=i;

for i:=1 to n do begin

for j:=1 to n+1 do

write (a[i, j],' ');writeln; end;

writeln;

for i:= 1 to n-1 do

for j:=i+1 TO n do begin

c:=a[i, i]/a[j, i];

for k:=i to n+1 do

a[j, k]: =a[i, k]-a[j, k]*c;

end;

for i:=1 to n do begin

for j:=1 to n+1 do

write (a[i, j],' ');writeln; end;

writeln;

X[n]: =a[n, n+1]/a[n, n];

write (X[n]:0:6,' ');

for i:= n-1 downto 1 do begin

x[i]: =a[i, n+1]/a[i, i];

for j:=i+1 to n do

x[i]: =x[i]-a[i, j]*x[j]/a[i, i];

write (x[i]:0:7,' ');

end;

close (t);

end.

Результаты

Корни уравнений данной системы: х1= 2; х2= -1; х3= 3; х4= 1;

4.2 Метод Ньютона

Теоретические сведения Чтобы численно решить уравнение методом простой итерации, его необходимо привести к следующей форме:, где — сжимающее отображение.

Для наилучшей сходимости метода в точке очередного приближения должно выполняться условие. Решение данного уравнения ищут в виде

тогда:

В предположении, что точка приближения «достаточно близка» к корню, и что заданная функция непрерывна, окончательная формула для такова:

С учётом этого функция определяется выражением:

Эта функция в окрестности корня осуществляет сжимающее отображение, и алгоритм нахождения численного решения уравнения сводится к итерационной процедуре вычисления:

Рис. Иллюстрация метода Ньютона

уравнение алгебраическое интеграл фурье Рис иллюстрирует работу метода Ньютона. В данном случае вторая производная функции положительна, поэтому в качестве начального приближения выбрана точка хо = b. Как видно из рисунка, метод имеет очень быструю сходимость среди всех методов решения нелинейных уравнений: обычно заданная точность достигается за 2−3 итерации

Блок схема Текст задачи

var a: array[1.100,1.100] of integer;

t:text;

i, j, m, n, k:byte;

e, ed, alfa: real;d:boolean; x, xx: array[1.100] of real;

begin

ed:=0.001;

assign (t,'H:Phoenix4.txt');

reset (t);

while not eof (t) do

begin

i:=i+1;

j:=0;

while not eoln (t) do

begin

j:=j+1;

read (t, a[i, j]);

end;

readln (t);

end;

close (t);

n:=i;

m:=j;

for i:=1 to n do begin

for j:=1 to m do

write (a[i, j],' ');

writeln;

end;

repeat

xx:=x;

k:=k+1;

writeln ('k=', k,' ');

for i:=1 to n do

begin

x[i]: =a[i, n+1]/a[i, i];

for j:=1 to n do

if j<>i then x[i]: =x[i]-(a[i, j]*xx[j])/a[i, i];

write ('x[', i,']=', x[i]:0:2,' ');

end;

d:=false;

for i:=1 to n do

begin

alfa:=0;

for j:=1 to n do

alfa:=alfa+(a[i, j]*x[j]);

e:=abs ((alfa-a[i, n+1])/(a[i, n+1]));

write ('e=', e:2:5,' ');

if e>=ed then d:=true;

end;

writeln;

until not d; end.

Результаты

Корни уравнений х1=-1; х2=1;

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

4.3. Метод Зейделя

Теоретические сведения Метод Зейделя (иногда называемый методом Гаусса-Зейделя) является модификацией метода простой итерации, заключающейся в том, что при вычислении очередного приближения x^(k+1) его уже полученные компоненты x₁^(k+1), …, x_{i — 1}^(k+1) сразу же используются для вычисления x_i^(k+1).

В координатной форме записи метод Зейделя имеет вид:

x₁^(k+1) = c₁₁x₁^(k) + c₁₂x₂^(k) + … + c_1n-1x_n-1^(k) + c_1nx_n^(k) + d₁ x₂^(k+1) = c₂₁x₁^(k+1) + c₂₂x₂^(k) + … + c_2n-1x_n-1^(k) + c_2nx_n^(k) + d₂ … x_n^(k+1) = c_n1x₁^(k+1) + c_n2x₂^(k+1) + … + c_nn-1x_n-1^(k+1) + c_nnx_n^(k) + d_n где x⁽⁰⁾ — некоторое начальное приближение к решению.

Таким образом i-тая компонента (k+1)-го приближения вычисляется по формуле

x_i^(k+1) =? _j=1^i-1 c_ijx_j^(k+1) +? ⁿ_j=i c_ijx_j^(k) + d_i, i = 1, …, n

Условие окончания итерационного процесса Зейделя при достижении точности е в упрощенной форме имеет вид:

|| x ^(k+1) — x ^(k) ||? е.

Блок схема

Текст задачи

var a: array[1.100,1.100] of integer;

t:text;

i, j, m, n, k:byte;

e, ed, alfa: real;d:boolean; x: array[1.100] of real;

begin

ed:=0.001;

assign (t,'H:Phoenix4.txt');

reset (t);

while not eof (t) do

begin

i:=i+1;

j:=0;

while not eoln (t) do

begin

j:=j+1;

read (t, a[i, j]);

end;

readln (t);

end;

close (t);

n:=i;

m:=j;

for i:=1 to n do begin

for j:=1 to m do

write (a[i, j],' ');

writeln;

end;

repeat

k:=k+1;

writeln ('k=', k,' ');

for i:=1 to n do

begin

x[i]: =a[i, n+1]/a[i, i];

for j:=1 to n do

if j<>i then x[i]: =x[i]-(a[i, j]*x[j])/a[i, i];

write ('x[', i,']=', x[i]:0:5,' ');

end;

d:=false;

for i:=1 to n-1 do

begin

alfa:=0;

for j:=1 to n do

alfa:=alfa+(a[i, j]*x[j]);

e:=abs ((alfa-a[i, n+1])/(a[i, n+1]));

write ('e=', e:2:9,' ');

if e>=ed then d:=true;

end;

writeln;

until not d;

end.

Результаты

Корни уравнений х1=2; х2=-1; х3=0; х4=4;

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

5. Решение дифференциальных уравнений

5.1 Метод Эйлера

Теоретические сведения. Решить дифференциальное уравнение у^/=f (x, y) численным методом — это значит для заданной последовательности аргументов х₀, х₁…, х_n и числа у₀, не определяя функцию у=F (x), найти такие значения у₁, у₂,…, у_n, что у_i=F (x_i)(i=1,2,…, n) и F (x₀)=y₀. Таким образом, численные методы позволяют вместо нахождения функции У=F (x) получить таблицу значений этой функции для заданной последовательности аргументов. Величина h=x_k-x_k-1 называется шагом интегрирования.

Текст задачи

const

x0=0;

y0=1;

n=5;

Xn=1;

deltax=(Xn-x0)/n;

var x, y: array[1.100] of real;

i, j: integer;

k:integer;

t:text;

function f (x, y: real):real;

begin

f:=2*y-2*x+1;

{f:=exp (2*x)+x;}

end;

begin

assign (t,'H:Phoenixeyler.txt');

rewrite (t);

x[1]: =x0; y[1]: =y0;

writeln (x[1]:2:5,' ', y[1]: 2:5);

for i:=2 to n+1 do begin

x[i]: =x[i-1]+deltax;

y[i]:=y[i-1]+f (x[i-1], y[i-1])*deltax;

writeln (x[i]:2:5,' ', y[i]: 2:5);

end;

for i:=1 to n+1 do

writeln (t, x[i]: 2:5,' ', y[i]: 2:5);

close (t);

end.

Результаты

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

5.2 Метод Рунге-Кутта

Теоретические сведения

1. Эти методы являются одноступенчатыми: чтобы найти уm+1, нужна информация о предыдущей точке xm, ym.

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

3. Они не требуют вычисления производных от f (x, y), а требуют вычисления самой функции.

Текст программы

const

x0=0;

y0=1;

n=5;

Xn=1;

deltax=(Xn-x0)/n;

var x, y: array[1.100] of real;

i, j: integer;

k:integer;

t:text;

k1,k2,k3,k4,deltay:real;

function f (x, y: real):real;

begin

f:=2*y-2*x+1;

{f:=exp (2*x)+x;}

end;

begin

assign (t,'H:Phoenixrunge-kut.txt');

rewrite (t);

x[1]: =x0; y[1]: =y0;

writeln (x[1]:2:5,' ', y[1]: 2:5);

for i:=2 to n+1 do begin

x[i]: =x[i-1]+deltax;

k1:=deltax*f (x[i-1], y[i-1]);

k2:=deltax*f (x[i-1]+deltax/2,y[i-1]+k½);

k3:=deltax*f (x[i-1]+deltax/2,y[i-1]+k2/2);

k4:=deltax*f (x[i-1]+deltax, y[i-1]+k3);

deltay:=(k1+2*k2+2*k3+k4)/6;

y[i]:=y[i-1]+deltay;

writeln (x[i]:2:5,' ', y[i]: 2:5);

end;

for i:=1 to n+1 do

writeln (t, x[i]: 2:5,' ', y[i]: 2:5);

close (t);

end.

Результаты

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

Общий вывод

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

6. Ряды Фурье

Теоретические сведения Ряд Фурье является частным случаем функциональных рядов. Функциональным рядом называется выражение вида где U₁(x), U₂(x), …, U_k(x), … — функции, зависящие от одной переменной х или от нескольких переменных (тогда через х обозначают точку в пространстве: х = (х₁, …, x_n)).

Текст задачи

const

ee=1e-3;

a=0;

b=3;

period=b-a;

n=20;

id=1e-3;

var k: integer;

ex:text;

a0,t, dt, st, Ampl, phi, i: real;

ak, bk: array [1.n] of real;

function f (t:real):real;

begin

if (t>=0) and (t<1) then f:=t;

if (t>=1) and (t<2) then f:=1;

if (t>=2) then f:=6−2*t;

end;

function si (c, k: integer):real;

var t, dt, e, ss, s: real; z: integer;

begin

dt:=(b-a)/50;

repeat

ss:=s;

dt:=dt/2;

case c of

1: s:=f (a)+f (b);

2: s:=f (a)*cos (2*Pi*k*a/period)+f (b)*cos (2*Pi*k*b/period);

3: s:=f (a)*sin (2*Pi*k*a/period)+f (b)*sin (2*Pi*k*b/period);

4: s:=sqr (f (a))+sqr (f (b));

end;

t:=a+dt;

z:=-1;

while t

begin

z:=-z;

case c of

1: s:=s+f (t)*(3+z);

2: s:=s+f (t)*cos (2*Pi*k*t/period)*(3+z);

3: s:=s+f (t)*sin (2*Pi*k*t/period)*(3+z);

4: s:=s+sqr (f (t))*(3+z);

end;

t:=t+dt;

end;

s:=2*(s*dt/3)/period;

if abs (s)

begin

s:=0; e:=0; end

else e:=abs ((ss-s)/s);

until e<=ee; si:=s;

end; begin

{assign (ex, 'H:PhoenixAfi.txt');

rewrite (ex);} a0:=si (1,0);

writeln ('a0=', a0:0:7);

{for k:=1 to n do}

while i<=id do

begin

ak[k]: =si (2,k);

bk[k]:=si (3,k);

Ampl:=sqrt (sqr (ak[k])+sqr (bk[k]));

if (ak[k]>0) and (bk[k]>0) then phi:=arctan (ak[k]/bk[k]);

if (ak[k]>0) and (bk[k]<0) then phi:=arctan (ak[k]/bk[k])+pi;

if (ak[k]>0) and (bk[k]=0) then phi:=pi/2;

if (ak[k]<0) and (bk[k]>0) then phi:=arctan (ak[k]/bk[k])+2*pi;

if (ak[k]<0) and (bk[k]<0) then phi:=arctan (ak[k]/bk[k])+pi;

if (ak[k]<0) and (bk[k]=0) then phi:=3*pi/2;

if (ak[k]=0) and (bk[k]>0) then phi:=0;

if (ak[k]=0) and (bk[k]<0) then phi:=pi;

if (ak[k]=0) and (bk[k]=0) then writeln ('net');

writeln ({ex,}k:2,ak[k]: 8:4,bk[k]:8:4,Ampl:6:2,phi*180/pi:6:0);

{writeln (k:2,ak[k]:8:4,bk[k]:8:4,Ampl:6:2,phi*180/pi:6:0);}

end;

dt:=(b-a)/(40*n);

t:=a;

while t<=b do begin

st:=a0/2;

for k:=1 to n do

st:=st+ak[k]*cos (2*pi*k*t/period)+bk[k]*sin (2*pi*k*t/period);

writeln (t:0:5,' ', f (t):0:5,' ', st:0:5);

t:=t+dt; end; {close (ex);} end.

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

1. К. Ланцош, Практические методы прикладного анализа, 1961.

2. И. С. Березин и Н. П. Жидков, Методы вычислений в 2 томах, 1959;1962.

3. Демидович Б. П. и Марон И. А. Основы вычислительной математики. — Наука, 1970. — С. 664.

4. Бахвалов, Жидков, Кобельков Численные методы. — Наука.

5. Каханер Д., Моулер К., Нэш С. Численные методы и программное обеспечение (пер. с англ.). М.: Мир, 2001, 575 c.

6. Самарский А. А., Гулин А. В. Численные методы: Учеб. пособие для вузов. — М.: Наука. Гл. ред. физ-мат. лит., 1989. — 432 с.

7. Пискунов Н. С. Дифференциальное и интегральное исчисления для вузов. — 13-е изд. — М.: Наука. Гл. ред. физ-мат. лит., 1985. — 432 с.

8. И. П. Натансон: Краткий курс высшей математики

9. И. И. Валуцэ, Г. Д. Дилигун: Математика для техникумов

10. И. А. Сахарников: Высшая математика

11. П. П. Коровнин: Математический анализ

12. Л. И. Лихтарников, А. Н. Поволоцкий: основы математического анализа

13. http://www.ph4s.ru/book_pc_chisl.html