octave> version
ans = 4.0.1
$> octave -q <<NB
> eval ("2 + 4.1")
> NB
ans = 6.1000
octave knows about arithmetic operations (+ ,- , *, /), exponentiation (^), natural logarithms/exponents (log, exp), and the trigonometric functions (sin, cos, tan)
octave calculations work on real or imaginary numbers (i,j)
you can save vars to file (save) and then you can load them from this file (load) S
octave:> load "filenameIn.m" octave:> t = 1:0.001:2; octave:> plot (t, sin(t)); octave:> saveas (1, "filenamePic1", "fmt"); octave:> function f (a) > 2 * a + 1 > end octave:> x = [1,2,3,4,5,6]; octave:> y = [0.16, 0.16, 0.16, 0.16,0.16, 0.16]; octave:> bar (x, y); octave:> saveas (2, "filenamePic2", "fmt"); octave:> save "filenameOut.m" octave:> quit
the base of the natural logarithm (e) and the ratio of a circle's circumference to its diameter (pi) are pre-defined
to display the value of a variable, simply type the name of the variable at the prompt
you can recall, edit, and reissue previous commands using Emacs/vi style editing commands
once you know the name of the function you wish to use, you can get more help on the function by simply including the name as an argument to help. for example,
help plotwill display the help text for the plot function.
ending a command with a semicolon tells Octave not to print the result of the command.
a = menu ("example", "say bye", "say good", "quit"); if (a == 3) disp ("quit") elseif (a == 1) disp ("bye") else disp ("good") endif; quit
for example,
myfile = fopen ("splat.dat", "r"); the possible values `mode' may have are `r' for reading `w' for writing. the previous contents are discared `a' open or create a file for writing at the end `r+' open an existing file for reading and writing `w+' open a file for reading or writing. the previous contents are discared `a+' open or create a file for reading or writing at the end
$> cat a.dat 1 2 3 4 5 6 7 8 9 $> octave -q octave:1> a = load("a.dat") a = 1 2 3 4 5 6 7 8 9 octave:2> a a = 1 2 3 4 5 6 7 8 9
if there's a text header, the file won't load. if you have more than two fields at any point in the file, it won't load. the input file has to have the same number of columns (and each column must have the same number of rows), otherwise octave will complain about the data
> fid = fopen ('data.dat', 'r') ; > data = fscanf (fid, '%f') ; > nRows = data (1) ; > data = reshape (data(2:end), nRows, nRows) ;
you can use the function dlmread():
so you can simply use the one-liner:
data = dlmread ('input_file', ' ', 1, 0);by calling the function with r0 set to 1, you're effectively skipping the first line, which contains the (now useless) number of rows
to plot line from (a1,b1) to (a2,b2) use
octave:> plot ([3,4],[7,8])
octave:> t = linspace(-1,2,100); octave:> plot (t, arrayfun(@ (x) -.5*(x^2 - x) , t))
to plot a vertical line
octave:> plot ([x; x], [0; 5])
to plot a horizontal line
octave:> plot ([0; 40], [y; y])
text (10, 80, strcat ("slope is ", num2str(a), "\nintersect is ", num2str(b)));
which corresponds to> p = [1 2 -1]
x² + 2*x - 1
> p = [1, -1, 1]; > polyval (p, 4) ans = 13
> p = [1, -1, 1]; > roots (p) ans = 0.50000 + 0.86603i 0.50000 - 0.86603i
> p p = 2 1 > q q = 1 2 > % (2*x + 1) * (x + 2) = 2*x^2 + 5*x + 2 > conv(p, q) ans = 2 5 2
the first calling form computes the partial fraction expansion for the quotient of the polynomials, B and A
K - vector is a polynomial of order N-1 representing the direct contribution,
E -vector specifies the multiplicity of the m-th residue's pole.
> b = [1, 1, 1] ; > a = [1, -5, 8, -4] ; > [r, p, k, m] = residue (b, a) ; > r r = -2.0000 7.0000 3.0000 > p p = 2.00000 2.00000 1.00000 > k k = [](0x0) > m m = 1 2 1
which represents the following partial fraction expansion
s^2 + s + 1 -2 7 3 ------------------- = ----- + ------- + ----- s^3 - 5s^2 + 8s - 4 (s-2) (s-2)^2 (s-1)
> f = [3:7]' f = 3 4 5 6 7
octave:> A = [ 1, 1, 2; 3, 5, 8; 13, 21, 34 ] ;
for example,
B = rand (3, 2);will create a 3 row, 2 column matrix with each element set to a uniformly distributed random value between zero and one
functions like exp, log, sin, cos, sqrt, etcetera, can also be applied to matrices. in this case, they are simply applied to every element of the matrix
the first and the second also have a linear algebraic counterpart: the matrix exponential and matrix logarithm. the names of the corresponding functions is formed with the suffix m, so the matrix exponential is denoted by expm, matrix logarithm by logm
function Outputs = function_name (Par1, Par2, ...)
...
end
return
exit
exit(Status)
function [z] = vect_mul (x , y) d = length (x) ; if (d != length (y)) printf ("vectors are not equal\n") ; z = [0] ; return ; else for i = 1 : d z(i) = x(i) * y(i) ; endfor endif endfunction
octave> function last_words () > disp ("bye bye kanzas"); > endfunction octave> atexit ("last_words");
octave:> f = @(x) x * x + 2 ; octave:> f(5) ans = 27
arrayfun (@ (x) 1/x, Param)
function r = computeit (f, g, v)
r = f (g (v));
endfunction
computeit (@exp, @sin, pi/3)
computeit (@log, @cos, pi/6)
or pass the string name of the function and use the feval primitive.
function r = computeit2 (f, g, v)
r = f (feval (g, v));
endfunction
computeit2(@exp, "sin", pi/3)
f - the name of the function, the function handle or the inline function.
a and b - limits of integration. either or both may be infinite
the optional argument tol is a vector that specifies the desired accuracy of the result. the first element of the vector is the desired absolute tolerance, and the second element is the desired relative tolerance
the optional argument sing is a vector of values at which the integranl is known to be singular
the result of the integration is returned in v
ier contains an integer error code (0 indicates a successful integration)
nfun indicates how many function evaluations were required
err contains an estimate of the error in the solution
the function over with to integrate is defined by f, and the interval for the integration is defined by [xa, xb, ya, yb]. the function f must accept a vector x and a scalar y, and return a vector of the same length as x. If defined, tol defines the absolute tolerance to which to which to integrate each sub-integral. additional arguments, are passed directly to f
to use the default value for tol one may pass an empty matrix
dx/dt = f(x, t) x(0) = x₀
octave has built-in functions for solving nonlinear differential equations of the form dx/dt = f (x, t)
x(t0) = x0
f (x, t)
it will create a file called foo.eps that contains a rendering of the current plot in Encapsulated PostScript formatoctave:> function [y] = f (x, t) > r = 0.25; k = 1.4; > a = 1.5; b = 0.16; > c = 0.9; d = 0.8; > > y(1) = r * x(1) * (1 - x(1) / k) - a * x(1) * x(2) / (1 + b * x(1)); > y(2) = c * a * x(1) * x(2) / (1 + b * x(1)) - d * x(2); > endfunction > x0 = [1; 2]; > t = linspace (0, 50, 200); > x = lsode ("f", x0, t); > plot (t, x); > print -deps foo.eps
octave can handle both 'stiff' and 'nonstiff' equations
lsode expects three arguments:
the function returns a matrix of x in times, specified by user
options can be specified by lsode_options function
after a computation, each row of the matrix A corresponds to an element of the vector t
after a successful computation, the value of istate will be 2 (consistent with the Fortran version of Lsode). if the computation is not successful, istate will be something other than 2 and msg will contain additional information
#!/home/guest/.nix-profile/bin/octave
% equation of moving in 3Dspace
% dr/dt = v(t) = [x/t, y/t, z/t]
1;
function v = some_func (x, t)
d = exp (-t) ;
v = (eye (3) * d) * x ;
end;
function diff3 (fdot, x0, t)
f = lsode (fdot, x0, t) ;
plot (t, f, 'linewidth', 2) ; grid ;
legend (" x " , " y ", " z ") ;
xlabel ("time") ; ylabel ("position") ;
end;
x0 = [1; 2; 3] ;
t = linspace (1, 10, 200) ;
diff3 (@some_func, x0, t) ; pause ;
y" = -sin(y) + sin(5*t) and the initial conditions: y (0) = 1 y'(0) = 0
now, let y₁=y and y₂=y'
this gives you the first order system:
y₁' = y₂ y₂' = -sin(y₁) + sin(5*t) with the initial conditions y₁(0) = 1 y₂(0) = 0
f (x, x', t) = 0 x (0) = x₀ x'(0) = x'₀
the function is called dassl, and it expects four arguments:
the first element of t should be t0 and correspond to the initial state of the system x0 and its derivative xdot0, so that the first row of the output A is x0 and the first row of the output B is xdot0
the first argument, fcn, is a string, that names the function f to call to compute the vector of residuals for the set of equations. it must have the form res = f (x, xdot, t) in which x, xdot, and res are vectors, and t is a scalar
#!/home/guest/.nix-profile/bin/octave
function yourdif (fdot, x0, d0, t)
[A, B, c, d] = dassl (fdot, x0, d0, t) ;
plot ( t , A , "r" , 'linewidth' , 2 , t , B , "b" , 'linewidth' , 2 ) ;
legend ("func" , "deriv") ; pause ;
end
function r = some_func_here (x, d, t) r(1) = 3 - d ; end
x0 = 5 ; d0 = 3 ; t = linspace (1, 10, 200) ;
yourdif ("some_func_here", x0, d0, t) ;
defining a sinusoidal signal:
freq1 = 100;
period1 = 1 / freq1;
w1 = 2 * pi * freq1;
num_tsteps = 1000;
num_periods = 2;
tstep = num_periods * period1 / num_tsteps;
t = 0:tstep:(num_periods * period1);
x1 = sin(w1 * t);
plot(t, x1);
modulate it with a 2.5kHz signal:
freq2 = 2500;
period2 = 1 / freq2;
w2 = 2 * pi * freq2;
x2 = sin(w2 * t);
x3 = x1 .* x2;
plot(t, x3);
the function
the list returned from the function fft isn't dependent on any particular timescale (provided the signal is sampled at a constant rate, called the sampling frequency). it's up to you to set a timescale that makes sense
you know the list from fft has n points - the list of associated frequencies will therefore contain n points. since the signal was sampled with a constant f[s], the spacing between the frequencies will also be equal
now I declare the frequency list as a range of n numbers, which will then be scaled by an appropriate amount:
f = (0:n-1) * scale;
what's the scaling factor? it's easiest to first normalize the frequencies so they run from 0 to 1 (non-inclusive) by dividing by n
the sampling frequency is the maximum frequency because one period of the Fourier transform represents an amount of time equal to the sampling period (1/f_s). I can therefore rewrite the last statement as:
f = f_s * (0 : n-1)/n;
in our case, f_s = num_tsteps / (period1 * num_periods)
. pulling all of this together, let's define list y as the Fourier transform of the AM signal, and f as its associated frequencies:
n = 2**16;
y = abs (fft (x3, n));
f_s = num_tsteps / (period1 * num_periods);
f = f_s * (0 : n-1) / n;
the list returned from fft contains complex numbers, but we are only interested in the magnitude of each frequency component of the signal
let's plot the FFT:
axis ([0, 5000]);
plot (f, y);
this symmetry, however, is not a reflection across f = 0 Hz, as it should be. to view the true spectrum it's necessary to use the function fftshift (v), which shifts the values in v, a list containing the results of a Fourier transform, such that the frequency 0 Hz is at the center of the list. you'll also have to shift the corresponding frequencies:
y = abs (fft (x3, n));
y = fftshift (y);
f = f_s * (-n/2 : n/2-1) / n;
axis ([-5000, 5000]);
plot (f, y);
there are more fft related functions
injecting some random noise by using the randn (rows, columns), which returns a matrix of random numbers with normal distribution
alternatively, we could use rand (rows, columns), but for this purpose we want values centered about zero
let's plot the new modulated signal:
x3 = x3 + randn (1, length(x3));
plot(t, x3);
%% defining a sinusoidal signal ...
freq1 = 100 ; period1 = 1 / freq1 ; w1 = 2 * pi * freq1 ;
num_steps = 1000 ; num_periods = 4 ;
dt = num_periods * period1 / num_steps ;
t = 0:dt:dt * num_steps ;
x1 = sin (w1 * t) ; plot (t, x1, 'linewidth', 2) ; pause
%% ... and modulate it with a 250 Hz signal
freq2 = 250 ; period2 = 1 / freq2 ; w2 = 2 * pi * freq2 ;
x2 = sin (w2 * t) ;
x3 = x1 .* x2 ;
plot (t, x3, 'linewidth', 2) ; pause
%% ... and its spectrum now will be
n = 2^8 ;
y = fft (x3, n) ; df = 1/dt ; f = df * (0:n-1) / n ;
plot (f, y, 'linewidth', 2) ; pause
octave can be used to solve a linear program (LP)
LP is finding minimax for functions like that:
f(x₁,x₂,...) = c₁*x₁ + c₂*x₂ + .. a₁₁*x₁ + a₁₂*x₂ + ... ≥|≤ b₁ a₂₁*x₁ + a₂₂*x₂ + ... ≥|≤ b₂ ... 0 < xₙ < ∞ or lower < xₙ < ∞ or lower < xₙ < upper or 0 < xₙ < upper
that's all about it. example:
> A = [2,1;1,1;1,0] A = 2 1 1 1 1 0 > c = [3,2] c = 3 2 > b = [100,80,40] b = 100 80 40 > [x_max, z_max] = glpk (c, A, b, [0,0], [], "UUU", "CC", -1) x_max = 20 60 z_max = 180
there is glpk function to do the job:
LP : minimax c⁺ * x subject to a_{1i} * x cond b_{1}, ... a_{mi} * x cond b_{m}, x ∈ ℝⁿ c ∈ ℝⁿ b ∈ ℝⁿ
octave can be used to solve a quadratic program (QP)
function qp solves the quadratic program
QP: minimize .5 * x' * H * x + x' * q subject to A * x = b lb ≤ x ≤ ub A_lb ≤ A_in * x ≤ A_ub
any bound (A, B, lb, ub, A_LB, A_UB) may be set to the empty matrix [] if not present. if the initial guess is feasible the algorithm is faster
#!/usr/bin/octave -q
% real integrator with damping
x = zeros (64);
y = zeros (64);
x (2) = 1;
for t = 2:64
y (t) = x (t) + .9 * y (t - 1);
end;
plot (1:64, y, 'linewidth', 2)
pause;
a script file is a file containing (almost) any sequence of Octave commands. it is read and evaluated just as if you had typed each command at the Octave prompt, and provides a convenient way to perform a sequence of commands that do not logically belong inside a function it is possible to define more than one function in a single script file and load (but not execute) all of them at once
to do this, the first token in the file (ignoring comments and other white space) must be something other than function. if you have no other statements to evaluate, you can use a statement that has no effect, like this:
# prevent Octave from thinking that this is a function file:
1;
# define function one:
function one ()
…
to run this script just input filename in octave prompt
if the first token in a file (ignoring comments) is function, Octave will compile the function and try to execute it
unlike a function file, a script file must not begin with the keyword function. if it does, Octave will assume that it is a function file, and that it defines a single function that should be evaluated as soon as it is defined
a script file also differs from a function file in that the variables named in a script file are not local variables, but are in the same scope as the other variables that are visible on the command line
%% octave --persist lprg_1.m
%% or
%% octave:1> eval lprg_1
t = 0:1:30;
function y = a (t)
y = 50 - 2 * t;
end
function y = b (t)
y = 40 - 4 * t / 3;
end
function y = c (t)
y = 50 - 6 * t / 5;
end
hold; grid; xlabel ("factory A, days"); ylabel ("factory B, days");
plot (t, a(t), 'r', 'linewidth', 2)
plot (t, b(t), 'g', 'linewidth', 2)
plot (t, c(t), 'm', 'linewidth', 2)
plot([24,24],[0,24], 'linewidth', 2)
plot([0,24],[24,24], 'linewidth', 2)
set (0, "defaultlinelinewidth", 2)
octave> pkg prefix /usr/local/lib/octave/packages/4.0.0 octave> pkg install some-pkg-name.tar.gz octave> pkg install --forge some-pkg-name.tar.gz octave> pkg list octave> pkg uninstall some-pkg-name