! 木高研の試験体(半断面)を、トラス部分をB32、木床版をC3D8で作ってみる。 ! ※このプログラムで書いたinpファイルでは節点番号が共有されないので、 ! inpファイルの木部材要素の節点番号を共有させるトラス部材の ! 節点番号に書き換える必要がある。 ! ※そのため、木部材の高さ(y軸)方向の分割数は偶数にする。 implicit real*8(a-h,o-z) dimension x(1000000),y(1000000),z(1000000) ! 節点座標:x(節点番号),y(節点番号),z(節点番号) print'("**")' print'("** Structure: Pratt truss")' print'("** Test objective: B32 and C3D8 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=1 ! 上弦材(下弦材)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'(i6,a,f9.4,a,f9.4,a,f9.4)',& & nset,',',x(nset),',',y(nset),',',z(nset) end do end do ! ! ※木部材について(上部:ue,下部:sita) print'("** 木部材の節点番号と座標")' print'("*NODE, NSET=Wood_ue")' ! 床版の幅(x方向,半断面で):bmu(m) ! 板厚(y方向):hmu(m) ! 支間長(z方向):ell(m) bmu=420.d-3 hmu=120.d-3 ! ell=4900.d-3 ! x,y,z方向の要素分割数(偶数で):mx,my,mz mx=6 my=8 mz=56 ! 1要素のx,y,z方向の長さ:xn,yn,zn xm=bmu/mx ym=hmu/my zm=ell/mz ! 節点番号:nset nset=(nk*nz*2+1)*ny do k=0,mz do j=0,my do i=0,mx nset=nset+1 x(nset)=xm*i y(nset)=ym*j-(hmu/2) z(nset)=zm*k print'(i6,a,f9.4,a,f9.4,a,f9.4)',& & nset,',',x(nset),',',y(nset),',',z(nset) end do end do end do ! print'("*NODE, NSET=Wood_sita")' ! 節点番号:nset nset=(nk*nz*2+1)*ny+((mx+1)*(my+1)*(mz+1)) do k=0,mz do j=0,my do i=0,mx nset=nset+1 x(nset)=xm*i y(nset)=ym*j-(hmu/2)+h z(nset)=zm*k print'(i6,a,f9.4,a,f9.4,a,f9.4)',& & nset,',',x(nset),',',y(nset),',',z(nset) end do 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'(i6,a,i6,a,i6,a,i6)',& & 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'(i6,a,i6,a,i6,a,i6)',& & 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'(i6,a,i6,a,i6,a,i6)',& & 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'(i6,a,i6,a,i6,a,i6)',& & 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'(i6,a,i6,a,i6,a,i6)',& & nyou,',',n1,',',n2,',',n3 end do ! ! ※木部材について(C3D8要素) print'("** 木部材の要素番号と構成する節点番号")' print'("** トラス部材と共有させる節点番号を書き換える必要あり")' print'("*ELEMENT, TYPE=C3D8, ELSET=Banmen_ue")' nyou=nk*nz+nk*nz+nk+1+nk do k=1,mz do j=1,my do i=1,mx nyou=nyou+1 n1=(mx+1)*(my+1)*(k-1)+(mx+1)*(j-1)+i+3*(nk*nz*2+1) n2=n1+1 n3=n1+(mx+1) n4=n3+1 n5=n1+(mx+1)*(my+1) n6=n2+(mx+1)*(my+1) n7=n3+(mx+1)*(my+1) n8=n4+(mx+1)*(my+1) print'(i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6)',& & nyou,',',n1,',',n3,',',n7,',',n5,',',n2,',',n4,',',n8,',',n6 end do end do end do ! !syouban.txtに上の床版の前半部分をコピー open(7,file='syouban.txt') write(7,*) mx write(7,*) my write(7,*) mz write(7,*) nk write(7,*) ny write(7,*) nz write(7,*) (mz*my*mx)/2 nyou=nk*nz+nk*nz+nk+1+nk do k=1,mz/2 do j=1,my do i=1,mx nyou=nyou+1 n1=(mx+1)*(my+1)*(k-1)+(mx+1)*(j-1)+i+3*(nk*nz*2+1) n2=n1+1 n3=n1+(mx+1) n4=n3+1 n5=n1+(mx+1)*(my+1) n6=n2+(mx+1)*(my+1) n7=n3+(mx+1)*(my+1) n8=n4+(mx+1)*(my+1) write(7,*) nyou,n1,n3,n7,n5,n2,n4,n8,n6 end do end do end do close(7) ! !syouban.txtに上の床版の後半部分をコピー open(7,file='syouban2.txt') write(7,*) mx write(7,*) my write(7,*) mz write(7,*) nk write(7,*) ny write(7,*) nz write(7,*) (mz*my*mx)/2 nyou=nk*nz+nk*nz+nk+1+nk+(mx*my*mz)/2 do k=mz/2+1,mz do j=1,my do i=1,mx nyou=nyou+1 n1=(mx+1)*(my+1)*(k-1)+(mx+1)*(j-1)+i+3*(nk*nz*2+1) n2=n1+1 n3=n1+(mx+1) n4=n3+1 n5=n1+(mx+1)*(my+1) n6=n2+(mx+1)*(my+1) n7=n3+(mx+1)*(my+1) n8=n4+(mx+1)*(my+1) write(7,*) nyou,n1,n3,n7,n5,n2,n4,n8,n6 end do end do end do close(7) ! ! print'("*ELEMENT, TYPE=C3D8, ELSET=Banmen_sita")' nyou=nk*nz+nk*nz+nk*nz+1+nk+(mx*my*mz) do k=1,mz do j=1,my do i=1,mx nyou=nyou+1 n1=(mx+1)*(my+1)*(k-1)+(mx+1)*(j-1)+i+3*(nk*nz*2+1)& & +((mx+1)*(my+1)*(mz+1)) n2=n1+1 n3=n1+(mx+1) n4=n3+1 n5=n1+(mx+1)*(my+1) n6=n2+(mx+1)*(my+1) n7=n3+(mx+1)*(my+1) n8=n4+(mx+1)*(my+1) print'(i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6)',& & nyou,',',n1,',',n3,',',n7,',',n5,',',n2,',',n4,',',n8,',',n6 end do end do end do ! !syouban3.txtに下の床版の前半をコピー ! open(7,file='syouban3.txt') write(7,*) mx write(7,*) my write(7,*) mz write(7,*) nk write(7,*) ny write(7,*) nz write(7,*) (mz*my*mx)/2 nyou=nk*nz+nk*nz+nk*nz+1+nk+(mx*my*mz) do k=1,mz/2 do j=1,my do i=1,mx nyou=nyou+1 n1=(mx+1)*(my+1)*(k-1)+(mx+1)*(j-1)+i+3*(nk*nz*2+1)& & +((mx+1)*(my+1)*(mz+1)) n2=n1+1 n3=n1+(mx+1) n4=n3+1 n5=n1+(mx+1)*(my+1) n6=n2+(mx+1)*(my+1) n7=n3+(mx+1)*(my+1) n8=n4+(mx+1)*(my+1) write(7,*) nyou,n1,n3,n7,n5,n2,n4,n8,n6 end do end do end do close(7) ! !syouban4.txtに下の床版の後半をコピー ! open(7,file='syouban4.txt') write(7,*) mx write(7,*) my write(7,*) mz write(7,*) nk write(7,*) ny write(7,*) nz write(7,*) (mz*my*mx)/2 nyou=nk*nz+nk*nz+nk*nz+1+nk+(mx*my*mz)+(mz*my*mx)/2 do k=mz/2+1,mz do j=1,my do i=1,mx nyou=nyou+1 n1=(mx+1)*(my+1)*(k-1)+(mx+1)*(j-1)+i+3*(nk*nz*2+1)& & +((mx+1)*(my+1)*(mz+1)) n2=n1+1 n3=n1+(mx+1) n4=n3+1 n5=n1+(mx+1)*(my+1) n6=n2+(mx+1)*(my+1) n7=n3+(mx+1)*(my+1) n8=n4+(mx+1)*(my+1) write(7,*) nyou,n1,n3,n7,n5,n2,n4,n8,n6 end do end do end do close(7) ! ! ! <境界条件> ! 拘束節点の入力は、 ! (節点番号,拘束する最初の自由度,拘束する最後の自由度) ! の順。 ! ※トラスについて print'("** トラス部材の境界条件")' ! ヒンジ支承 print'("*NSET,NSET=HingeT")' ! 拘束節点の節点番号 print'(i6,a)',ny,',' ! ローラー支承 print'("*NSET,NSET=RollerT")' print'(i6,a)',(nk*nz*2+1)*ny,',' ! ! ※木部材について(下の部材のみ) print'("** 木部材の境界条件")' ! ヒンジ支承 print'("*NSET,NSET=HingeW")' ! 拘束節点の節点番号 HW=(nk*nz*2+1)*ny+((mx+1)*(my+1)*(mz+1))& & +(mx+1)*my+1 do i=HW,HW+mx print'(i6,a)',i,',' end do ! ローラー支承 print'("*NSET,NSET=RollerW")' RW=(nk*nz*2+1)*ny+((mx+1)*(my+1)*(mz+1))& & +(mx+1)*(my+1)*mz+(mx+1)*my+1 do i=RW,RW+mx print'(i6,a)',i,',' 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'("*NSET,NSET=Danmen_ue")' do k=0,mz do j=1,my+1 print'(i6,a)',(nk*nz*2+1)*ny+(mx+1)*j+(mx+1)*(my+1)*k,',' end do end do ! print'("*NSET,NSET=Danmen_sita")' do k=0,mz do j=1,my+1 print'(i6,a)',(nk*nz*2+1)*ny+(mx+1)*(my+1)*(mz+1)+(mx+1)*j& & +(mx+1)*(my+1)*k,',' end do end do ! ! 上記(NSET)の拘束節点の拘束変位 print'("** それぞれの拘束変位")' print'("*BOUNDARY")' print'(" HingeT, 1,2")' print'(" HingeT, 5,6")' print'(" RollerT, 1,2")' print'(" RollerT, 5,6")' print'(" HingeW, 1,2")' print'(" RollerW, 1,2")' ! print'(" SaikaZ, 3")' print'(" Danmen_ue, 1")' print'(" Danmen_sita, 1")' ! ! ! <材料の定義> ! 3行目にヤング率(MPa),ポアソン比 ! ※トラスについて print'("** トラス部材に用いる材料")' print'("*MATERIAL,NAME=METAL")' print'("*ELASTIC")' print'(" 206.d3, 0.3")' ! ※木部材について(C3D8要素) print'("** 木部材に用いる材料")' print'("*MATERIAL,NAME=WOOD")' ! *ELASTICは等方性の意味 print'("*ELASTIC")' print'(" 7.11d3, 0.4")' ! ! ! <材料と要素の関連づけ> ! ※プラットトラスについて ! (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. ! ※木部材について print'("*SOLID SECTION,ELSET=Banmen_ue,MATERIAL=WOOD")' print'("*SOLID SECTION,ELSET=Banmen_sita,MATERIAL=WOOD")' ! ! <載荷条件> ! 節点ごとに異なる載荷条件を与える場合をふまえて、 ! 載荷条件の種類に対して適当な名前(LASTとか)を付ける ! ! ※プラットトラスについて(中央載荷) print'("** トラス部材の載荷節点")' print'("*NSET,NSET=SaikaT")' print*,nk*nz*ny+1,"," ! ※木部材について(中央載荷) print'("** 木部材の載荷節点")' print'("*NSET,NSET=SaikaWhasi")' print*,(nk*nz*2+1)*ny+(mx+1)*(my+1)*mz/2+1,"," print*,(nk*nz*2+1)*ny+(mx+1)*(my+1)*mz/2+mx+1,"," ! print'("*NSET,NSET=SaikaWnaka")' do i=2,mx print*,(nk*nz*2+1)*ny+(mx+1)*(my+1)*mz/2+i,"," end do ! ! print'("*STEP")' print'("*STATIC,SOLVER=SPOOLES")' print'("*CLOAD")' ! 載荷節点のy(つまり2)方向に荷重を分配して載荷する。 ! P/(mx+2)載荷 write(*,*)'SaikaT, 2,',27.5d-3/(mx+2) ! P/(2*(mx+2))載荷 write(*,*)'SaikaWhasi,2,',27.5d-3/((mx+2)*2) ! P/(4*(mx+2))載荷 write(*,*)'SaikaWnaka,2,',27.5d-3/((mx+2)*4) ! ! print'("**BUCKLE")' ! print'("**10,0.1")' print'("*NODE PRINT,NSET=Truss")' print'("U")' print'("*NODE PRINT,NSET=Wood_ue")' print'("U")' print'("*NODE PRINT,NSET=Wood_sita")' 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'("*EL PRINT,ELSET=Banmen_ue")' print'("S")' print'("*EL PRINT,ELSET=Banmen_sita")' print'("S")' print'("*NODE FILE,output=3D")' print'("U")' print'("*EL FILE,POSITION=AVERAGED AT NODES")' print'("S")' print'("*END STEP")' ! ! ! end