Week 5 Challenge- Calculating Transformation Matrix and Fixed end Force/Moments
Objective:-
- Calculating Transformation Matrix and Fixed end Force/Moments
Solution procedure:-
1) Calculate the transformation matrix for a frame shown below Develop a matrix that will transform forces from their local axis x’-y’-z’ to x-y-z axis as shown in the figure below.


`q1' = q1cos(thetax) -q2cos(thetay)`
`q2' = q1cos(thetay) +q2cos(thetax)`

`lamdax = cos(thetax)= (x(i)-xj)/L`
`lamday = cos(thetay)= (y(i)-yj)/L`
`[T] = [q1' q2'q3'] = [lamdax lamday0;-lamdaylamdax0;001]*[q1q2q3]`
________________________________________________________________________________________________________________
2) Develop a transformation matrix for the portal frame shown below. Using the transformation matrix assemble global stiffness matrix for the structure given below

Code for the Transformation Matrix: -
function[k_element] = Frame_Stiffness_Matrix1(node_xy,A,I,length_element)
%Step 1: To construct element stiffness matrix
E =2*10^11;
frame_stiffness=zeros(6,6);
k11 =A*E/length_element;
k22 =12*E*I/length_element^3;
k23= 6*E*I/length_element^2;
k33 = 4*E*I/length_element;
k36 = 2*E*I/length_element;
frame_stiffness(:,:) = [k11 0 0 -k11 0 0;0 k22 k23 0 -k22 k23;0 k23 k33 0 -k23 k36;-k11 0 0 k11 0 0;0 -k22 -k23 0 k22 -k23;0 k23 k36 0 -k23 -k33]
%To compute Transformation matrix
lamda_x=(node_xy(2,1)-node_xy(1,1))/length_element;
lamda_y=(node_xy(2,2)-node_xy(1,2))/length_element;
T_rot=[lamda_x -lamda_y 0;lamda_y lamda_x 0;0 0 1];
T=[T_rot zeros(3);zeros(3) T_rot]
k_element=T'*frame_stiffness*T
end
Explanation of the above code:-
function[k_element] = Frame_Stiffness_Matrix1(node_xy,A,I,length_element) - function to perform task
%Step 1: To construct element stiffness matrix - comments
E =2*10^11; - assing E yogs modulus vaue
frame_stiffness=zeros(6,6); - creating zeros matrix
k11 =A*E/length_element;- calculating k11
k22 =12*E*I/length_element^3;
k23= 6*E*I/length_element^2;
k33 = 4*E*I/length_element;
k36 = 2*E*I/length_element;
frame_stiffness(:,:) = [k11 0 0 -k11 0 0;0 k22 k23 0 -k22 k23;0 k23 k33 0 -k23 k36;-k11 0 0 k11 0 0;0 -k22 -k23 0 k22 -k23;0 k23 k36 0 -k23 -k33] - assing the matrix y rows and columns
%To compute Transformation matrix
lamda_x=(node_xy(2,1)-node_xy(1,1))/length_element; - calculating lamda for x
lamda_y=(node_xy(2,2)-node_xy(1,2))/length_element; - clculating lamda for y
T_rot=[lamda_x -lamda_y 0;lamda_y lamda_x 0;0 0 1]; - calcuating in rows and columns
T=[T_rot zeros(3);zeros(3) T_rot] - transformation matrix calculation
k_element=T'*frame_stiffness*T
end
Output:

Global stiffness matrix for above code:-

_______________________________________________________________________________________________________________
3) Using a direct stiffness method calculates displacements at A, B. Determine the global stiffness matrix fixed end forces and moments and calculate the displacements at A and B

