Ac� te mando los scrips que hice de Scilab... no me acuerdo
cuales te mand� ya, pero te mando algunos. Son de este a�o.
Espero te sirvan.
Salu2
JANUS
P/D: el formato es P.E.[p].sce donde P=pr�ctica, E=ejercicio,
p=punto
--
El computador naci� para resolver problemas que antes no
exist�an.
-- Www.frases.com.
// Pr�ctica 1 - Ejercicio 4iii
function [n]=epsilon
c = 1;
while (1+c) > 1
n=c;
c = c/2;
end
// Pr�ctica 1 - Ejercicio 4ii
function [n]=maxpos
c=1;
while 1/c
n=c;
c=c*2;
end
// Pr�ctica 1 - Ejercicio 4i
function [n]=pos
c=1;
while c
n=c;
c=c/2;
end
// Ejercicio 17 - Pr�ctica 2
// X es el vector vertical inicial
function [X]=newt(X)
function J=Jinv(Y)
J=[2*Y(1,1)+Y(2,1)^3 3*Y(1,1)*Y(2,1)^2; 6*Y(1,1)*Y(2,1) 3*Y(1,1)^2-3*Y(2,1)^2];
J=inv(J);
endfunction
function FX=F(Y)
FX=[Y(1,1)^2+Y(1,1)*Y(2,1)^3-9; 3*Y(1,1)^2*Y(2,1)-4-Y(2,1)^3];
endfunction
K=F(X)
Tolerancia=10E-9;
while abs(K(1,1))>Tolerancia & abs(K(2,1))>Tolerancia
X=X-Jinv(X)*F(X);
K=F(X);
end
// Ejercicio 10b - Pr�ctica 2
function x=newt
x=0.5;
while abs(-log(x)-x)>10^(-9)
x=x-(-log(x)-x)/(-1/x-1);
end
// Pr�ctica 2 - Ejercicio 14
// (a,b) es el intervalo donde comenzar con el m�todo de bisecci�n
function r=raiz(a,b)
// Definici�n de la funci�n que estamos trabajando
// Para esta funci�n las ra�ces est�n en 0 y m�s o menos en 1.934, as� que usar un intervalo
// apropiado
function y=f(x)
y=x^2/4-sin(x);
endfunction
if ~f(a)
r=a;
elseif ~f(b)
r=b;
elseif f(a)*f(b)>0
disp ("Error en el intervalo elegido, f(a)*f(b) debe ser <= 0")
r=0;
else
// Bisecci�n hasta una tolerancia de 10E-3
tolerancia=10E-3
r=(a+b)/2;
while abs(f(r))>tolerancia
if f(r)*f(a)<0
b=r;
else
a=r;
end
r=(a+b)/2;
end
// Secante
while a~=b
r=a-f(a)*(a-b)/(f(a)-f(b));
a=b;
b=r;
end
end
// Ejercicio 10a - Pr�ctica 2
function x=punto_fijo
x=0.5;
while abs(-log(x)-x)>10^(-9)
x=exp(-x);
end
// Ejercicio 9 - Pr�ctica 2
function [rs,rn,e]=raiz
// rs: Ra�z calculada por el m�todo de la Secante
// rn: Ra�z calculada por el m�todo de Newton
// e: Diferencia entre ambas
// x0: Valor inicial
// M�todo de la Secante
x0=1;
x1=2;
rs=x1;
while tan(rs)-0.5*rs
rs=x1-(tan(x1)-0.5*x1)*((x1-x0)/((tan(x1)-0.5*x1)-(tan(x0)-0.5*x0)));
x0=x1;
x1=rs;
end
// M�todo de Newton
rn=1;
while tan(rn)-0.5*rn
rn=rn-(tan(rn)-0.5*rn)/(0.5+tan(rn)^2);
end
// Diferencia
e=abs(rs-rn);
// Ejercicio 20 - Pr�ctica 2
function i=veinte()
x=0.25;
i=1;
izq=(2^(i-1)*x^(2^(i-1)-1)-2^(i)*x^(2^(i)-1))/(1-x^(2^(i-1))+x^(2^(i)));
while abs(izq-((1+2*x)/(1+x+x^2)))>10E-6
i=i+1;
izq=izq+(2^(i-1)*x^(2^(i-1)-1)-2^(i)*x^(2^(i)-1))/(1-x^(2^(i-1))+x^(2^(i)));
end
// Pr�ctica 3 - Ejercicio 1
// Sisntaxis: X = gauss(T)
// Donde T es una matriz de n * n+1 ya que contiene al vector
// de t�rminos independientes.
///////////////////////
// M�todo de Gauss //
///////////////////////
function [X] = gauss(T)
M = T;
N = size(M);
n=N(1,1);
m=N(1,2);
if n~=m-1
disp('Error en las dimensiones de la matriz (debe ser de n x n+1)')
else
for k=1:n-1
j=k+1;
if M(k,k)==0
while j<=n & M(j,k)==0
j=j+1;
end
if j<=n
aux=M(j,:);
M(j,:)=M(k,:);
M(k,:)=aux;
else
j=0;
end
else
if j
for j=k+1:N(1,1)
M(j,:) = M(j,:)-(M(j,k)/M(k,k))*M(k,:);
end
end
end
end
j=1;
while j & j<=n
if ~M(j,j)
j=0;
disp ('Sistema con m�ltiples soluciones o sin soluci�n')
else
j=j+1;
end
end
if j
X=zeros(n,1);
j=n;
while j
h=0;
for k=j:n
h=h+M(j,k)*X(k);
end
X(j)=(M(j,m)-h)/M(j,j);
j=j-1;
end
end
end
// Pr�ctica 4 - Ejercicio 12
function X=gc(A,b)
N=size(A);
n=N(1,1);
p=zeros(n,1);
X=zeros(n,1);
r=zeros(n,1);
rant=zeros(n,1);
k=0;
while r
k=k+1;
if k==1
p(:,1)=r(:,1);
else
B=(r'*r)/(rant'*rant);
p=r+B*p;
end
t=(r'*r)/(p'*A*p);
X=X+t*p;
rant=r;
function p=fact(k)
p=prod(1:k);
endfunction
H=zeros(n,n);
Hi=zeros(n,n);
for i=1:n
for j=1:n
H(i,j)=1/(i+j-1);
Hi(i,j)=((-1)^(i+j)*fact(n+i-1)*fact(n+j-1))/((i+j-1)*(fact(i-1)*fact(j-1))^2*fact(n-i)*fact(n-j));
end
end
// K(H) para norma 1
norma1H=0;
norma1Hi=0
for j=1:n
sumH=0;
sumHi=0;
for i=1:n
sumH=sumH+abs(H(i,j));
sumHi=sumHi+abs(Hi(i,j));
end
if (sumH>norma1H)
norma1H=sumH;
end
if (sumHi>norma1Hi)
norma1Hi=sumHi;
end
end
K1=norma1H*norma1Hi;
// K(H) para norma 2
// Tomo como valor aprox de norma 2 a:
// raiz del radio espectral de H*H' donde H' es H tranpuesta
// Pero, como H es sim�trica, H=H' (lo mismo con Hi)
autovaloresH=spec(H*H);
autovaloresHi=spec(Hi*Hi);
radespH=0;
radespHi=0;
for i=1:n
if (abs(autovaloresH(i))>radespH)
radespH=abs(autovaloresH(i));
end
if (abs(autovaloresHi(i))>radespHi)
radespHi=abs(autovaloresHi(i));
end
end
norma2H=sqrt(radespH);
norma2Hi=sqrt(radespHi);
K2=norma2H*norma2Hi;
// K(H) para norma infinito
norma_infH=0;
norma_infHi=0
for i=1:n
sumH=0;
sumHi=0;
for j=1:n
sumH=sumH+abs(H(i,j));
sumHi=sumHi+abs(Hi(i,j));
end
if (sumH>norma_infH)
norma_infH=sumH;
end
if (sumHi>norma_infHi)
norma_infHi=sumHi;
end
end
Kinf=norma_infH*norma_infHi;
// Pr�ctica 4 - Ejercicio 12
function X=gc(A,b)
N=size(A);
n=N(1,1);
X=zeros(n,1);
p=zeros(n,1);
r=b;
rant=r;
k=0;
while r
k=k+1;
if k==1
p=r;
else
B=r'*r/(rant'*rant);
p=r+B*p;
end
t=r'*r/(p'*A*p);
X=X+t*p;
rant=r;
r=r-t*A*p;
end
// A = [x1 y1;x2 y2; ... ;xn yn]
// x = valor a inter o extra polar
function y=Lagrange(A,x)
N=size(A);
n=N(1,1);
ma=A(1,1);
mi=A(1,1);
for k=2:n
if (A(k,1)>ma)
ma=A(k,1);
end
if (A(k,1)<mi)
mi=A(k,1);
end
end
if (x>ma | x<mi)
disp("Esta es una extrapolaci�n")
end
y=0;
for k=1:n
produ=1;
for i=1:n
if(i~=k)
produ=produ*(x-A(i,1))/(A(k,1)-A(i,1));
end
end
y=y+produ*A(k,2);
end
