Energy method with base is based by piecewise linear functions
4 views (last 30 days)
Show older comments
The differential problem is:
-(x^2 y')'+6y=8x+1/x
y(-2)-2y'(-2)=-4
y(-1/2)=-5/4
The output isn't little different in te point zero but i don't see the motivation. Can you help me?
function [x,u] = ex_PWLener2(m,n)
%metodo energia con base delle funzioni lineari a tratti con condizione
%essenziale in b
B=11/14;
C=-6/7;
h=1.5/m;
xC=-2:h:-0.5; %x di reazione della base
A=zeros(m);
b=zeros(m,1);
A(1,1)=2+integral(@(x)fA(x,0,0,xC),xC(1),xC(2));
A(1,2)=integral(@(x)fA(x,0,1,xC),xC(1),xC(2));
b(1)=integral(@(x)fb(x,0,xC),xC(1),xC(2));
for k=1:m-2
A(k+1,k+1)=integral(@(x)fA(x,k,k,xC),xC(k),xC(k+2));
A(k+1,k+2)=integral(@(x)fA(x,k,k+1,xC),xC(k+1),xC(k+2));
A(k+1,k)=integral(@(x)fA(x,k,k-1,xC),xC(k),xC(k+1));
b(k+1)=integral(@(x)fb(x,k,xC),xC(k),xC(k+2));
end
A(m,m)=integral(@(x)fA(x,m-1,m-1,xC),xC(m-1),xC(m+1));
A(m,m-1)=integral(@(x)fA(x,m-1,m-2,xC),xC(m-1),xC(m));
b(m)=integral(@(x)fb(x,m-1,xC),xC(m-1),xC(m+1));
delta=A\b;
h=1.5/n;
x=-2:h:-0.5;
u=zeros(1,n+1);
for k=1:m-1
u=u+delta(k)*Cfun(x,k,xC);
end
u=u+B.*x+C;% u=u+l(x) dove l(x) va calcolata
end
function y = Cfun(x,k,xC)
m=length (xC)-1;
y=zeros(size(x));
for j=1:length (x)
if k==0
if xC(1)<=x(j) && x(j)<xC(2)
y(j)=(xC(2)-x(j))/(xC(2)-xC(1));
end
elseif k==m
if xC(m)<x(j) && x(j)<=xC(m+1)
y(j)=(x(j)-xC(m))/(xC(m+1)-xC(m));
end
else

if xC(k)<x(j) && x(j)<=xC(k+1)
y(j)=(x(j)-xC(k))/(xC(k+1)-xC(k));
elseif xC(k+1)<=x(j) && x(j)<xC(k+2)
y(j)=(xC(k+2)-x(j))/(xC(k+2)-xC(k+1));
end
end
end
end
function y=Cder(x,k,xC)
%mettere derivata rispetto x(j)
m=length (xC)-1;
y=zeros(size(x));
for j=1:length (x)
if k==0
if xC(1)<=x(j) && x(j)<xC(2)
y(j)=-1/(xC(2)-xC(1));
end
elseif k==m
if xC(m)<x(j) && x(j)<=xC(m+1)
y(j)=1/(xC(m+1)-xC(m));
end
else
if xC(k)<x(j) && x(j)<=xC(k+1)
y(j)=1/(xC(k+1)-xC(k));
elseif xC(k+1)<=x(j) && x(j)<xC(k+2)
y(j)=-1/(xC(k+2)-xC(k+1));
end
end
end
end
function y=fA(x,k,l,xC)
y=(x.^2).*Cder(x,k,xC).*Cder(x,l,xC)+6.*Cfun(x,k,xC).*Cfun(x,l,xC);
end
function y=fb(x,k,xC)
B=11/14;
C=-6/7;
y=Cfun(x,k,xC).*(8.*x+(1./x)-4*B*x-6*C);
end
2 Comments
Answers (0)
See Also
Categories
Find more on Other Formats in Help Center and File Exchange
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!