Code covered by the BSD License  

Highlights from
ABD matrix

ABD matrix

by

 

to calculate ABD matrix in composite materials

ABD_mat_revised.m
%%Program for Calculation of A, B and D matrices for a Uniform laminate.
clc;
close all;
clear all;
digits(32);
format short e;
disp('Please Select the Material')
disp('1  NCT 301')
disp('2  New Material')
comm=input('Enter the number correspond to your choice>> ');
n=input('Enter Number of Layers (Total Number)>> ');
disp('----------------------------------------------------------')
disp('please enter layer stacking, one example is provided below');
disp('If the Layer Stacking(LS) is [+-30 0]s enter " [30 -30 0 0 -30 30] " Do Not Forget to Use Brackets [] !!!')
disp('---------------------------------------------------------------------------------------------------------')
l=input('Enter LS >> ');
clc

%Mechanical Properties of NCT 301
if comm==1
    E1=155E+9;
    E2=12.1E+9;
    E3=E2; 
    n21=0.248; 
    n31=n21;
    n23=0.458; 
    G12=4.4E+09;
    G13=G12; 
    G23=3.20E+09;
    
else
    disp('Please enter the material properties')
    disp('Enter the values for extensional and shear modulus in GPa')
    mm=1e+9;
    E1=input('E1 = ');
    E1=E1*mm;
    E2= input('E2 = ');
    E2=E2*mm;
    E3= input('E3 = ');
    E3=E3*mm;
    n21= input('n21 = ');
    n31= input('n31 = ');
    n23= input('n23 = ');
    G12= input('G12 = ');
    G12=G12*mm;
    G13= input('G13 = ');
    G13=G13*mm;
    G23= input('G23 = ');
    G23=G23*mm;
end

t   = 125e-6;

all_ang = l;
k = {};

N=length(all_ang);
h=zeros(N+1);
h(N+1)=N*t/2;
lstf=zeros(3*N,3);

for k = 1 : N
    

        %%Compliance Matrix elements
        S = [1/E1     -n21/E1    -n31/E1    0         0       0;
            -n21/E1   1/E2       -n23/E2    0         0       0;
            -n31/E1   -n23/E2     1/E3      0         0       0;
            0           0           0       1/G23     0       0;
            0           0           0        0       1/G13    0;
            0            0          0        0       0       1/G12];

        %%Stiffness Matrix elements
        C = inv(S);
        
        %Reduced Stiffness
        Q = zeros(3,3);
        Q(1,1) = C(1,1) - (C(1,3))^2 / C(3,3);
        Q(2,2) = C(2,2) - (C(2,3))^2 / C(3,3);
        Q(1,2) = C(1,2) - C(1,3)*C(2,3) / C(3,3);
        Q(2,1) = Q(1,2);
        Q(3,3) = C(6,6);

%Transformatioin Matrix        

m = cosd (all_ang(k));
n = sind (all_ang(k));
       
       
        T = zeros(3,3);
        T(1,1)=m^2;
        T(1,2)=n^2;
        T(1,3)=2*m*n;
        T(2,1)=n^2;
        T(2,2)=m^2;
        T(2,3)=-2*m*n;
        T(3,1)=-m*n;
        T(3,2)=m*n;
        T(3,3)=m^2-n^2;
     
%Transformed Stiffness
q=zeros(3,3);
q([1:2],[1:2])=Q([1:2],[1:2]);
q(3,3)=2*Q(3,3);

q=inv(T)*q*T;
    for i = 1:3
        q(i,3)=q(i,3)/2;
    end

lstf([3*k-2:3*k],[1:3])=q([1:3],[1:3]);

h(k)=(k-N/2-1)*t;
end

% To Calculate A B and D Matrices for Uniform Laminate

A=zeros(3,3);
B=zeros(3,3);
D=zeros(3,3);
    
    for i=1:3
        for j=1:3
            q([1:3],[1:3])=lstf([1:3],[1:3]);
            A(i,j) = q(i,j) * (h(2) - h(1));
            B(i,j) = 1/2*(q(i,j) * (h(2)^2 - h(1)^2));
            D(i,j) = 1/3*(q(i,j) * (h(2)^3 - h(1)^3));
            
            for k = 2 : N
            q([1:3],[1:3]) = lstf( [3*k-2:3*k] , [1:3] );
            A(i,j) = q(i,j) * (h(k+1) - h(k)) + A(i,j);
            B(i,j) = 1/2*(q(i,j) * (h(k+1)^2 - h(k)^2)) + B(i,j);
            D(i,j) = 1/3*(q(i,j) * (h(k+1)^3 - h(k)^3)) + D(i,j);    
        end
    end
    end

LamnStf=zeros(6,6);
a=zeros(3,3);
b=zeros(3,3);
d=zeros(3,3);
LamnStf([1:3],[1:3])=A([1:3],[1:3]);
LamnStf([4:6],[4:6])=D([1:3],[1:3]);
LamnStf([1:3],[4:6])=B([1:3],[1:3]);
LamnStf([4:6],[1:3])=B([1:3],[1:3]);
LamnCmp=inv(LamnStf);
a([1:3],[1:3])=LamnCmp([1:3],[1:3]);
b([1:3],[1:3])=LamnCmp([1:3],[4:6]);
d([1:3],[1:3])=LamnCmp([4:6],[4:6]);

if rem(length(all_ang),2)==0
A= A.*[1 1 0 ; 1 1 0 ; 0 0 1];
B= B.*[0 0 0 ; 0 0 0 ; 0 0 0];
D= D.*[1 1 1 ; 1 1 1 ; 1 1 1];

a= a.*[1 1 0 ; 1 1 0 ; 0 0 1];
b= b.*[0 0 0 ; 0 0 0 ; 0 0 0];
d= d.*[1 1 1 ; 1 1 1 ; 1 1 1];

end;

A
B
D

a
b
d

% To Calculate Force and Moment Resultants or Reference Strains or
% reference curvatures.

Nx= A(

Contact us