Elastic beam solver (engineers only ;))

Post Reply
Henko
Posts: 814
Joined: Tue Apr 09, 2013 12:23 pm
My devices: iPhone,iPad
Windows
Location: Groningen, Netherlands
Flag: Netherlands

Elastic beam solver (engineers only ;))

Post by Henko »

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.

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$)

Operator
Posts: 138
Joined: Mon May 06, 2013 5:52 am

Re: Elastic beam solver (engineers only ;))

Post by Operator »

Very, very interesting concept and nice
compact code.

Would you mind posting the main loop in the
"expanded" version, without "!". Every time I
try to rearrange the code without "!" I get it
broken...
Thx.


iPhone 4 users, here is a micro mod to get all
fields on our tiny display : )

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
'
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+30*nequ+20 ! old_equ=0
DRAW FONT SIZE 30 ! DRAW TEXT system_title$ AT 50,1 ! 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
'
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+35*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+35*i SIZE 150,30
  FIELD "var" & i TEXT .VALUE(.equ(nr,i)) AT 170,.y_var+35*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$)

Henko
Posts: 814
Joined: Tue Apr 09, 2013 12:23 pm
My devices: iPhone,iPad
Windows
Location: Groningen, Netherlands
Flag: Netherlands

Re: Elastic beam solver (engineers only ;))

Post by Henko »

Hi,
Here the program with the main program without the ! tokens.

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$)

Henko
Posts: 814
Joined: Tue Apr 09, 2013 12:23 pm
My devices: iPhone,iPad
Windows
Location: Groningen, Netherlands
Flag: Netherlands

Re: Elastic beam solver (engineers only ;))

Post by Henko »

Operator wrote:Very, very interesting concept and nice
compact code.

iPhone 4 users, here is a micro mod to get all
fields on our tiny display : )
The compactness illustrates how powerfull this Basic has become. Mr. K. deserves a big compliment for it, but we don't want to make him vain, do we? ;)

As for the output on the iPhone screens, even in your mod, the "keyboard" may cover the inputfield. The solution would be to use the (scrollable!) "list" object in stead of "buttons" and "fields".

User avatar
Mr. Kibernetik
Site Admin
Posts: 4786
Joined: Mon Nov 19, 2012 10:16 pm
My devices: iPhone, iPad, MacBook
Location: Russia
Flag: Russia

Re: Elastic beam solver (engineers only ;))

Post by Mr. Kibernetik »

Henko wrote:The compactness illustrates how powerfull this Basic has become. Mr. K. deserves a big compliment for it, but we don't want to make him vain, do we? ;)
Thanks! :D

Post Reply