%=========================================
%PARAMETRIC MODEL ANALYSIS RESULTS CONTROL
%=========================================
% _________________________________________________________________________
%|____10NODE_20DOF ISOPARAMETRIC TRIANGULAR PLANE STRESS FEM___[A.]_(C)(R)|
%| |
%| Shape function:Homogen |
%| Stifness :Topology+Accumulate method |
%|_______________________________Solve Equation:Cholesky-[L][D][u]_________|
%| Maxima 5.9.0.9beta2 http://maxima.sourceforge.net |
%| Using Lisp Kyoto Common Lisp GCL 2.6.3 (aka GCL) |
%| Distributed under the GNU Public License. See the file COPYING. |
%| Dedicated to the memory of William Schelter. |
%| This is a development version of Maxima. The function bug_report() |
%| provides bug reporting information. |
%| |
%|_____PARAMETRIC FEM ANALYSIS SUBPROGRMAM_______Maxima5.9.0.9beta2_[GPL]__|
%========================================This procedure runing MAXIMA
%PROCEDURE(1):Jacobian matrix
%PROCEDURE(2):Element Stiffness matrix depend the Jacobian transform
%_____________________________________[PROCEDURE(1)]
%N1(x,y):=1/2*e1*(3*e1-1)*(3*e1-2);
%N2(x,y):=1/2*e2*(3*e2-1)*(3*e2-2);
%N3(x,y):=1/2*e3*(3*e3-1)*(3*e3-2);
%N4(x,y):=9/2*e1*e2*(3*e1-1);
%N5(x,y):=9/2*e1*e2*(3*e2-1);
%N6(x,y):=9/2*e2*e3*(3*e2-1);
%N7(x,y):=9/2*e2*e3*(3*e3-1);
%N8(x,y):=9/2*e3*e1*(3*e3-1);
%N9(x,y):=9/2*e3*e1*(3*e1-1);
%N10(x,y):=27*e1*e2*e3;
%H1e(e,n):=diff(N1(e,n),e1)-diff(N1(e,n),e3);
%H2e(e,n):=diff(N2(e,n),e1)-diff(N2(e,n),e3);
%H3e(e,n):=diff(N3(e,n),e1)-diff(N3(e,n),e3);
%H4e(e,n):=diff(N4(e,n),e1)-diff(N4(e,n),e3);
%H5e(e,n):=diff(N5(e,n),e1)-diff(N5(e,n),e3);
%H6e(e,n):=diff(N6(e,n),e1)-diff(N6(e,n),e3);
%H7e(e,n):=diff(N7(e,n),e1)-diff(N7(e,n),e3);
%H8e(e,n):=diff(N8(e,n),e1)-diff(N8(e,n),e3);
%H9e(e,n):=diff(N9(e,n),e1)-diff(N9(e,n),e3);
%H10e(e,n):=diff(N10(e,n),e1)-diff(N10(e,n),e3);
%H1n(e,n):=diff(N1(e,n),e2)-diff(N1(e,n),e3);
%H2n(e,n):=diff(N2(e,n),e2)-diff(N2(e,n),e3);
%H3n(e,n):=diff(N3(e,n),e2)-diff(N3(e,n),e3);
%H4n(e,n):=diff(N4(e,n),e2)-diff(N4(e,n),e3);
%H5n(e,n):=diff(N5(e,n),e2)-diff(N5(e,n),e3);
%H6n(e,n):=diff(N6(e,n),e2)-diff(N6(e,n),e3);
%H7n(e,n):=diff(N7(e,n),e2)-diff(N7(e,n),e3);
%H8n(e,n):=diff(N8(e,n),e2)-diff(N8(e,n),e3);
%H9n(e,n):=diff(N9(e,n),e2)-diff(N9(e,n),e3);
%H10n(e,n):=diff(N10(e,n),e2)-diff(N10(e,n),e3);
%Pos(x,y):=matrix([x1,y1],[x2,y2],[x3,y3],[x4,y4],[x5,y5],[x6,y6],[x7,y7],[x8,y8],[x9,y9],[x10,y10]);
%N(e,n):=matrix([H1e(e,n),H2e(e,n),H3e(e,n),H4e(e,n),H5e(e,n),H6e(e,n),H7e(e,n),H8e(e,n),H9e(e,n),H10e(e,n)],
% [H1n(e,n),H2n(e,n),H3n(e,n),H4n(e,n),H5n(e,n),H6n(e,n),H7n(e,n),H8n(e,n),H9n(e,n),H10n(e,n)]);
%J(e,n):=N(e,n).Pos(x,y);
%J(e,n);
%Fortmx (J,J(e,n));
%_____________________________________[PROCEDURE(2)]
%Ju(e,n):=matrix([g11,g12],[g21,g22]);
%Ju(e,n);
%Nx(e,n):=Ju(e,n).N(e,n);
%Nx(e,n);
%fortmx(H,Nx(e,n));
%=====================================================[END MAXIMA]
% _________________________________________________________________________
%|_____10NODE_20DOF ISOPARAMETRIC TRIANGULAR PLANE STRESS FEM___[A.]_(C)(R)|
%| |
%| Shape function:Homogen |
%| Stifness :Topology+Accumulate method |
%|_______________________________Solve Equation:Cholesky-[L][D][u]_________|
%| |
%|_____PARAMETRIC FEM ANALYSIS SUBPROGRMAM_______Maxima5.9.0.9beta2_[GPL]__|
clc
clear
%===============================================INPUT DATA
%===============MATERIALS PROPERTIES
E=2E6; %Element elasticity constant;
v=0.3; %Element material poission ratio;
th=0.10; %Element thickness
%===============
%==================POSITION MATRIX
%Pos(Element No,:)=[Element Node Number]
Pos(1,:)= [1 34 31 12 23 33 32 21 11 22];
Pos(2,:)= [1 4 34 2 3 14 24 23 12 13];
Pos(3,:)= [4 37 34 15 26 36 35 24 14 25];
Pos(4,:)= [4 7 37 5 6 17 27 26 15 16];
Pos(5,:)= [7 40 37 18 29 39 38 27 17 28];
Pos(6,:)= [7 10 40 8 9 20 30 29 18 19];
Pos(7,:)= [31 34 51 32 33 43 52 51 41 42];
Pos(8,:)= [61 34 64 52 43 44 54 63 62 53];
Pos(9,:)= [34 37 64 35 36 46 55 54 44 45];
Pos(10,:)=[64 37 67 55 46 47 57 66 65 56];
Pos(11,:)=[37 40 67 38 39 49 58 57 47 48];
Pos(12,:)=[67 40 70 58 48 50 60 69 68 59];
%==================
Number=size(Pos);
No=Number(1);
%_____________Automatic coordinate function
%Cor=Element Node Castesian Coordinate for Global System Axis
% Cor=[ xi yi ]
nomin=0;
for sut=1:7;
for sat=1:10;
nomin=nomin+1;
Cor(nomin,1)=0.25*(sat-1);
Cor(nomin,2)=0.20*(sut-1);
end
end
clear sut sat nomin;
Number=size(Cor);
Node=Number(1); %System Node Number
for i=1:Node;
Re(i,:)=[1 1];
end
%===============SYSTEM SUPPORT
%Re(Node number,:)=[ux vy]
Re(1,:) =[0 0];
Re(10,:) =[1 0];
Re(61,:) =[0 1];
%===============
Number=size(Re);
Nom=Number(2); %Plane Element node d.o.f Number
%_______________Topology and accumulate method
sayman=0;
for i=1:Node;
for j=1:Nom;
if Re(i,j)==1 ;
sayman = sayman +1;
Re(i,j) = sayman;
end
end
end
Item=sayman; %Plane system sum deplacement value
%_________________________________Modal deplacement parameter
for i=1:No;
for j=1:Nom;
R(i,j+0*Nom)= Re(Pos(i,1),j);
R(i,j+1*Nom)= Re(Pos(i,2),j);
R(i,j+2*Nom)= Re(Pos(i,3),j);
R(i,j+3*Nom)= Re(Pos(i,4),j);
R(i,j+4*Nom)= Re(Pos(i,5),j);
R(i,j+5*Nom)= Re(Pos(i,6),j);
R(i,j+6*Nom)= Re(Pos(i,7),j);
R(i,j+7*Nom)= Re(Pos(i,8),j);
R(i,j+8*Nom)= Re(Pos(i,9),j);
R(i,j+9*Nom)= Re(Pos(i,10),j);
end
end
%============SYSTEM LOAD
%Re(Node number, Freedoom number)[u v]
P(Item)=0;
P(Re(40,1)) =-5;
P(Re(64,2))=-10;
P(Re(67,2))=-10;
%=============
%___________________________________________________END_INPUT_DATA
%===================================================SYSTEM_ANALYSIS
%open dimension local element stiffness matrix
K(20,20,No)=0;
for s=1:No;%Element number
X1=Cor(Pos(s,1),1); Y1=Cor(Pos(s,1),2);
X2=Cor(Pos(s,2),1); Y2=Cor(Pos(s,2),2);
X3=Cor(Pos(s,3),1); Y3=Cor(Pos(s,3),2);
X4=Cor(Pos(s,4),1); Y4=Cor(Pos(s,4),2);
X5=Cor(Pos(s,5),1); Y5=Cor(Pos(s,5),2);
X6=Cor(Pos(s,6),1); Y6=Cor(Pos(s,6),2);
X7=Cor(Pos(s,7),1); Y7=Cor(Pos(s,7),2);
X8=Cor(Pos(s,8),1); Y8=Cor(Pos(s,8),2);
X9=Cor(Pos(s,9),1); Y9=Cor(Pos(s,9),2);
X10=Cor(Pos(s,10),1); Y10=Cor(Pos(s,10),2);
Area=0.5*(-X1*(Y3-Y2)+X2*Y3-X3*Y2+(X3-X2)*Y1);
%12 Point Gauss_Legendre_Numerical_Integraton_Coefficients
%Multiply value
ha=0.050844906370207;hb=0.116786275726379;hc=0.082851075618374;
%Weight function value
wa=0.872831971016996;
wb=0.063089014491502;
wc=0.501426509659179;
wd=0.249286745170910;
we=0.636502499121399;
wf=0.310352451033785;
wg=0.053145049844816;
W=[wa wb wb; %ha
wb wa wb;
wb wb wa; %hb
wc wd wd;
wd wc wd;
wd wd wc; %hc
we wf wg;
we wg wf;
wf wg we;
wf we wg;
wg we wf;
wg wf we];
for j=1:12
if j==1 ;h=ha;end
if j==4 ;h=hb;end
if j==7 ;h=hc;end
e1=W(j,1);
e2=W(j,2);
e3=W(j,3);
%====================================[PROCEDURE(1)]Maxima
%Jacobian matrix terms
J11 = (9.0*(3*e1-1)*e3/2.0+27.0*e1*e3/2.0+(-9.0)*e1*(3*e1-1)/2.0)*X9+(9.0*e3*(3*e3-1)/2.0+(-9.0)*e1*(3*e3-1)/2.0+(-27.0)*e1*e3/2.0)*X8+((-9.0)*e2*(3*e3-1)/2.0+(-27.0)*e2*e3/2.0)*X7+(-9.0)*e2*(3*e2-1)*X6/2.0+9.0*e2*(3*e2-1)*X5/2.0+(9.0*(3*e1-1)*e2/2.0+27.0*e1*e2/2.0)*X4+(-(3*e3-2)*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-2)/2.0)*X3+(27*e2*e3-27*e1*e2)*X10+((3*e1-2)*(3*e1-1)/2.0+3.0*e1*(3*e1-1)/2.0+3.0*e1*(3*e1-2)/2.0)*X1;
J12 = (9.0*(3*e1-1)*e3/2.0+27.0*e1*e3/2.0+(-9.0)*e1*(3*e1-1)/2.0)*Y9+(9.0*e3*(3*e3-1)/2.0+(-9.0)*e1*(3*e3-1)/2.0+(-27.0)*e1*e3/2.0)*Y8+((-9.0)*e2*(3*e3-1)/2.0+(-27.0)*e2*e3/2.0)*Y7+(-9.0)*e2*(3*e2-1)*Y6/2.0+9.0*e2*(3*e2-1)*Y5/2.0+(9.0*(3*e1-1)*e2/2.0+27.0*e1*e2/2.0)*Y4+(-(3*e3-2)*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-2)/2.0)*Y3+(27*e2*e3-27*e1*e2)*Y10+((3*e1-2)*(3*e1-1)/2.0+3.0*e1*(3*e1-1)/2.0+3.0*e1*(3*e1-2)/2.0)*Y1;
J21 = (-9.0)*e1*(3*e1-1)*X9/2.0+((-9.0)*e1*(3*e3-1)/2.0+(-27.0)*e1*e3/2.0)*X8+(9.0*e3*(3*e3-1)/2.0+(-9.0)*e2*(3*e3-1)/2.0+(-27.0)*e2*e3/2.0)*X7+(9.0*(3*e2-1)*e3/2.0+27.0*e2*e3/2.0+(-9.0)*e2*(3*e2-1)/2.0)*X6+(9.0*e1*(3*e2-1)/2.0+27.0*e1*e2/2.0)*X5+9.0*e1*(3*e1-1)*X4/2.0+(-(3*e3-2)*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-2)/2.0)*X3+((3*e2-2)*(3*e2-1)/2.0+3.0*e2*(3*e2-1)/2.0+3.0*e2*(3*e2-2)/2.0)*X2+(27*e1*e3-27*e1*e2)*X10;
J22 = (-9.0)*e1*(3*e1-1)*Y9/2.0+((-9.0)*e1*(3*e3-1)/2.0+(-27.0)*e1*e3/2.0)*Y8+(9.0*e3*(3*e3-1)/2.0+(-9.0)*e2*(3*e3-1)/2.0+(-27.0)*e2*e3/2.0)*Y7+(9.0*(3*e2-1)*e3/2.0+27.0*e2*e3/2.0+(-9.0)*e2*(3*e2-1)/2.0)*Y6+(9.0*e1*(3*e2-1)/2.0+27.0*e1*e2/2.0)*Y5+9.0*e1*(3*e1-1)*Y4/2.0+(-(3*e3-2)*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-2)/2.0)*Y3+((3*e2-2)*(3*e2-1)/2.0+3.0*e2*(3*e2-1)/2.0+3.0*e2*(3*e2-2)/2.0)*Y2+(27*e1*e3-27*e1*e2)*Y10;
%Jacobian matrix
Jacobi=[J11 J12;
J21 J22];
%Invers to Jacobian matrix
InvJacobi=Jacobi^-1;
g11=InvJacobi(1,1); g12=InvJacobi(1,2);
g21=InvJacobi(2,1); g22=InvJacobi(2,2);
%====================================[PROCEDURE(2)]Maxima
Hx1 = ((3*e1-2)*(3*e1-1)/2.0+3.0*e1*(3*e1-1)/2.0+3.0*e1*(3*e1-2)/2.0)*g11;
Hx2 = ((3*e2-2)*(3*e2-1)/2.0+3.0*e2*(3*e2-1)/2.0+3.0*e2*(3*e2-2)/2.0)*g12;
Hx3 = (-(3*e3-2)*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-2)/2.0)*g12+(-(3*e3-2)*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-2)/2.0)*g11;
Hx4 = 9.0*e1*(3*e1-1)*g12/2.0+(9.0*(3*e1-1)*e2/2.0+27.0*e1*e2/2.0)*g11;
Hx5 = (9.0*e1*(3*e2-1)/2.0+27.0*e1*e2/2.0)*g12+9.0*e2*(3*e2-1)*g11/2.0;
Hx6 = (9.0*(3*e2-1)*e3/2.0+27.0*e2*e3/2.0+(-9.0)*e2*(3*e2-1)/2.0)*g12+(-9.0)*e2*(3*e2-1)*g11/2.0;
Hx7 = (9.0*e3*(3*e3-1)/2.0+(-9.0)*e2*(3*e3-1)/2.0+(-27.0)*e2*e3/2.0)*g12+((-9.0)*e2*(3*e3-1)/2.0+(-27.0)*e2*e3/2.0)*g11;
Hx8 = ((-9.0)*e1*(3*e3-1)/2.0+(-27.0)*e1*e3/2.0)*g12+(9.0*e3*(3*e3-1)/2.0+(-9.0)*e1*(3*e3-1)/2.0+(-27.0)*e1*e3/2.0)*g11;
Hx9 = (-9.0)*e1*(3*e1-1)*g12/2.0+(9.0*(3*e1-1)*e3/2.0+27.0*e1*e3/2.0+(-9.0)*e1*(3*e1-1)/2.0)*g11;
Hx10 = (27*e1*e3-27*e1*e2)*g12+(27*e2*e3-27*e1*e2)*g11;
Hy1 = ((3*e1-2)*(3*e1-1)/2.0+3.0*e1*(3*e1-1)/2.0+3.0*e1*(3*e1-2)/2.0)*g21;
Hy2 = ((3*e2-2)*(3*e2-1)/2.0+3.0*e2*(3*e2-1)/2.0+3.0*e2*(3*e2-2)/2.0)*g22;
Hy3 = (-(3*e3-2)*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-2)/2.0)*g22+(-(3*e3-2)*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-1)/2.0+(-3.0)*e3*(3*e3-2)/2.0)*g21;
Hy4 = 9.0*e1*(3*e1-1)*g22/2.0+(9.0*(3*e1-1)*e2/2.0+27.0*e1*e2/2.0)*g21;
Hy5 = (9.0*e1*(3*e2-1)/2.0+27.0*e1*e2/2.0)*g22+9.0*e2*(3*e2-1)*g21/2.0;
Hy6 = (9.0*(3*e2-1)*e3/2.0+27.0*e2*e3/2.0+(-9.0)*e2*(3*e2-1)/2.0)*g22+(-9.0)*e2*(3*e2-1)*g21/2.0;
Hy7 = (9.0*e3*(3*e3-1)/2.0+(-9.0)*e2*(3*e3-1)/2.0+(-27.0)*e2*e3/2.0)*g22+((-9.0)*e2*(3*e3-1)/2.0+(-27.0)*e2*e3/2.0)*g21;
Hy8 = ((-9.0)*e1*(3*e3-1)/2.0+(-27.0)*e1*e3/2.0)*g22+(9.0*e3*(3*e3-1)/2.0+(-9.0)*e1*(3*e3-1)/2.0+(-27.0)*e1*e3/2.0)*g21;
Hy9 = (-9.0)*e1*(3*e1-1)*g22/2.0+(9.0*(3*e1-1)*e3/2.0+27.0*e1*e3/2.0+(-9.0)*e1*(3*e1-1)/2.0)*g21;
Hy10 = (27*e1*e3-27*e1*e2)*g22+(27*e2*e3-27*e1*e2)*g21;
%Connection Matrix
A(:,:,s)= [Hx1 0 Hx2 0 Hx3 0 Hx4 0 Hx5 0 Hx6 0 Hx7 0 Hx8 0 Hx9 0 Hx10 0 ;
0 Hy1 0 Hy2 0 Hy3 0 Hy4 0 Hy5 0 Hy6 0 Hy7 0 Hy8 0 Hy9 0 Hy10 ;
Hy1 Hx1 Hy2 Hx2 Hy3 Hx3 Hy4 Hx4 Hy5 Hx5 Hy6 Hx6 Hy7 Hx7 Hy8 Hx8 Hy9 Hx9 Hy10 Hx10 ];
%Elasticity Matrix
C=E/(1-v^2)*[1 v 0 ;
v 1 0 ;
0 0 0.5*(1-v)];
%Gauss-Legendre Numerical Integration
K(:,:,s)=K(:,:,s)+0.5*th*det(Jacobi)*h*A(:,:,s)'*C*A(:,:,s);
end
end
%open dimension system stiffness matrix dimension
Ksis(Item,Item)=0;
for n=1:No;
for sat=1:20;
for sut=1:20;
if (R(n,sat)~=0)
if (R(n,sut)~=0);
Ksis(R(n,sut),R(n,sat))=Ksis(R(n,sut),R(n,sat)) + K(sat,sut,n);
end
end
end
end
end
%System global stiffness matris singularity check
equation=size(Ksis);
if equation(1)~=rank(Ksis)
display('This system stiffness matrix is badly scaled')
R
error('Control system support boundary conditions')
else
Ku = inv(Ksis);
D = Ku *P';
end
clear equation
%Element per unit node topolog matrix
for v = 1 : No;
for m = 1 :20;
u = R(v, m);
if u ~=0
Hu(v, m) = D(u) ;
else
Hu(v,m)=0;
end
end
end
%Global system node displacement is moving local element nodes.
for s=1:No;
for i=1:10*Nom;
if R(s,i)~=0 ;
Dep(s,i)=D(R(s,i));
else
Dep(s,i)=0;
end
end
end
display('Global displacement [u(x) v(y)]T');
Dep'
for v=1 :No;
Pg(:,v) = K(:,:,v)*Hu(v,:)';
end
display('Global system node reactions');
Pg
contourf(Ksis)
title('Global system stiffness band matrix');
%=================================Node Patch Testing Code
for sutun=1:No;
for sat=1:10;
Pgx(sat,sutun)=Pg(2*sat-1,sutun);
Pgy(sat,sutun)=Pg(2*sat,sutun);
end
end
Pdx(Node)=0;
Pdy(Node)=0;
for Elemanno=1:No;
for sutun=1:10;
Pdx(Pos(Elemanno,sutun))=Pdx(Pos(Elemanno,sutun))+Pgx(sutun,Elemanno);
Pdy(Pos(Elemanno,sutun))=Pdy(Pos(Elemanno,sutun))+Pgy(sutun,Elemanno);
end
end
%Information:This Pdx,Pdy different the zero terms Support Reactions
('Manual Patch Test (node kin. eq. met. [Pgx Pgy])')
[Pdx' Pdy']
%================================================INFORMATION
%____________________PARAMETRIC Element stiffness matrix value (1)
%1.0e+006 *
%
% Columns 1 through 9
%
% 0.0409 0 0 -0.0067 -0.0084 0.0067 0.0036 -0.0548 0.0036
% 0 0.1168 -0.0058 0 0.0058 -0.0240 -0.0470 0.0103 0.0198
% 0 -0.0058 0.0747 0 -0.0154 0.0058 0.0066 0.0198 0.0066
% -0.0067 0 0 0.0262 0.0067 -0.0054 0.0231 0.0023 -0.0548
% -0.0084 0.0058 -0.0154 0.0067 0.1156 -0.0607 -0.0102 0.0054 -0.0102
% 0.0067 -0.0240 0.0058 -0.0054 -0.0607 0.1429 0.0054 -0.0126 0.0054
% 0.0036 -0.0470 0.0066 0.0231 -0.0102 0.0054 0.4590 -0.1205 -0.0918
% -0.0548 0.0103 0.0198 0.0023 0.0054 -0.0126 -0.1205 0.5674 -0.0427
% 0.0036 0.0198 0.0066 -0.0548 -0.0102 0.0054 -0.0918 -0.0427 0.4590
% 0.0231 0.0103 -0.0470 0.0023 0.0054 -0.0126 -0.0538 -0.1135 -0.1205
% -0.0036 0 -0.1187 0.0548 0.0629 -0.0251 0.0325 -0.0241 -0.1623
% 0 -0.0103 0.0470 -0.0415 -0.0284 0.0311 -0.0241 0.0927 0.1205
% -0.0036 0 0.0593 -0.0231 -0.1151 0.0416 0.0325 -0.0241 0.0325
% 0 -0.0103 -0.0198 0.0208 0.0495 -0.0312 -0.0241 0.0927 -0.0241
% 0.0325 -0.0198 -0.0066 0 -0.0583 0.0495 0.0593 -0.0241 0.0593
% -0.0231 0.0927 0 -0.0023 0.0416 -0.1831 -0.0241 0.0208 -0.0241
% -0.0649 0.0470 -0.0066 0 0.0390 -0.0284 -0.2967 0.1205 0.0593
% 0.0548 -0.1854 0 -0.0023 -0.0251 0.0950 0.1205 -0.1038 -0.0241
% 0 0 0 0 0 0 -0.1947 0.1446 -0.3560
% 0 0 0 0 0 0 0.1446 -0.5563 0.1446
% Columns 10 through 18
% 0.0231 -0.0036 0 -0.0036 0 0.0325 -0.0231 -0.0649 0.0548
% 0.0103 0 -0.0103 0 -0.0103 -0.0198 0.0927 0.0470 -0.1854
% -0.0470 -0.1187 0.0470 0.0593 -0.0198 -0.0066 0 -0.0066 0
% 0.0023 0.0548 -0.0415 -0.0231 0.0208 0 -0.0023 0 -0.0023
% 0.0054 0.0629 -0.0284 -0.1151 0.0495 -0.0583 0.0416 0.0390 -0.0251
% -0.0126 -0.0251 0.0311 0.0416 -0.0312 0.0495 -0.1831 -0.0284 0.0950
% -0.0538 0.0325 -0.0241 0.0325 -0.0241 0.0593 -0.0241 -0.2967 0.1205
% -0.1135 -0.0241 0.0927 -0.0241 0.0927 -0.0241 0.0208 0.1205 -0.1038
% -0.1205 -0.1623 0.1205 0.0325 -0.0241 0.0593 -0.0241 0.0593 -0.0241
% 0.5674 0.1205 -0.4636 -0.0241 0.0927 -0.0241 0.0208 -0.0241 0.0208
% 0.1205 0.4590 -0.1205 -0.2698 0.0909 0 0.0241 0 0.0241
% -0.4636 -0.1205 0.5674 0.1020 -0.1758 0.0241 0 0.0241 0
% -0.0241 -0.2698 0.1020 0.4590 -0.1205 0 -0.1205 0 0.0241
% 0.0927 0.0909 -0.1758 -0.1205 0.5674 -0.1205 0 0.0241 0
% -0.0241 0 0.0241 0 -0.1205 0.4590 -0.1205 -0.1891 0.0909
% 0.0208 0.0241 0 -0.1205 0 -0.1205 0.5674 0.1020 -0.3916
% -0.0241 0 0.0241 0 0.0241 -0.1891 0.1020 0.4590 -0.1205
% 0.0208 0.0241 0 0.0241 0 0.0909 -0.3916 -0.1205 0.5674
% 0.1446 0 -0.1446 -0.1947 0.1446 -0.3560 0.1446 0 -0.1446
% -0.1246 -0.1446 0 0.1446 -0.5563 0.1446 -0.1246 -0.1446 0
% Columns 19 through 20
% 0 0
% 0 0
% 0 0
% 0 0
% 0 0
% 0 0
% -0.1947 0.1446
% 0.1446 -0.5563
% -0.3560 0.1446
% 0.1446 -0.1246
% 0 -0.1446
% -0.1446 0
% -0.1947 0.1446
% 0.1446 -0.5563
% -0.3560 0.1446
% 0.1446 -0.1246
% 0 -0.1446
% -0.1446 0
% 1.1015 -0.2893
% -0.2893 1.3619
%__________________________ISOPARAMETRIC element stiffness matrix value(1)
%K(:,:,1) =
%
% 1.0e+006 *
%
% Columns 1 through 9
%
% 0.0406 0.0000 -0.0000 -0.0067 -0.0084 0.0067 0.0036 -0.0545 0.0036
% 0.0000 0.1161 -0.0057 -0.0000 0.0057 -0.0239 -0.0467 0.0101 0.0197
% -0.0000 -0.0057 0.0743 0.0000 -0.0153 0.0057 0.0066 0.0197 0.0065
% -0.0067 -0.0000 0.0000 0.0260 0.0067 -0.0054 0.0230 0.0023 -0.0545
% -0.0084 0.0057 -0.0153 0.0067 0.1149 -0.0604 -0.0102 0.0054 -0.0102
% 0.0067 -0.0239 0.0057 -0.0054 -0.0604 0.1421 0.0054 -0.0126 0.0054
% 0.0036 -0.0467 0.0066 0.0230 -0.0102 0.0054 0.4578 -0.1204 -0.0913
% -0.0545 0.0101 0.0197 0.0023 0.0054 -0.0126 -0.1204 0.5669 -0.0427
% 0.0036 0.0197 0.0065 -0.0545 -0.0102 0.0054 -0.0913 -0.0427 0.4583
% 0.0230 0.0103 -0.0467 0.0023 0.0054 -0.0126 -0.0538 -0.1129 -0.1204
% -0.0036 -0.0000 -0.1181 0.0545 0.0627 -0.0250 0.0323 -0.0240 -0.1615
% -0.0000 -0.0103 0.0467 -0.0413 -0.0283 0.0310 -0.0240 0.0923 0.1201
% -0.0036 0.0000 0.0591 -0.0230 -0.1145 0.0414 0.0325 -0.0241 0.0323
% 0.0000 -0.0103 -0.0197 0.0207 0.0492 -0.0312 -0.0241 0.0928 -0.0240
% 0.0323 -0.0197 -0.0066 0.0000 -0.0581 0.0492 0.0590 -0.0240 0.0594
% -0.0230 0.0924 0.0000 -0.0023 0.0414 -0.1822 -0.0240 0.0206 -0.0241
% -0.0646 0.0467 -0.0066 -0.0000 0.0389 -0.0283 -0.2956 0.1201 0.0590
% 0.0545 -0.1845 -0.0000 -0.0023 -0.0250 0.0947 0.1201 -0.1033 -0.0240
% 0.0001 -0.0000 0.0001 0.0000 0.0002 -0.0001 -0.1947 0.1446 -0.3560
% -0.0000 0.0002 0.0000 0.0000 -0.0001 0.0002 0.1446 -0.5563 0.1446
% Columns 10 through 18
% 0.0230 -0.0036 -0.0000 -0.0036 0.0000 0.0323 -0.0230 -0.0646 0.0545
% 0.0103 -0.0000 -0.0103 0.0000 -0.0103 -0.0197 0.0924 0.0467 -0.1845
% -0.0467 -0.1181 0.0467 0.0591 -0.0197 -0.0066 0.0000 -0.0066 -0.0000
% 0.0023 0.0545 -0.0413 -0.0230 0.0207 0.0000 -0.0023 -0.0000 -0.0023
% 0.0054 0.0627 -0.0283 -0.1145 0.0492 -0.0581 0.0414 0.0389 -0.0250
% -0.0126 -0.0250 0.0310 0.0414 -0.0312 0.0492 -0.1822 -0.0283 0.0947
% -0.0538 0.0323 -0.0240 0.0325 -0.0241 0.0590 -0.0240 -0.2956 0.1201
% -0.1129 -0.0240 0.0923 -0.0241 0.0928 -0.0240 0.0206 0.1201 -0.1033
% -0.1204 -0.1615 0.1201 0.0323 -0.0240 0.0594 -0.0241 0.0590 -0.0240
% 0.5657 0.1201 -0.4619 -0.0240 0.0923 -0.0241 0.0208 -0.0240 0.0206
% 0.1201 0.4575 -0.1198 -0.2692 0.0907 -0.0000 0.0240 0.0000 0.0241
% -0.4619 -0.1198 0.5655 0.1018 -0.1752 0.0240 -0.0000 0.0241 0.0000
% -0.0240 -0.2692 0.1018 0.4581 -0.1206 0.0002 -0.1202 -0.0000 0.0240
% 0.0923 0.0907 -0.1752 -0.1206 0.5670 -0.1202 0.0002 0.0240 -0.0000
% -0.0241 -0.0000 0.0240 0.0002 -0.1202 0.4584 -0.1206 -0.1886 0.0907
% 0.0208 0.0240 -0.0000 -0.1202 0.0002 -0.1206 0.5661 0.1018 -0.3908
% -0.0240 0.0000 0.0241 -0.0000 0.0240 -0.1886 0.1018 0.4574 -0.1198
% 0.0206 0.0241 0.0000 0.0240 -0.0000 0.0907 -0.3908 -0.1198 0.5657
% 0.1446 -0.0001 -0.1446 -0.1948 0.1447 -0.3561 0.1447 -0.0001 -0.1446
% -0.1246 -0.1446 -0.0000 0.1447 -0.5563 0.1447 -0.1248 -0.1446 -0.0002
% Columns 19 through 20
% 0.0001 -0.0000
% -0.0000 0.0002
% 0.0001 0.0000
% 0.0000 0.0000
% 0.0002 -0.0001
% -0.0001 0.0002
% -0.1947 0.1446
% 0.1446 -0.5563
% -0.3560 0.1446
% 0.1446 -0.1246
% -0.0001 -0.1446
% -0.1446 -0.0000
% -0.1948 0.1447
% 0.1447 -0.5563
% -0.3561 0.1447
% 0.1447 -0.1248
% -0.0001 -0.1446
% -0.1446 -0.0002
% 1.1014 -0.2893
% -0.2893 1.3618
