simple frame analysis by stiffness matrix approach 1

[civil2013]

hi, i am analysing a simple frame structure - composed of a horizontal beam and a vertical column (the frame looks like '7' with 6 dofs restrained at the bottom of column)... i have prepared the stiffness matrix for the strcture, now need to prepare the force matrix..... i need to analyse the structure only for the dead weight and wish to know how to consider the weight of column in the force matrix. understand that the beam weight can be applied as momement and vertical shear at the beam ends, how to consider the column weight - is is just applied as equal point loads at the top and bottom ? thanks

 

[JAE]

What sort of problem is this? Why would you not use other loads besides dead load? Is this some kind of problem given in a course of study?

Check out Eng-Tips Forum's Policies here:
faq731-376

 

[civil2013]

hi, its not part of course...i am trying to develop a fortran code for analysing frame structures, just worried how to handle the self weight in the analysis.. all the other external loads are resolved as joint loads.. doubt is regarding how to take care of the column weight - if the column weight is 'W', is it ok to consider W/2 on the top node and W/2 at bottom node ? in my post, i wanted to simplify the problem, so mentioned that i am not using external loads

 

[awa5114]

You could use equivalent nodal loads to convert distributed loads to joint loads, I believe this lends itself well to programming.

I have a question for you though, as I am also trying to develop a similar program using another language (not Fortran). Are you trying to make your program be able to output complete shear, moment and deflection diagrams for all your elements? Or just nodal forces and displacements?

 

[civil2013]

your comment noted.


regarding my pgrm,i have not decided what all to extract as output. if your worry is about plotting the outputs, u cud use libraries from softwares like plplot/gnplot etc... they sync with fortran/python etc

 

[BAretired]

hi, its not part of course...i am trying to develop a fortran code for analysing frame structures, just worried how to handle the self weight in the analysis.. all the other external loads are resolved as joint loads.. doubt is regarding how to take care of the column weight - if the column weight is 'W', is it ok to consider W/2 on the top node and W/2 at bottom node ? in my post, i wanted to simplify the problem, so mentioned that i am not using external loads

Place W at the bottom node and find the maximum axial load in the column. To find the axial load elsewhere on the column, simply subtract the weight of column between the point in question and the bottom node.



BA

 

[medeek]

When I wrote my matrix analysis program for the truss designer I applied the distributed loads as you indicated (equivalent shear and moments loads at the nodes):

Just as the other distributed loads are applied as a equivalent shear and moments I would think one would do the same for a vertical member. If you then look at your final output (axial, shear, moment diagrams) you should have an axial load that increases linearly from top to bottom, due to the distributed dead weight along its length.

To check you work create an equivalent model in STAAD or RISA and compare with your own program, that is the only way I was able to work the bugs out of my own programming and to finally arrive at the correct output.

I'm sure you already have it but a copy of Structural Analysis by Hibbeler was very helpful.

A confused student is a good student.
Nathaniel P. Wilkerson, PE

http://www.medeek.com

 

[BAretired]

If W/2 is placed at the top and bottom nodes of a vertical column, the axial force at the upper node will be too high by W/2.

BA

 

[KootK]

As others have intimated, the usual way to go about this programmatically is to convert member self weight into statically equivalent fixed end forces that pay homage to the orientation of the member relative to the direction of gravitational acceleration. So, for vertical members in a static run, that does become w/2 top and bottom. If accurate, mid-height, self weight axial demand is important, you pretty much have to divide the columns up into multiple members.

I like to debate structural engineering theory -- a lot. If I challenge you on something, know that I'm doing so because I respect your opinion enough to either change it or adopt it.

 

[KootK]

Are you trying to make your program be able to output complete shear, moment and deflection diagrams for all your elements? Or just nodal forces and displacements?

I've done this and found the post-processing to be quite a bit more difficult than the matrix stuff. Let me know if there's anything that I can do to help with your project.

I like to debate structural engineering theory -- a lot. If I challenge you on something, know that I'm doing so because I respect your opinion enough to either change it or adopt it.

 

[medeek]

