Page 3 - Demo

Planar TrussExamplefor ComrelAdd-onRCPConsult, 2019-2025Page 3## Global flags to manage plot facilities of engineStrurelPlot=TrueStrurelPlotName='PlanarTruss_';StrurelPlotType='.png';StrurelPlotMode=3## Code based on Max Ehre work https://github.com/maxehre/polynomial_surrogates## Used with his kind permission as version from 2017 that adapted to Strurel.importnumpyasnp # https://numpy.orgifStrurelPlot:importmatplotlib.pyplotasplt # https://matplotlib.orgimportmatplotlib.cmascmximportmatplotlibasmplimportctypesdefPlanar_Truss(XP):r'''============================================================================Planar Truss Example for Python============================================================================Inputs:X = [A1, A2, E1, E2, P1, P2, P3, P4, P5, P6, u_y]E1,E2: Youngs modulus of horizontal and inclined bars resp. lognormally distributedA1,A2: cross section of horizontal and inclined bars resp. lognormally distributedP1P6: negative vertical forces attacking in nodes 813. Gumbel distributedu_y: maximal allowed vertical deflection in node 4. constantOutput:uout: vertical truss deflection at bottom center (N04)============================================================================Truss description taken from Lee, S.H., Kwak, B.M. (2006).Response surface augmented moment method for efficient reliability analysis.Structural Safety 28(3), 261 272.https://www.sciencedirect.com/science/article/abs/pii/S0167473005000421Based on Diego Andr%u00e9s Alvarez Mar%u00edn work https://github.com/diegoandresalvarez/elementosfinitos/tree/master/codigo/repaso_matricial/cercha_2dN: NodesR: Rods============================================================================N13__R4___N12__R8___N11__R12__N10__R16__N09__R20__N08/\\/\\/\\/\\/\\/\\/ \\/ \\/ \\/ \\/ \\/ \\R1 R3 R5 R7 R9 R11 R13 R15 R17 R19 R21 R23/ \\/ \\/ \\/ \\/ \\/ \\/___R2___\\/___R6___\\/__R10___\\/__R14___\\/__R18___\\/__R22___\\N01 N02 N03 N04 N05 N06 N07||\\/uout => u_y @ N04============================================================================'''#### Planar truss model#### Vector inputs order from Stochastic ModelA1=XP[0]# cross section of horizontal barsA2=XP[1]# cross section of inclined barsE1=XP[2]# Youngs modulus of horizontal barsE2=XP[3]# Youngs modulus of inclined barsP=XP[4:10]# negative vertical forces attacking in nodes 813. (1x6 vector)## Element, nodes and dofs association## IEN: connectivity matrix, nfe x nodes -bar 1 has nodes 1 and 13, bar 2 has nodes 1 and 2 ...IEN =np.array([[1,13],[1,2],[13,2],[13,12],[2,12],[2,3],[12,3],[12,11],[3,11],[3,4],[11,4],[11,10],[4,10],[4,5],[10,5],[10,9],[5,9],[5,6],[9,6],[9,8],[6,8],[6,7],[8,7]])-1## Deterministic rod propertiesgrid =4# size of grid side in [m]ng =12# count of grid cellsang =np.arctan2(grid,grid)# inclination angle of the truss [rad]theta =np.tile([ang,0,-ang,0],int(ng/2))# inclination angle [rad]theta =np.delete(theta,-1)leng =np.tile([grid/np.sqrt(2),grid],ng)# bar length [m]leng =np.delete(leng,-1)## FEM constantsned =2# number of dof per nodennp =13# number of nodal pointsnfe =np.shape(IEN)[0]# number of barsndof =ned*nnp # number of degrees of freedom (dof)## Dof: degrees of freedom (rows are the nodes, cols are dofs)dof =np.array(range(ndof)).reshape(nnp,ned)## Stochastic rod propertiesarea =np.tile([A2,A1],ng)# bar cross sectional area [m2]area =np.delete(area,-1)E =np.tile([E2,E1],ng)# Young's modulus [N/m^2]E =np.delete(E,-1)## Material propertiesk =E*area/leng # stiffness of each bar## Boundary conditions (supports)cc =np.array([1,2,14])1# dof fixed/supportsdd =np.array(list(set(range(ndof))-set(cc)))# dof free## Boundary conditions (applied forces)fc =np.zeros(ndof)fc[15::2]=np.array(P)# negative vertical forces [N]fc =np.delete(fc,cc)fc =np.expand_dims(fc,axis=1)## Global stiffness matrixK =np.zeros((ndof,ndof))# stiffness matrix