identification division.
       program-id. elastic.
      *  aus D.Nowottny p. 77
      * Mathematik am Computer
       author. "JD".
       date-written. 24.9.2006.
       date-compiled.
       data division.
       working-storage section.
         78 anz pic 9(4) value 100.
         78 n   pic 9(4) value 100.
         78 n2  pic 9(4) value 200.
         78 m   pic 9(4) value 100.
         78 m2   pic 9(4) value 200.
         01 occurs n.
           05 mm   pic 9(4)V99 occurs n.
         77 i   pic 9(4).
         77 j   pic 9(4).
         77 k   pic 9(4).
         77 te  pic 9(4).
         77 c   pic 9(4).
         77 hi  pic 9(4).
         77 l   pic 9(4).
         77 h   pic 9(4).
         77 y   pic 9(4).
         77 nenner   pic 9(4).
         01 diff     pic 9(4) occurs n2.
       linkage section.
         01 p   pic 9(4)V99 occurs n2.
         01 f   pic 9(4)V99 occurs m2.
         01 a   pic 9(4)V99 occurs n2.
         01 x   pic 9(4)V99 occurs m2.
       procedure division using p f a x.
      *  Tragwerk auf 2 Stützen
      * p= koordinaten f=Kräfte, a=Topologie, x=Verschiebungen
         perform varying i from 1 until i > n
           perform varying j from 1 until i > n
             move 0 to mm(i,j)
           end-perform
         end-perform
         perform varying k from 1 until k > anz
           move a(1,k) to i
           move a(2,k) to j
           move a(3,k) to c
           perform varying l from 1 until l > n
             compute diff(l) = p(j)-p(i)
           end-perform
           move 0 to nenner
           perform varying l from 1 until l > n
             compute nenner=nenner+diff(l)*diff(l)
           end-perform
      *    dyadisches produkt
           perform varying i from 1 until i > n
             perform varying j from 1 until i > n
               compute hi(i,j) = c/nenner * x(i)*y(j)
             end-perform
           end-perform
      *    [ Hi -Hi; -Hi Hi]t
           perform varying i from 1 until i > n
             perform varying j from 1 until i > n
               compute h(i,j)     = Hi(i,j)
               compute h(i+n,j)   = - Hi(i,j)
               compute h(i,j+n)   = Hi(i,j)
               compute h(i+n,j+n) = Hi(i,j)
             end-perform
           end-perform
           compute te=2*i-1
           compute m(te,te)=M(te,te)+H 
         end-perform.
         perform varying i from 1 until i > n
           perform varying j from 1 until i > n
             if (i>n) or (j>n) then
               move 0 to m(i,j)
        move 0 to Hi(i,j)
             end-perform
           end-perform
         call agev using X M F
      *  X=Reshape(X,2,m)
         call reshape using X 2 M
      *  nun sind die Verschiebungen berechnet
        end-program elastic.