1 -- Glasow Haskell 0.403 : FINITE ELEMENT PROGRAM V2
    2 -- **********************************************************************
    3 -- *                                                                    *
    4 -- * FILE NAME : elemstif.hs            DATE : 13-3-1991                *
    5 -- *                                                                    *
    6 -- * CONTENTS : Computes element stiffness matrix.                      *
    7 -- *                                                                    *
    8 -- **********************************************************************
    9 
   10 -- Two dimensional Beam Element:
   11 --
   12 --     Left Node nl (xl,yl)
   13 --     Right Node nr (xr,yr)
   14 --     Material Property Parameter (ea,ei)
   15 --
   16 --
   17 --     Ke = [ T * K11 * T~  T * K12 * T~ ]
   18 --          [ T * K21 * T~  T * K22 * T~ ]
   19 --
   20 --          [ cos a   -sin a    0 ]
   21 --     T  = [ sin a    cos a    0 ]
   22 --          [   0        0      1 ]
   23 --
   24 --          [ ea/l       0       0       ]
   25 --     K11= [  0     12ei/l"3  6ei/l"2   ]
   26 --          [  0      6ei/l"2   4ei/l    ]
   27 --
   28 --          [ -ea/l       0       0      ]
   29 --     K12= [  0    -12ei/l"3  6ei/l"2   ]
   30 --          [  0      -6ei/l"2  2ei/l    ]
   31 --
   32 --     K21= K12~
   33 --
   34 --          [ ea/l       0       0       ]
   35 --     K22= [  0     12ei/l"3  -6ei/l"2  ]
   36 --          [  0      -6ei/l"2   4ei/l   ]
   37 --
   38 --     where
   39 --          l = length of the beam element
   40 --          sin a = (yr-yl) / l
   41 --          cos a = (xr-xl) / l
   42 
   43 module Elemstif( beam2d ) where         
   44 
   45 import Basics
   46 import Vector
   47 import Matrix
   48 import DB_interface
   49 
   50 beam2d :: (Array Int Int, Array Int Float) -> Int -> Mat Float
   51 
   52 beam2d str element =
   53         makemat (6,6) 
   54 
   55               ( \ (i,j) -> 
   56                 if      ( i==1 && j==1 ) || ( i==4 && j==4 ) then b1
   57                 else if ( i==2 && j==1 ) || ( i==5 && j==4 ) ||
   58                         ( i==1 && j==2 ) || ( i==4 && j==5 ) then b2
   59                 else if ( i==2 && j==2 ) || ( i==5 && j==5 ) then b3
   60                 else if ( i==1 && j==4 ) || ( i==4 && j==1 ) then -b1
   61                 else if ( i==2 && j==4 ) || ( i==5 && j==1 ) ||
   62                         ( i==1 && j==5 ) || ( i==4 && j==2 ) then -b2
   63                 else if ( i==2 && j==5 ) || ( i==5 && j==2 ) then -b3
   64                 else 0.0  )
   65         where
   66         a1      = ea / l
   67         b1      = a1 * c * c
   68         b2      = a1 * c * s
   69         b3      = a1 * s * s
   70         dx      = xr - xl
   71         dy      = yr - yl
   72         l       = sqrt (dx*dx + dy*dy)
   73         c       = dx / l
   74         s       = dy / l
   75         (ea,ei) = getmpro str (getemat str element)
   76         (nl,nr) = getenlr str element
   77         (xl,yl) = getnxy str nl
   78         (xr,yr) = getnxy str nr
   79