function [kg]=frameg(x1,y1,x2,y2,el,xi,area,p);

%--------------------------------------------------------------
%  Synopsis:
%     [kg]=frame(x1,y1,x2,y2,el,xi,area,p)
%
%  Variable Description:
%     kg  - element global geometric stiffness matrix (size of 6x6)   
%     x1 - the X coordinate of the node 1
%     x2 - the X coordinate of the node 2
%     y1 - the Y coordinate of the node 1
%     y2 - the Y coordinate of the node 2
%     el - elastic modulus 
%     xi - second moment of inertia of cross-section
%     area - area of beam cross-section
%     p - the normalized axial force in the member
%
%  Internal Variables:
%     kgl   - local element geomteric stiffness matrix (size of 6x6)
%     leng - element length
%     t    - orientation transformation matrix (size of 6x6)
%     beta - angle between the local and global axes  ipt = 1: consistent mass matrix
%            is positive if the local axis is in the ccw direction from
%            the global axis
%--------------------------------------------------------------------------
   
% calculate element length and rotation angle

 leng = sqrt((x2-x1)^2 + (y2-y1)^2);
 
 dx = x2-x1;
 dy = y2-y1;
 if dx == 0 
     if dy > 0
         beta = pi/2.;
     else
         beta = -pi/2.;
     end
 else
     beta = atan (dy/dx);
 end
 
% stiffness matrix at the local axis

 a=el*area/leng;
 c=el*xi/(leng^3);
 kgl=(p/leng)*[0   0           0            0		0          0;...
     		0   6/5        -leng/10      0    -6/5	-leng/10;...
     		0  -leng/10	2/15*leng^2     0     leng/10   -leng^2/30;...
    	 	0   0           0             0     0          0;...
     		0  -6/5        leng/10       0     6/5        leng/10;...
     		0  -leng/10    -leng^2/30    0     leng/10   2*leng^2/15];


% transformation matrix

 t=[ cos(beta)  sin(beta)  0   0          0         0;...
    -sin(beta)  cos(beta)  0   0          0         0;...
     0          0          1   0          0         0;...
     0          0          0   cos(beta)  sin(beta) 0;...
     0          0          0  -sin(beta)  cos(beta) 0;...
     0          0          0   0          0         1];

 
% stiffness matrix in the global coordinate system

 kg=t'*kgl*t;


