Численные методы в инженерных расчетах

Размещено на http://www.allbest.ru/

СОДЕРЖАНИЕ

ВВЕДЕНИЕ

1. РЕШЕНИЕ НЕЛИНЕЙНЫХ УРАВНЕНИЙ

1.1 Задание

1.2 Решение уравнения методом половинного деления

1.3 Метод секущих

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

1.5 Метод Matlab

2. РЕШЕНИЕ СИСТЕМ ЛИНЕЙНЫХ АЛГЕБРАИЧЕСКИХ УРАВНЕНИЙ

2.1 Задание

2.2 Решение СЛАУ

3. ИНТЕРПОЛЯЦИЯ ФУНКЦИЙ

3.1 Задание

3.2 Интерполяция заданной функции

4. РЕШЕНИЕ СИСТЕМЫ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ

4.1 Задание

4.2 Решение системы О.Д.У.

СПИСОК ЛИТЕРАТУРЫ

Размещено на http://www.allbest.ru/

ВВЕДЕНИЕ

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

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



1. РЕШЕНИЕ НЕЛИНЕЙНЫХ УРАВНЕНИЙ

1.1 Задание

нелинейный алгебраический дифференциальный уравнение

Определить корни графически.

Уточнить один из корней с точностью до 10^(-6) методами:

а) методом половинного деления;

б) методом секущих;

в) методом Ньютона;

г) используя стандартную функцию Matlab.

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

Визуализировать итерационную последовательность.

Заданные уравнения:

2ex-2x-3=0

(x-2)cos(x)=1

-2П<=x<=2П

3x4+8x3+6x2-10=0

1.2 Решение уравнения методом половинного деления

В методе половинного деления предполагается, что функция f(x) прерывна на отрезке [a,b] и на концах отрезка имеет разные знаки: f(a)f(b)<0.

Определяется середина отрезка: x0=(a+b)/2.

Если f(a)f(x0)<0, то [a1,b1]=[a,x0]; если f(a)f(x0)>0, то [a1,b1]=[x0,b].

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

Метод половинного деления медленный, но всегда сходится. На n шаге отрезок уменьшится в раз:

,

где - точность вычисления.

Программа в Matlab для уравнения 2ex-2x-3=0

Во - первых, вводим функцию в Matlab следующим образом:

function y=myf1(x)

y=2*exp(x)-2*x-3;

Теперь строим график данной функции:

x=0:0.001:3;

y=myf1(x);

plot (x,y);grid



Рис.1. График функции 2ex-2x-3=0

Программа для определения корня уравнения

clc;

a=0.5;

b=1;

eps=10^(-6);

l=b-a;

m=1;

format long

while l>=eps;

x=(a+b)/2;

c(m)=x;

if myf1(a)*myf1(x)<0

b=x;

else

a=x;

end

m=m+1;

l=b-a;

end

disp('число итераций длЯ уточнениЯ корнЯ');

disp(m);

disp('последовательное приближение')

disp(c);

plot(c,'*');

Получен результат: число итераций длЯ уточнениЯ корнЯ

последовательное приближение

Columns 1 through 4

0.75000000000000 0.87500000000000 0.81250000000000 0.84375000000000

Columns 5 through 8

0.85937500000000 0.85156250000000 0.85546875000000 0.85742187500000

Columns 9 through 12

0.85839843750000 0.85791015625000 0.85766601562500 0.85778808593750

Columns 13 through 16

0.85772705078125 0.85769653320313 0.85768127441406 0.85767364501953

Columns 17 through 19

0.85767745971680 0.85767555236816 0.85767650604248



Рис.2. Итерационная последовательность

Программа в Matlab для уравнения (x-2)cos(x)=1.

Во - первых, вводим функцию в Matlab следующим образом:

function y=myf2(x)

y=(x-2).*cos(x)-1;

Теперь строим график данной функции:

x=-2*pi:0.001:2*pi;

y=myf2(x);

plot (x,y);grid



Рис.3. График функции y= (x-2)cos(x)-1

Программа для определения корня уравнения

clc;

a=-6;

b=-4;

eps=10^(-6);

l=b-a;

m=1;

format long

while l>=eps;

x=(a+b)/2;

c(m)=x;

if myf2(a)*myf2(x)<0

b=x;

else

a=x;

end

m=m+1;

l=b-a;

end

