MATLAB PROGRAM FOR THE CALCULATION OF EFFECTIVE RIGIDITIES

OF A COUPLED BENDING TORSION BEAM USING CLASSICAL LAMINATE THEORY

MATLAB CODE STARTS

% This is the user friendly way of inputing the beam characteristics
prompt={'1. Ef (Modulus of Fibre)','2. Em Modulus of Matrix ','3. Gf Shear Modulus of Fibre','4. Gm Shear Modulus of Matrix','5. vf Poisson of fibre','6 vm Poisson of Matrix','7. V Volume fraction','8. Df Density of fibre','9. Dm Density of Matrix'}
def={'275.6e9','2.76e9','114.8e9','1.036e9','0.2','0.33','0.5','1900','1600'};
TITLE='Define Composite Properties'
line=1
ANSWER=inputdlg(prompt,TITLE,line,def)
convertc=char(ANSWER);
prec=str2num(convertc);
double(prec);
Ef=prec(1);                %%%%Elastic Modulus of Fiber
Em=prec(2);              %%%%Elastic Modulus of Matrix
Gf=prec(3);                %%%%Shear Modulus of Fiber
Gm=prec(4);              %%%%Elastic Modulus of Matrix
vf=prec(5);                 %%%%Poison's Ratio of Fiber
vm=prec(6);               %%%%Poison's Ratio of Matrix
V=prec(7);                 %%%%Volume Fiber Fraction
Df=prec(8);                 %%%%Density of Fiber
Dm=prec(9)                 %%%%Density Ratio of Matrix
 

%%%%%%%%%%Calculation of Macromechanical Properties%%%%%%%%
 

    E11=Ef*V+Em*(1-V);
E22=Em*(Ef+Em+(Ef-Em)*V)/(Ef+Em-(Ef-Em)*V);
v12=vf*V+vm*(1-V);
v23=vf*V+vm*(1-V)*(1+vm-v12*Em/E11)/(1-vm^2+vm*v12*Em/E11)
G12=Gm*(Gf+Gm+(Gf-Gm)*V)/(Gf+Gm-(Gf-Gm)*V);
G23=E22/(2*(1+v23));
Den=Df*V+Dm*(1-V);

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

NL=1                 %%%%Number of Layers
t=0.005                 %%%%Thickness of each Layer
b=0.1                 %%%%Base of the beam
L=.5                 %%%%Length of the Beam
E1=E11
E2=E22
v=v12
G=G12
pause

    TT=30
TTT=deg2rad(TT);
TH(1)=TTT;                 %%%%This is an array of angles default set for one

                                       %%%%Generally TT=[angle1 angle2 angle3] from top

                                                    to bottom of Laminate

j=0;
for j=1:NL+1
if j==1
range(j)=NL*t/2;
end
range(j)=NL*t/2-t*(j-1);
end

Q11=E1/(1-E2/E1*v^2);
Q22=E2/E1*Q11;
Q12=v*Q22;
Q66=G;
D=0;
for n=1:NL
if TH(n)==0
Qp=[Q11 Q12 0
Q12 Q22 0
0 0 Q66];

else

Q11p(n)=Q11*(cos(TH(n)))^4+Q22*(sin(TH(n)))^4+2*(Q12+2*Q66)*(sin(TH(n)))^2*(cos(TH(n)))^2;
Q12p(n)=(Q11+Q22-4*Q66)*(sin(TH(n)))^2*(cos(TH(n)))^2+Q12*((cos(TH(n)))^4+(sin(TH(n)))^4);
Q16p(n)=(Q11-Q12-2*Q66)*(sin(TH(n)))*(cos(TH(n)))^3+(Q12-Q22+2*Q66)*(sin(TH(n)))^3*(cos(TH(n)));
Q22p(n)=Q11*(sin(TH(n)))^4+Q22*(cos(TH(n)))^4+2*(Q12+2*Q66)*(sin(TH(n)))^2*(cos(TH(n)))^2;
Q26p(n)=(Q11-Q12-2*Q66)*(sin(TH(n)))^3*(cos(TH(n)))+(Q12-Q22+2*Q66)*(sin(TH(n)))*(cos(TH(n)))^3;
Q66p(n)=(Q11+Q22-2*(Q12+Q66))*(sin(TH(n)))^2*(cos(TH(n)))^2+Q66*((cos(TH(n)))^4+(sin(TH(n)))^4);



Qp=[Q11p(n) Q12p(n) Q16p(n)
Q12p(n) Q22p(n) Q26p(n)
Q16p(n) Q26p(n) Q66p(n)];

end
h1=range(n);
h2=range(n+1);

Ddummy=1/3*Qp*(h1^3-h2^3);

D=D+Ddummy;
end


EI=b*(D(2,2)-D(1,2)^2/D(1,1))                 %%%%Effective Bending Rigidity
GJ=4*b*(D(3,3)-(D(1,3))^2/D(1,1))        %%%%Effective Torsional Rigidity
K=2*b*(D(2,3)-D(1,2)*D(1,3)/D(1,1))     %%%%Effective Coupled Rigidity