Hi ,i have problem whit the function quad
10 visualizzazioni (ultimi 30 giorni)
Mostra commenti meno recenti
% -u''=(3/2)*t^(1/2)+(2/3)*t^(3/2)+(2/3)-(2/3)*t, int(u'v')=int(fv)
% u(0)=u(1)=0
clear;figure;
n=10;% Nombre des noeuds
a=0;b=1;
x=linspace(a,b,n);h=(b-a)/(n-1);% Noeuds
syms t ;
f=-(3/2)*t^(1/2)+(2/3)*t^(3/2)+(2/3)-(2/3)*t;
for i=1:n-2
H1=@(t)(x(i+2)-t)/h;
H2=@(t)(t-x(i))/h;
k(i,i)=2/h;
f1=@(t)f.*(x(i+2)-t)/h;
f2=@(t)f.*(t-x(i))/h;
F(i)=quad(f2,x(i),x(i+1))+quad(f1,x(i+1),x(i+2));
end
% Calcul des éléments non diagonaux de la matrice de régidité
for i=2:n-2
k(i-1,i)=-1/h;
k(i,i-1)=k(i-1,i);
end
% Matrice globale de régidité
A=zeros(n-2);
for i=1:n-2
A(i,i)=k(i,i);
end
for i=2:n-2
A(i-1,i)=k(i-1,i);
A(i,i-1)=A(i-1,i);
end
%Vecteur de second memebre du système global
B(:,1)=F(1:n-2);
Uh=zeros(n,1);
Uh(2:n-1)=A\B;% Solution approchée
t=a:.001:b;
u=@(t)2*t.^(3/2)/3-2*t/3 ;% Solution exacte
y=feval(u,t);
plot(t,y),hold on
%m=length(x);
u=@(x)2*x.^(3/2)/3-2*x/3;
U(1:n,1)=feval(u,x(1:n));
Uh(1)=0;Uh(n)=0;
for i=1:n-1
t=x(i):.001:x(i+1);
H2=(t-x(i))/h;
H1=(x(i+1)-t)/h;
uh=Uh(i)*H1+Uh(i+1)*H2;%Expression analytique de la solution approchée
plot(t,uh,'-r')
end
%Affichage des résultats
er=abs(Uh-U);
T(:,1)=x';T(:,2)=U;T(:,3)=Uh;T(:,4)=er;
disp([' x ',' u ',' uh ',' Erreur'])
disp(T)
2 Commenti
Geoff Hayes
il 20 Apr 2020
Elkalai - please copy and paste the full error message to this question (rather than creating tags of the error message).
Risposte (1)
Vedere anche
Categorie
Scopri di più su Object Containers in Help Center e File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
