Hi, my loop function is not working and giving me index error in my member stiffens matrix. any idea what is the issue? thank you

23 views (last 30 days)
%define coordinates of each nodes%
Coord_nodes = [0 0; 10 0; 10 10; 0 10];
% define connection of each element%
Connect_elements = [1 2; 1 3; 1 4; 4 3; 4 2; 2 3];
%define degree of Freedom of each nodes
Deg_free=[1 2; 6 5; 7 8; 3 4];
% define number of reactions at supports
Nr= 3;
%define number of members
Nm= 6;
%define number of nodes
Nn= 4;
%define EA
EA= 1;
%define Lenght of each member
L=zeros(Nm,1);
L= [10; 14.14 ;10 ;10 ; 14.14 ;10]
L = 6×1
10.0000 14.1400 10.0000 10.0000 14.1400 10.0000
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
% creation of ID Array
ID = zeros(Nm,4);
for k = 1:Nm
i = Connect_elements(k,1);
j = Connect_elements(k,2);
ID(k,1:2)=Deg_free(i,1:2);
ID(k,3:4)=Deg_free(j,1:2);
end
% assemble Global Stiffness matrix
Kg=zeros(8,8);
for k = 1:Nm
i = Connect_elements(k,1);
j = Connect_elements(k,2);
Lamda_x= (Coord_nodes(j,1)- Coord_nodes(i,1))/L(k)
Lamda_y= (Coord_nodes(j,2)- Coord_nodes(i,2))/L(k)
end
Lamda_x = 1
Lamda_y = 0
Lamda_x = 0.7072
Lamda_y = 0.7072
Lamda_x = 0
Lamda_y = 1
Lamda_x = 1
Lamda_y = 0
Lamda_x = 0.7072
Lamda_y = -0.7072
Lamda_x = 0
Lamda_y = 1
%member Stiffness matrix
for k = 1:Nm
k11(k)= (Lamda_x(k)^2/L(k));
k12(k)= (Lamda_x(k)*Lamda_y(k))/L(k);
k13(k)= -k11(k);
k14(k)= (-Lamda_x(k)*Lamda_y(k))/L(k);
k21(k)= k12(k);
k22(k)= (Lamda_y(k))^2/L(k);
k23(k)= -k12(k);
k24(k)= -k22(k);
k31(k)= -k11(k);
k32(k)= k23(k);
k33(k)= k11(k);
k34(k)= k12(k);
k41(k)= k14(k);
k42(k)= k24(k);
k43(k)= k12(k);
k44(k)= k22(k);
k_elem(:,:,k)= [k11(k) k21(k) k31(k) k41(k);
k12(k) k22(k) k32(k) k42(k);
k13(k) k23(k) k33(k) k43(k);
k14(k) k24(k) k34(k) k44(k)]
end
k_elem = 4×4
0 0 0 0 0 0.1000 0 -0.1000 0 0 0 0 0 -0.1000 0 0.1000
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
Index exceeds the number of array elements. Index must not exceed 1.

Accepted Answer

Stephen23
Stephen23 on 29 Mar 2025
Edited: Stephen23 on 29 Mar 2025
% define coordinates of each nodes
Coord_nodes = [0,0; 10,0; 10,10; 0,10];
% define connection of each element
Connect_elements = [1,2; 1,3; 1,4; 4,3; 4,2; 2,3];
% define degree of Freedom of each nodes
Deg_free = [1,2; 6,5; 7,8; 3,4]; % [x-direction DOF, y-direction DOF]
% define number of reactions at supports
Nr = 3;
% define number of members
Nm = size(Connect_elements,1);
% define number of nodes
Nn = 4;
% define EA
EA = 1;
% creation of ID Array
ID = zeros(Nm,4);
for k = 1:Nm
i = Connect_elements(k,1);
j = Connect_elements(k,2);
ID(k,1:2) = Deg_free(i,1:2);
ID(k,3:4) = Deg_free(j,1:2);
end
% assemble Global Stiffness matrix
k_elem = nan(4,4,Nm);
for k = 1:Nm
i = Connect_elements(k,1);
j = Connect_elements(k,2);
dx = Coord_nodes(j,1) - Coord_nodes(i,1);
dy = Coord_nodes(j,2) - Coord_nodes(i,2);
L = norm([dx,dy]);
Lambda_x = dx/L;
Lambda_y = dy/L;
%
k11 = Lambda_x^2/L;
k12 = (Lambda_x*Lambda_y)/L;
k13 = -k11;
k14 = -k12; % not (-Lambda_x*Lambda_y)/L;
k21 = k12;
k22 = Lambda_y^2/L;
k23 = -k12;
k24 = -k22;
k31 = -k11;
k32 = k23;
k33 = k11;
k34 = k12;
k41 = k14;
k42 = k24;
k43 = k12;
k44 = k22;
k_elem(:,:,k) = EA * [...
k11, k21, k31, k41;...
k12, k22, k32, k42;...
k13, k23, k33, k43;...
k14, k24, k34, k44];
end
%
% Get total number of DOFs
Nd = max(Deg_free(:));
% Initialize global stiffness matrix
GSM = zeros(Nd,Nd);
% Assemble global stiffness matrix
for k = 1:Nm
for i = 1:4
for j = 1:4
GSM(ID(k,i), ID(k,j)) = GSM(ID(k,i), ID(k,j)) + k_elem(i,j,k);
end
end
end
Giving
display(GSM)
GSM = 8×8
0.1354 0.0354 0 0 0 -0.1000 -0.0354 -0.0354 0.0354 0.1354 0 -0.1000 0 0 -0.0354 -0.0354 0 0 0.1354 -0.0354 0.0354 -0.0354 -0.1000 0 0 -0.1000 -0.0354 0.1354 -0.0354 0.0354 0 0 0 0 0.0354 -0.0354 0.1354 -0.0354 0 -0.1000 -0.1000 0 -0.0354 0.0354 -0.0354 0.1354 0 0 -0.0354 -0.0354 -0.1000 0 0 0 0.1354 0.0354 -0.0354 -0.0354 0 0 -0.1000 0 0.0354 0.1354
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>

More Answers (0)

Categories

Find more on Numerical Integration and Differential Equations in Help Center and File Exchange

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!