!  木高研の試験体(半断面)のトラス部分。
!  上弦材、下弦材の断面換算はなし。
!  上弦材、下弦材のPC鋼棒の位置で、x軸方向を拘束。
         implicit real*8(a-h,o-z)
         dimension x(10000),y(10000),z(10000)
!  節点座標：x(節点番号),y(節点番号),z(節点番号)
        print'("**")'
        print'("** Structure: Pratt truss")'
        print'("** Test objective: B32 elements.")'
        print'("** 木高研試験体(トラス部分)")'
        print'("**")'
! <節点番号ごとの節点座標>
!  ※プラットトラスについて
        print'("** トラス部材の節点番号と座標")'
        print'("*NODE, NSET=Truss")'
!  主構高(y方向)：h(m)
!  支間長(z方向)：ell(m)
         h=380.d-3
         ell=4900.d-3
!  部材の断面(長方形へ換算した断面)
!  (1) 上弦材(x方向:bj、y方向:hj)
         bj=6.d-3
         hj=163.d-3
!  (2) 下弦材(x方向:bk、y方向:hk)
         bk=6.d-3
         hk=163.d-3
!  (3) 垂直材(x方向:bs、y方向:hs)
         bs=12.d-3
         hs=68.d-3
!  (3) 左の斜材(x方向:bsL、y方向:hsL)
         bsL=12.d-3
         hsL=68.d-3
!  (4) 右の斜材(x方向:bsR、y方向:hsR)
         bsR=12.d-3
         hsR=68.d-3
!  y方向の節点数：ny
         ny=3
!  上弦材(下弦材)の数：nk
         nk=14
!  上弦材(下弦材)の1本あたりの要素分割数：nz
         nz=10
!  上弦材(下弦材)1本あたりの長さ：zn
         zn=ell/nk
!
!  節点番号：nset
           nset=0
        do k=0,nk*nz*2
        do j=0,ny-1
           nset=nset+1
           x(nset)=0
           y(nset)=(h*j)/2
           z(nset)=(zn/(2*nz))*k
       print'(i4,a,f9.4,a,f9.4,a,f9.4)',
     &         nset,',',x(nset),',',y(nset),',',z(nset)
        end do
        end do
!
!
!
! <要素番号とその要素を構成する節点番号>
!  ※プラットトラスについて(B32要素1分割)
!  わかりやすいように各部材に分別して表示させる
       print'("** トラス部材の要素番号と構成する節点番号")'
!  (1) 上弦材
       print'("*ELEMENT, TYPE=B32, ELSET=jougen")'
          nyou=0
        do k=1,nk*nz
          nyou=nyou+1
         n1=6*k-5
         n2=6*k-2
         n3=6*k+1
       print'(i4,a,i4,a,i4,a,i4)',
     &   nyou,',',n1,',',n2,',',n3
        end do
!  (2) 下弦材
       print'("*ELEMENT, TYPE=B32, ELSET=kagen")'
          nyou=nk*nz
        do k=1,nk*nz
          nyou=nyou+1
         n1=6*k-3
         n2=6*k
         n3=6*k+3
       print'(i4,a,i4,a,i4,a,i4)',
     &   nyou,',',n1,',',n2,',',n3
        end do
!  (3) 垂直材
       print'("*ELEMENT, TYPE=B32, ELSET=suityoku")'
          nyou=nk*nz+nk*nz
        do k=1,nk+1
          nyou=nyou+1
         n1=6*nz*k-(6*nz-1)
         n2=6*nz*k-(6*nz-2)
         n3=6*nz*k-(6*nz-3)
       print'(i4,a,i4,a,i4,a,i4)',
     &   nyou,',',n1,',',n2,',',n3
        end do
!  (4) 左の斜材
       print'("*ELEMENT, TYPE=B32, ELSET=syazaiL")'
          nyou=nk*nz+nk*nz+nk+1
        do k=1,nk/2
          nyou=nyou+1
         n1=6*nz*k-(6*nz-1)
         n2=6*nz*k-(3*nz-2)
         n3=6*nz*k+3
       print'(i4,a,i4,a,i4,a,i4)',
     &   nyou,',',n1,',',n2,',',n3
        end do
!  (5) 右の斜材
       print'("*ELEMENT, TYPE=B32, ELSET=syazaiR")'
          nyou=nk*nz+nk*nz+nk+1+nk/2
        do k=nk/2+1,nk
          nyou=nyou+1
         n1=6*nz*k-(6*nz-3)
         n2=6*nz*k-(3*nz-2)
         n3=6*nz*k+1
       print'(i4,a,i4,a,i4,a,i4)',
     &   nyou,',',n1,',',n2,',',n3
        end do
!
!
! <境界条件>
!  拘束節点の入力は、
!  (節点番号,拘束する最初の自由度,拘束する最後の自由度)
!  の順。
       print'("** トラス部材の境界条件")'
!  ヒンジ支承
       print'("*NSET,NSET=HingeT")'
!  拘束節点の節点番号
       print'(i4,a)',ny,','
