function tappered_bar(At,Ab,l,fb,E,n,filename)
y0=At;
K=zeros(n);k=linspace(0,0,n);f=linspace(0,0,n)';f(n)=fb; )%At=area at fixed end Ab=area at free end
op=fopen(filename,'wt'); %fb=axial
force at free end and n=no. of element
fprintf(op,'=======================================================================\n');
fprintf('=======================================================================\n');
fprintf(op,'\t\t\tTappered
bar solution by FEM\n'); fprintf('\t\t\t\t\tTappered bar solution by FEM\n');
fprintf(op,'-----------------------------------------------------------------------\n');
fprintf('-----------------------------------------------------------------------\n');
fprintf(op,'Maximum area =
%f\t\t\t',At); fprintf(op,'Minimum area = %f\n',Ab);
fprintf('Maximum area =
%f\t\t\t\t',At); fprintf('Minimum area = %f\n',Ab);
fprintf(op,'Force at
minimum area = %f\t',fb); fprintf(op,'Lenght = %f\n',l);
fprintf('Force at
minimum area = %f\t',fb); fprintf('Lenght = %f\n',l);
fprintf(op,'Modulus of
elasticity = %f\t',E); fprintf(op,'No. of element = %f\n',n);
fprintf('Modulus of
elasticity = %f\t',E); fprintf('No. of element = %f\n',n);
fprintf(op,'-----------------------------------------------------------------------\n');
fprintf('-----------------------------------------------------------------------\n');
for i = 1:n
yi = (Ab-At)*i/n+At;
Ai= (y0 +yi)/2;
k(i)=n*E*Ai/l;
y0=yi;
fprintf(op,'A%d = %f\t',i,Ai);fprintf('A%d = %f\t',i,Ai);
fprintf(op,'k%d = %f\n',i,k(i));fprintf('k%d = %f\n',i,k(i));
end
fprintf(op,'-----------------------------------------------------------------------\n');
fprintf('-----------------------------------------------------------------------\n');
for i=1:n
for j=1:n
if abs(i-j)==1
if i>j
K(i,j)=-k(i);
elseif j>i
K(i,j)=-k(j);
end
elseif i==n && j==n
K(i,j)=k(i);
elseif i~=n && j~=n && i==j
K(i,j)=k(i)+k(i+1);
end
end
end
u = inv(K)*f;
disp('K = ');fprintf(op,'K = \n');
disp(K);
for i=1:n
for j=1:n
fprintf(op,'%f\t',K(i,j));
end
fprintf(op,'\n');
end
fprintf(op,'-----------------------------------------------------------------------\n');
fprintf('-----------------------------------------------------------------------\n');
for i=1:n
fprintf(op,'u%d = %f\n',i,u(i));fprintf('u%d = %f\n',i,u(i));
%u_approx=ui;
end
fprintf(op,'=======================================================================\n');
fprintf('=======================================================================\n');
fprintf(op,'\t\t\tExact
solution of Tappered bar\n');fprintf('\t\t\t\t\tExact solution of Tappered
bar\n');
fprintf(op,'-----------------------------------------------------------------------\n');
fprintf('-----------------------------------------------------------------------\n');
syms x;
Area_at_x = (Ab-At)*x/l+At;
del= (fb/E)*int(1/Area_at_x,x,0,l);
fprintf(op,'\t\t\t\tdel = %f\n',del);fprintf('\t\t\t\t\t\t\tdel
= %f\n',del);
fprintf(op,'=======================================================================\n');
fprintf('=======================================================================\n');
plot(n,u(n),'r*',n,del,'gx');
xlabel('no. of Elements for FEM solution');ylabel('Displacement');title('Comapring FEM
and Exact solution');
legend('FEM solution','Exact solution');
fclose(op);
end
OUTPUT
when n= 3
=======================================================================
Tappered
bar solution by FEM
--------------------------------------------------------------------------------------------------------------------
Maximum area = 0.031416 Minimum area = 0.007854
Force at minimum area = 10.000000 Lenght = 2.000000
Modulus of elasticity = 1.000000 No. of element = 3.000000
-------------------------------------------------------------------------------------------------------------------
A1 = 0.027489 k1 = 0.041233
A2 = 0.019635 k2 = 0.029452
A3 = 0.011781 k3 = 0.017671
-------------------------------------------------------------------------------------------------------------------
K =
0.070686 -0.029452 0.000000
-0.029452 0.047124 -0.017671
0.000000 -0.017671 0.017671
-------------------------------------------------------------------------------------------------------------------
u1 = 242.521818
u2 = 582.052363
u3 = 1147.936605
=======================================================================
Exact
solution of Tappered bar
-------------------------------------------------------------------------------------------------------------------
del
= 1176.723201
=======================================================================
OUTPUT
when n=4
=======================================================================
Tappered
bar solution by FEM
--------------------------------------------------------------------------------------------------------------------
Maximum area = 0.031416 Minimum area = 0.007854
Force at minimum area = 10.000000 Lenght = 2.000000
Modulus of elasticity = 1.000000 No. of element = 4.000000
--------------------------------------------------------------------------------------------------------------------
A1 = 0.028471 k1 = 0.056941
A2 = 0.022580 k2 = 0.045160
A3 = 0.016690 k3 = 0.033379
A4 = 0.010799 k4 = 0.021598
--------------------------------------------------------------------------------------------------------------------
K =
0.102102 -0.045160 0.000000 0.000000
-0.045160 0.078540 -0.033379 0.000000
0.000000 -0.033379 0.054978 -0.021598
0.000000 0.000000 -0.021598 0.021598
--------------------------------------------------------------------------------------------------------------------
u1 = 175.619248
u2 = 397.052212
u3 = 696.637987
u4 = 1159.634185
=======================================================================
Exact
solution of Tappered bar
-------------------------------------------------------------------------------------------------------------------
del
= 1176.723201
=======================================================================
Solutions has been obtained by taking E=1. To get
solution the values of displacement obtained should be divided by E in KN/m2.
No comments:
Post a Comment