Its been awhile (Sept. 15, 2013) since I wrote this so looking through my own perl code (matrix analysis module) is a bit like trying to read latin, its ugly. The thing to remember is that the final values for axial, shear and moments given by your output will be a result of the matrix analysis and the superposition of the loads applied at the nodes for the analysis, see page 622 Structural Analysis Hibbeler 7th Edition. So it should all work out in the end as expected.

###############################################
#
# Sub Routines for Truss Designer (FINK TRUSS MATRIX)
#
###############################################




########################
#
# Rigid Peak Joint
#
########################

sub Matrix {



	##############################
	#
	# Overhang Loading and Graphs
	#
	##############################

	if($Overhang > 0)
	{
		$Moment_overhang = ($Tctl_pli*$Overhang*$Overhang)/2;
		$Moment_overhang_out = sprintf("%.1f", $Moment_overhang);

		$Lgtoverhang = $Overhang/(cos($Phi));
		
		$coefaover = -1*($Tctl_pli/2) * cos($Phi) * cos($Phi);
		$coefaoverv = -1*($Tctl_pli) * cos($Phi) * cos($Phi);
		$coefaovera = $Tctl_pli * sin($Phi) * cos($Phi);

		$coefaover_out = sprintf("%.1f", $coefaover);
		$coefaoverv_out = sprintf("%.1f", $coefaoverv);
		$coefaovera_out = sprintf("%.1f", $coefaovera);

		@momentsoverlist = ();
		@xoverlist = ();
		$Y_axisover = "";
		$X_axisover = "";
		$Y_axisoverv = "";
		$Y_axisovera = "";

		$Eqn_momentover = 'Overhang:  M(x) = ' . "$coefaover_out" . 'x^2';
		$Eqn_shearover = 'Overhang:  V(x) = ' . "$coefaoverv_out" . 'x';
		$Eqn_axialover = 'Overhang:  A(x) = ' . "$coefaovera_out" . 'x';

		for ($xover = 0; $xover < $Lgtoverhang; $xover +=.5)
		{

			$Momentoverx = $coefaover * $xover**2;
			$Shearoverx = $coefaoverv * $xover;
			$Axialoverx = $coefaovera * $xover;
			push (@momentsoverlist, $Momentoverx);
			push (@xoverlist, $xover);

			$Y_axisover .= "$Momentoverx,";
			$Y_axisoverv .= "$Shearoverx,";
			$Y_axisovera .= "$Axialoverx,";
			$X_axisover .= "\"$xover\",";

		}

		# Add point at panel point

		$xover = $Lgtoverhang;
		$Momentoverx = $coefaover *  $xover**2;
		$Shearoverx = $coefaoverv *  $xover;
		$Axialoverx = $coefaovera *  $xover;
		push (@momentsoverlist, $Momentoverx);
		push (@xoverlist, $xover);
	
		$Maxm_over_out = sprintf("%.1f", $Momentoverx);
		$Maxv_over_out = sprintf("%.1f", $Shearoverx);
		$Maxa_over_out = sprintf("%.1f", $Axialoverx);

		$Y_axisover .= "$Momentoverx";
		$Y_axisoverv .= "$Shearoverx";
		$Y_axisovera .= "$Axialoverx";
		$X_axisover .= "\"$xover\"";
	


		
		
	}
	else
	{

		$Moment_overhang = 0;
		$Overhang_graphs_html = "";
	}



	####################################
	#
	# DISTRIBUTED LOADS MATRIX ANALYSIS
	#
	####################################


	############################
	#
	# Known External Loads (Qk)
	#
	############################
	

	
	$Shear1 = $Shear4 = $Tctl_pli*$L12*.5;
	$Shear2 = $Shear3 = $Tctl_pli*$Ltcp*.5;

	$Moment1 = $Moment4 = ($Tctl_pli*$L12*$L12)/12;
	$Moment2 = $Moment3 = ($Tctl_pli*$Ltcp*$Ltcp)/12;
	

	$Shear5 = $Shear7 = $Bctl_pli*$L17*.5;
	$Shear6 = $Bctl_pli*$Lbcp*.5;	
	
	

	$Moment5 = $Moment7 = ($Bctl_pli*$L17*$L17)/12;
	$Moment6 = ($Bctl_pli*$Lbcp*$Lbcp)/12;

	$Moment1_neg = $Moment4_neg = -1*$Moment1;
	$Moment2_neg = $Moment3_neg = -1*$Moment2;
	$Moment5_neg = $Moment7_neg = -1*$Moment5;
	$Moment6_neg = -1*$Moment6;
	

	$Q1 = 0;
	$Q2 = -1*($Shear1 + $Shear2);
	$Q3 = $Moment1 - $Moment2;
	$Q4 = 0;
	$Q5 = -1*($Shear2 + $Shear3);
	$Q6 = $Moment2 - $Moment3;
	$Q7 = 0;
	$Q8 = -1*($Shear3 + $Shear4);
	$Q9 = $Moment3 - $Moment4;
	$Q10 = 0;
	$Q11 = -1*($Shear5 + $Shear6);
	$Q12 = $Moment6 - $Moment5;
	$Q13 = 0;
	$Q14 = -1*($Shear6 + $Shear7);
	$Q15 = $Moment7 - $Moment6;
	$Q16 = 0;
	$Q17 = ($Moment4 + $Moment5) - $Moment_overhang;
	$Q18 = -1*($Moment1 + $Moment7) + $Moment_overhang;


  	$Qk_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ $Q1 ]
      	[ $Q2 ]
      	[ $Q3 ]
      	[ $Q4 ]
      	[ $Q5 ]
      	[ $Q6 ]
      	[ $Q7 ]
      	[ $Q8 ]
      	[ $Q9 ]
      	[ $Q10 ]
      	[ $Q11 ]
      	[ $Q12 ]
      	[ $Q13 ]
      	[ $Q14 ]
      	[ $Q15 ]
      	[ $Q16 ]
      	[ $Q17 ]
      	[ $Q18 ]  	
