Basis of Mathematical Modeling

Basis of Mathematical Modeling LECTURE 3 Numerical Analysis with MATLAB Dr. N.K. Sakhnenko, PhD, Professor Associate

Outline • Finding roots with MATLAB • Interpolation • Numerical Integration • Optimization Problems: Minimizing Functions of One Variable; Minimizing Functions of Several Variables; Linear Programming

Root finding Computing roots of the polynomials Polynomials are represented in MATLAB by their coefficients in the descending order of powers. For instance, the cubic polynomial p(x) = 3x3 + 2x2+ 1 is represented as p = [3 2 0 1]. Its roots can be found using function roots. >> format short >> r = roots(p) r = -1.0000 0.1667 + 0.5528i 0.1667 - 0.5528i

Finding zeros using MATLAB functionfzero Let now f be an one variable function f=f(x). MATLAB function fzero computes a zero ofthe function f using user supplied initial guess. In the following example let f(x) = cos(x) – x. First we define a function y = f1(x) function y = f1(x) y = cos(x) - x; To compute its zero we use MATLAB function fzero >> format long >>r = fzero('f1', 0.5) r =0.73908513321516 Name of the function whose zero is computed is entered as a string. Second argument of functionfzerois the initial approximation of r. f=0

Functionfzero A user can enter a two-element vector that designates a starting interval. In our example we choose [ 0 1] as a starting interval to obtain r = fzero('f1', [0 1]) r = 0.73908513321516 One can check last result using function feval err = feval('f1', r) err =0 By adding the third input parameter tol you can force MATLAB to compute the zero of afunction with the relative error tolerance tol. In our example we let tol = 10-3 to obtain rt = fzero('f1', .5, 1e-3) rt = 0.73886572291538

Functionfzero functionf = myf(x) f = sin(x) - х.^2.*cos(x); >> x1=fzero ('myf', -5) x1 = -4.7566 >> x2=fzero ('myf', 0.1) x2 = 7.9161e-019 >> x3=fzero ('myf', 5) x3 = 4.6665 f=0 f=0 f=0

Systems of nonlinear equations MATLAB function fsolvecomputes solution of this system as well. For more information help fsolve

Systems of nonlinear equations function F=mysys(x) F(1)=x(1)*(2-x(2))-cos(x(1))*exp(x(2)); F(2)=2+x(1)-x(2)-cos(x(1))-exp(x(2)); >> [x,f]=fsolve('mysys', [0 0]) x = 0.7391 0.4429 f = 1.0e-011 * -0.4702 -0.6404

Interpolation • Given set of n+1 points xk , yk, k=(0…n), with x0 < x1 < … < xn, find a function f(x) whosegraph interpolates the data points, i.e., f(xk) = yk, for k = 0, 1, …, n. • The general form of the function interp1 is yi = interp1(x, y, xi, method), where the vectors xand y are the vectors holding the x- and the y- coordinates of points to be interpolated,respectively, xi is a vector holding points of evaluation, i.e., yi = f(xi) and method is an optionalstring specifying an interpolation method. The following methods work with the function interp1: • Nearest neighbor interpolation, method = 'nearest'. Produces a locally piecewise constantinterpolant. • Linear interpolation method = 'linear'. Produces a piecewise linear interpolant. • Cubic spline interpolation, method = 'spline'. Produces a cubic spline interpolant. • Cubic interpolation, method = 'cubic'. Produces a piecewise cubic polynomial.

Interpolation In this example, the following points (xk, yk) = (k/5, sin(2k)), k = 0, 1, … , 5, x = 0:pi/5:pi; y = sin(2.*x); are interpolated using the 'cubic' method of interpolation. The interpolant is evaluated at the following points xi = 0:pi/100:pi; yi = interp1(x, y, xi, 'cubic'); Points of interpolation together with the resulting interpolant are displayed below plot(x, y, 'o', xi, yi), title('Cubic interpolant of y = sin(2x)')

Numeric Integration Two MATLAB functions quad('f ', a, b, tol) andquad8('f ', a, b, tol) are designed for numerical integration of the univariatefunctions. The input parameter 'f' is a string containing the name of the function to be integratedfrom a to b. The fourth input parameter tol is optional and specifies user's chosen relative error inthe computed integral. >> quad('sin',0,pi/2) ans = 0.99999999787312 >> quad('sin',0,pi/2,0.0000001) ans = 0.99999999946896 The default accuracy is 0.001 function f=fint(x) f=exp(-x.^2); >> quad('fint',0,1) ans = 0.74682418072642