Calculating the Fixed end moments by using
`(-WL^2)/12` and `(WL^2)/12` - for udl loads comes out to be 6.25kn/m and -6.25kn/m`
`(-Pab^2)/l^2` and `(Pa^2b)/l^2`- for point loads be 20kn/m and -20kn/m
then we calculate fixed end forces
- Then we calculate the v1,theata1 and v3
v1 : `(12EI)/l^3` `(-12EI)/l^3` `(-6EI)/l^2`
`theta1`:`(-12EI)/l^3` `(12EI)/l^3` `(6EI)/l^2`
v3: `(-6EI)/l^2` `(6EI)/l^2` `(4EI)/l`
v3 : `(12EI)/l^3` `(-12EI)/l^3` `(-6EI)/l^2`
v2: `(-12EI)/l^3` `(12EI)/l^3` `(6EI)/l^2`
`theta2`: `(-6EI)/l^2` `(6EI)/l^2` `(4EI)/l`
On solving the above we get:
v1 `theta1` v3 v2 `theta2`
v1 0.096 -0.24 -0.096 0 0
`theta1` -0.24 0.8 0.24 0 0
v3 -0.096 0.24 0.096,0.02 -0.02 -0.09
v2 0 0 -0.02 0.02 0.09
`theta2` 0 0 -0.09 -0.09 0.5
____________________________________________________________
K1+K2 is Kg
On calculating theta1 and theta2 is `-1.95*10^-4` and `10*10^-4`
4) Calculate fixed end forces and moments for loading

a)Fixed-End Forces
P/2
Moment
`-(PL)/8`
`(PL)/8`
b)Fixed-End Forces
p/2
Moment
`(-Pab^2)/L^2`
`(Pa^2b)/L^2`
c)Fixed-End Forces
`(WL)/2`
Moment
`(-WL^2)/12`
`(WL^2)/12`
d)
Moment
`(-Pab^2)/L^2`
`(Pa^2b)/L^2`
e)Fixed end force
0.5*l*p
Moment
`(-WL^2)/30`
`(WL^2)/30`
_______________________________________________________________________________________________________________
5) Develop a Matlab module to compute fixed end forces and moments when the load is applied between nodal points. Compute displacements and forces for the following beam

function [global_stiffness_matrix]=ComputeGlobalStiffnessMatrix()
%Get Element Info
[node_coord,element_no,element_node,length_element,element_geometry,element_dof]=GetElementInfo();
%Get Boundary Conditions Info
boundary_conditions=GetBoundary_Conditions();
%Get Nodal loads Info
node_loads = GetNodalLoads();
%Calculating Number Of Elements And Nodes
no_elements = size(element_no,1);
no_nodes = size(node_coord,1);
%Intializing Matrix with Zeros
global_stiffness_matrix= zeros(3*no_nodes,3*no_nodes);
displacement_matrix= zeros(3*no_nodes,1);
force_matrix = zeros(3*no_nodes,1);
%Boundary Conditions
all_dof= 1:3*no_nodes;
dof_specified =[]
for i=1:no_nodes
n_dof=3*(boundary_conditions(i,1)-1)+boundary_conditions(i,2);
dof_specified=[dof_specified n_dof];
displacement_matrix(n_dof)=boundary_conditions(i,3);
end
free_dof=all_dof;
free_dof(dof_specified)=[]
%Compute Nodal Loads
for j =1:size(node_loads,1)
force_matrix(3*(node_loads(j,1)-1)+node_loads(j,2))=node_loads(j,3);
end
%Global Stiffness Matrix
for k=1:no_elements
ele_nodes=element_node(k,:);
node_xy=node_coord(ele_nodes,:);
A = element_geometry(k,1)*element_geometry(k,2);
I=element_geometry(k,1)*element_geometry(k,2).^3/12;
length=length_element(k);
k_element=Frame_Stiffness_Matrix1(node_xy,A,I,length);
eldofs=element_dof(k,:);
global_stiffness_matrix(eldofs,eldofs)=global_stiffness_matrix(eldofs,eldofs)+k_element;
end
displacement_matrix(free_dof)=global_stiffness_matrix(free_dof,free_dof)\(force_matrix(free_dof)-global_stiffness_matrix(free_dof,dof_specified)*displacement_matrix(dof_specified));
force_matrix(dof_specified)=global_stiffness_matrix(dof_specified,:)*displacement_matrix;
disp('Displacement Vector');
disp(displacement_matrix);
disp('Force Vector');
disp(force_matrix);
end
Explanation of the above code:-
function [global_stiffness_matrix]=ComputeGlobalStiffnessMatrix()
%Get Element Info
[node_coord,element_no,element_node,length_element,element_geometry,element_dof]=GetElementInfo(); - taking data from the functions
%Get Boundary Conditions Info
boundary_conditions=GetBoundary_Conditions();
%Get Nodal loads Info
node_loads = GetNodalLoads();
%Calculating Number Of Elements And Nodes
no_elements = size(element_no,1);
no_nodes = size(node_coord,1);
%Intializing Matrix with Zeros
global_stiffness_matrix= zeros(3*no_nodes,3*no_nodes);
displacement_matrix= zeros(3*no_nodes,1);
force_matrix = zeros(3*no_nodes,1);
%Boundary Conditions
all_dof= 1:3*no_nodes;
dof_specified =[]
for i=1:no_nodes
n_dof=3*(boundary_conditions(i,1)-1)+boundary_conditions(i,2);
dof_specified=[dof_specified n_dof];
displacement_matrix(n_dof)=boundary_conditions(i,3);
end
free_dof=all_dof;
free_dof(dof_specified)=[]
%Compute Nodal Loads
for j =1:size(node_loads,1)
force_matrix(3*(node_loads(j,1)-1)+node_loads(j,2))=node_loads(j,3);
end
%Global Stiffness Matrix
for k=1:no_elements - for loop to iteration
ele_nodes=element_node(k,:); - assigning the elements nodes
node_xy=node_coord(ele_nodes,:); - assigning the coordinates for node x and node y
A = element_geometry(k,1)*element_geometry(k,2);
I=element_geometry(k,1)*element_geometry(k,2).^3/12;
length=length_element(k);
k_element=Frame_Stiffness_Matrix1(node_xy,A,I,length);
eldofs=element_dof(k,:);
global_stiffness_matrix(eldofs,eldofs)=global_stiffness_matrix(eldofs,eldofs)+k_element;
end
displacement_matrix(free_dof)=global_stiffness_matrix(free_dof,free_dof)\(force_matrix(free_dof)-global_stiffness_matrix(free_dof,dof_specified)*displacement_matrix(dof_specified));
force_matrix(dof_specified)=global_stiffness_matrix(dof_specified,:)*displacement_matrix;
disp('Displacement Vector');
disp(displacement_matrix);- display the dispacement amtrix
disp('Force Vector');- display the force vector name
disp(force_matrix);- display the force matrix
end
The output of the above code:-




________________________________________________________________________________________________________________
6) Calculate deflections for the beam under the loading specified

a) Point load in Cantilever Beam
This deflection is calculated as:
`deltamax =(FL^3)/(3EI)`
Where:
F = force per unit length (N/m)
L = unsupported length (m)
E = modulus of elasticity (N/m2)
I = planar moment of inertia (m4)
b)
Udl in Cantilever Beam
This deflection is calculated as:
`deltamax =(qL^4)/(8EI)`
Where:
q = force per unit length (N/m)
L = unsupported length (m)
E = modulus of elasticity (N/m2)
I = planar moment of inertia (m4)
c) Point load at center in the simply supported beam
`deltamax =(FL^3)/(48EI)`
Where:
F = force per unit length (N/m)
L = unsupported length (m)
E = modulus of elasticity (N/m2)
I = planar moment of inertia (m4)
d) Udl in the simply supported beam
`deltamax =(5qL^4)/(384EI)`
Where:
q = force per unit length (N/m)
L = unsupported length (m)
E = modulus of elasticity (N/m2)
I = planar moment of inertia (m4)