MATRIX

	$Qk_matrix_out = sprintf("%s", $Qk_matrix);

	

	#################################
	#
	# Superimposed Loading TC and BC
	#
	#################################

	$q1super_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ 0 ]
      	[ $Shear1 ]
      	[ $Moment1 ]
      	[ 0 ]
      	[ $Shear1 ]
      	[ $Moment1_neg ]
MATRIX



	$q2super_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ 0 ]
      	[ $Shear2 ]
      	[ $Moment2 ]
      	[ 0 ]
      	[ $Shear2 ]
      	[ $Moment2_neg ]
MATRIX


###################################
#
# Corrected Superposition Matrices
#
###################################
	
	$Shear1new = $Tctl_pli*$L12*.5*cos($Phi);
	$Axial1new = $Tctl_pli*$L12*.5*sin($Phi);
	$Moment1new = $Moment1;
	$Moment1_negnew = $Moment1_neg;
	$Shear2new = $Tctl_pli*$Ltcp*.5*cos($Phi);
	$Axial2new = $Tctl_pli*$Ltcp*.5*sin($Phi);
	$Moment2new = $Moment2;
	$Moment2_negnew = $Moment2_neg;
	
	$q1super_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ $Axial1new ]
      	[ $Shear1new ]
      	[ $Moment1new ]
      	[ $Axial1new ]
      	[ $Shear1new ]
      	[ $Moment1_negnew ]
MATRIX


	
	$q2super_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ $Axial2new ]
      	[ $Shear2new ]
      	[ $Moment2new ]
      	[ $Axial2new ]
      	[ $Shear2new ]
      	[ $Moment2_negnew ]
MATRIX

	
###################################
	

	$q3super_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ 0 ]
      	[ $Shear3 ]
      	[ $Moment3 ]
      	[ 0 ]
      	[ $Shear3 ]
      	[ $Moment3_neg ]
MATRIX

	$q4super_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ 0 ]
      	[ $Shear4 ]
      	[ $Moment4 ]
      	[ 0 ]
      	[ $Shear4 ]
      	[ $Moment4_neg ]
