Elastic beam solver (engineers only ;))
Posted: Fri Jun 12, 2015 8:56 am
For an elastic (prismatic) beam of certain dimensions and for various loads, results such as bending stress, deflection and displacement can be calculated using a formula for each result.
In this example, a beam of steel is fixed on one end and loaded with a moment, a force and an uniform load on the entire beam.
Three formula are available to calculate the bending stress at the root of the beam, the deflection at the free end of the beam, and the axial force to cause buckling of the beam.
The system contains 10 variables (Young's modulus of the chosen material, 3 dimensions of the beam, 3 types of external loads, and 3 primary results).
However, if you want to be able to calculate any of these variables, the other ones being given, those formulas have to be rewritten and programmed for each variable, hence the program will contain 3x10=30 formulas (and of course a selection mechanism to choose the correct formula and the variable to be calculated).
This little program does the trick with only the 3 basic formulas. No need to rewrite each formula for every variable in it.
By using indirect adressing, the user interface need not be adapted to a specific problem.
First choose one of the 3 calculation groups. The relevant variables will then be displayed, together with their initial values (some of wich may be best guesses).
Enter values of the variables in the input fields on the right side (don't forget the "return" to confirm the entry), and press a variable name on the left side to have that variable being calculated and displayed in the field on the right side.
In this example I've used old-fashioned units kgf (= 9,81 N) and cm. Other units can be used when entering the (initial) values, but consistancy is needed of course.
In this example, a beam of steel is fixed on one end and loaded with a moment, a force and an uniform load on the entire beam.
Three formula are available to calculate the bending stress at the root of the beam, the deflection at the free end of the beam, and the axial force to cause buckling of the beam.
The system contains 10 variables (Young's modulus of the chosen material, 3 dimensions of the beam, 3 types of external loads, and 3 primary results).
However, if you want to be able to calculate any of these variables, the other ones being given, those formulas have to be rewritten and programmed for each variable, hence the program will contain 3x10=30 formulas (and of course a selection mechanism to choose the correct formula and the variable to be calculated).
This little program does the trick with only the 3 basic formulas. No need to rewrite each formula for every variable in it.
By using indirect adressing, the user interface need not be adapted to a specific problem.
First choose one of the 3 calculation groups. The relevant variables will then be displayed, together with their initial values (some of wich may be best guesses).
Enter values of the variables in the input fields on the right side (don't forget the "return" to confirm the entry), and press a variable name on the left side to have that variable being calculated and displayed in the field on the right side.
In this example I've used old-fashioned units kgf (= 9,81 N) and cm. Other units can be used when entering the (initial) values, but consistancy is needed of course.
Code: Select all
' equations solver / system analyzer
' programmed by Henko in smartbasic, june 2015
' in yellow: specific data and coding for "elastic beam" example
'
option base 0
dim varname$(16),value(16),status(16)
dim equname$(6),equ(6,11)
graphics ! graphics clear ! draw color 0,0,0
read system_title$
read nequ,nvar
for i=1 to nequ
read equname$(i),equ(i,0)
for j=1 to equ(i,0) ! read equ(i,j) ! next j
next i
for i=1 to nvar ! read void,varname$(i),value(i) ! next i
'y'
data "elastic beam" ' name of problem
data 3,10 ' number of formulas and variables
data "stress",8,1,2,3,4,5,6,7,8 ' variables in the stress formula
data "displacement",8,1,2,3,4,5,6,7,9 ' same in the displacement
data "buckling",5,1,2,3,4,10 ' same in the buckling formula
data 1,"elast. modulus",2100000
data 2,"beam width",1
data 3,"beam height",10
data 4,"beam length",100
data 5,"bending moment",0
data 6,"force",1000
data 7,"load",0
data 8,"stress",10
data 9,"displacement",1
data 10,"axial force",100
''
y_eq=10 ! y_var=y_eq+50*nequ+20 ! old_equ=0
draw font size 30 ! draw text system_title$ at 50,10 ! draw font size 20
equ_disp()
loop: slowdown
for i=1 to nequ ! if bp("eq" & i) then cur_equ=i ! next i
if cur_equ <> old_equ then
if old_equ>0 then var_del(old_equ)
var_disp(cur_equ) ! old_equ=cur_equ
end if
for i=1 to equ(cur_equ,0)
fn$="var" & i
if field_changed(fn$) then
value(equ(cur_equ,i))=field_text$(fn$) ! break
end if
if bp(fn$) then field fn$ text solve(value,equ,cur_equ,i)
next i
goto loop
end
def equi(value(),equ(,),nre) ' coding the equations
'y'
imom=value(2)*value(3)^3/12 ' moment of inertia
on nre goto eq1,eq2,eq3
eq1: ' calculate bending stress
res=value(equ(nre,5))+value(equ(nre,6))*value(equ(nre,4))
res+=.5*value(equ(nre,7))*value(equ(nre,4))^2
return .5*res*value(equ(nre,3))/imom-value(equ(nre,8))
eq2: ' calculate vertical displacement of beam's end
res=value(equ(nre,5))/2+value(equ(nre,6))*value(equ(nre,4))/3
res+=value(equ(nre,7))*value(equ(nre,4))^2/8
res=res*value(equ(nre,4))^2/value(equ(nre,1))/imom
return res-value(equ(nre,8))
eq3: ' calculate buckling force (Euler)
if value(2)<value(3) then imom=value(3)*value(2)^3/12
res=2.4674*value(equ(nre,1))*imom/value(equ(nre,4))^2
return res-value(equ(nre,5))
''
end def
def equ_disp ' display equation buttons
for i=1 to .nequ
button "eq" & i text .equname$(i) at 10,.y_eq+50*i size 300,30
next i
end def
def var_disp(nr) ' display variable name-buttons and value-fields
for i=1 to .equ(nr,0)
button "var" & i text .varname$(.equ(nr,i)) at 10,.y_var+50*i size 150,30
field "var" & i text .value(.equ(nr,i)) at 170,.y_var+50*i size 140,30
next i
end def
def var_del(nr) ' remove variable buttons and fields
for i=1 to .equ(nr,0)
button "var" & i delete ! field "var" & i delete
next i
end def
def solve(value(),equ(,),nre,nrv) ' newton-raphson iteration
imax=10 ! eps=.001
for i=1 to imax
delta=equi(value,equ,nre)/deriv(value,equ,nre,nrv)
value(equ(nre,nrv))-=delta
if abs(delta)<eps then return value(equ(nre,nrv))
next i
return -999
end def
def deriv(value(),equ(,),nre,nrv) ' partial derivate of nre w.r.t. nrv
fx=equi(value,equ,nre)
dx=.001*value(equ(nre,nrv)) ! if abs(dx)<.001 then dx=.001
value(equ(nre,nrv))+=dx
fxplus=equi(value,equ,nre)
value(equ(nre,nrv))-=dx
return (fxplus-fx)/dx
end def
def bp(b$) = button_pressed(b$)
