From b8e4ff72039227a57ee773962685e13011b66e79 Mon Sep 17 00:00:00 2001 From: Giovanni Di Grezia Date: Wed, 17 Dec 2014 19:33:04 +0100 Subject: [PATCH] aggiunto note, piccoli fix a eigenvalues, aggiunto polinomio di chebyshev --- functions/{eigenvalues.m => eigenvalues_qr.m} | 4 +- functions/external/ChebyshevRoots.m | 154 ++++++++++++++++++ functions/{ => external}/float2bin.m | 0 samples/autovalori,eig ,qr e power_method.m | 24 +++ samples/chebyshev.m | 2 + samples/plot,functions,polyval.m | 15 +- samples/quadratura e integrali.m | 4 + samples/spline cubic.m | 5 + samples/spline.m | 15 -- 9 files changed, 201 insertions(+), 22 deletions(-) rename functions/{eigenvalues.m => eigenvalues_qr.m} (69%) create mode 100644 functions/external/ChebyshevRoots.m rename functions/{ => external}/float2bin.m (100%) create mode 100644 samples/autovalori,eig ,qr e power_method.m create mode 100644 samples/chebyshev.m create mode 100644 samples/quadratura e integrali.m create mode 100644 samples/spline cubic.m delete mode 100644 samples/spline.m diff --git a/functions/eigenvalues.m b/functions/eigenvalues_qr.m similarity index 69% rename from functions/eigenvalues.m rename to functions/eigenvalues_qr.m index 0819f12..cbc246d 100644 --- a/functions/eigenvalues.m +++ b/functions/eigenvalues_qr.m @@ -1,10 +1,10 @@ -function [ values, err ] = eigenvalues( A, ite, toll ) +function [ values, err, A ] = eigenvalues_qr( A, ite, toll ) i=0; err = inf; while (i < ite) && (err > toll) [q,r] = qr(A); - A1= r * q + A1= r * q; err = norm(diag(A1) - diag(A),inf) / norm (diag(A1),inf); A = A1; i=i+1; diff --git a/functions/external/ChebyshevRoots.m b/functions/external/ChebyshevRoots.m new file mode 100644 index 0000000..3e6ccaf --- /dev/null +++ b/functions/external/ChebyshevRoots.m @@ -0,0 +1,154 @@ +%{ +Copyright (c) 2009, Russell Francis +All rights reserved. + +Redistribution and use in source and binary forms, with or without +modification, are permitted provided that the following conditions are +met: + + * Redistributions of source code must retain the above copyright + notice, this list of conditions and the following disclaimer. + * Redistributions in binary form must reproduce the above copyright + notice, this list of conditions and the following disclaimer in + the documentation and/or other materials provided with the distribution + +THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS" +AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE +IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE +ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT OWNER OR CONTRIBUTORS BE +LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR +CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF +SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS +INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN +CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) +ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE +POSSIBILITY OF SUCH DAMAGE. +%} + +function x = ChebyshevRoots( n, type, range ) +% USAGE: x = ChebyshevRoots( n [,type [, range]] ) +% +% This method returns the values of the roots of the Chebyshev polynomial, +% either type one or type two of degree n. In the literature, these are +% often referred to as T_n(x) for type 1 or U_n(x) for type 2. If the +% optional parameter type is omitted, it is assumed to be type 1. This +% function also supports scaling and translating the roots of the +% polynomial to lie within a range specified. The polynomials T_n(x) and +% U_n(x) have n roots which lie within [1,-1]. It is often useful to scale +% and translate these roots to have the same relative distance but lie over +% a different range so that they may be used as the nodes for interpolation +% of a function which does not lie within this [-1, 1] range. +% +% PARAMETERS: +% n - The degree of the Chebyshev polynomial whose roots we wish to +% find. +% +% type [optional] - Either 'Tn' for type 1 or 'Un' for type 2 +% depending on whether you wish to generate the roots of the type 1 +% or type 2 polynomial. 'Tn' is the default if this parameter is +% omitted. +% +% range [optional] - The Chebyshev polynomial is defined over [-1, 1] +% this parameter allows the roots of the polynomial to be translated +% to be within the range specified. ie. The relative distance of the +% Chebyshev nodes to each other will be the same but their values +% will span the provided range instead of the range [-1, 1]. +% +% RETURNS: +% A vector with the n roots of the Chebyshev polynomial T_n(x) or +% U_n(x) of degree n, optionally scaled to lie within the range +% specified. +% +% AUTHOR: +% Russell Francis +% +% THANKS: +% John D'Errico - Provided reference to the Abramowitz and Stegun book +% which provides a rigorous definition of the two types of Chebyshev +% polynomials and is suggested reading particularly chapter 22 for +% interested parties. +% +% REFERENCES: +% +% [1] Abramowitz, M.; Stegun, I.A. Handbook of Mathematical Functions with +% Formulas, Graphs, and Mathematical Tables. U.S. Department of Commerce. +% Online version available at: +% http://www.knovel.com/knovel2/Toc.jsp?BookID=528&VerticalID=0 +% +% [2] http://en.wikipedia.org/wiki/Chebyshev_polynomials +% +% [3] Burden, Richard L.; Faires, J. Douglas Numerical Analysis 8th ed. +% pages 503-512. +% + +%% +% Verify our parameters. +%% + +% The degree must be specified and must be greater than or equal to 1. +if( nargin() < 1 ) + error( 'You must provide the parameter n which is the degree of the polynomial you wish to calculate the roots of.' ); +else + if( n < 1 ) + error( 'The parameter n must be greater than or equal to 1.' ); + end +end + +% The type of the Chebyshev polynomial to calculate the roots of, optional +% and defaults to T_n(x) +if( nargin() < 2) + type = 'Tn'; +else + if( (strcmp( type, 'Un' ) ~= 1 ) && (strcmp( type, 'Tn' ) ~= 1) ) + error( 'The type parameter which was specified is not valid!, Please specify either "Tn" or "Un"' ); + end +end + +% The range which we wish to scale and translate the result to, optional +% and the default is to not scale or translate the result. +if( nargin() < 3 ) + range = [-1 1]; +else + if( length( range ) ~= 2 ) + error( 'The parameter range must contain two values.' ); + end + + if( range(1) == range(2) ) + error( 'The parameter range must contain two distinct values.' ); + end + + range = sort(range); +end + +%% +% Begin to compute our Chebyshev node values. +%% +if( n == 1 ) + x = [0]; +else + if( strcmp( type, 'Tn' ) == 1 ) + x = [(pi/(2*n)):(pi/n):pi]; + else + x = [(pi/(n+1)):(pi/(n+1)):((n*pi)/(n+1))]; + end + x = sort( cos(x) ); +end + +%% +% x now contains the roots of the nth degree Chebyshev polynomial, +% we need to scale and translate the result if necessary. +%% +if( (range(1) ~= -1) || (range(2) ~= 1) ) + M = eye(n+1); + % Calculate the scaling factor to apply to the nodes. + sf = (range(2) - range(1)) / 2; + % Calculate the translation to apply to the nodes. + tl = range(1) + sf; + % Generate our transformation matrix. + M(1:n,1:n) = M(1:n,1:n) * sf; + M(n+1,1:n) = tl; + % Apply to our earlier result. + x = [x 1] * M; + x = x(1:n); +end +return; % The x values of the Chebyshev nodes. diff --git a/functions/float2bin.m b/functions/external/float2bin.m similarity index 100% rename from functions/float2bin.m rename to functions/external/float2bin.m diff --git a/samples/autovalori,eig ,qr e power_method.m b/samples/autovalori,eig ,qr e power_method.m new file mode 100644 index 0000000..ab7e590 --- /dev/null +++ b/samples/autovalori,eig ,qr e power_method.m @@ -0,0 +1,24 @@ +eig(a) autovalori della matrice a + +qr(a) fattorizzazione matric a in q*r con q ortogonale(inversa e trasposta uguale) e r triangolare superiore + 2/3 n 3 operazioni. + +eningenvalues_qr(a, iterazioni,tolleranza) ccalcola tuttgli gli autovalori usando qr + se a ha autovalori reali la matrice converge a triangolare superiore + e si richiede che abs(a1) > abs(a2) > abs(a3) ... + la velocita con cui gli elementi sotto la diagonale convergono a 0 e' bassa se gli autovalori sono molto vicini tra loro + se a ha autovalori complessi la matrice converge a quasi triangolare superiore con elmenti non zero sotto la diagonale + + extra: + Teorema: A reale simmetrica, esiste Q ortogonale t.c. + D = QTA Q diagonale + (la trasformazione non avviene in un numero finito di passi) + + Teorema: A reale, esiste Q ortogonale t.c. QTA Q è di + 2 + Hessenberg (tridiagonale se A è simmetrica) + (la trasformazione avviene in un numero finito di passi) + + +power_method funziona solo se il massimo in valore assoluto degli autovalori di una matrice e' unico + e' lento se gli autovalori sono molto vicini tra loro \ No newline at end of file diff --git a/samples/chebyshev.m b/samples/chebyshev.m new file mode 100644 index 0000000..9c02c1b --- /dev/null +++ b/samples/chebyshev.m @@ -0,0 +1,2 @@ +ChebyshevRoots(6,'Tn',[0.1 0.6]) +input: grado, tipo, range \ No newline at end of file diff --git a/samples/plot,functions,polyval.m b/samples/plot,functions,polyval.m index a153998..5a32f43 100644 --- a/samples/plot,functions,polyval.m +++ b/samples/plot,functions,polyval.m @@ -1,10 +1,15 @@ -fplot( @(x) x^2 * cos(x), [-4 4],'r*' ) +f=@(x) x^2 * cos(x) usa ./ .* .^ per divisioni moltiplicazioni elevazioni elemento per elemento +fplot(f, [-4 4],'r*' ) hold on -p= [1,2,3,4,5] //coefficinti del polinomio x=1:0.1:10 oppure x=linspace(-5,5,100) 100 punti tra -5,5 -y=polyval(p,x) -plot(x,y) +y=f(x) + +p= polyfit(x,y,n) // n è il grado, 1 in meno dei punti +p= [1,2,3,4,5] //coefficinti del polinomio + +xx=1:0.01:10 +yy=polyval(p,xx) +plot(xx,yy) -p= polyfit(x,y,n) // n è il grado, 1 in meno dei punti \ No newline at end of file diff --git a/samples/quadratura e integrali.m b/samples/quadratura e integrali.m new file mode 100644 index 0000000..9e90821 --- /dev/null +++ b/samples/quadratura e integrali.m @@ -0,0 +1,4 @@ +f = @(x) x.^2 + 2; +integral(f,-6,6) integrale matlab +quad(f,-6,6) o quad(f,-6,6,1e-5) quadratura integrale matlab + diff --git a/samples/spline cubic.m b/samples/spline cubic.m new file mode 100644 index 0000000..389139c --- /dev/null +++ b/samples/spline cubic.m @@ -0,0 +1,5 @@ +x=[-2 0 2 3 4 5]; +y=[4 0 -4 -30 -40 -30]; +xx=-2:0.1:5; +yy=spline(x,y,xx); +plot(xx,yy,'g') diff --git a/samples/spline.m b/samples/spline.m deleted file mode 100644 index d943046..0000000 --- a/samples/spline.m +++ /dev/null @@ -1,15 +0,0 @@ -x=[-2 0 2 3 4 5]; -y=[4 0 -4 -30 -40 -30]; -plot(x,y,'sr') -hold on -p=polyfit(x,y,5) - -p = - - -0.2286 2.5333 -6.2286 -10.1333 26.5714 0.0000 - -xx=-2:0.1:5; -yy=polyval(p,xx); -plot(xx,yy,'b') -ys=spline(x,y,xx); -plot(xx,ys,'g')