disp('число итераций длЯ уточнениЯ корнЯ');

disp(m);

disp('последовательное приближение')

disp(c);

plot(c,'*');

Получен результат:

число итераций длЯ уточнениЯ корнЯ

22

последовательное приближение

Columns 1 through 4

-5.00000000000000 -4.50000000000000 -4.75000000000000 -4.62500000000000

Columns 5 through 8

-4.56250000000000 -4.53125000000000 -4.54687500000000 -4.55468750000000

Columns 9 through 12

-4.55859375000000 -4.56054687500000 -4.55957031250000 -4.55908203125000

Columns 13 through 16

-4.55932617187500 -4.55944824218750 -4.55938720703125 -4.55935668945313

Columns 17 through 20

-4.55934143066406 -4.55933380126953 -4.55933761596680 -4.55933952331543

Column 21

-4.55933856964111



Рис.4. Итерационная последовательность

Программа в Matlab для уравнения 3x4+8x3+6x2-10=0

Во - первых, вводим функцию в Matlab следующим образом:

function y=myf3(x)

y=3*x.^4+8*x.^3+6*x.^2-10;

Теперь строим график данной функции:

x=1:0.001:3;

y=myf3(x);

plot (x,y);grid



Рис.5. График функции y=3x4+8x3+6x2-10

Программа для определения корня уравнения

clc;

a=0.5;

b=1;

eps=10^(-6);

l=b-a;

m=1;

format long

while l>=eps;

x=(a+b)/2;

c(m)=x;

if myf3(a)*myf3(x)<0

b=x;

else

a=x;

end

m=m+1;

l=b-a;

end

disp('число итераций длЯ уточнениЯ корнЯ');

disp(m);

disp('последовательное приближение')

disp(c);

plot(c,'*');

Получен результат:

число итераций длЯ уточнениЯ корнЯ

20

последовательное приближение

Columns 1 through 4

0.75000000000000 0.87500000000000 0.81250000000000 0.84375000000000

Columns 5 through 8

0.82812500000000 0.82031250000000 0.82421875000000 0.82617187500000

Columns 9 through 12

0.82519531250000 0.82568359375000 0.82592773437500 0.82604980468750

Columns 13 through 16

0.82598876953125 0.82601928710938 0.82603454589844 0.82604217529297

Columns 17 through 19

0.82603836059570 0.82604026794434 0.82604122161865



Рис.6. Итерационная последовательность

1.3 Метод секущих

В методе секущих предполагается, что функция f(x) прерывна на отрезке [a,b] и на концах отрезка имеет разные знаки: f(a)f(b)<0.Предположим, что

f(a)>0, f(b)<0.

Уравнение секущей:

.

Тогда

.



Получен новый отрезок [a1,b1]=[a,x0]. Вычисления производятся до тех пор, пока значение x0 не приблизится к искомому корню с заданной точностью.

В ряде случаев алгоритм имеет быструю сходимость.

Программа в Matlab для уравнения 2ex-2x-3=0

clc;

a=0.5;

b=1;

eps=10^(-6);

x=a-((b-a)/((myf1(b)-myf1(a))))*myf1(a);

m=1;

format long

while abs(myf1(x))>eps;

c(m)=x;

if myf1(x)*myf1(a)<0;

b=x;

else

a=x;

end

m=m+1;

x=a-((b-a)/((myf1(b)-myf1(a))))*myf1(a);

end

disp('число итераций длЯ уточнениЯ корнЯ');

disp(m);

disp('последовательное приближение')

disp(c);

plot(c,'*');

Получен результат:

число итераций длЯ уточнениЯ корнЯ

7

последовательное приближение

Columns 1 through 4

0.80837697986152 0.85188251839544 0.85700982100534 0.85760011315751

Columns 5 through 6

0.85766788654576 0.85767566538909

Рис.7. Итерационная последовательность

Программа в Matlab для уравнения (x-2)cos(x)=1

clc;

a=-6;

b=-4;

eps=10^(-6);

x=a-((b-a)/((myf2(b)-myf2(a))))*myf2(a);

m=1;

format long

while abs(myf2(x))>eps;

c(m)=x;

if myf2(x)*myf2(a)<0;

b=x;

else

a=x;

end

m=m+1;

x=a-((b-a)/((myf2(b)-myf2(a))))*myf2(a);

end

