Page 4 - Demo

Planar Truss Examplefor ComrelAdd-onRCPConsult, 2021-2025Page 4enddoforn fromupperbound(cc)by1to1dofc :=Remove(fc,cc[n]):enddo## Global stiffness matrixK :=Matrix(ndof,ndof):# stiffness matrixT:=Array(1..nfe):# transformation matrixforn tonfe do# sin and cos of the inclinationc :=cos(theta[n]):s :=sin(theta[n]):# coordinate transformation matrix for the bar nT[n]:=Matrix([[c,s,0,0],[s,c,0,0],[0,0,c,s],[0,0,s,c]]):# local stiffness matrix for the bar nKloc :=k[n]*Matrix([[1,0,-1,0],[0,0,0,0],[-1,0,1,0],[0,0,0,0]]):# global stiffness matrix assemblyidx:=[entries(LM[n,1..],'nolist')]:K[idx,idx]:=K[idx,idx]+Array(evalm(transpose(T[n]).Kloc.T[n])):enddo## Solve the system of equationsec :=[entries(cc,'nolist')]:ed :=[entries(dd,'nolist')]:ac :=Vector([0,0,0]):Kcc :=SubMatrix(K,ec,ec):Kcd :=SubMatrix(K,ec,ed):Kdc :=SubMatrix(K,ed,ec):Kdd :=SubMatrix(K,ed,ed):ad :=LinearSolve(Kdd,fc Multiply(Kdc,ac)):qd :=Multiply(Kcc,ac)+Multiply(Kcd,ad):## Final displacementsa :=Vector(ndof):q :=Vector(ndof):a[ec]:=ac:q[ec]:=qd:a[ed]:=ad:N :=Vector(nfe):## Axial loads in barsforn tonfe doN[n]:=Multiply(Multiply(Multiply(k[n],Array([1,0,1,0])),T[n]),a(LM[n])):enddo## VisualisationifStrurelPlotthen# begin of plot blockifStrurelMode=0orStrurelMode=2thenglobalStrurelPlotCount:ifStrurelMode=0thentit:=\Deterministic Solution\# 0 -initial deterministic solutionelifStrurelMode=2thentit:=\Stochastic Solution\# 2 -final stochastic solutionendif# Get a screen resolutiongetmetric:=define_external(\t::integer[4],RETURN::integer[4],LIB=\## Plot Function to be used in MapletTrussPlot:=proc(np)# X and Y coordinates in mXY :=Matrix([[0,0],[4,0],[8,0],[12,0],[16,0],[20,0],[24,0],[22,2],[18,2],[14,2],[10,2],[6,2],[2,2]]):ifnp =1thenPlotT:=textplot([12.5,6.,\BarL :=Array(1..nfe):# rods as linesBarR :=Array(1..nfe):# rod numbersNodN :=Array(1..nnp):# nodal numbers# Plot initial state of trussfori tonfe dox1:=XY[IEN[i,1],1]:y1:=XY[IEN[i,1],2]:x2:=XY[IEN[i,2],1]:y2:=XY[IEN[i,2],2]:BarL[i]:=line([x1,y1],[x2,y2],color=blue):# blue color for initial stateBarR[i]:=textplot([(x1+x2)/2,(y1+y2)/2,sprintf(\# rod number in the middle of rod lengthenddofori tonnp doNodN[i]:=textplot([XY[i,1],XY[i,2],sprintf(\# nodal numbersenddoBarD :=Array(1..nfe):# deformed rods as lines# Plot deformed state of truss -scale factor 2fori tonfe don1:=IEN[i,1]:n2:=IEN[i,2]:x1:=XY[n1,1]+a[n1*21]*2:y1:=XY[n1,2]+a[n1*2]*2:x2:=XY[n2,1]+a[n2*21]*2:y2:=XY[n2,2]+a[n2*2]*2:BarD[i]:=line([x1,y1],[x2,y2],color=red):# red color for deformed state