MATRIX

	$q5super_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ 0 ]
      	[ $Shear5 ]
      	[ $Moment5 ]
      	[ 0 ]
      	[ $Shear5 ]
      	[ $Moment5_neg ]
MATRIX

	$q6super_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ 0 ]
      	[ $Shear6 ]
      	[ $Moment6 ]
      	[ 0 ]
      	[ $Shear6 ]
      	[ $Moment6_neg ]
MATRIX

	$q7super_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ 0 ]
      	[ $Shear7 ]
      	[ $Moment7 ]
      	[ 0 ]
      	[ $Shear7 ]
      	[ $Moment7_neg ]
MATRIX



	###########################
	#
	# Known Displacements (Dk)
	#
	###########################

	$Dk_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
      	[ 0 ]
      	[ 0 ]
      	[ 0 ]
MATRIX
		
	$Dk_matrix_out = sprintf("%s", $Dk_matrix);


$Matrix_email_output2 = qq{
Qk = 
$Qk_matrix_out

Dk = 
$Dk_matrix_out

	};
	
	###########################################
	#
	# Structure Stiffness Matrices for Members
	#
	###########################################

	%k1map = (1 => 20, 2 => 21, 3 => 18, 4 => 1, 5 => 2, 6 => 3);
	%k2map = (1 => 1, 2 => 2, 3 => 3, 4 => 4, 5 => 5, 6 => 6);
	%k3map = (1 => 4, 2 => 5, 3 => 6, 4 => 7, 5 => 8, 6 => 9);
	%k4map = (1 => 7, 2 => 8, 3 => 9, 4 => 16, 5 => 19, 6 => 17);
	%k5map = (1 => 10, 2 => 11, 3 => 12, 4 => 16, 5 => 19, 6 => 17);
	%k6map = (1 => 13, 2 => 14, 3 => 15, 4 => 10, 5 => 11, 6 => 12);
	%k7map = (1 => 20, 2 => 21, 3 => 18, 4 => 13, 5 => 14, 6 => 15);
	%k8map = (1 => 13, 2 => 14, 3 => 15, 4 => 1, 5 => 2, 6 => 3);
	%k9map = (1 => 13, 2 => 14, 3 => 15, 4 => 4, 5 => 5, 6 => 6);
	%k10map = (1 => 10, 2 => 11, 3 => 12, 4 => 4, 5 => 5, 6 => 6);
	%k11map = (1 => 10, 2 => 11, 3 => 12, 4 => 7, 5 => 8, 6 => 9);



	%Thetax = (1 => $Phi,
	2 => $Phi,
	3 => -$Phi,
	4 => -$Phi,
	5 => 0,
	6 => 0,
	7 => 0,
	8 => 3.14159265358 - $Alpha,
	9 => $Beta,
	10 => 3.14159265358 - $Beta,
	11 => $Alpha);

	
	%Area_hash = (1 => $A12,
	2 => $A23,
	3 => $A34,
	4 => $A45,
	5 => $A56,
	6 => $A67,
	7 => $A17,
	8 => $A27,
	9 => $A37,
	10 => $A36,
	11 => $A46);

	
	%Inertia_hash = (1 => $Ix12,
	2 => $Ix23,
	3 => $Ix34,
	4 => $Ix45,
	5 => $Ix56,
	6 => $Ix67,
	7 => $Ix17,
	8 => $Ix27,
	9 => $Ix37,
	10 => $Ix36,
	11 => $Ix46);


	%Length_hash = (1 => $Lgt12,
	2 => $Lgt23,
	3 => $Lgt34,
	4 => $Lgt45,
	5 => $Lgt56,
	6 => $Lgt67,
	7 => $Lgt17,
	8 => $Lgt27,
	9 => $Lgt37,
	10 => $Lgt36,
	11 => $Lgt46);


	if (($E12) && ($E17) && ($E27) && ($E37))
	{

	
	%E_hash = (1 => $E12,
	2 => $E12,
	3 => $E12,
	4 => $E12,
	5 => $E17,
	6 => $E17,
	7 => $E17,
	8 => $E27,
	9 => $E37,
	10 => $E37,
	11 => $E27);
	}
	else
	{
	%E_hash = (1 => 1500000,
	2 => 1500000,
	3 => 1500000,
	4 => 1500000,
	5 => 1500000,
	6 => 1500000,
	7 => 1500000,
	8 => 1500000,
	9 => 1500000,
	10 => 1500000,
	11 => 1500000);
	}

	%Lambdax = ();
	%Lambday = ();
	%T_matrix_hash = ();
	%TT_matrix_hash = ();
	%kprime_matrix_hash = ();
	%k_matrix_hash = ();
	%k_assembly_matrix_hash = ();

	$K_matrix = new Math::MatrixReal(21,21); 


	for ($Membernum = 1; $Membernum <= 11; $Membernum +=1)
	{

		$lmdx = $Lambdax{"$Membernum"} = cos($Thetax{"$Membernum"});
		$lmdy = $Lambday{"$Membernum"} = sin($Thetax{"$Membernum"});

		$lmdy_neg = $lmdy * -1;

		#########################################
		#
		# Displacement Transformation Matrix (T)
		#
		#########################################

		$T_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
		[ $lmdx     $lmdy 0  0         0    0 ]
      		[ $lmdy_neg $lmdx 0  0         0    0 ]
      		[  0         0    1  0         0    0 ]
      		[  0         0    0 $lmdx     $lmdy 0 ]
		[  0         0    0 $lmdy_neg $lmdx 0 ]
		[  0         0    0  0         0    1 ]
MATRIX
		
		$T_matrix_out = sprintf("%s", $T_matrix);
		$T_matrix_hash{"$Membernum"} = $T_matrix;

		$Matrix_email_output2 .= qq{
T$Membernum = 
$T_matrix_out
		};

		#####################################
		#
		# Force Transformation Matrix (TT)
		#
		#####################################

		$TT_matrix = new Math::MatrixReal(6,6); 
		$TT_matrix->transpose($T_matrix);
		$TT_matrix_out = sprintf("%s", $TT_matrix);
		$TT_matrix_hash{"$Membernum"} = $TT_matrix;

		$Matrix_email_output2 .= qq{
TT$Membernum = 
$TT_matrix_out
		};


		#####################################
		#
		# Member Stiffness Matrix (k')
		#
		#####################################

		$AE_L = ($Area_hash{"$Membernum"} * $E_hash{"$Membernum"})/$Length_hash{"$Membernum"};
		$EI12_LLL = (12 * $Inertia_hash{"$Membernum"} * $E_hash{"$Membernum"})/($Length_hash{"$Membernum"})**3;
		$EI6_LL = (6 * $Inertia_hash{"$Membernum"} * $E_hash{"$Membernum"})/($Length_hash{"$Membernum"})**2;
		$EI4_L = (4 * $Inertia_hash{"$Membernum"} * $E_hash{"$Membernum"})/($Length_hash{"$Membernum"});
		$EI2_L = (2 * $Inertia_hash{"$Membernum"} * $E_hash{"$Membernum"})/($Length_hash{"$Membernum"});
		
		$AE_L_neg = $AE_L * -1;
		$EI12_LLL_neg = $EI12_LLL * -1;
		$EI6_LL_neg = $EI6_LL * -1;



		if ($Membernum >= 8)
		{

			$kprime_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
			[ $AE_L     0  0  $AE_L_neg 0  0 ]
      			[  0        0  0   0        0  0 ]
      			[  0        0  0   0        0  0 ]
      			[ $AE_L_neg 0  0  $AE_L     0  0 ]
			[  0        0  0   0        0  0 ]
      			[  0        0  0   0        0  0 ]
MATRIX

		}
		else
		{

			
			
			$kprime_matrix = Math::MatrixReal->new_from_string(<<"MATRIX");
			[ $AE_L       0              0          $AE_L_neg  0              0          ]
      			[  0         $EI12_LLL      $EI6_LL      0        $EI12_LLL_neg  $EI6_LL     ]
      			[  0         $EI6_LL        $EI4_L       0        $EI6_LL_neg    $EI2_L      ]
      			[ $AE_L_neg   0              0          $AE_L      0              0          ]
			[  0         $EI12_LLL_neg  $EI6_LL_neg  0        $EI12_LLL      $EI6_LL_neg ]
      			[  0         $EI6_LL        $EI2_L       0        $EI6_LL_neg    $EI4_L      ]
MATRIX

			
		}
		
		$kprime_matrix_out = sprintf("%s", $kprime_matrix);
		$kprime_matrix_hash{"$Membernum"} = $kprime_matrix;

		$Matrix_email_output2 .= qq{
k'$Membernum = 
$kprime_matrix_out
		};


		#####################################
		#
		# Global Member Stiffness Matrix (k)
		#
		#####################################

		$k_matrix = ($TT_matrix * $kprime_matrix) * $T_matrix;

		$k_matrix_out = sprintf("%s", $k_matrix);
		$k_matrix_hash{"$Membernum"} = $k_matrix;

		$Matrix_email_output2 .= qq{
k$Membernum = 
$k_matrix_out
		};


		#######################################################
		#
		# Global Member Stiffness Matrix (Assembly) (k padded)
		#
		#######################################################


		$kmap_name = 'k'."$Membernum".'map';
		%kmap = %{$kmap_name};
		
		$k_assembly_matrix = new Math::MatrixReal(21,21); 

		for ($ith = 1; $ith <= 6; $ith +=1)
		{
			$rowmap = $kmap{"$ith"};
			for ($jth = 1; $jth <= 6; $jth +=1)
			{
				$colmap = $kmap{"$jth"};
				$kvalue = $k_matrix->element($ith,$jth); 
				$k_assembly_matrix->assign($rowmap,$colmap,$kvalue); 
			}

		}

		$k_assembly_matrix_out = sprintf("%s", $k_assembly_matrix);
		$k_assembly_matrix_hash{"$Membernum"} = $k_assembly_matrix;

		$Matrix_email_output2 .= qq{
k$Membernum Assembly = 
$k_assembly_matrix_out
		};


		################################
		#
		# Structure Stiffness Matrix (K)
		#
		################################



    		$K_matrix = $K_matrix + $k_assembly_matrix;

	

	}



	$K_matrix_out = sprintf("%s", $K_matrix);
		
	$Matrix_email_output2 .= qq{
K Structure Stiffness Matrix = 
$K_matrix_out
	};
	

	#############################
	#
	# K11 Matrix (K11)
	#
	#############################

	$K11_matrix = new Math::MatrixReal(18,18); 

	for ($ith = 1; $ith <= 18; $ith +=1)
	{	
		for ($jth = 1; $jth <= 18; $jth +=1)
		{
				$K11value = $K_matrix->element($ith,$jth); 
				$K11_matrix->assign($ith,$jth,$K11value); 
		}
	}

	$K11_matrix_out = sprintf("%s", $K11_matrix);
	#print qq{<br>$K11_matrix_out};
	
	$Matrix_email_output2 .= qq{
K11 Matrix = 
$K11_matrix_out
	};


	#############################
	#
	# K21 Matrix (K21)
	#
	#############################

	$K21_matrix = new Math::MatrixReal(3,18); 

	for ($ith = 19; $ith <= 21; $ith +=1)
	{
		$ith21 = $ith - 18; 
					
		for ($jth = 1; $jth <= 18; $jth +=1)
		{
				$K21value = $K_matrix->element($ith,$jth);
				$K21_matrix->assign($ith21,$jth,$K21value); 
		}
	}

	$K21_matrix_out = sprintf("%s", $K21_matrix);
		
	$Matrix_email_output2 .= qq{
K21 Matrix = 
$K21_matrix_out
	};


	#############################
	#
	# K11INV Matrix (K11 Inverse)
	#
	#############################

	$K11INV_matrix = $K11_matrix->decompose_LR->invert_LR;

	$K11INV_matrix_out = sprintf("%s", $K11INV_matrix);
	
	$I_matrix = $K11_matrix * $K11INV_matrix;

	$I_matrix_out = sprintf("%s", $I_matrix);
	
	$Matrix_email_output2 .= qq{
K11 Inverse Matrix = 
$K11INV_matrix_out

Indentity Matrix Check =
$I_matrix_out
	};


	#############################
	#
	# Du Matrix (Du)
	#
	#############################


	$Du_matrix = $K11INV_matrix * $Qk_matrix;

	$Du_matrix_out = sprintf("%s", $Du_matrix);
		
	$Matrix_email_output2 .= qq{
Du Matrix = 
$Du_matrix_out
	};


	#############################
	#
	# Qu Matrix (Qu)
	#
	#############################


	$Qu_matrix = $K21_matrix * $Du_matrix;

	$Qu_matrix_out = sprintf("%s", $Qu_matrix);
		
	$Matrix_email_output2 .= qq{
Qu Matrix = 
$Qu_matrix_out
	};



	####################################
	#
	# Displacement Matrix (D = Du + Dk)
	#
	####################################

	$DuT_matrix = new Math::MatrixReal(1,18); 
	$DuT_matrix->transpose($Du_matrix);
	
	$DkT_matrix = new Math::MatrixReal(1,3); 
	$DkT_matrix->transpose($Dk_matrix);

	$DT_matrix = $DuT_matrix . $DkT_matrix;

	$D_matrix = new Math::MatrixReal(21,1); 
	$D_matrix->transpose($DT_matrix);

	$D_matrix_out = sprintf("%s", $D_matrix);
		
	$Matrix_email_output2 .= qq{
D Matrix = 
$D_matrix_out
	};


	
	####################################
	#
	# Internal Loading of Members (q)
	#
	####################################

	%D_member_matrix_hash = ();
	%q_matrix_hash = ();
	%qmod_matrix_hash = ();
	
	for ($Membernum = 1; $Membernum <= 11; $Membernum +=1)
	{

		$kmap_name = 'k'."$Membernum".'map';
		%kmap = %{$kmap_name};
		
		$D_member_matrix = new Math::MatrixReal(6,1); 

		for ($ith = 1; $ith <= 6; $ith +=1)
		{
			$rowmap = $kmap{"$ith"};
			$kvalue = $D_matrix->element($rowmap,1); 
			$D_member_matrix->assign($ith,1,$kvalue); 
		}

		$D_member_matrix_out = sprintf("%s", $D_member_matrix);
		$D_member_matrix_hash{"$Membernum"} = $D_member_matrix;

		
		$q_matrix = ($kprime_matrix_hash{"$Membernum"} * $T_matrix_hash{"$Membernum"}) * $D_member_matrix;
		
		$q_matrix_out = sprintf("%s", $q_matrix);
		$q_matrix_hash{"$Membernum"} = $q_matrix;		


		if ($Membernum =~ /^[1234567]$/)
		{
			# $q1super_matrix
			
			$super_matrix_name = 'q'."$Membernum".'super_matrix';
			$qsuper_matrix = ${$super_matrix_name};	
			
			$qmod_matrix = $qsuper_matrix + $q_matrix;
			
			$qmod_matrix_out = sprintf("%s", $qmod_matrix);
			$qmod_matrix_hash{"$Membernum"} = $qmod_matrix;

			$Matrix_email_output_full .= qq{
D$Membernum = 
$D_member_matrix_out

q$Membernum =
$q_matrix_out

q$Membernum (mod) =
$qmod_matrix_out
			};


$Matrix_email_output .= qq{
q$Membernum (mod) =
$qmod_matrix_out
			};

		}
		else
		{

			$Matrix_email_output_full .= qq{
D$Membernum = 
$D_member_matrix_out

q$Membernum =
$q_matrix_out
			};

$Matrix_email_output .= qq{
q$Membernum =
$q_matrix_out
			};

		}	

	}


	# TEST AND DEBUG CODE

	#$Matrix_email_output .= join("\n", %Thetax);
	#$Matrix_email_output .= join("\n", %Lambdax);
	#$Matrix_email_output .= join("\n", %Lambday);
	$Matrix_email_output .= join("\n", %E_hash);