disp('число итераций длЯ уточнениЯ корнЯ');

disp(m);

disp('последовательное приближение')

disp(c);

plot(c,'*');

Получен результат:

число итераций длЯ уточнениЯ корнЯ

4

последовательное приближение

-4.50362928795517 -4.56127779222470 -4.55931106364103



Рис.8. Итерационная последовательность

Программа в Matlab для уравнения 3x4+8x3+6x2-10=0.

clc;

a=0.5;

b=1;

eps=10^(-6);

x=a-((b-a)/((myf3(b)-myf3(a))))*myf3(a);

m=1;

format long

while abs(myf3(x))>eps;

c(m)=x;

if myf3(x)*myf3(a)<0;

b=x;

else

a=x;

end

m=m+1;

x=a-((b-a)/((myf3(b)-myf3(a))))*myf3(a);

end

disp('число итераций длЯ уточнениЯ корнЯ');

disp(m);

disp('последовательное приближение')

disp(c);

plot(c,'*');

Получен результат:

число итераций длЯ уточнениЯ корнЯ

10

последовательное приближение

Columns 1 through 4

0.75545851528384 0.81290500229492 0.82367706378330 0.82561839681742

Columns 5 through 8

0.82596568870814 0.82602773448238 0.82603881669747 0.82604079604705

Column 9

0.82604114956799



Рис.9. Итерационная последовательность

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

Требуется, чтобы функция f(x) была дифференцируема в некоторой окрестности точки с (искомого корня уравнения). Производная не должна менять свой знак на промежутке [c, x0].

Уравнение касательной

.

Тогда



.

Алгоритм продолжается до тех пор, пока не будет достигнута заданная точность

.

Программа в Matlab для уравнения 2ex-2x-3=0

clc;

format long

eps=0.00000001;

m=0.5;

n=1;

while abs(myf1(m))>eps;

m=m-myf1(m)/(2*exp(m)-2);

c(n)=m;

n=n+1;

end

disp('число итерации длЯ уточнениЯ корнЯ');

disp (n);

disp('последовательное приблежение');

disp(c);

plot(c,'*');

Получен результат:

число итерации длЯ уточнениЯ корнЯ

6

последовательное приблежение

Columns 1 through 4

1.04149408253680 0.88225679128705 0.85818785708263 0.85767690070957

Column 5

0.85767667394594

Рис.10. Итерационная последовательность

Программа в Matlab для уравнения (x-2)cos(x)=1

clc;

format long

eps=0.00000001;

m=-6;

n=1;

while abs(myf2(m))>eps;

m=m-myf2(m)/(cos(m)-(m-2)*sin(m));

c(n)=m;

n=n+1;

end

disp('число итерации длЯ уточнениЯ корнЯ');

disp (n);

disp('последовательное приблежение');

disp(c);

plot(c,'*');

Получен результат:

число итерации длЯ уточнениЯ корнЯ

5

последовательное приблежение

-4.75166784825219 -4.56522808749884 -4.55934580524512 -4.559337734920

Рис.11. Итерационная последовательность



Программа в Matlab для уравнения 3x4+8x3+6x2-10=0

clc;

format long

eps=0.00000001;

m=0.5;

n=1;

while abs(myf3(m))>eps;

m=m-myf3(m)/(12*m.^3+24*m.^2+12*m);

c(n)=m;

n=n+1;

end

disp('число итерации длЯ уточнениЯ корнЯ');

disp (n);

disp('последовательное приблежение');

disp(c);

plot(c,'*');

Получен результат:

число итерации длЯ уточнениЯ корнЯ

6

последовательное приблежение

Columns 1 through 4

1.04166666666667 0.86731136332084 0.82789983686378 0.82604519906920

Column 5

0.82604122645581



Рис.12. Итерационная последовательность

1.5 Метод Matlab

для уравнения 2ex-2x-3=0:

fzero('myf1',0.5)

ans =

0.85767667394590

для уравнения (x-2)cos(x)=1

fzero('myf2',-4)

ans =

-4.55933773490524

для уравнения 3x4+8x3+6x2-10=0

fzero('myf3',1)

ans =

0.82604122643762135



2. РЕШЕНИЕ СИСТЕМ ЛИНЕЙНЫХ АЛГЕБРАИЧЕСКИХ УРАВНЕНИЙ

2.1 Задание

Записать систему линейных уравнений 5х5.

