Хлыбов А.А., Хлыбова О.Н. MathCAD и MATLAB Начальный курс Учебно-методическое пособие

Подождите немного. Документ загружается.

Глава 2. MATLAB

d=b*b-4*a*c;

% проверим знак дискриминанта

if d<0

disp('действительных корней нет');

else

x1=(-b-sqrt(d))/(2*a)

x2=(-b+sqrt(d))/(2*a)

end

2. Определим рабочий или выходной день недели по его номеру

% вводим номер дня недели

n=input('n=');

switch n

case 1

disp('рабочий')

case 2

disp('рабочий')

case 3

disp('рабочий')

case 4

disp('рабочий')

case 5

disp('рабочий')

otherwise

disp('выходной')

end

3. Ввод вектора осуществляется с помощью функции strcat,

предназначенной для объединения строк в одну строку, int2str – преобразует

число в строку символов.

N=input('N=')

for i=1:N

Глава 2. MATLAB

x(i)=input(strcat('x(',int2str(i),')='));

end

4. Ввод матрицы

n=input('n=');

m=input('m=');

for i=1:n

for j=1:m

a(i,j)=input(strcat('a(',int2str(i),',',int2str(j),')='));

end

5. Найти максимальный элемент массива

% положим первый элемент массива максимальным, его номер = 1,

затем все элементы, начиная со второго сравниваем с первым в цикле.

N=input('N=')

for i=1:N

x(i)=input(strcat('x(',int2str(i),')='));

max=x(1);

Nmax=1;

for i=2:N

if x(i)>max

max=x(i);

Nmax=i;

end

max

6. Найти минимальный элемент последовательности

N=input('N=')

for i=1:N

x(i)=input(strcat('x(',int2str(i),')='));

Глава 2. MATLAB

min=x(1);

Nmin=1;

for i=2:N

if min>x(i)

min=x(i);

Nmin=i;

end

min

Nmin

7. Вычисление факториала

n=input('n=')

y=1;

for i=1:n

y=y*i;

end

8. Пример кода, являющегося частью функции magic.m,

генерирующего

магический квадрат M для нечётных значений размера

стороны n:

[J,I] = meshgrid(1:n);

A = mod(I+J-(n+3)/2,n);

B = mod(I+2*J-2,n);

M = n*A + B + 1;

9. Пример кода, рисующего график функции sinc-функции

sin

[X,Y] = meshgrid(-8:.5:8);

R = sqrt(X.^2 + Y.^2) + eps;

Z = sin(R)./R;

mesh(X,Y,Z)

10. Вычисление факториала с помощью М-функции

Глава 2. MATLAB

function y=f(n)

n=input('n=');

for i=1:n

if i<=1

y=1;

else

y=i*y;

end

Можно просто воспользоваться встроенной функцией factorial(n).

2.5. Решение уравнений и систем уравнений средствами

MATLAB

2.5.1. Вычисление нулей функции одной переменной

Ряд функций системы MATLAB предназначен для работы с

функциями. Под функциями понимаются как встроенные функции, например

sin(x) или ехр(х), так и функции пользователя, например f(x), задаваемые как

m-файлы-функции. По аналогии с дескрипторами графических объектов

могут использоваться объекты класса дескрипторов функций, задаваемых с

помощью символа @, например: fe=@exp.

Численные значения таких функций, заданных дескрипторами,

вычисляются с помощью функции feval:

fe=@exp;

feval(fe,1.0)

ans =

2.7183

Довольно часто возникает задача решения нелинейного уравнения вида

f(x) = 0

Соответствующая функция MATLAB, решающая данную задачу,

приведена ниже:

Глава 2. MATLAB

• fzero(@fun,x) — возвращает уточненное значение х, при котором

достигается нуль функции fun, представленной в символьном виде, при

начальном значении аргумента х. Возвращенное значение близко к

точке, где функция меняет знак, или равно NaN, если такая точка не

найдена;

• fzero(@fun,[xl x2]) — возвращает значение х, при котором fun(x)=0 с

заданием интервала поиска с помощью вектора x=[xl х2], такого, что

знак fun(x(D)) отличается от знака fun(x(2)). Если это не так, выдается

