% -d^2f(x)/dx^2 = 3+2*x^2;
clc;clear all;close all;syms k(x) ;k(x) = 3 +2*x^2; % known fn
nb = 4; % selection of heighest no of basis fn
%% Analytical solution
syms g(x) eqn = -diff(g(x),2) == 3 +2*x^2; % diff eqn
cond = [g(0)==0, g(1)==0]; % boundary cond
a(x) = dsolve(eqn,cond); % diff soln
%% calculating z_matrix and v_vector % Numerical solution %
for N = 2:nb z = zeros(2,2); v = zeros(1,2); for n = 1:N b(n) = (x-(x^(n+1))); % basis fn
w(n) = b(n); % galerkin's MOM
l(n) = -diff(b(n),x,2); % operator L
kr(n) = w(n).*k; % inner product
v(n) = int(kr(n),[0 1]); % v_vector using inner product b/w x->0 to 1
end for i = 1:N for j = 1:N pr(i,j) = w(i).*l(j); z(i,j) = int(pr(i,j),[0 1]); % z_matrix
end end%% LU_decomposition using gauss_elimination
zv = [z v'];L = eye(N,N); for i = 2:1:N alpha(i,1) = zv(i,1)/zv(1,1); L(i,1) = alpha(i,1); % L_matrix elements
zv(i,:) = zv(i,:)- (alpha(i,1)*zv(1,:)); % making zeros in col 1 using pivot 1
end for r = 3:1:N alpha(r,2) = zv(r,2)/zv(2,2); L(r,2) = alpha(r,2); zv(r,:) = zv(r,:)-(alpha(r,2)*zv(2,:)); % making zero in col 2 using pivot 2
endU = zv(1:N,1:N);%% using dummy variable_y finding_I
zb = [L v']; % Augmented matrix
y = zeros(N,1);y(1) = zb(1,end)/zb(1,1); for p = 2:1:N y(p) = (zb(p,end)-(zb(p,1:p-1)*y(1:p-1)))/zb(p,p); % forward substitution
endI = zeros(N,1);I(N) = y(N)/U(N,N); % I_vector calculation
for q = N-1:-1:1 I(q) = (y(q)-(U(q,2:N)*I(2:N)))/U(q,q); end f(N) = b*I ; %finding f(x)
endfor h = 2:nb figure(h) d(x) = f(h); x = 0:0.01:1; plot(x,d(x),'b-','linewidth',2); hold on grid on plot(x,a(x),'r-','linewidth',1); xlabel('x') ylabel('f(x)') legend('Numerical method','Analytical method')end
I have ' f ' as a vector of different functions and I want to plot is simultaneously in different figures but I'm not able to do it, error is showing ''Invalid indexing or function definition. Indexing must follow MATLAB indexing. Function arguments must be symbolic variablesand function body must be sym expression''. But when I ommitting the loop in plot and by separate plotting I'm getting the plot for each h = 2/3/4. So please help me how to fix this problem.
Thank you for your help!
Best Answer