Проверить обусловленность систем:

а) вычислить определитель;

б) найти обратную матрицу;

в) вычислить числа обусловленности системы в трёх нормах.

Решить систему методом Зейделя.

Решить систему методом Matlab.

2.2 Решение СЛАУ

При большом числе уравнений (100 и более) прямые методы решения СЛАУ (за исключением метода прогонки) становятся труднореализуемыми на ЭВМ из-за сложности хранения и обработки матриц большой размерности. Здесь разумнее применять итерационные методы.

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

Метод Зейделя состоит в том, что итерации производятся по формуле:

,

где -произвольны, i=1, 2, …, n; k=1, 2,…

Итерации по методу Зейделя отличаются от простых итераций тем, что при нахождении i-й компоненты k-го приближения сразу используются уже найденные компоненты k-го приближения с меньшими номерами.

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

Исходная система линейных уравнений 5x5

-+4+8+6+4=6

-5+2-8+5+4=-5

4+-5++5+10=3

-4+3+-7+5=5

5+10+7-4-3=3

Программа в Matlab:

clc

a=[-1 4 8 6 4; -5 2 -8 5 4;4 -5 1 5 10; -4 3 1 -7 5;5 10 7 -4 -3]

b=[6; -5; 3; 5; 3]

w=det(a); %определитель матрицы

a1=inv(a); %обратная матрица

x=a\b%решение системы методом Matlab

masa=abs(a);