сообщение об ошибке. Вызов функции fzero с интервалом гарантирует,

что fzero возвратит значение, близкое к точке, где fun изменяет знак;

• fzero(@fun,x,tol) — возвращает результат с заданной погрешностью

tol;

• fzero(@fun,x,tol ,trace) — выдает на экран информацию о каждой

итерации;

• fzero(@fun,х,tol, trace,Р1,Р2,...) — предусматривает дополнительные

аргументы, передаваемые в функцию fun(x,Pl,P2,...). При задании

пустой матрицы для tol или trace используются значения по

умолчанию. Пример: fzero(fun,x,[ ],[ ],Р1).

Для функции fzero ноль рассматривается как точка, где график

функции fun пересекает ось х, а не касается ее. В зависимости от формы

задания функции fzero реализуются следующие хорошо известные

численные методы поиска нуля функции: деления отрезка пополам, секущей

и обратной квадратичной интерполяции. Приведенный ниже пример

показывает приближенное вычисление π/2 из решения уравнения cos(x)=0 с

представлением косинуса дескриптором:

x= fzero(@cos,[1 3])

x =

1.5708

Глава 2. MATLAB

В более сложных случаях рекомендуется строить график функции f(x)

для приближенного определения корней и интервалов, в пределах которых

они находятся.

Пример:

%Функция, корни которой

ищутся

function f=funl(x)

x=0:0.1:10;

f=0.25*x+sin(x)-1;

plot(x,f) grid

end

Из рисунка нетрудно заметить, что значения корней заключены в

интервалах [0.5 1], [2 3] и [5 6]. Найдем их, используя функцию fzero:

xl=fzero(@funl,[0.5 1]) xl = 0.8905

x2=fzero(@funl,[2 3]) x2 = 2.8500

x3=fzero(@funl,[5 6]) x3 = 5.8128

К сожалению, сразу найти все корни функция fzero не в состоянии.

Решим эту же систему при помощи функции fsolve из пакета Optimization

Toolbox, которая решает систему нелинейных уравнений вида f(x)=0 методом

наименьших квадратов, ищет не только точки пересечения, но и точки

касания.

fsolve(@funl,0:10 )

ans =

Columns 1 through 7

0.8905 0.8905 2.8500 2.8500 2.8500 5.8128 5.8128

Columns 8 through 11

5.8128 2.8500 2.8500 10.7429

Глава 2. MATLAB

2.5.2. Минимизация функции одной переменной

Еще одна важная задача численных методов — поиск минимума

функции f(x) в некотором интервале изменения х — от х

до х

. Если нужно

найти максимум такой функции, то достаточно поставить знак «минус»

перед функцией.

Для решения этой задачи используется функция [X,fval,exitflag,output]

= fminbnd(@fun.x1,x2,options, p1,p2,...)

• fminbnd(@fun,xl,x2) — возвращает значение х, которое является

локальным минимумом функции fun(x) на интервале xl<x<x2;

• fminbnd(@fun,xl,x2,options) — сходна с описанной выше

формой функции, но использует параметры tolX, maxfuneval, maxiter, display

из вектора options, предварительно установленные при помощи команды

optimset (смотрите описание lsqnonneg);

• fminbnd(@fun,xl,x2,options,P1,P2...) — сходна с описанной

выше, но передает в целевую функцию дополнительные аргументы: Р1, Р2.....

Если требуется использовать параметры вычислений по умолчанию, то

вместо options перед P1, Р2 необходимо ввести [ ] (пустой массив);

• [x,fval] = fminbnd(...) — дополнительно возвращает значение

целевой функции fval в точке минимума;

• [x,fval ,exitflag] = fminbndL) —дополнительно возвращает

параметр exitflag, равный 1, если функция сошлась с использованием

options.tolX, и 0, если достигнуто максимальное число итераций

options.maxiter.

Используются следующие обозначения: xl:х2 — интервал, на котором

ищется минимум функции; Р1,Р2... — дополнительные, помимо х,

передаваемые в функцию аргументы; fun — строка, содержащая название

функции, которая будет минимизирована; options — вектор параметров

вычислений. В зависимости от формы задания функции fminbnd вычисление

минимума выполняется известными методами золотого сечения или