Numeric Integration of the bivariate functions Functiondblquadcomputes a numerical approximation of the double integral where D = {(x, y): a < x< b, c< y <d} is the domain of integration. Syntax of the functionis dblquad (fun, a, b, c, d, tol), where the parameter tol has the same meaning as in thefunction quad.

Numeric Integration of the bivariate functions The M-file esin is used to evaluate function f: function z = esin(x,y); z = exp(-x.*y).*sin(x.*y); >>result = dblquad('esin', -1, 1, 0, 1) result = -0.22176646183245

Optimization The problem of optimality criterion finding is in a base of most managing protocols in modern telecom technologies. The optimization problem is the problem of finding the best solution from all feasible solutions (from Wikipedia, the free encyclopedia) This may correspond to maximizing profit in a company, or minimizing loss in a conflict. For Telecom, for example, it is necessary to minimize middle delivery of packages in a network; to maximize the productivity of telecom system; to minimize the cost of the using of resources of networks; to maximize a profit of the sale of telecommunications services.

Minimizing Functions of One Variable For a given mathematical function of a single variable coded in an M-file, you can use the fminbnd function to find a local minimum of the function in a given interval. If one needs a maximum he should to put minus before the function. [x,fval] = fminbnd(@myfun,x1,x2) returns value of x in which the function myfun(x) has a minimum on the interval x1<x<x2 and returns fval, that is the value of the function in the minimum. function f = myfun(x) f = ... % the function under consideration For minimizing one can use the following functions [x,fval] = fminsearch(@myfun,x0) [x,fval] = fminunc(@myfun, x0) where x0 is initial guess.

Minimizing Functions of One Variable Find the minimum of the functionf(x)=x3-2x-5on the interval from 0 to 2 clear all; fplot('fun1',[0,2]); grid on; title('fun1'); xlabel('x'); ylabel('f(x)'); [x,fval]=fminbnd(@fun1,0,2) % fun1.m function f = fun1(x) f=x^3-2*x-5 Solution: x=0.8165; fval=-6.0887

Minimizing Functions of One Variable For minimizing problem one can use the following functions [x,fval] = fminsearch(@myfun,x0) [x,fval] = fminunc(@myfun, x0) where x0 is initial guess. E.g., Find the minimum of the function f(x)=arctg(x3-2x-5)on the interval from -2to 2 clear all; fplot('fun2',[-2,2]); grid on; title('fun2'); xlabel('x'); ylabel('f(x)'); [x,fval]=fminunc(@fun2,0) Solution: x=0.8165; fval=-1.4080 [x,fval]=fminunc(@fun2,-1) Solution: x =-2; fval=-1.5708 % fun2.m function f=fun2(x) f=atan(x^3-2*x-5);

Minimizing Functions of Several Variables The fminsearch function minimizes a function of several variables. Syntax x = fminsearch(fun,x0) x = fminsearch(fun,x0,options) [x,fval] = fminsearch(...) Description fminsearch finds the minimum of a scalar function of several variables, starting at an initial estimate. This is generally referred to as unconstrained nonlinear optimization. x = fminsearch(fun,x0) starts at the point x0 and finds a local minimum x of the function described in fun. x0 can be a scalar, vector, or matrix.

Minimizing Functions of Several Variables Find the minimum of the function f(x)= 3x2 + 2xy + y2+5 near the point [20, -20] clear all; x0 = [20,-20]; [x,fval] = fminunc(@fun3,x0) % fun3.m function f = fun3(x) f = 3*x(1)^2 + 2*x(1)*x(2) + x(2)^2+5;

Linear Programming A linear programming problem (or LP in brief) is any decision problem of the form

Linear Programming First, enter the coefficients f = [-5; -4; -6] A = [1 -1 1 3 2 4 3 2 0]; b = [20; 42; 30]; lb = zeros(3,1); Next, call a linear programming routine: [x,fval] = linprog(f,A,b,[],[],lb) Optimization terminated successfully. x = 0.0000 15.0000 3.0000 fval = -78.0000

Basis of Mathematical Modeling