!  ローラー支承
       print'("*NSET,NSET=RollerT")'
       print'(i4,a)',(nk*nz*2+1)*ny,','
!  PC鋼棒の位置(左)
       print'("*NSET,NSET=PCL")'
        do i=0,3
         nlt=3*nz+1+12*nz*i
         nlb=3*nz+3+12*nz*i
         print'(i4,a)',nlt,','
         print'(i4,a)',nlb,','
        end do
!  PC鋼棒の位置(右)
       print'("*NSET,NSET=PCR")'
        do i=4,7
         nrt=3*nz+1+12*nz*i-6*nz
         nrb=3*nz+3+12*nz*i-6*nz
         print'(i4,a)',nrt,','
         print'(i4,a)',nrb,','
        end do
!  格点の位置
!       print'("*NSET,NSET=Kakuten")'
!        do i=0,14
!         nkakut=6*nz*i+1
!         nkakub=6*nz*i+3
!         print'(i4,a)',nkakut,','
!         print'(i4,a)',nkakub,','
!        end do
!  すべての点
!       print'("*NSET,NSET=Nall")'
!         nall=0
!        do i=1,(nz*2*nk+1)*ny
!         nall=nall+1
!         print'(i4,a)',nall,','
!        end do
!       
!  載荷部のz方向変位
!       print'("*NSET,NSET=SaikaZ")'
!       print*,nz*ny+1,","
!        do i=1,mx+1
!         print*,(nz*2+1)*ny+(mx+1)*(my+1)*mz/2+i,","
!        end do
!
       print'("** それぞれの拘束変位")'
       print'("*BOUNDARY")'
       print'(" HingeT,      1,3")'
       print'(" HingeT,      5,6")'
       print'(" RollerT,     1,2")'
       print'(" RollerT,     5,6")'
       print'(" PCL,         1")'
       print'(" PCR,         1")'
!       print'(" Kakuten,     5,6")'
!       print'(" Nall,        5,6")'
!       print'(" SaikaZ,      3")'
!
!
! <材料の定義>
!  3行目にヤング率(MPa),ポアソン比
! ※トラスについて
       print'("** トラス部材に用いる材料")'
       print'("*MATERIAL,NAME=METAL")'
       print'("*ELASTIC")'
       print'(" 206.d3, 0.3")'
!
!
! <材料と要素の関連づけ>
!  ※プラットトラスについて
!  (1) 上弦材
       print'("*BEAM SECTION,ELSET=jougen,MATERIAL=METAL,SECTION=RECT")'
!  断面の第1軸方向の長さ,第2軸方向の長さ
!  断面の第1軸を定義するための単位ベクトル
!        print'(i6,a)',i,','
       print*,bj,',',hj
       print*,1.,',',0.,',',0.
!  (2) 下弦材
       print'("*BEAM SECTION,ELSET=kagen,MATERIAL=METAL,SECTION=RECT")'
       print*,bk,',',hk
       print*,1.,',',0.,',',0.
!  (3) 垂直材
       print'("*BEAM SECTION,ELSET=suityoku,MATERIAL=METAL,SECTION=RECT"
     &  )'
       print*,bs,',',hs
       print*,1.,',',0.,',',0.
!  (4) 左の斜材
      print'("*BEAM SECTION,ELSET=syazaiL,MATERIAL=METAL,SECTION=RECT")'
       print*,bsL,',',hsL
       print*,1.,',',0.,',',0.
!  (5) 右の斜材
      print'("*BEAM SECTION,ELSET=syazaiR,MATERIAL=METAL,SECTION=RECT")'
       print*,bsR,',',hsR
       print*,1.,',',0.,',',0.
!
! <載荷条件>
!  節点ごとに異なる載荷条件を与える場合をふまえて、
!  載荷条件の種類に対して適当な名前(LASTとか)を付ける
!
!  ※プラットトラスについて(中央載荷)
       print'("** トラス部材の載荷節点")'
       print'("*NSET,NSET=SaikaT")'
       print*,nk*nz*ny+1,","
!
!
       print'("*STEP")'
       print'("*STATIC,SOLVER=SPOOLES")'
       print'("*CLOAD")'
!  載荷節点のy(つまり2)方向に荷重を分配して載荷する。
!  座屈荷重を出すので、倍率bは1
       write(*,*)'SaikaT,    2,',1.
!
       print'("*BUCKLE")'
       print'("5,0.1")'
       print'("*NODE PRINT,NSET=Truss")'
       print'("U")'
       print'("*EL PRINT,ELSET=jougen")'
       print'("S")'
       print'("*EL PRINT,ELSET=kagen")'
       print'("S")'
       print'("*EL PRINT,ELSET=suityoku")'
       print'("S")'
       print'("*EL PRINT,ELSET=syazaiL")'
       print'("S")'
       print'("*EL PRINT,ELSET=syazaiR")'
       print'("S")'
       print'("*NODE FILE,output=")'
       print'("U")'
       print'("*EL FILE,POSITION=AVERAGED AT NODES")'
       print'("S")'
       print'("*END STEP")'
!
!
!
      end