параболической интерполяции.

Глава 2. MATLAB

Пример:

options=optimset('tolX',1.e-10);

[x]=fminbnd(@cos,3,4,options)

x = 3.1416

Значительно сложнее задача минимизации функций нескольких

переменных f(х

,...). При этом значения переменных представляются

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

Для

минимизации функций ряда переменных MATLAB обычно использует

разновидности симплекс-метода Нелдера-Мида. Мы не будем рассматривать

эту задачу.

2.5.3. Графический способ решения уравнений

1. Найти корни уравнения на отрезке [0,3].

0)sin(

=⋅ xx

Программа:

1;x=0:.01:3; f=x.*sin(x.^2); plot(x,[f;0*f]), grid

2;ginput (% в команде ginput точка снимается нажатием левой клавиши

мыши, Enter – выход из ginput).

Полученные корни и график

функции:

1.7869 -0.0439

2.4919 -0.0088

2. Построим график неявной функции f(x,y)≡x

y-2xy

+y-0.2=0,

x,y=[0,1].

Глава 2. MATLAB

1;h=.02; x=0:h:1; [X,Y]=meshgrid(x);

f=X.^3.*Y-2*X.*Y.^2+Y-.2;

2;v=[0,0];

Выясним, какой знак имеет f в

области G, для чего выполним:

mesh(x,x,f.*(f>0))

Это пример трехмерной графики. В ней цвет используется для

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

через голубой, зеленый и желтый до темно-красного.

Вычислим площадь S этой области:

4;S=h^2*sum(f(:)>=0) (S=0.7296)

5;V=h^2*sum(f(f>=0)) (V=0.1268)

2.5.4. Системы линейных алгебраических уравнений

Для реализации различных алгоритмов решения СЛАУ и связанных с

ними матричных операций применяются следующие операторы: +,-,*,/, \, *, '

Как отмечалось ранее, MATLAB имеет два различных типа арифметических

операций — поэлементные и для массивов (векторов и матриц) в целом.

Матричные арифметические операции определяются правилами линейной

алгебры.

Арифметические операции сложения и вычитания над массивами

выполняются поэлементно. Знак точки «.» отличает операции над

элементами массивов от матричных операций. Однако, поскольку операции

сложения и вычитания одинаковы для матрицы и элементов массива, знаки

«.+» и «.-» не используются. Рассмотрим другие операторы и выполняемые

ими операции:

• * — матричное умножение;

• С = А*В — линейное алгебраическое произведение матриц А и В:

() ()(

∑

⋅=

jkBkiAjiC

,,,

)

Глава 2. MATLAB

Для случая нескалярных А и В число столбцов матрицы А должно

равняться числу строк матрицы В. Скаляр может умножаться на матрицу

любого размера.

• / — правое деление. Выражение Х=В/А дает решение ряда

систем линейных уравнений АХ=В, где А — матрица размера тхп и В —

матрица размера nxk;

• \ — левое деление. Выражение Х=В\А дает решение ряда систем

линейных уравнений ХА=В, где А — матрица размера тхп и В — матрица

размера nxk. Если А — квадратная матрица, то А\В — примерно то же самое,

что и inv(A)*B, в остальных случаях возможны варианты, отмеченные ниже.

Если А — матрица размера пхп, а В — вектор-столбец с п

компонентами или матрица с несколькими подобными столбцами, тогда

Х=А\В — решение уравнения АХ=В, которое находится хорошо известным

методом исключения Гаусса.

Пример: Для невырожденной квадратной матрицы A и вектора-столбца

b вектор-стобец x=b/А есть решение системы линейных уравнений Ax=b.

Система Ax=b решается одним из трех способов.

>> A=[2 1 0 1;1 -3 2 4; -5 0 -1 -7; 1 -6 2 6];

>> b=[8 9 -5 0];

1 способ:

>> x=b/A

x = 3.0000 -4.0000 -1.0000 1.0000

2 способ:

>> x=b*A^-1

x = 3.0000 -4.0000 -1.0000 1.0000

3 способ:

>> x=b*inv(A)

x = 3.0000 -4.0000 -1.0000 1.0000

Как и следовало ожидать, результат решения одинаков.