A confused student is a good student.
Nathaniel P. Wilkerson, PE

http://www.medeek.com

 

[medeek]

This block of code coincides with the diagram of the fink truss above, and does not include any point loads, just the distributed loads shown. Not sure if this is very helpful or not but the same procedure to analyze a member horizontal or at some angle to the horizon should be used to analyze a member that is vertical in my opinion. If the system of local and general coordinate systems is consistent the loads should work out as expected.

A confused student is a good student.
Nathaniel P. Wilkerson, PE

http://www.medeek.com

 

[IDS]

Place W at the bottom node and find the maximum axial load in the column. To find the axial load elsewhere on the column, simply subtract the weight of column between the point in question and the bottom node.

That will work for finding the forces in a determinant structure, but it won't give the right answer for deflections or for forces in an indeterminant structure.

The procedure for a stiffness analysis is:
[ol 1]
[li]For each member find the fixed end actions for all applied loads[/li]
[li]Calculate the node deflections due to these end actions[/li]
[li]Calculate the reactions at the supports[/li]
[li]Knowing the support reactions, calculate the actions along each member, taking account of distributed loads, including self weight.[/li]
[/ol]

The axial fixed end forces for the self weight of a vertical member are +-W/2.

Doug Jenkins
Interactive Design Services

http://newtonexcelbach.wordpress.com/

 

