{VERSION 6 0 "IBM INTEL NT" "6.0" } {USTYLETAB {CSTYLE "Maple Input" -1 0 "Courier" 0 1 255 0 0 1 0 1 0 0 1 0 0 0 0 1 }{CSTYLE "2D Math" -1 2 "Times" 0 1 0 0 0 0 0 0 2 0 0 0 0 0 0 1 }{CSTYLE "2D Output" 2 20 "" 0 1 0 0 255 1 0 0 0 0 0 0 0 0 0 1 } {PSTYLE "Normal" -1 0 1 {CSTYLE "" -1 -1 "Times" 1 12 0 0 0 1 2 2 2 2 2 2 1 1 1 1 }1 1 0 0 0 0 1 0 1 0 2 2 0 1 }{PSTYLE "Maple Output" -1 11 1 {CSTYLE "" -1 -1 "Times" 1 12 0 0 0 1 2 2 2 2 2 2 1 1 1 1 }3 3 0 0 0 0 1 0 1 0 2 2 0 1 }} {SECT 0 {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 3113 "# Finte Element Tr uss Analysis Program. Text-only version\n# Written for ENGN1300 at Bro wn Univ. \n# Adapted from AF Bower's 2-D FEM code by JA Blume\n# resta rt:\n#\n#********** GIVE THE NAME OF THE INPUTFILE IN THE LINE BELOW * **********\n# Use frontslash/ in place of backslash \\ \n# Enclose the name in single backquotes`\n# \+ End the line with a colon:\n#\ninfilename:=`D:/Do cuments and Settings/blume/Desktop/tutorial.inp.txt`:\n#\n#*********** **************** THEN HIT RETURN *************************** \n# \n#** ********************************************************************* \n# load linear algebra maple package\n#\nwith(linalg):\n#\n# defin e an output file name\n#\nlifn:=length(infilename):\noutfilename:=cat( substring(infilename,1..lifn-8),`.out.txt`):\n#\n# ================= R ead data from the input file ==================\n#\n infile:=fopen (infilename,READ):\n#\n# Title:\n#\nruntitle := readline(infile):\n #\n# No. nodes and nodal coordinates\n#\nnnode := fscanf(infile,`%d `)[1]:\n coord := array(1..nnode,1..2):\n for i from 1 to nnode d o\n coord[i,1] := fscanf(infile,`%f`)[1]:\n coord[i,2] := \+ fscanf(infile,`%f`)[1]:\n od:\n#\n# No. elements connectivity, \+ and material data.\n# connect[i,1] and connect[i,2] hold the tw o joint numbers for the\n# ith member. connect[i,3] holds EA fo r that member \n#\n nelem := fscanf(infile,`%d`)[1]:\n connect := \+ array(1..nelem,1..3):\n for i from 1 to nelem do\n for j fro m 1 to 2 do\n connect[i,j] := fscanf(infile,`%d`)[1]:\n \+ od:\n connect[i,3] := fscanf(infile,`%e`)[1]:\n od:\n\n#\n# No. degrees of freedom with prescribed displacements; prescribed \+ \n# displacements\n# \n nfix := fscanf(infile,`%d`)[1]:\n fixno des := array(1..nfix,1..3):\n for i from 1 to nfix do\n#\n# fixno des[i,1] holds the number of the ith fixed node\n# fixnodes[i,2] holds the direction (1 or 2) \n# fixnodes[i,3] holds the valus of the presc ribed displacement\n#\n fixnodes[i,1] := fscanf(infile,`%d`)[1]: \n fixnodes[i,2] := fscanf(infile,`%d`)[1]:\n fixnodes[i,3 ] := fscanf(infile,`%e`)[1]:\n od:\n#\n# No. loaded nodes, wit h the loads\n#\n nload := fscanf(infile,`%d`)[1]:\n loads := arr ay(1..nload,1..3):\n#\n# Node number of ith loaded node is in loads [i,1]\n# F1 and F2 for the ith loaded node are in loads[i,2] and lo ads[i,3]\n#\n for i from 1 to nload do\n loads[i,1] := f scanf(infile,`%d`)[1]:\n loads[i,2] := fscanf(infile,`%e`)[1]: \n loads[i,3] := fscanf(infile,`%e`)[1]:\n od:\n fclo se(infile):\n# Print input data to the screen:\nprint(`number of nodes and elements:`,nnode, nelem);\nnodecoords:=array(1..nnode,1..3,\n \+ [seq([nn,coord[nn,1],coord[nn,2]],nn=1..nnode)]);\nprint(`coor ds: [node #, x,y]:`,nodecoords);\nelcon:=array(1..nelem,1..4, \+ [seq([nn,connect[nn,1],connect[nn,2],connect[nn,3]],nn=1..nelem)]); \nprint(`Connectivity: [el.#, node 1, node 2, EA]:`,elcon);\nprint(`[f ixed node #, dof, value]:`,fixnodes);\nprint(`[loaded node #, Px, Py]: `,loads);\n#unassign(`nodecoords`,`elcon`):" }}{PARA 11 "" 1 "" {XPPMATH 20 "6%%>number~of~nodes~and~elements:G\"\"$F$" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%+nodecoordsGK%'matrixG6#7%7%\"\"\"$\"\"!F,F+7%\" \"#$\"$2(!\"$F/7%\"\"$$\"$T\"!\"#F+Q(pprint06\"" }}{PARA 11 "" 1 "" {XPPMATH 20 "6$%7coords:~[node~#,~x,y]:GK%'matrixG6#7%7%\"\"\"$\"\"!F+ F*7%\"\"#$\"$2(!\"$F.7%\"\"$$\"$T\"!\"#F*Q(pprint06\"" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#>%&elconGK%'matrixG6#7%7&\"\"\"F*\"\"#$\"#5!\"\"7& F+F+\"\"$F,7&F0F0F*F,Q(pprint06\"" }}{PARA 11 "" 1 "" {XPPMATH 20 "6$% JConnectivity:~[el.#,~node~1,~node~2,~EA]:GK%'matrixG6#7%7&\"\"\"F)\" \"#$\"#5!\"\"7&F*F*\"\"$F+7&F/F/F)F+Q(pprint06\"" }}{PARA 11 "" 1 "" {XPPMATH 20 "6$%<[fixed~node~#,~dof,~value]:GK%'matrixG6#7%7%\"\"\"F)$ \"\"!F+7%F)\"\"#F*7%\"\"$F-F*Q(pprint06\"" }}{PARA 11 "" 1 "" {XPPMATH 20 "6$%9[loaded~node~#,~Px,~Py]:GK%'matrixG6#7#7%\"\"#$\"\"!F +$\"#5!\"\"Q(pprint06\"" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 4793 "#\n#================ ELEMENT STIFFNESS MATRIX ================== =\n#\nelstif := proc (xa,ya,xb,yb,EA) \n local elstif, nx, ny, leng th;\n#\n# Define the element length and direction vector \n#\n l ength := sqrt((xa-xb)*(xa-xb)+(ya-yb)*(ya-yb));\n# \n nx := (xb-xa )/length;\n ny := (yb-ya)/length;\n \n#\n# Define the element \+ stiffness\n#\n elstif := evalm((EA/length)*array(1..4,1..4,\n \+ [[nx*nx,nx*ny,-nx*nx,-nx*ny],\n [ny*nx,ny *ny,-ny*nx,-ny*ny],\n [-nx*nx,-nx*ny,nx*nx,nx*ny],\n \+ [-ny*nx,-ny*ny,ny*nx,ny*ny]]));\nend:\n#\n#========== =Assemble the global stiffness matrix==============\n#\nStif := array( 1..2*nnode,1..2*nnode):\n for i from 1 to 2*nnode do\n for j fro m 1 to 2*nnode do\n Stif[i,j] := 0.0:\n od:\n od:\n#\n# \+ Loop over all the elements\n#\nfor lmn from 1 to nelem do\n#\n# Set up the stiffness for the current element\n#\n a := connect[lmn,1]: \n b := connect[lmn,2]:\n EA := connect[lmn,3]:\n k := array( 1..4,1..4):\n k:=elstif(coord[a,1],coord[a,2],coord[b,1],coord[b,2] ,EA):\n#\n# Add the current element stiffness to the global stiffnes s\n#\n for i from 1 to 2 do\n for ii from 1 to 2 do\n f or j from 1 to 2 do\n for jj from 1 to 2 do\n rw : = 2*(connect[lmn,i]-1)+ii:\n cl := 2*(connect[lmn,j]-1)+jj: \n Stif[rw,cl] := Stif[rw,cl] + k[2*(i-1)+ii,2*(j-1 )+jj]:\n od:\n od:\n od:\n od:\n od:\n#\n # ==================== Assemble global residual vector ============\n# \n# Define the residual\n#\n resid := array(1..2*nnode):\n for i from 1 to 2*nnode do\n resid[i] := 0.:\n od:\n#\n# Loop l oaded noads\n#\n for i from 1 to nload do\n a := loads[i,1]:\n f1 := loads[i,2]:\n f2 := loads[i,3]:\n resid[2*a-1] : = resid[2*a-1] + f1:\n resid[2*a] := resid[2*a] + f2:\n od: \n#\n# Modify the global stiffness and residual to include constrain ts\n#\n for i from 1 to nfix do\n rw := 2*fixnodes[i,1]+fixno des[i,2]-2:\n for j from 1 to 2*nnode do\n resid[j] := r esid[j] - Stif[j,rw]*fixnodes[i,3]:\n Stif[rw,j] := 0:\n \+ Stif[j,rw] := 0:\n od;\n Stif[rw,rw] := 1.0:\n \+ resid[rw] := fixnodes[i,3]:\n od:\n#\n# ================== So lve the FEM equations ===================\n#\nu := linsolve(Stif,resid ):\n#\n# Define a procedure to calculate the element strains\n#\ne lstrain := proc (xa,ya,xb,yb,uax,uay,ubx,uby) \n local elstrain, nx , ny, length, uel;\n#\n# length and direction vector for the element \n#\nlength := sqrt((xa-xb)*(xa-xb)+(ya-yb)*(ya-yb));\n# \n nx := \+ (xb-xa)/length;\n ny := (yb-ya)/length;\n#\n# Element displaceme nt vector\n#\n uel := array(1..4,[uax,uay,ubx,uby]);\n#\n# Eleme nt strains\n#\n elstrain := ((ubx-uax)*nx+(uby-uay)*ny)/length;\nend :\n# \n# compute element forces and strains:\n#\nfor lmn from 1 \+ to nelem do\n a := connect[lmn,1]:\n b := connect[lmn,2]:\n \+ EA := connect[lmn,3]:\n xa := coord[a,1]:\n ya := coord [a,2]:\n xb := coord[b,1]:\n yb := coord[b,2]:\n uxa := u[2*a-1]:\n uya := u[2*a]:\n uxb := u[2*b-1]:\n uyb := u[2*b]:\n Strain[lmn] := elstrain(xa,ya,xb,yb,uxa,uya,uxb,uyb): \n Force[lmn]:=EA*Strain[lmn]; \n od: \n#\n# Print nodal displacements, element strains and forces to a file\n#\n\n\noutfile:= fopen(outfilename,WRITE):\nfprintf(outfile, `%s\\n`,runtitle):\nfprint f(outfile,`%s %3d %3d\\n` ,`Number of Nodes and Elements:`,nnode, nele m):\nfprintf(outfile,`%s\\n`,`Node# x1 x2`):\n for i from 1 t o nnode do\n fprintf(outfile,`%5d %9.4f %9.4f\\n`,\n \+ i,coord[i,1],coord[i,2]):\n od:\nfprintf(outfile,`%s\\n`,`El# \+ Node 1 Node 2 EA`):\n for i from 1 to nelem do\n fprintf (outfile,`%4d %10d %10d %9.4f\\n`,\n i,connect[i,1], connect[i,2],connect[i,3]):\n od:\nfprintf(outfile,`%s\\n`,`Fixed No de# DOF Value`):\n for i from 1 to nfix do\n fprintf(outfile, `%9d %12d %10.4f\\n`,\n fixnodes[i,1],fixnodes[i,2],f ixnodes[i,3]):\n od:\nfprintf(outfile,`%s\\n`,`Loaded Node# P1 \+ P2`):\n for i from 1 to nload do\n fprintf(outfile,`%14d %1 3.4f %9.4f\\n`,\n loads[i,1],loads[i,2],loads[i,3]): \n od:\nfprintf(outfile,`%s\\n`,` `):\nfprintf(outfile,`%s\\n`,`**** *****************************`):\nfprintf(outfile,`%s\\n`,` `):\n#\nfp rintf(outfile,`%s\\n`,`Nodal Displacements:`):\nfprintf(outfile,`%s\\n `,` Node u1 u2`):\n for i from 1 to nnode do\n fprintf (outfile,`%3d %8.4f %8.4f\\n`,i,u[2*i-1],u[2*i]):\n od:\nfprintf(out file,`%s\\n`,`Strains and Forces:`):\nfprintf(outfile,`%s\\n`,`Element Strain Force`):\n for i from 1 to nelem do\n fprintf(out file,`%3d %9.4f %9.4f\\n`,i,Strain[i],Force[i]):\n od:\nfclose(outfi le):\nprint('finished!');\n" }}{PARA 11 "" 1 "" {XPPMATH 20 "6#-%*fact orialG6#%)finishedG" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}} {EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 " " {MPLTEXT 1 0 0 "" }}}{EXCHG {PARA 0 "> " 0 "" {MPLTEXT 1 0 0 "" }}}} {MARK "0 0 0" 674 }{VIEWOPTS 1 1 0 1 1 1803 1 1 1 1 }{PAGENUMBERS 0 1 2 33 1 1 }