|
Основы редактирования и отладки sci-файлов
Для набора, редактирования, отладки и запуска sci-файлов служит специальный редактор, который можно вызвать из командной строки (команда scipad();), либо командой Editor (кнопкаSciNotesв пятой версии). Для запуска файла его необходимо предварительно записать на диск, используя команду Save As в меню File редактора. После записи файла на диск его надо загрузить в среду Scilab командой редактора Execute/Load into Scilabили из главного меню Scilab вызвать команду Execи указать имя файла-сценария. С помощью редактора, в частности, можно устанавливать в тексте файла специальные метки – точки прерывания. Операторы и функции Некоторые специальные символы: a(: , j) – j-й столбец матрицы a; a(i , :) – i-я строка матрицы a; a(j : k) – элементы aj, aj+1, …, ak; a( : ) – записывает все элементы массива a в виде столбца; a(m , :) =[ ] – удаляет из матрицы строку m; a’ – транспонированная матрица а; a.’ – транспонирование массива; prod(a) – произведение элементов массива; prod(a, dim) – произведение элементов столбцов (dim=1) или строк (dim=2) ; sum(a) – сумма элементов массива; sum(a, dim) - сумма элементов столбцов (dim=1) или строк (dim=2). Матричные операции линейной алгебры · det(a)– возвращает определитель квадратной матрицы а; · rank(a)– возвращает ранг матрицы; · norm(a)– возвращает норму матрицы а; · b=orth(a)– возвращает ортонормальный базис матрицы а; · inv(a)– возвращает матрицу, обратную матрице а; · spec(a)– возвращает вектор собственных значений матрицы а. Примеры: Вычислить: определитель матрицы d; ранг матрицы c; обратную матрицу к матрице x. Решение: В командном окне в строках ввода набираем: --> d=[2 9 9 4; 2 -3 12 8; 4 8 3 -5; 1 2 6 4]; --> c=[1, -2, -3, 0; 2, 3, 8, 7; -1, 1, 1, -1]; -->x=[3 2 1; 1 0 2; 4 1 3]; --> x1=det(d), x2=rank(c), x3=inv(x) В строках вывода получаем: x1 = 147. x2 = 2. x3 = -0.4 -1. 0.8 1. 1. -1. 0.2 1. -0.4 Численные методы и обработка данных Решение систем линейных уравнений Пример: решить систему линейных уравнений Решение возможно одним из способов (s1, s2, s3 или s4 – см. приведенную ниже программу). --> a=[2, 1, 0, 1; 1, -3, 2, 4; -5, 0, -1, -7; 1, -6, 2, 6]; b=[8 9 -5 0]; --> s1=b/a', s2=a\b', s3=b*a'^(-1), s4=b*inv(a') s1 = 8.1481481 - 1.5185185 11.703704 - 6.7777778 s2 = 8.1481481 - 1.5185185 11.703704 - 6.7777778 s3 = 8.1481481 - 1.5185185 11.703704 - 6.7777778 s4 = 8.1481481 - 1.5185185 11.703704 - 6.7777778 Решить систему линейных уравнений вида можно с помощью функции linsolve : linsolve(a,b). В нашем примере: a=[2, 1, 0, 1; 1, -3, 2, 4; -5, 0, -1, -7; 1, -6, 2, 6]; b=[8; 9 ;-5 ;0]; x=linsolve(a,-b) x = 8.1481481 - 1.5185185 11.703704 - 6.7777778 Если система уравнений имеет бесчисленное множество решений, то выводится одно из них: a=[2, 1; 4, 2]; b=[-7; -14 ];x=linsolve(a,b) x = 2.8 1.4 Если система не имеет решений: a=[2, 1; 4, 2]; b=[-7; -13 ];x=linsolve(a,b) WARNING:Conflicting linear constraints! x = [] Вычисление корней полинома Функция roots(c) возвращает вектор-столбец из корней полинома с. Пример: решить уравнение --> x=[7, 0, 12, 23]; d=roots(x) d = 0.5564046 + 1.6257442i 0.5564046 - 1.6257442i - 1.1128093 Примечание: Коэффициенты полинома следует вводить в порядке убывания степеней переменной x. Если в уравнении отсутствует слагаемое, содержащее, например, x2, то в векторе коэффициентов на соответствующем месте надо ввести 0. Решение нелинейных уравнений вида f(x)=0 Уравнения бывают алгебраическими и трансцендентными. Алгебраическим называют уравнение вида . Если уравнение нельзя свести к алгебраическому заменой переменных, то его называют трансцендентным. Пример: Для решения уравнений, в том числе трансцендентных, в Scilab применяют функцию fsolve(x0,f) где x0 - начальное приближение, f - функция, описывающая левую часть уравнения f(x)=0. Пример: решить уравнение Набираем в окне редактора файл: function y=f(x) y=7*x.^3+45*x.^2+12*x+23; endfunction и сохраняем его под именем f.sci.Загружаем его в Scilab(Execute/Load into Scilab). Для нахождения отрезка [a, b], на котором отделен корень данного уравнения, построим график функции .
-->x=-8:0.1:-5; plot(x, f(x)); xgrid() Из графика видно, что корень отделен на отрезке [-6.5 , -6]. Найдем его, используя функциюfsolve: -->x0=-6.5;x1= fsolve(x0,f) Получаем: x1 = - 6.2381997 Систему нелинейных уравнений также можно решить, используя функцию fsolve.
clc function [y]=ff(x) y(1)=x(1)^2+x(2)^2-1; y(2)=x(1)^3-x(2); endfunction t=fsolve([-.5,-.5],ff) t = - 0.8260314 - 0.5636242 Поиск минимума функции y=f(x) на интервале [a, b] Функция [f,xopt]=optim(costf,x0) возвращает локальный минимум функции costf. Функция возвращает минимум функции (значение f) и точку, в которой этот минимум достигается (xopt). Главной особенностью функции optim является структура функции costf, которая должна быть такой: function [f,g,ind]=costf(x,ind) f=gg(x);//функция, минимум которой мы ищем g=numdiff(gg,x); //градиент функции f Endfunction Если возвращаемое сформированной функцией costf значение ind равно 2, 3 или 4, то функция costf обеспечивает поиск минимума, т.е. в качестве результата функции optim возвращается f и xopt. Если ind=1, то в функции optim ничего не считается, условие ind<0 означает, что минимум f(x) не может быть оценен, а ind=0 прерывает оптимизацию. Вообще говоря, значение параметра ind является внутренним параметром для связи между optim и costf, для использования optim необходимо помнить, что параметр ind должен быть определен в функции costf. Пример: найти минимум функции . Решение. Построим график функции для определения интервалов [a, b], на которых находятся экстремумы этой функции. -->x=-3:.1:3; y=x.^4-3*x.^2-5*x-4; plot(x, y); xgrid() Из графика видно, что это отрезок [1, 2]. Набираем в окне редактора и отправляем на выполнение файл function [f,g,r]=z(x,r) f=x.^4-3*x.^2-5*x-4 g=4*x.^3-6*x-5 endfunction x0=1; [fmin,xmin]=optim(z,x0) Получаем --> xmin = 1.5233402 fmin = - 13.193373 Возможен другой вариант, без ручного вычисления производной: function y=gg(x) y=x.^4-3*x.^2-5*x-4; endfunction function [f,g,r]=z(x,r) f=gg(x) g=numdiff(gg,x) endfunction x0=1; [fmin,xmin]=optim(z,x0)
xmin = 1.5233402 fmin = - 13.193373 В случае функции двух переменных: (Поиск минимума функции Розенброка ) clc x0=[-2;2]; function y=gg(x) y=100*(x(2)-x(1)^2)^2+(1-x(1))^2; endfunction function [f,g,r]=z(x,r) f=gg(x) g=numdiff(gg,x) endfunction [fmin,xmin]=optim(z,x0) xmin = 0.9999955 0.9999910 fmin = 2.010D-11
Численное интегрирование При вычислении определенных интегралов от функций, заданных в виде таблицы или в явном виде ( ), одним из численных методов используют функции intsplin, inttrap, integrate, intg. Способ 1. С помощью команды intsplin. Это интегрирование экспериментальных данных с помощью сплайн-интерполяции. Известно значение интегрируемой функции в дискретных точках (узлах). Пример: вычислить интеграл от таблично заданной функции
-->x=1:.4:5; -->y=exp((x-3).^2/8)//значения у в таблице получены табулированием этой функции -->v=intsplin(x,y) Получаем: v = 4.7799684 Или: -->x=[1 1.4 1.8 2.2 2.6 3 3.4 3.8 4.2 4.6 5]; -->y=[1.6487 1.3771 1.1972 1.0833 1.0202 1 1.0202 1.0833 1.1972 1.3771 1.6487]; -->v=intsplin(x,y) v = 4.7799328 Способ 2. С помощью команды inttrap. Интегрирование экспериментальных данных с помощью трапецеидальной интерполяции (метод трапеций). При вычислении интеграла между соседними узлами функция интерполируется линейно. Этот метод вычислений называется методом трапеций. Вычислим интеграл от той же самой функции -->x=1:.4:5; -->y=exp((x-3).^2/8) -->v=inttrap(x,y) Получаем:
v = 4.8017553 Способ 3. С помощью команды integrate. Это интегрирование по квадратуре. Может задаваться требуемая точность вычислений. Синтаксис [x]=integrate(expr,v,x0,x1 [,ea [,er]]) Параметры expr : подынтегральная функция; v : переменная интегрирования; x0, x1 : пределы интегрирования; ea, er : действительные числа. Первое из них – еа - предельная абсолютная ошибка. По умолчанию принимает значение 0. Второе - er -, - предельная относительная ошибка. По умолчанию принимает значение 1.d-8. Пример: вычислить . -->integrate('exp((x-3)^2/8)','x',1,5) ans = 4.7798306 Можно это сделать и так. Набираем и сохраняем в окне редактора под именем a.sci файл function g=a(x) g=exp((x-3).^2/8); endfunction Загружаем этот файл в среду Scilab (Load into Scilab). Далее в строке ввода набираем: -->integrate('a','x',1,5) Получаем: ans = 4.7798306 Способ 4. С помощью команды intg. Интегрируемая функция задана извне: либо в виде набора дискретных точек, либо вычисляется с помощью внешней подпрограммы. Может задаваться требуемая точность вычислений. Синтаксис [v,err]=intg(a,b,f [,ea [,er]) Параметры a, b : действительные числа; f : внешняя функция или список строк; ea, er : - параметры, описанные выше; err : оценка абсолютной ошибки результата. Вычисляется определенный интеграл от a до b от f(t)dt.
Вычислим тот же интеграл: function g=a(x) g=exp((x-3).^2/8); endfunction intg(1,5,a) Вызвав этот файл на выполнение (Execute/Load into Scilab), получим ответ: ans = 4.7798306 |
|||||||||||||||||||||||||
|