n1=max(sum(masa')); % норма 1

masa1=abs(a1);

n11=max(sum(masa1')); %норма 1 для обратной матрицы

m1=n1*n11

n2=max(sum(masa)) % норма 2

n22=max(sum(masa1)); %норма 2для обратной матрицы

m2=n2*n22;

r1=a*a';

la1=eig(r1);

A1=max(abs(la1));

n3=sqrt(A1); % норма 3

r2=a1*a1';

la2=eig(r2);

A2=max(abs(la2));

n33=sqrt(A2); %норма 3 для обратной матрицы

m3=n3*n33

Er=rand(5,5);

aa=a+Er;

Erb=rand(5,1);

bb=b+Erb;

X=aa\bb

absp=max(abs(x-X))

fds=absp/(max(abs(x)))

Получаем результат:

a =

-1 4 8 6 4

-5 2 -8 5 4

4 -5 1 5 10

-4 3 1 -7 5

5 10 7 -4 -3

b =

6

-5

3

5

3

x =

-0.22786051834609

-0.11319065451120

0.77419480408262

-0.24929348291513

0.38177614113214

m1 =

8.45727031382594

n2 =

27

m3 =

3.26291163081868

X =

-0.24060047867289

-0.04637852582081

0.79476339451552

-0.26209671096621

0.39532812354543

absp =

0.06681212869039

fds =

0.08629885958684

Функция Зейделя:

function [z1,z2]=zeide1(a,b,eps) % векторное присвоение

a=[-1 4 8 6 4; -5 2 -8 5 4;4 -5 1 5 10; -4 3 1 -7 5;5 10 7 -4 -3]

b=[6; -5; 3; 5; 3]

eps=10^(-6);

N=size(a,1) % возвращает число строк массива А

C=a'*a

D=a'*b % приведение системы к нормальному виду

for i=1:N % размерность задачи

D1(i)=D(i)/C(i,i);

end

% D11=D1'; % D11-начальное приближение (столбец правой части)

d1=D1;

for i=1:N

for j=1:N

if i==j % проверка на равенство

C1(i,j)=0;

else

C1(i,j)=-C(i,j)/C(i,i);

end

end

end

% решение системы методом Зейделя

R1=d1;

k=0;

while k==0

for i=1:N

v=C1(i,1:N); % вырезка

A=dot(v',d1); % скалярное произведение

d1(i)=A+D1(i);

R2=d1;

S=max(abs(R2-R1));

if S<eps

z1=d1; z2=S;

k=k+1;

end

R1=R2;

end

end

Получен результат:

a =

-1 4 8 6 4

-5 2 -8 5 4

4 -5 1 5 10

-4 3 1 -7 5

5 10 7 -4 -3

b =

6

-5

3

5

3

N =

5

C =

83 4 67 -3 -19

4 154 84 -52 -41

67 84 179 -22 -6

-3 -52 -22 151 71

-19 -41 -6 71 166

D =

26

44

117

-21

50

ans =

Columns 1 through 4

-0.22784165101666 -0.11317732900043 0.77418152802107 -0.24929480082327

Column 5

0.38178538703810



3. ИНТЕРПОЛЯЦИЯ ФУНКЦИЙ

3.1 Задание

Из первой части курсовой работы взять функцию, которая имеет больше изгибов. На промежутке (-4р;4р) взять равномерную сетку из десяти значений(x). Вычислить значения функции в этих точках. Сформировать сетку c шагом 1/100.

Выполнить интерполяцию 3 методами.

3.2 Интерполяция заданной функции

Часто для восстановления зависимостей; сравнения значений двух или более функций, заданных на разной сетке; подготовки качественных графических материалов требуется построить интерполирующую функцию f(xi)=yi, которую называют интерполянтой. Интерполянту ищут в виде многочлена:

где -неопределённые коэффициенты.

Более удобным является представление полинома в виде полинома Лагранжа степени n:

Формула Лагранжа имеет вид:

.



Формат функции в Matlab:

yi=interp(x,y,xi,'метод').

Метод принимает значения:

nearest;

linear;

spline;

При заданной сеточной функции в узлах x принимаются значения, определённые компонентами y. Функция interp осуществляет интерполяцию функции в точке с координатой xi.

Принципы этих методов можно проследить при решении конкретной задачи. Для интерполяции по формуле Лагранжа используются функции polyfit и polyval.

Программа в Matlab:

xi=-4*pi:0.01:4*pi;

yi=(xi-2).*cos(xi)-1;

plot(xi,yi,'r'); grid;

v=[-15 15 -20 20];

axis(v);

hold on;

x=-4*pi:8*pi/10:4*pi;

y=(x-2).*cos(x)-1;

yi=interp1(x,y,xi,'nearest');

plot(x,y,'ob',xi,yi,'y')

pause

yi=interp1(x,y,xi,'linear')

plot(x,y,'ob',xi,yi,'m')

pause

yi=interp1(x,y,xi,'spline')

plot(x,y,'ob',xi,yi,'g')

pause

p=polyfit(x,y,11)

f=polyval(p,xi)

plot(x,y,'ob',xi,f,'b')

hold off

Рис.13. Сетка из десяти значений и интерполяция методом `nearest'



Рис.14. Интерполяция методом `linear'

Рис.15. Интерполяция методом `spline'



Рис.16 Полином Лагранжа



4. РЕШЕНИЕ СИСТЕМЫ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ

4.1 Задание

Привести уравнения к нормальной системе.

Решить систему методом Рунге-Кутты.

Y/V- Y///(Y//+1)-sin (t)Y/+YY/=3

Y(0)=0

Y/0)=0

Y//(0)=1

Y///(0)=-1

4.2 Решение системы О.Д.У.

Обыкновенным дифференциальным уравнением называется уравнение, связывающее одну независимую переменную x, искомую функцию y(x) и её производные до n порядка включительно:

.

Система О.Д.У. относительно n неизвестных функций y1(x), y2(x),… , yn(x) yназывается нормальной, если каждое уравнение, входящее в систему, содержит в левой части производную первого порядка от соответствующей искомой функции, а правая часть является функцией от независимой переменной x и n искомых функций

Одним из методов решения системы О.Д.У. является метод Рунге-Кутты. Вычислительная схема МРК для 4-го порядка имеет вид:



Приведение системы к нормальному виду

Y/V- Y///(Y//+1)-sin (t)Y/+YY/=3

Y(0)=0

Y/0)=0

Y//(0)=1

Y///(0)=-1

;

x4=y///=x3/

x1/ = x2

x2/=x3

x3/= x4

x4/=-x4(x3+1)-sin(t)x2+x1x2-3

Программа в Matlab:



function u=difur(t,x)

u=[x(2);x(3);x(4);-x(4)*(x(3)+1)-sin(t)*x(2)+x(2)*x(3)-3];

clc;

[t,X]=ode45('difur',[0,pi],[0,0,1,-1]);

y=X(:,1);

subplot(311)

plot(t,y)

y1=X(:,2);

subplot(312)

plot(t,y1)

y2=X(:,3);

subplot(313)

plot(t,y2)

Получен результат:

Рис.17 Решение уравнения

Размещено на Allbest.ru