문제는
-xu"-u'+u*x^2=cosx
u(0)=u(1)=0
0<x<1
일 때 u 의 해에 대해 plotting하는 건데요...
local 행렬을 구해서 u에 대한 global 행렬로 assembling해서 푸는 형식입니다.
저는 이 문제에서 전체 term을 x로 나눠서
-u"-(1/x)*u'+x*u=cosx/x
이렇게 해서 boundary condition을 더 간단하게 만들어서 푸는 중인데요
도저히 그 중의 한 term인 (cosx/x)((-x+b)/(b-a))부분과 (cosx/x)((x-a)/(b-a))부분의 적분결과를 설정할 수 가 없어서요..
cosx/x를 적분하면 element로 표현하기 어려워서 function으로 해보았는데,,,아직 초보라 제가 function을 제대로 설정하였는지도 모르겠어요..
여기서 a와 b는 0부터 1까지를 n등분한 것 중의 한 개의 element part이구요
일단 첫번째는
function lF1=int(N1,N2);
lF1=(cos(x)/x)*((-x+N2)/(N2-N1)); 이걸 엠파일로 만들어서
lF(1,1)= int(sym('int'),N1,N2); 이렇게 실행을 시켰는데.......
실행은 되요...맞는지는 모르겠지만..ㅠㅠ 이런식으로 두번째 적분도 했구요...
그다음엔..정확한 정답하고 제가 코딩한 거하고 비교해야하는데
여기서도 좀 잘못 되지 않았을까 싶어요.
x_1=eCoord;
% compare the solution with that obtained using Matlab built-in functions.
solinit = bvpinit(linspace(0.00000001,1,100), [0 1]);
sol = bvp4c(@problem_,@bc,solinit);
alpha=linspace(0.00000001,1,100);
beta=deval(sol,alpha);
title=sprintf('comparison of FEM solution=%3d',Nelem);
figure('Name',title);
plot(x_1,sU, alpha, beta(1,:),'-ro');
xlabel('x')
ylabel('y')
여기서
problem_.m은
function dydx = problem_(x,y)
dydx = [y(2) ; x*y(1)-(1/x)*y(2)-cos(x)/x];
bc.m은
function res = bc(ya,yb)
res = [ya(1)-0 yb(1)-0];로 했구요...
뭐가 문제일까요....
ㅠㅠ 그래프는 점점 산으로 가네요..도와주세요
