NUMERICS USER GUIDE Contents 1. Introduction 1 2. Installing
Contents
1. Numerics provides the functions explicit_euler_fs implicit_euler_fs crank_nicolson_fs rk_fs rkf_fs sdirk_fs nystroem_fs adams_bashforth_fs adams_moulton_fs bdf_fs ode fs for the solution of initial value problems for ordinary differential equa tions in R o t fGe0 to lt t lt T z to Lo with fixed step size methods All functions plot the computed solution according to their argument pc The function explicit_euler_fs implements the explicit Euler scheme see 6 Algorithmus 1 2 1 8 IV 1 10 1 2 2 It has the following arguments the first four being mandatory f force function f x0 initial value xo T final time T nt number of time steps the step size is h tO initial time to default is t 0 pc integer vector of length 2 fixing the components of x that will be plotted default settings are 0 1 for scalar odes i e plot x t versus t 1 2 for systems of odes i e plot the curve x t xa t T to nt 7 Note that the user defined function f must have the two arguments t and x in this order The function implicit_euler_fs implements the implicit Euler scheme NUMERICS USER GUIDE 13 see 6 Algorithmus 1 2 2 8 IV 1 10 1 2 2 It has the following arguments the first five being mandatory f force function f df derivative D f w r t x of the force function f x0 initial value xo T final time T nt number of time steps the step size is h2. best value of f at every iteration 3 print current best value of f and corresponding point x at every iteration itm maximal number N of iterations default is N 100 tol tolerance default is e 0 001 iteration terminates if the variation of the function values is less than e or if the number of iterations exceeds N REFERENCES 1 Micha l Baudin Introduction to Scilab http wiki scilab org Tutorials Programming in Scilab http wiki scilab org Tutorials 3 Arieh Iserles A First Course in the Numerical Analysis of Differential Equa tions Cambridge Texts in Applied Mathematics Cambridge University Press 1997 4 Eike Rietsch An Introduction to Scilab from a Matlab User s Point of View http wiki scilab org Tutorials 26 R VERFURTH R diger Verf rth Einf hrung in die Numerische Mathematik http www rub de num1 files lectures EinfNumerik pdf Numerik I Gew hnliche Differentialgleichungen und Differenzenver fahren f r partielle Differentialgleichungen http www rub de num1 files lectures NumDgl1l pdf Optimierung http www rub de num1 files lectures Optimierung pdf Numerische Mathemattk fiir Maschinenbauer Bauingenieure und Um welttechniker http www rub de num1 files lectures NumIng pdf Vertiefung Numerische Mathematik f r den Masterstudiengang UTRM http www rub de num1 files lectures VertNum pdf Numerical Methods and Stochastics Part I Numerical Met
3. equations The function poisson_FD implements the standard 5 point difference discretization of the two dimensional Poisson equation with homoge neous Dirichlet boundary conditions Au f inQ u 0 onT see 6 111 3 9 SIIL 1 The domain Q is the intersection of the rectangle x r1 X 12 amp r2 with bottom left corner x and top right corner xz with the set x R w x lt 0 the square 1 1 being the default setting The function poisson_FD plots the level lines of the discrete solution and its maximal error if the exact solution is known It has the following arguments the first two being mandatory f force function f 18 R VERFURTH np number of interior grid points in each coordinate direction y di a Trinti Tr2 T12 the mesh size in z and y direction is rem and ap TESp uex exact solution if known default setting is signifying an unknown exact solution xbl bottom left corner x default is x 1 1 xtr top right corner z default is x 1 1 fom function w describing the domain Q default is a les 1 nll number of level lines for plotting the solution default is 10 The function heat_FD_1D solves the one dimensional heat equation u Ou a ge ee 1 u t 0 0 n 0 7 u t 1 0 n 0 T u 0 uo in 0 1 using a second order symmetric difference discretization for the spatial derivative and the d scheme for the time discretizatio
4. explicit Euler scheme The function nystroem_fs has the same arguments as explicit_euler_fs The function adams_bashforth fs implements the Adams Bashforth two and three step schemes see 6 Beispiel 1 5 1 The first two or three resp steps are computed with the explicit Euler scheme The function adams_bashforth_fs has the same arguments as the function explicit_euler_fs except that pc now is the seventh argument and that a new sixth argument steps determines the method This argu ment may take the values 2 or 3 with 2 being its default value The function adams_moulton_fs implements the Adams Moulton two and three step schemes see 6 Beispiel 1 5 1 The first two or three resp steps are computed with the implicit Euler scheme The function 14 R VERFURTH adams_moulton fs has the same arguments as implicit_euler fs ex cept that pc now is the eighth argument and that a new seventh ar gument steps determines the method This argument may take the values 2 or 3 with 2 being its default value The function bdf_fs implements the backward difference two and three step schemes see 6 Beispiel 1 5 2 8 VI 8 The first two or three resp steps are computed with the implicit Euler scheme The function bdf_fs has the same arguments as implicit_euler fs except that pc now is the eighth argument and that a new seventh argument steps determines the method This argument may take the values 2 or 3 with 2 being its default value T
5. gradient of fctc gradients of the component functions of fctc should be column vectors NUMERICS USER GUIDE 25 nic number p of inequality constraints the values p 0 no inequality constraints and p m no equality constraints are allowed the first nic components of the vector yp must be positive pri print level default is 2 0 print optimal value of f 1 print optimal value of f and corresponding point z 2 print current value of f at every iteration 3 print current value of f and x at every iteration tol tolerance default is 1 0D 4 itmax maximal number N of iterations default is N 100 iteration terminates if the KKT conditions are satisfied up to tolerance e or if the number of iterations exceeds N The function nelder mead implements the simplex method of Nelder and Mead for unconstrained nonlinear optimization problems minimize f x in R see 7 Algorithmus III 8 1 It prints results according to its argument prl and has the following arguments the first two being mandatory fct objective function f inx either the dimension n of the problem or an n x n 1 matrix with the initial guesses xo 2 as column vectors in the first case the initial guesses are automatically set to the origin zo 0 and to the unit vectors 7 e 1 lt i lt n prl print level default is 2 0 print optimal value of f 1 print optimal value of f and corresponding point z 2 print current
6. the maximal interpolation error and plots the inter polation polynomial the original function and the error by drawing a piecewise linear interpolation of these functions at k equidistant points on the interval a 3 The function hermite performs Hermite interpolation using divided differences and Newton s scheme see 8 8IIL 4 When interpolating a given function the user must provide the following arguments the first three being mandatory top p n Isis NUMERICS USER GUIDE 7 f function that should be interpolated df derivative of the function that should be interpolated n number n of interpolation points at least 2 a first interpolation point a default is 1 b last interpolation point b default is 1 t spacing of interpolation points e equidistant points a b a 1 lt i lt c Cebysev points a b a 1 cos n i default is e ap left endpoint a of plot interval default is a bp right endpoint p of plot interval default is b np number k of plot points default is 101 n 1 lt i lt n n 1 When interpolating discrete data the role of the first three arguments changes f and df now are the vectors of the function values and corre sponding derivatives resp that should be interpolated and n is replaced by the vector x of the corresponding abscissa These three vectors must have the same dimensions The function hermite prints the maximal interpolation error and
7. 6 Algorithmus 1 2 3 8 IV 1 10 81 2 2 It has the same arguments as implicit_euler_vs The function rk_vs implements the classical Runge Kutta scheme see 6 81 2 8 1V 6 4 10 81 2 3 It has the same arguments as the function explicit_euler_vs The function rkf_vs implements Runge Kutta Fehlberg schemes of or der 2 with 3 stages of order 3 with 4 stages and of order 4 with 6 stages see 6 Beispiel 1 4 7 8 Beispiel VI 7 It has the same arguments as the function explicit_euler_vs except that pc now is the eleventh ar gument and that a new tenth argument order determines the method This argument may take the values 2 3 and 4 with 3 being its default value The function sdirk vs implements a strongly diagonally implicit Runge Kutta scheme of order 4 see 8 VI 4 10 1 2 4 It has the same arguments as implicit_euler_vs The function ode_vs allows a comparison of the above methods Its first argument is an array sm of strings which identify the methods that should be compared and which may take the following values ee for the explicit Euler scheme ie for the implicit Euler scheme cn for the Crank Nicolson scheme rk for the classical Runge Kutta scheme rkf2 for the Runge Kutta Fehlberg method of order 2 rkf3 for the Runge Kutta Fehlberg method of order 3 rkf4 for the Runge Kutta Fehlberg method of order 4 sdirk for the strongly diagonal implicit Runge Kutta method o
8. NUMERICS USER GUIDE R VERFURTH CONTENTS 1 Introduction 2 Installing Scilab 3 Installing Numerics 4 Using Numerics 4 1 Data structures 4 2 Structure of Numerics functions 4 3 Start 4 4 Help 4 5 Interpolation 4 6 Numerical integration 4 7 Nonlinear equations 4 8 Iterative solvers for linear systems of equations 4 9 Eigenvalue problems 4 10 Initial value problems for odes fixed step size methods 4 11 Initial value problems for odes variable step size methods 4 12 Boundary value problems for odes 4 13 Finite difference methods for pdes 4 14 Linear optimization 4 15 Discrete optimization 4 16 Nonlinear optimization References Index 1 INTRODUCTION Or omooww rnwnnrnndsksr mm 17 17 19 20 23 25 27 Numerics is a library of Scilab functions implementing many of the algorithms presented in our courses on numerical analysis for math ematicians and engineers 5 6 7 8 9 10 It covers the following areas e interpolation e numerical integration e nonlinear equations Date January 5 2012 2 R VERFURTH iterative solvers for linear systems of equations eigenvalue problems initial value problems for odes boundary value problems for odes finite difference methods for pdes linear optimization integer optimization e nonlinear optimization A detailed description of Numerics is given in section 4 below In sections 2 and 3 we first explain how to in
9. RTH xp initial guess xo fct objective function f gfct gradient Df of f gfct should be a column vector Hfct Hessian D f of f if Hfct is provided the line search is an exact one otherwise it is an Armijo one default is Hfct i e Armijo line search prl print level default is 2 0 print optimal value of f 1 print optimal value of f and corresponding point z 2 print current value of f at every iteration 3 print current value of f and x at every iteration tol tolerance default is 1 0D 4 itmax maximal number N of iterations default is N 100 iteration terminates if Df x lt or if the number of itera tions exceeds N The function augmented_lagrangian implements an augmented La grangian algorithm for constrained nonlinear optimization problems minimize f x subject to fiz lt 0 1 lt i lt p fi z 0 p 1 lt j lt m with augmented Lagrange function NOU eric yon ro a i 1 m 1 y 2 gt an Tj j p 1 Sui ve 2 Tk see 7 Algorithmus II 6 9 It prints results according to its argument prl and has the following arguments the first eight being mandatory xp initial guess for the vector x yp initial guess for the vector y rp augmentation factor the vector r is given by r rp l 1 fct objective function f gfct gradient Df of f gfct should be a column vector fctc vector valued constraint function f1 fm fctc should be a column vector gfctc
10. a and b first number of starting node s default is 1 last number of terminal node t default is size Adjc c pri print level default is 2 0 no print output 1 print the value of the maximal flow z 2 print the value of the maximal flow x and the corresponding matrix where the element i j gives the flow through the arc connecting nodes i and j 3 print in every iteration the value of the current flow z 4 print in every iteration the value of the current flow x and the corresponding matrix where the element i j gives the flow through the arc connecting nodes i and j The function minimal_cost_flow implements algorithm 11 4 16 of 7 for finding an s t flow x with prescribed value and minimal cost in NUMERICS USER GUIDE 23 a network It prints results according to its argument prl and has the following arguments the first three being mandatory Adjc one of the following two objects e a square matrix such that Adjc i j gives the capacity of the arc connecting nodes i and j of the graph and equals inf if nodes i and j are not connected by an arc or if i j e a list of vectors a b c describing the properties of the arcs such that c is the capacity of the arc connecting nodes a and b Adjp one of the following two objects e a square matrix such that Adjpli j gives the cost of the arc connecting nodes i and j of the graph and equals inf if nodes i and j are not connected by an arc or if i j e a list o
11. ctors D V and M see 7 Algorithmus 11 3 3 22 R VERFURTH The function floyd_warshall implements the Floyd Warshall algo rithm for finding the shortest paths between all pairs of nodes in a graph see 7 Algorithmus II 3 11 It prints results according to its argument prl and has the following arguments the first one be ing mandatory Adjl one of the following two objects e a square matrix such that Adjl i j gives the length of the arc connecting nodes i and j of the graph and equals hinf if nodes i and j are not connected by an arc or if i j e a list of vectors a b c describing the properties of the arcs such that c is the length of the arc connecting nodes a and b prl print level default is 2 0 no print output 1 print matrix of path lengths 2 print matrices of path lengths and predecessors on shortest paths The function ford_fulkerson implements the Ford Fulkerson algo rithm for finding a maximal s t flow x in a network see 7 Algo rithmus I1 4 6 It prints results according to its argument prl and has the following arguments the first one being mandatory Adjc one of the following two objects e a square matrix such that Adjc i j gives the capacity of the arc connecting nodes i and j of the graph and equals hinf if nodes i and j are not connected by an arc or if i J e a list of vectors a b c describing the properties of the arcs such that c is the capacity of the arc connecting nodes
12. f order 4 The remaining arguments of ode_vs are the same as for the function implicit_euler_vs Note that the argument df is also present for NUMERICS USER GUIDE 17 explicit methods as the explicit Euler scheme When comparing explicit schemes exclusively you may enter the void argument for df 4 12 Boundary value problems for odes Numerics provides the function sturm_FD for solving the Sturm Liouville problem pu qu f in a b u a a u b 6 with a finite difference discretization on a uniform mesh with mesh size h see 6 I1 4 9 81 6 10 1 5 It has the following arguments the first two being mandatory nx number of interior points the mesh size is h f force function f p diffusion function p default sp 1 q reaction function q default is q 0 uex exact solution u if known default is signifying an unknown exact solution a left end point a default is a 0 b right end point b default is b 1 al left boundary value a default is a 0 be right boundary value default is 6 0 eit nx 1 The function sturm_FD plots the discrete solution and additionally the exact solution and the error if uex is provided In this case sturm_FD also prints the L and H norms of the error 4 13 Finite difference methods for pdes Numerics provides the functions poisson FD heat_FD_1D lax_friedrichs lax_wendroff for finite difference discretizations of partial differential
13. f vectors a b c describing the properties of the arcs such that c is the cost of the arc connecting nodes a and b val requested value of the flow first number of starting node s default is 1 last number of terminal node t default is size Adjc c prl print level default is 2 0 no print output 1 print the value of the final flow x 2 print the value of the final flow x and the corresponding matrix where the element i j gives the flow through the arc connecting nodes i and j 3 print in every iteration the value of the current flow z 4 print in every iteration the value of the current flow x and the corresponding matrix where the element i j gives the flow through the arc connecting nodes i and j 4 16 Nonlinear optimization Numerics provides the functions descent augmented_lagrangian nelder_mead for nonlinear optimization The function descent realizes a descent algorithm for unconstrained nonlinear optimization problems minimize f x in R see 7 Algorithmus II 1 5 The search direction is the negative gradient of the objective function f The line search is either an exact one or an Armijo one with parameters c ca 0 5 see 7 Korollar I11 1 10 In the first case the user must provide the Hessian matrix D f of the objective function The function descent prints results according to its argument prl and has the following arguments the first three being mandatory 24 R VERFU
14. he function ode_fs allows a comparison of the above methods Its first argument is an array sm of strings which identify the methods that should be compared and which may take the following values ee for the explicit Euler scheme ie for the implicit Euler scheme cn for the Crank Nicolson scheme rk for the classical Runge Kutta scheme rkf2 for the Runge Kutta Fehlberg method of order 2 rkf3 for the Runge Kutta Fehlberg method of order 3 rkf4 for the Runge Kutta Fehlberg method of order 4 sdirk for the strongly diagonal implicit Runge Kutta method of order 4 ny for the Nystroem method ab2 for the Adams Bashforth two step method ab3 for the Adams Bashforth three step method am2 for the Adams Moulton two step method am3 for the Adams Moulton three step method bdf2 for the backward difference two step method bdf3 for the backward difference three step method The remaining arguments of ode_fs are the same as for the function implicit_euler_fs Note that the argument df is also present for explicit methods as the explicit Euler or Adams Bashforth schemes When comparing explicit schemes exclusively you may enter the void argument for df 4 11 Initial value problems for odes variable step size meth ods Numerics provides the functions explicit_euler_vs implicit_euler_vs crank_nicolson_vs rk_vs rkf_vs sdirk_vs ode_vs for the solution of initial value problems for o
15. he secant rule for scalar nonlinear equations f x 0 in R see 5 Algorithmus III 3 1 8 Algorithmus 11 3 It has the following arguments the first three being mandatory x1 first initial value xo x2 second initial value 7x1 f function f pr if t results are printed in each iteration default is tol error tolerance e iteration is stopped once f x lt default is 1 0D 8 nmax maximal number of iterations default is 10 The function secant_rule prints the required number of iterations the final residual f x and the approximate solution x The function regula_falsi implements the regula falsi for scalar non linear equations f x 0 in R see 5 Algorithmus III 3 2 It has the same arguments as secant_rule and produces the same output but the function values f xo and f x of the two initial values must now be of different sign 4 8 Iterative solvers for linear systems of equations Numerics provides the functions richardson jacobi gauss_seidel ssor cg pcg _ssor bicg_stab pbicg_stab_ssor compare_solvers for the iterative solution of linear systems of equations Ar b With the exception of compare_solvers all functions print intermediate and final results depending on their argument prl and plot the relative Aan b Azn b residuals on I convergence rates Aen Ol and mean convergence Azo 6 Atn 1 b Axn b 5 rates over the number n of itera
16. hods http www rub de num1 files lectures NMS pdf adams_bashforth_fs 13 adams_moulton_fs 13 augmented_lagrangian 24 autostart_simplex 19 bdf_fs 14 bicg_stab 10 boolean 2 branch_and_bound 21 cg 10 compare_solvers 10 crank nicolson_fs 13 crank _nicolson_vs 16 cubic_spline 7 cutting_planes 21 descent 23 dijkstra 21 dual_simplex 20 explicit_euler fs 12 explicit_euler_vs 15 floyd_warshall 22 ford_fulkerson 22 function 4 gauss_rule 8 gauss_seidel 9 goertzel 7 heat_FD_1D 18 hermite 6 implicit_euler_fs 12 implicit_euler_vs 15 inner_point 20 inverse _power_iteration 11 inverse_rayleigh 12 jacobi 9 lagrange 6 lax_friedrichs 18 lax_wendroff 18 matrix 3 midpoint_rule 8 minimal_cost_flow 22 nelder_mead 25 neville 5 newton 8 INDEX number 3 numerics_help 5 nystroem_fs 13 ode fs 14 ode_vs 16 pbicg_stab_ssor 10 pcg_ssor 10 poisson_FD 17 power_iteration 11 rayleigh 12 regula_falsi 9 richardson 9 rkf_fs 13 rk fs 13 rkf_vs 16 rk_vs 16 romberg 8 sdirk fs 13 sdirk_vs 16 secant_rule 9 simplex 19 simpson_rule 8 sparse matrix 3 ssor 9 string 3 trapezoidal_rule 8 vector 3 27
17. idistant points on the interval 1 1 Then you simply enter lagrange exp 10 in the Scilab console This invokes the Numerics function lagrange with the built in exponential function exp and the number 10 as argu ments Suppose next that you want to interpolate the function gt in 10 equidistant points on the interval 1 1 Then you first have to create the corresponding function by entering function y myfunction x y x 2 x x endfunction in the Scilab console Then you may start the interpolation by typing without the full stop at the end lagrange myfunction 10 4 2 Structure of Numerics functions Numerics functions never have a return value they only print the required time for the calcula tion In addition depending on their arguments they may print and eventually plot results All functions have a certain number of manda tory and optional arguments Mandatory arguments must be specified by the user Optional arguments may be omitted corresponding de fault values will be used then Note that arguments are identified by order Consequently you must enter a void argument for all those op tional arguments that you don t want to change and that come prior to an optional argument which you want to change As an example con sider the function lagrange It requires a function as mandatory first argument the number of interpolation points as mandatory second ar gument the left endpoint of the interpolati
18. lean variable myboolean by typing myboolean t or myboolean f in the Scilab console and then NUMERICS USER GUIDE 3 calling the Numerics function with myboolean as argument String Strings are enclosed by quotation marks They are transferred to a Numerics function by either directly entering name for the string value name or by creating an own string variable mystring by typ ing mystring name in the Scilab console and then calling the Numerics function with mystring as argument Number Numbers may be integers like 2 and 3 or fixed point reals like 3 14 and 2 87 or floating point reals like 0 578D4 0 689D5 0 203D 7 and 0 1763D 5 They are transferred to a Numerics function by ei ther directly entering the numerical value or by creating an own variable mynumber by typing mynumber 2 or similar in the Scilab console and then calling the Numerics function with mynumber as argument Vector Matrix A vector is either an mx 1 matrix for a column vector or a 1 x n matrix for a row vector Matrices are created by typing the appropriate sequence of numbers in lexicographic order from top left to bottom right in brackets with columns separated by blanks or commas and rows separated by semicolons Thus the matrices 1 2 6 _ and 3 4 5 6 are represented by 1 23 456 and 1 2 3 4 5 6 A is the transpose of the matrix A Ax B is the product of the matrices A and B ones m n and zeros m n create m x n matrices
19. lt n y scalar vector or matrix of points in which the interpolation polynomial should be evaluated f vector of function values at the interpolation points at least 2 x vector of interpolation points at least 2 Note that both vectors f and x must have the same dimensions In both cases neville prints the vectors of the evaluation points and of the corresponding values of the interpolation polynomial The function lagrange performs Lagrange interpolation using divided differences and Newton s scheme see 5 Algorithmus 1 2 11 1 2 13 8 Algorithmus II 2 IIL 3 When interpolating a given function the user must provide the following arguments the first two being mandatory f function that should be interpolated number n of interpolation points at least 2 first interpolation point a default is 1 last interpolation point b default is 1 spacing of interpolation points e equidistant points a b a 1 lt i lt c Cebysev points a b a 1 cos n default is e ap left endpoint a of plot interval default is a bp right endpoint 6 of plot interval default is b np number k of plot points default is 101 When interpolating discrete data the role of the first two arguments changes f now is the vector of the function values that should be interpolated and n is replaced by the vector x of the corresponding abscissa Both vectors must have the same dimensions The function lagrange prints
20. m mesh refinements each refinement step doubling the number of grid points in each direction de fault is 5 rtol relative error tolerance for iterative solvers iteration stops omer lt e default is 0 01 om relaxation parameter for ssor iteration default is 1 5 Note that compare_solvers may be called with an empty argument list and that the strings identifying methods are case sensitive The function compare_solvers prints and plots the computing times of the chosen methods and their complexity defined as the ratio of the logarithm of the computing time over the logarithm of the number of unknowns 4 9 Eigenvalue problems Numerics provides the functions power_iteration rayleigh inverse_power_iteration inverse_rayleigh for computing the maximal and minimal eigenvalues of a general qua dratic matrix A or of a symmetric positive definite matrix The functions power_iteration and inverse_power_iteration im plement the power iteration see 5 Algorithmus V 3 1 8 V 2 and the inverse power iteration see 5 Algorithmus V 3 9 8 V 4 for computing the in absolute value largest and smallest resp eigenvalue of a general quadratic matrix A and a corresponding eigenvector They both have the following arguments the first one being mandatory A matrix A x initial guess for an eigenvector default is a random vector prl print level default is 0 0 print number of required iterations and final approximati
21. me arguments as richardson with the default value 1 5 for the relaxation parameter om Notice that for both algorithms the matrix A must by symmetric positive definite The functions bicg_stab and pbicg_stab_ssor implement the stabi lized bi conjugate gradient algorithm see 9 Algorithmus VI 7 1 10 Algorithm IJI 5 1 and the preconditioned stabilized bi conjugate gra dient algorithm with SSOR preconditioning They have the same ar guments as their symmetric counterparts cg and pcg_ssor up to an additional optional final argument nmaxrestart which limits the max imal number of restarts with the default value being 10 The function compare_solvers compares various direct and iterative solvers for the solution of the linear system of equations arising from the five point difference discretization of the Poisson equation with ho mogeneous Dirichlet boundary conditions on the unit square see 6 SIII 3 9 III 1 It has the following arguments none being manda tory methods a vector of strings specifying the chosen methods available methods are LU sparse LU factorization CG the conjugate gradient algorithm PCG the preconditioned conjugate gradient algorithm with ssor preconditioning GS gauss seidel iteration NUMERICS USER GUIDE 11 SSOR ssor iteration default setting is to use all methods nfirst number of interior grid points in each direction on the coarsest grid default is 3 maxref number of unifor
22. n see 6 SIIL 4 9 8IIL 2 It plots the discrete solution at each intermediate time step and has the following arguments the first three being mandatory T final time T nt number of time steps the fixed time step size is T Z nx number of interior grid points in x direction the fixed spatial mesh size is h a theta parameter 0 for the temporal discretization 0 0 is the explicit Eule scheme 0 0 5 is the Crank Nicolson scheme 0 1 is the implicit Euler scheme default is 0 0 5 u0 initial function ug default is ug 0 f force function f default is f 0 The functions lax_friedrichs and lax_wendroff implement the Lax Friedrichs 3 Exercise 14 4 and Lax Wendroff 3 Equation 14 29 finite difference schemes for the numerical solution of first order partial differential equations of the form Ou Ou ra 0 in 0 7 x 0 1 42 0 u t 1 in 0 7 u 0 uo in 0 1 They both plot the numerical solution at each time step and have the following arguments the first three being mandatory T final time T NUMERICS USER GUIDE 19 nx number of interior grid points in x direction the fixed spatial mesh size is h 23 the number of time steps is automatically set to nt T gt ax 1 which result in a constant time step T Sh u0 initial function uo a advection function a default is a 1 rp clear plot window every rp time steps default is 20 4 14 Linear optimization Numerics provide
23. of arguments including an empty ar gument list Typing numerics_help yields a list of all user relevant functions for which help is available When specifying arguments these must be strings with the names of Numerics functions Thus typing numerics help newton secant_rule regula falsi provides information on the functions newton secant_rule and regula_falsi 4 5 Interpolation Numerics provides the functions neville lagrange hermite cubic_spline goertzel for interpolation These functions can be used in two different ways which are determined by the type of their arguments e interpolation of a function e interpolation of discrete data The function neville realizes interpolation by Neville s scheme see 5 Algorithmus 1 2 7 When interpolating a given function the user must provide the following arguments the first three being mandatory y scalar vector or matrix of points in which the interpolation polynomial should be evaluated f function that should be interpolated 6 R VERFURTH number n of interpolation points at least 2 first interpolation point a default is 1 last interpolation point b default is 1 spacing of interpolation points e equidistant points a b a 1 1 lt i lt 1 top p n 1 c Cebysev points a 3 b a 1 cos r default is e When interpolating discrete data the user must provide the following three mandatory arguments i n 1 n Isi
24. on for the largest resp smallest eigenvalue 1 print number of actual iteration and actual approximation for the largest resp smallest eigenvalue 2 same as 0 additionally print the final approximation of the eigenvector 3 same as 1 additionally print the actual approximation of the eigenvector tol error tolerance e default is 1D 8 nmax maximal number N of iterations default is 100 iteration is stopped if two consecutive eigenvalue approxima tions differ less than and after N iterations latest 12 R VERFURTH The functions rayleigh and inverse_rayleigh realize the Rayleigh quotient iteration see 5 Algorithmus V 3 5 8 V 3 and the inverse Rayleigh quotient iteration see 5 Algorithmus V 3 11 8 V 5 for computing the largest and smallest resp eigenvalue of a symmetric positive definite matrix A They both have the following arguments the first one being mandatory A matrix A x initial guess for an eigenvector default is a random vector prl if t print iteration number and actual approximation for the largest resp smallest eigenvalue otherwise print this infor mation only at the end of the algorithm default is tol error tolerance default is 1D 8 nmax maximal number N of iterations default is 100 iteration is stopped if two consecutive eigenvalue approxima tions differ less than e and after N iterations latest 4 10 Initial value problems for odes fixed step size methods
25. on interval as optional third argument with default value 1 the right endpoint of the interpolation interval as optional fourth argument with default value 1 and further optional arguments which we do not need for this example Then the inputs NUMERICS USER GUIDE 5 lagrange exp 10 lagrange exp 10 2 lagrange exp 10 2 lagrange exp 10 4 3 yield the Lagrange interpolation of the built in function exp in 10 equidistant points on the intervals 1 1 2 1 1 2 and 4 3 respectively 4 3 Start Before using any Numerics function you must make it known to Scilab by first entering mylib lib path to the console Here mylib may be any other name different from the name of any Numerics function any built in function and any of your user defined functions path is a string giving the full path to the Numerics library Suppose for example that you are using Mac 0sX that Gargantua is your user name and that you have put the Numerics folder into the subfolder Pantagruel of your documents folder Then the above command must read without the full stop at the end mylib lib users Gargantua documents Pantagruel Numerics If everything is correctly installed Scilab returns the names of all functions in the Numerics library This includes auxiliary functions which are not user relevant and which will not be described below 4 4 Help Numerics comes with a help function numerics_help which may be used with any number
26. ot x t versus t 1 2 for systems of odes i e plot the curve 21 t z2 t Note that the user defined function f must have the two arguments t and x in this order The function implicit_euler_vs implements the implicit Euler scheme see 6 Algorithmus 1 2 2 8 IV 1 10 81 2 2 It has the following arguments the first five being mandatory f force function f df derivative D f w r t x of the force function f x0 initial value xp T final time T dt initial step size tO initial time to default is t 0 tol required tolerance e the function strives at obtaining an approximation n t for the solution x t such that n t x t lt e x t holds for all t default value is 0 001 16 R VERFURTH ntmax maximal number of time steps default is 10000 sfdt safety factor for step size variation reduce the step size at most by the factor at most by the factor sfdt default value is 10 maxr maximal number of step size reductions in a single time step default is 10 pe integer vector of length 2 fixing the components of x that will be plotted default settings are 0 1 for scalar odes i e plot z t versus t 1 2 for systems of odes i e plot the curve x1 t r2 t 1 ET and increase it Note that the user defined functions f and df must have the two argu ments t and x in this order The function crank_nicolson_vs implements the Crank Nicolson me thod see
27. plots the interpolation polynomial the original function and the error in the same way as the function lagrange The function cubic_spline realizes interpolation by natural cubic splines see 5 Satz 1 3 4 8 Algorithmus III 8 When interpolat ing a given function the user must provide the following arguments the first two being mandatory f function that should be interpolated number n of interpolation points at least 2 first interpolation point a default is 1 last interpolation point b default is 1 spacing of interpolation points e equidistant points a b a 1 lt i lt n c Cebysev points a b a 1 cos m 4 1 lt i lt n default is e np number k of plot points per subinterval default is 10 top p When interpolating discrete data the role of the first two arguments changes f now is the vector of the function values that should be interpolated and n is replaced by the vector x of the corresponding abscissa Both vectors must have the same dimensions The function cubic_spline prints the maximal interpolation error and plots the in terpolation polynomial the original function and the error by drawing a piecewise linear interpolation of these functions at k equidistant points in each subinterval of a b The function goertzel performs trigonometric interpolation using the Goertzel scheme see 5 Algorithmus I 5 3 When interpolating a given function the user must provide
28. ptimal value c x and vector x at the end 3 print value c x at every iteration 4 5 print value x and vector x at every iteration print value ctz vector x and simplex tableau at every iter ation 20 R VERFURTH The function dual_simplex implements the dual simplex method see 7 Algorithmus 1 5 5 It has the following arguments the first three being mandatory clp cost vector c R clp should be a row vector Alp constraint matrix A R blp constraint vector b R blp should be a column vector Jlp first basis J default setting is J n m 1 n prl print level determines the print output default is 2 0 no print output print optimal value c x and vector x at the end print value c x at every iteration print value x and vector x at every iteration print value ctz vector x and simplex tableau at every iter ation The function inner_point implements an inner point method for lin ear optimization see 7 Algorithmus I 7 7 It has the following ar guments the first three being mandatory clp cost vector c R clp should be a row vector Alp constraint matrix A R blp constraint vector b R blp should be a column vector xys initial column vector of the form z y s satisfying Ax b Ay s c z gt 0 s gt 0 default setting is xys inner_point then trie
29. rature formula with m interior Gau points and order 2m 1 see 5 811 3 8 1V 4 The value of m is passed as optional sixth argument ng The value of ng may be 2 3 or 4 with 4 being the default value The function romberg realizes the Romberg scheme see 5 Definition 11 4 6 8 IV 6 It has six arguments the first five being the same as those of midpoint_rule The optional sixth argument ne determines the number of columns of the Romberg scheme Its default value is 5 and it should be at least 2 and at most min 5 nr 4 7 Nonlinear equations Numerics provides the functions newton secant_rule regula_falsi for the solution of nonlinear systems of equations The function newton implements Newton s method with damping for systems of nonlinear equations f z 0 in R see 5 Definition 111 2 1 8 88I1 3 I1 5 It has the following arguments the first three being mandatory x initial value f function f NUMERICS USER GUIDE 9 df derivative Df of f pr if t results are printed in each iteration default is nod if t damping is switched off default is i e damping is performed tol error tolerance g iteration is stopped once f x f x lt default is 1 0D 8 nmax maximal number of iterations default is 10 The function newton prints the required number of iterations the final residual f x f x and the approximate solution x The function secant_rule realizes t
30. rdinary differential equa tions in R a t ft elt to lt t lt T x to To NUMERICS USER GUIDE 15 with variable step size methods The function rkf_vs controls the step size by comparing methods of different order see 6 Algorithmus 1 4 6 8 Algorithmus VI 5 the other functions perform the step size con trol by comparing results for different step sizes see 6 Algorithmus 1 4 3 8 Algorithmus VI 6 All functions plot the computed solution according to their argument pc The function explicit_euler_vs implements the explicit Euler scheme see 6 Algorithmus 1 2 1 8 IV 1 10 1 2 2 It has the following arguments the first four being mandatory f force function f x0 initial value xo T final time T dt initial step size tO initial time to default is t 0 tol required tolerance e the function strives at obtaining an approximation n t for the solution x t such that n t x t lt e x t holds for all t default value is 0 001 ntmax maximal number of time steps default is 10000 sfdt safety factor for step size variation reduce the step size at most by the factor X and increase it at most by the factor sfdt default value is 10 maxr maximal number of step size reductions in a single time step default is 10 pc integer vector of length 2 fixing the components of x that will be plotted default settings are 0 1 for scalar odes i e pl
31. row vector Alp constraint matrix A R blp constraint vector b R blp should be a column vector prl print level determines the print output default is 2 0 no print output 1 print optimal value c x at the end 2 print optimal value c x and vector x at the end 3 print value ctz at every branch and bound or cutting planes iteration 4 print value cx and vector x at every branch and bound or cutting planes iteration 5 print value c x at every simplex iteration 6 print value c x and vector x at every simplex iteration The function dijkstra implements the Dijkstra algorithm for finding a shortest s t path in a graph see 7 Algorithmus II 3 3 It prints results according to its argument prl and has the following arguments the first one being mandatory Adjl one of the following two objects e a square matrix such that Adjl i j gives the length of the arc connecting nodes i and j of the graph and equals hinf if nodes i and j are not connected by an arc or if i j e a list of vectors a b c describing the properties of the arcs such that c is the length of the arc connecting nodes a and b first number of starting node s default is 1 last number of terminal node t default is size Adj1 c pri print level default is 2 0 no print output 1 print length of shortest path at end 2 print length of shortest path and path itself at end 3 print all intermediate values of ve
32. s the functions simplex autostart_simplex dual_simplex inner_point for linear optimization problems minimize c x subject to Ax b x gt 0 All functions print values of the objective function ctz and of the solu tion vector x depending on their argument prl The function simplex implements the basic simplex method see 7 Al gorithmus 1 4 5 It has the following arguments the first three being mandatory clp cost vector c R clp should be a row vector Alp constraint matrix A R blp constraint vector b R blp should be a column vector Jlp first basis J default setting is J n m 1 n prl print level determines the print output default is 2 0 1 2 print optimal value c x and vector x at the end 3 print value c x at every iteration 4 print value c x and vector x at every iteration 5 print value cz vector x and simplex tableau at every iter ation The function autostart_simplex implements the simplex method with automatic initialization see 7 Algorithmus I 4 14 It has the follow ing arguments the first three being mandatory clp cost vector c R clp should be a row vector Alp constraint matrix A R blp constraint vector b R blp should be a column vector prl print level determines the print output default is 2 0 no print output 1 print optimal value c x at the end 2 print o
33. s to automat ically determine vectors x y and s with the above properties prl print level determines the print output default is 2 0 no print output 1 print optimal value c x at the end 2 print optimal value c x and vector x at the end 3 print value c x at every iteration 4 print value c x and vector x at every iteration exact default is f if t inner point tries to compute the exact solution of the optimization problem by guessing an admissible basis from the computed approximate solution and starting a simplex algo rithm with this basis itm maximal number of iterations default is 1000 tol tolerance default is 0 001 the iteration is stopped if the relative change of two consecu tive x vectors in the maximum norm is less than e or if the parameter u of inner point method is less than e 4 15 Discrete optimization Numerics provides the functions branch_and_bound cutting planes NUMERICS USER GUIDE 21 dijkstra floyd_warshall ford_fulkerson minimal_cost_flow for discrete optimization The functions branch_and_bound and cutting planes implement branch and bound and cutting planes methods resp for integer op timization problems minimize cz subject to x Z Ar b x gt 0 see 7 Algorithmen 11 1 8 II 1 12 They print results according to their argument prl and have the following arguments the first three being mandatory clp cost vector c R clp should be a
34. stall Scilab and Numerics respectively 2 INSTALLING SCILAB Scilab is a numerical programming environment comparable to Matlab with a similar syntax It was developed by INRIA France and is currently maintained by the Digiteo foundation in collaboration with INRIA Contrary to Matlab Scilab is completely free of charge There are versions for the operating systems Mac OsX Unix and Windows To use Numerics you first have to get the appropriate version of Scilab from http www scilab org and install it on your computer following the instructions of that web site Useful introductions to Scilab are e g 1 2 4 may be of par ticular interest for those familiar with Matlab 3 INSTALLING NUMERICS To install Numerics click the link Numerics zip on http www rub de num1 softwareE html or http www rub de num1 software html This installs a zip archive Numerics zip on your computer Unpack the zip archive and put the created folder Numerics at any place on your computer that suits you Now you can start Numerics as described in section 4 3 4 USING NUMERICS 4 1 Data structures We refer to 1 2 4 for a detailed introduction to Scilab Here we briefly present those data structures of Scilab that are used for input arguments of Numerics functions Boolean Boolean arguments are transferred to a Numerics function by either directly entering t for the value true or f for the value false or by creating an own boo
35. tO initial time to default is t 0 pe integer vector of length 2 fixing the components of x that will be plotted default settings are 0 1 for scalar odes i e plot z t versus t 1 2 for systems of odes i e plot the curve 21 t z2 t T to nt Note that the user defined functions f and df must have the two argu ments t and x in this order The function crank_nicolson_fs implements the Crank Nicolson me thod see 6 Algorithmus 1 2 3 8 IV 1 10 1 2 2 It has the same arguments as implicit_euler_fs The function rk_fs implements the classical Runge Kutta scheme see 6 81 2 8 1V 6 4 10 81 2 3 It has the same arguments as the function explicit_euler fs The function rkf_fs implements Runge Kutta Fehlberg schemes of or der 2 with 3 stages of order 3 with 4 stages and of order 4 with 6 stages see 6 Beispiel 1 4 7 8 Beispiel VI 7 It has the same arguments as the function explicit_euler_fs except that pc now is the seventh ar gument and that a new sixth argument order determines the method This argument may take the values 2 3 and 4 with 3 being its default value The function sdirk_fs implements a strongly diagonally implicit Runge Kutta scheme of order 4 see 8 VI 4 10 81 2 4 It has the same arguments as implicit_euler fs The function nystroem_fs implements the Nystr m two step method see 6 Beispiel 1 5 1 The first two approximations are computed with the
36. the following arguments the first one being mandatory 8 R VERFURTH f function f that should be interpolated n interpolate f at the points oe Oa 20 np number k of plot points default is 101 When interpolating discrete data the role of the first argument changes f now is the vector of the function values that should be interpolated The function goertzel prints the maximal interpolation error and plots the interpolation polynomial the original function and the error by drawing a piecewise linear interpolation of these functions at k equidis tant points on the interval 0 27 4 6 Numerical integration Numerics provides the functions midpoint _rule trapezoidal_rule simpson_rule gauss _rule romberg for numerical integration The functions midpoint_rule trapezoidal_rule simpson_rule im plement the corresponding composite quadrature formulae see 5 Bei spiel 11 2 2 3 8 IV 2 They all have the following arguments the first one being mandatory and print the computed approximate value of the integral f function that should be integrated a left endpoint of integration interval default is 0 b right endpoint of integration interval default is 1 nc initial number n of subintervals default is 1 nr number k of refinement steps default is 0 in the th refinement step 0 lt lt k the integration interval is split into n2 subintervals The function gauss_rule implements a composite Gau quad
37. tions The functions richardson jacobi gauss_seidel and ssor imple ment the corresponding classical stationary iterative methods see 5 Beispiel IV6 3 Algorithmus IV 7 14 10 8811 2 2 11 2 5 They all have the following arguments the first two being mandatory A matrix A 10 R VERFURTH b right hand side vector b x initial guess xo default is zo 0 om damping parameter default is ae A for richardson 1 for jacobi and gauss_seidel and 1 5 for ssor pri print level default is 0 0 only print required number of iterations n final residual Axn b Az b and mean convergence rate 1 print actual iteration number n actual residual Ar Agn b De 2 same as 0 additionally print the solution vector n 3 same as 1 additionally print the current vector n tol tolerance default is 1 0D 8 nmax maximal number N of iterations default is 100 Azn bll lt rn N Aroa and actual convergence rate iteration is stopped if The functions cg and pcg_ssor realize the conjugate gradient algorithm see 5 Algorithmus IV 7 7 9 Algorithmus VI 4 1 10 III 3 2 and the preconditioned conjugate gradient algorithm with SSOR pre conditioning see 5 Algorithmen IV 7 10 IV 7 14 9 Algorithmen V1 4 2 V1 4 3 10 88I11 3 3 III 3 4 cg has the same arguments as richardson except the missing relaxation parameter om pcg_ssor has the sa
38. with all elements equal to 1 and 0 respectively eye n n creates the n x n identity matrix Matrices are transferred to a Numerics function by either directly entering the matrix or by creating an own matrix valued variable mymatrix and then calling the Numerics function with mymatrix as argument Sparse matrix Scilab has a built in data structure for sparse matrices which differs from the corresponding structure of Matlab A sparse matrix with n non zero entries is stored in an n x 1 column vector containing the non zero entries and an n x 2 matrix containing the row and column indices of the non zero entries Suppose you want to store the matrix 10 0 2 00 3 0 40 10 in sparse format Then you must enter in the Scilab console without the full stop at the end myentries 1 2 3 4 1 myindices 1 1 1 4 2 3 4 1 4 3 mysparsematrix sparse myindices myentries 4 R VERFURTH Note that the input myoentries 1 1 2 3 4 myoindices 4 3 1 1 1 4 2 3 4 1 myosparsematrix sparse myoindices myoentries represents the same matrix This matrix in sparse format is transferred to a Numerics function by using mysparsematrix or myosparsematrix or an equivalent construct as argument Function Functions user defined or built in ones can act as argu ments of Numerics functions To clarify this consider two examples Suppose first that you want to interpolate the built in exponential func tion in 10 equ
Download Pdf Manuals
Related Search