Hi ,i have problem whit the function quad
조회 수: 13 (최근 30일)
이전 댓글 표시
% -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
Geoff Hayes
2020년 4월 20일
Elkalai - please copy and paste the full error message to this question (rather than creating tags of the error message).
답변 (1개)
참고 항목
카테고리
Help Center 및 File Exchange에서 Object Containers에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!
