Here is my code along with the frame structure with Edof....
I am having a value for the moment (supposed to be null) at the join ball located at the top center of the structure (x=4.5m ; y=10m)
%----------------------------------------------------------------
echo on
%----- Topology -------------------------------------------------
Edof=[1 1 2 3 7 8 9;
2 4 5 6 16 17 18;
3 7 8 9 10 11 12;
4 10 11 12 13 14 15;
5 13 14 15 16 17 18;
6 7 8 9 19 20 21;
7 16 17 18 25 26 27;
8 19 20 21 22 23 24;
9 22 23 24 25 26 27]
%----- Stiffness matrix K and load vector f ---------------------
K=zeros(27); f=zeros(27,1);
f(11)=-20e3;
f(14)=-20e3;
%----- Element stiffness and element load matrices -------------
Eh=36e6;
Ah=0.16;
Ih=2.133e-3;
Ev=70e6;
Av=0.1;
Iv=1.333e-3;
eph=[Eh Ah Ih];
epv=[Ev Av Iv];
ex1=[0 0]; ex2=[9 9]; ex3=[0 3];
ey1=[0 5]; ey2=[0 5]; ey3=[5 5];
ex4=[3 6]; ex5=[6 9]; ex6=[0 0];
ey4=[5 5]; ey5=[5 5]; ey6=[5 10];
ex7=[9 9]; ex8=[0 4.5]; ex9=[4.5 9];
ey7=[5 10]; ey8=[10 10]; ey9=[10 10];
eq1=[0 0];
eq2=[0 0];
eq3=[0 0];
eq4=[0 0];
eq5=[0 0];
eq6=[0 0];
eq7=[0 0];
eq8=[0 -20e3];
eq9=[0 -20e3];
Ke1=beam2e(ex1,ey1,epv);
Ke2=beam2e(ex2,ey2,epv);
Ke3=beam2e(ex3,ey3,eph);
Ke4=beam2e(ex4,ey4,eph);
Ke5=beam2e(ex5,ey5,eph);
Ke6=beam2e(ex6,ey6,epv);
Ke7=beam2e(ex7,ey7,epv);
[Ke8,fe8]=beam2e(ex8,ey8,eph,eq8);
[Ke9,fe9]=beam2e(ex9,ey9,eph,eq9);
%----- Assemble Ke into K ---------------------------------------
K=assem(Edof(1,:),K,Ke1);
K=assem(Edof(2,:),K,Ke2);
K=assem(Edof(3,:),K,Ke3);
K=assem(Edof(4,:),K,Ke4);
K=assem(Edof(5,:),K,Ke5);
K=assem(Edof(6,:),K,Ke6);
K=assem(Edof(7,:),K,Ke7);
[K,f]=assem(Edof(8,:),K,Ke8,f,fe8);
[K,f]=assem(Edof(9,:),K,Ke9,f,fe9);
%----- Solve the system of equations and compute reactions ------
bc=[1 0; 2 0; 3 0; 4 0; 5 0; 24 0; 6 0];
[a,r]=solveq(K,f,bc)
%----- Section forces -------------------------------------------
Ed=extract(Edof,a);
es1=beam2s(ex1,ey1,epv,Ed(1,:),eq1,21)
es2=beam2s(ex2,ey2,epv,Ed(2,:),eq2,21)
es3=beam2s(ex3,ey3,eph,Ed(3,:),eq3,21)
es4=beam2s(ex4,ey4,eph,Ed(4,:),eq4,21)
es5=beam2s(ex5,ey5,eph,Ed(5,:),eq5,21)
es6=beam2s(ex6,ey6,epv,Ed(6,:),eq6,21)
es7=beam2s(ex7,ey7,epv,Ed(7,:),eq7,21)
es8=beam2s(ex8,ey8,eph,Ed(8,:),eq8,21)
es9=beam2s(ex9,ey9,eph,Ed(9,:),eq9,21)
%----- Draw deformed frame ---------------------------------------
figure(1)
plotpar=[2 1 0];
eldraw2(ex1,ey1,plotpar);
eldraw2(ex2,ey2,plotpar);
eldraw2(ex3,ey3,plotpar);
eldraw2(ex4,ey4,plotpar);
eldraw2(ex5,ey5,plotpar);
eldraw2(ex6,ey6,plotpar);
eldraw2(ex7,ey7,plotpar);
eldraw2(ex8,ey8,plotpar);
eldraw2(ex9,ey9,plotpar);
sfac=scalfact2(ex7,ey7,Ed(7,:),0.1);
plotpar=[1 2 1];
eldisp2(ex1,ey1,Ed(1,:),plotpar,sfac);
eldisp2(ex2,ey2,Ed(2,:),plotpar,sfac);
eldisp2(ex3,ey3,Ed(3,:),plotpar,sfac);
eldisp2(ex4,ey4,Ed(4,:),plotpar,sfac);
eldisp2(ex5,ey5,Ed(5,:),plotpar,sfac);
eldisp2(ex6,ey6,Ed(6,:),plotpar,sfac);
eldisp2(ex7,ey7,Ed(7,:),plotpar,sfac);
eldisp2(ex8,ey8,Ed(8,:),plotpar,sfac);
eldisp2(ex9,ey9,Ed(9,:),plotpar,sfac);
axis([-5 15 -1 15]);
pltscalb2(sfac,[1e-2 5 0]);
axis([-5 15 -1 15]);
title('displacements')
%----- Draw normal force diagram --------------------------------
figure(2)
plotpar=[2 1];
sfac=scalfact2(ex3,ey3,es3(:,1),0.2);
eldia2(ex1,ey1,es1(:,1),plotpar,sfac);
eldia2(ex2,ey2,es2(:,1),plotpar,sfac);
eldia2(ex3,ey3,es3(:,1),plotpar,sfac);
eldia2(ex4,ey4,es4(:,1),plotpar,sfac);
eldia2(ex5,ey5,es5(:,1),plotpar,sfac);
eldia2(ex6,ey6,es6(:,1),plotpar,sfac);
eldia2(ex7,ey7,es7(:,1),plotpar,sfac);
eldia2(ex8,ey8,es8(:,1),plotpar,sfac);
eldia2(ex9,ey9,es9(:,1),plotpar,sfac);
axis([-5 15 -1 15]);
pltscalb2(sfac,[3e4 5 0]);
title('normal force')
%----- Draw shear force diagram ---------------------------------
figure(3)
plotpar=[2 1];
sfac=scalfact2(ex5,ey5,es5(:,2),0.2);
eldia2(ex1,ey1,es1(:,2),plotpar,sfac);
eldia2(ex2,ey2,es2(:,2),plotpar,sfac);
eldia2(ex3,ey3,es3(:,2),plotpar,sfac);
eldia2(ex4,ey4,es4(:,2),plotpar,sfac);
eldia2(ex5,ey5,es5(:,2),plotpar,sfac);
eldia2(ex6,ey6,es6(:,2),plotpar,sfac);
eldia2(ex7,ey7,es7(:,2),plotpar,sfac);
eldia2(ex8,ey8,es8(:,2),plotpar,sfac);
eldia2(ex9,ey9,es9(:,2),plotpar,sfac);
axis([-5 15 -1 15]);
pltscalb2(sfac,[3e4 5 0]);
title('shear force')
%----- Draw moment diagram --------------------------------------
figure(4)
plotpar=[2 1];
sfac=scalfact2(ex4,ey4,es1(:,3),0.1);
eldia2(ex1,ey1,es1(:,3),plotpar,sfac);
eldia2(ex2,ey2,es2(:,3),plotpar,sfac);
eldia2(ex3,ey3,es3(:,3),plotpar,sfac);
eldia2(ex4,ey4,es4(:,3),plotpar,sfac);
eldia2(ex5,ey5,es5(:,3),plotpar,sfac);
eldia2(ex6,ey6,es6(:,3),plotpar,sfac);
eldia2(ex7,ey7,es7(:,3),plotpar,sfac);
eldia2(ex8,ey8,es8(:,3),plotpar,sfac);
eldia2(ex9,ey9,es9(:,3),plotpar,sfac);
axis([-5 15 -1 15]);
pltscalb2(sfac,[10e4 5 0]);
title('moment')
%------------------------ end -----------------------------------
echo off
React
Typescript
Django
Laravel
Facebook
Microsoft
Google
Alibaba
D3
Tencent