[BAretired]

Upon further reflection, I agree with a fixed end force of W/2 at each node due to the column weight, W.

BA

 

[medeek]

Set up your global displacement and loading (Dk and Qk) matrices carefully, things get a bit interesting at the reaction points. Example 16-1 and 16-2 in Hibbeler's text are very close to your problem, they would be worth taking a look at.

Where the column meets the ground (column bottom) this appears to be zero degrees of freedom at this node. Since this is the case no external forces or moments are placed at this node. At the other reaction point it appears your have a pin joint so only rotational degree of freedom is present. Based on your diagram your global displacement matrix should have 5 rows and your global loading matrix 4 rows.

There is no external lateral force so that entry in your loading matrix will be zero. The vertical load at node 2 will be Wb/2 + Wc/2, the other two entries in the loading matrix will be the moments at each end of the beam (with opposite signs).

A confused student is a good student.
Nathaniel P. Wilkerson, PE

http://www.medeek.com

 

[BAretired]

medeek,

If you want the program to output the correct reactions, why would you not place external forces and moments at a node having zero degrees of freedom?

BA

 

[IDS]

BA - Step 2 (as listed in my previous post above) is to find the deflections and rotations at each node, and from that you get the distribution of the applied loads to the members meeting at a node. At any fixed node you already know the deflection and/or rotation (zero), and the distribution of the applied nodal loads (100% to the fixed freedom), so you don't need to calculate these things. In fact in the standard stiffness method including fixed freedoms effectively results in dividing by zero in the solution process, so it won't work.

When you work out the support reactions you of course need to add in the fixed end actions at the supports.

Medeek - I think the horizontal member in the original post is a free cantilever, so there are 6 degrees of freedom, but if it was a pinned connection there would be 4 degrees of freedom, so 4 rows in the load matrix as well.

Doug Jenkins
Interactive Design Services

http://newtonexcelbach.wordpress.com/

 

[medeek]

To be perfectly honest I don't have a real intuitive feel for setting up the matrix analysis, I find it very abstract but I've reviewed numerous examples from Hibbeler's text and I fully understand the process and the algorithm that gets you there. The really cool thing is being able to take the output of the matrice calculations and then using these values as coefficients to develop actual equations for the axial, shear and moment diagrams. This post processing is not too bad when the loads are all distributed loads, however it gets really interesting when you add point loads into the mix. After numerous trial and error runs and head scratching I was also able to handle point loads by running two separate matrix analysis runs (one for point loads and one for distributed loads) and then combining the results using the principle of superposition to generate the equations and the graphs.

A confused student is a good student.
Nathaniel P. Wilkerson, PE

http://www.medeek.com