! EUC-jp
! 文字コードはEUC-jp
! 自由形式のfortranで書いています。
! g77でコンパイルする場合は、
!g77 -ffree-form -o smv smv.f
! のようにコンパイルしてください。
! gfortran でコンパイルする場合は文字コードをutf-8にしてください。
!program daiya
implicit real*8 (a-h, o-z)
pi=2.d0*asin(1.d0)
r=1.d0/2.d0/pi !半径
r=1.
h=0.4d0 !高さ
zure=0.5 !高さ方向のずれ
m=14 !周方向分割数
ny=4!高さ方向分割数
n=2
nset=3*m+3*m/2*(ny-1) !節点数
nyou=3*m*ny !要素数
       open(7,file='yousosuu.txt')
       write(7,*) nyou
       write(7,*) m,ny
       close(7)
th=2.d0*pi/m !ダイヤ1つぶんの中心角
dh=h/n !ダイヤ三角1段ぶんの高さ(斜めになってるから三角形の高さではない)
!
open(8,file='tubure.obj')
open(9,file='tubureccx.inp')
write(8,'("3DG1")')
write(8,*) nset
!do j=0,n   !高さ方向に
ph=(pi-th)/2. !多角形の1辺と中心を結ぶ三角形の底角
uegawa=1. !同じ段で上側を向いているか下側を向いているかを正負で区別
do i=0,m-1 !円周方向に1回りずつ
!xyz座標を出力
!上側右三角形の左下頂点、上、右下頂点を出して、
!thずつ座標変換で回転させる
thi=th*real(i)
j=1 !1段め
if(uegawa>0) then !上側を向いている場合、右三角形の上頂点だけz座標をdh下げる
write(8,*) sin(thi)*r*sin(ph), cos(thi)*r*sin(ph), (j-1)*dh
write(8,*) cos(thi)*(r*cos(ph)+dh/2./sin(th))+sin(thi)*r*sin(ph), &
&         -sin(thi)*(r*cos(ph)+dh/2./sin(th))+cos(thi)*r*sin(ph), (j-1)*dh
write(8,*) cos(thi)*(r*cos(ph)-dh/2./tan(th))+sin(thi)*(r*sin(ph)+dh/2.), &
&         -sin(thi)*(r*cos(ph)-dh/2./tan(th))+cos(thi)*(r*sin(ph)+dh/2.), j*dh
else !下側を向いている場合、右三角形の左右下頂点をdh下げる。上頂点は下げない。
write(8,*) sin(thi)*r*sin(ph), cos(thi)*r*sin(ph), j*dh
write(8,*) cos(thi)*(r*cos(ph)+dh/2./sin(th))+sin(thi)*r*sin(ph), &
&         -sin(thi)*(r*cos(ph)+dh/2./sin(th))+cos(thi)*r*sin(ph), j*dh
write(8,*) cos(thi)*(r*cos(ph)-dh/2./tan(th))+sin(thi)*(r*sin(ph)+dh/2.), &
&         -sin(thi)*(r*cos(ph)-dh/2./tan(th))+cos(thi)*(r*sin(ph)+dh/2.), (j-1)*dh
endif
uegawa=uegawa*(-1.)
end do
!end do
!
!
!
! 2段め
j=2
do i=0,m-2,2 !円周方向に1回りずつだけど1段めと節点を共有しないのは、
!下側右三角形の左下頂点とそれからth回転した上側中三角形の左右上頂点。
!この3頂点をthずつm/2回回転させればいいんでステップが2になってる。
thi=th*real(i)
th1=th*real(i+1) !下側右三角形からth回転した上側中三角形の頂点用
write(8,*) sin(thi)*r*sin(ph), cos(thi)*r*sin(ph), j*dh
write(8,*) cos(thi)*(r*cos(ph)+dh/2./sin(th))+sin(thi)*r*sin(ph), &
&         -sin(thi)*(r*cos(ph)+dh/2./sin(th))+cos(thi)*r*sin(ph), j*dh
write(8,*) cos(th1)*(r*cos(ph)-dh/2./tan(th))+sin(th1)*(r*sin(ph)+dh/2.), &
&         -sin(th1)*(r*cos(ph)-dh/2./tan(th))+cos(th1)*(r*sin(ph)+dh/2.), j*dh
end do
!
!
!
! 3段め
j=3
do i=0,m-2,2 !円周方向に1回りずつだけど2段めと節点を共有しないのは、
thi=th*real(i)
th1=th*real(i+1) !
write(8,*) cos(thi)*(r*cos(ph)-dh/2./tan(th))+sin(thi)*(r*sin(ph)+dh/2.), &
&         -sin(thi)*(r*cos(ph)-dh/2./tan(th))+cos(thi)*(r*sin(ph)+dh/2.), j*dh
write(8,*) sin(th1)*r*sin(ph), cos(th1)*r*sin(ph), j*dh
write(8,*) cos(th1)*(r*cos(ph)+dh/2./sin(th))+sin(th1)*r*sin(ph), &
&         -sin(th1)*(r*cos(ph)+dh/2./sin(th))+cos(th1)*r*sin(ph), j*dh
end do
!
!
!
! 4段め
j=4
do i=0,m-2,2 !円周方向に1回りずつだけど3段めと節点を共有しないのは、
thi=th*real(i)
th1=th*real(i+1) !
write(8,*) sin(thi)*r*sin(ph), cos(thi)*r*sin(ph), j*dh
write(8,*) cos(thi)*(r*cos(ph)+dh/2./sin(th))+sin(thi)*r*sin(ph), &
&         -sin(thi)*(r*cos(ph)+dh/2./sin(th))+cos(thi)*r*sin(ph), j*dh
write(8,*) cos(th1)*(r*cos(ph)-dh/2./tan(th))+sin(th1)*(r*sin(ph)+dh/2.), &
&         -sin(th1)*(r*cos(ph)-dh/2./tan(th))+cos(th1)*(r*sin(ph)+dh/2.), j*dh
end do
!
!
!
!1周ずつ上側三形と下側三角形の頂点の節点を外から見て左回りに出力
!do j=0,n-1 !高さ方向
!1段め
write(9,'("** 日本語")') !ccx用の載荷節点と拘束節点の出力
write(9,'("** 上端の載荷節点")')
do i=0,m-2,2 !周方向に
 if(i/=m-2) then !周の開始から一つ前まで
  !上側右三角形の左下頂点、右下頂点、上頂点
  write(8,*) 3,i*3,i*3+1,i*3+2,' ','0xcccccc'
  write(9,*) i*3,','
  write(9,*) i*3+1,','
  !下側左三角形の左下頂点、上頂点、右下頂点(下から見て左回り)
  write(8,*) 3,i*3+2,i*3+1,i*3+3,' ','0xaaaccc'
  !下側中三角形の左上頂点、右上頂点、下(下から見て左回り)
  write(8,*) 3,i*3+1,i*3+5,i*3+3,' ','0xaaaccc'
  !下側右三角形の左下頂点、上、右下頂点(下から見て左回り)
  write(8,*) 3,i*3+3,i*3+5,i*3+4,' ','0xaaaccc'
  !上側左三角形の左下頂点、右下頂点、上
  write(8,*) 3,i*3+5,i*3+6,i*3+4,' ','0xcccccc'
  write(9,*) i*3+5,','
  !上側中三角形の左上頂点、下、右上頂点
  write(8,*) 3,i*3+4,i*3+6,i*3+8,' ','0xcccccc'
 else !周の一番最後は、開始点の頂点とくっつくように
  !上側右三角形の左下頂点、右下頂点、上頂点
  write(8,*) 3,i*3,i*3+1,i*3+2,' ','0xcccccc'
  write(9,*) i*3,','
  write(9,*) i*3+1,','
  !下側左三角形の左下頂点、上頂点、右下頂点(下から見て左回り)
  write(8,*) 3,i*3+2,i*3+1,i*3+3,' ','0xaaaccc'
  !下側中三角形の左上頂点、右上頂点、下(下から見て左回り)
  write(8,*) 3,i*3+1,i*3+5,i*3+3,' ','0xaaaccc'
  !下側右三角形の左下頂点、上、右下頂点(下から見て左回り)
  write(8,*) 3,i*3+3,i*3+5,i*3+4,' ','0xaaaccc'
  !上側左三角形の左下頂点、右下頂点、上
  write(8,*) 3,i*3+5,0,i*3+4,' ','0xcccccc'
  write(9,*) i*3+5,','
  !上側中三角形の左上頂点、下、右上頂点
  write(8,*) 3,i*3+4,0,2,' ','0xcccccc'
 end if
end do
!end do
!
!
!2段め
do i=0,m-2,2 !周方向に
 i2=i/2
 if(i/=m-2) then !周の開始から一つ前まで
  !上側左三角形の左下頂点、右下頂点、上
  write(8,*) 3,i*3+2,i*3+3,m*3+i2*3+1,' ','0xcccccc'
  !上側中三角形の左上頂点、下、右上頂点
  write(8,*) 3,m*3+i2*3+1,i*3+3,m*3+i2*3+2,' ','0xcccccc'
  !上側右三角形の左下頂点、右下頂点、上頂点
  write(8,*) 3,i*3+3,i*3+4,m*3+i2*3+2,' ','0xcccccc'
  !下側左三角形の左下頂点、上頂点、右下頂点(下から見て左回り)
  write(8,*) 3,m*3+i2*3+2,i*3+4,m*3+i2*3+3,' ','0xaaaccc'
  !下側中三角形の左上頂点、右上頂点、下(下から見て左回り)
  write(8,*) 3,i*3+4,i*3+8,m*3+i2*3+3,' ','0xaaaccc'
  !下側右三角形の左下頂点、上、右下頂点(下から見て左回り)
  write(8,*) 3,m*3+i2*3+3,i*3+8,m*3+i2*3+4,' ','0xaaaccc'
 else !周の一番最後は、開始点の頂点とくっつくように
  !上側左三角形の左下頂点、右下頂点、上
  write(8,*) 3,i*3+2,i*3+3,m*3+i2*3+1,' ','0xcccccc'
  !上側中三角形の左上頂点、下、右上頂点
  write(8,*) 3,m*3+i2*3+1,i*3+3,m*3+i2*3+2,' ','0xcccccc'
  !上側右三角形の左下頂点、右下頂点、上頂点
  write(8,*) 3,i*3+3,i*3+4,m*3+i2*3+2,' ','0xcccccc'
  !下側左三角形の左下頂点、上頂点、右下頂点(下から見て左回り)
  write(8,*) 3,m*3+i2*3+2,i*3+4,m*3,' ','0xaaaccc'
  !下側中三角形の左上頂点、右上頂点、下(下から見て左回り)
  write(8,*) 3,i*3+4,2,m*3,' ','0xaaaccc'
  !下側右三角形の左下頂点、上、右下頂点(下から見て左回り)
  write(8,*) 3,m*3,2,m*3+1,' ','0xaaaccc'
 end if
end do
!
!
!
!3段め
write(9,'("** 日本語")') !ccx用の載荷節点と拘束節点の出力
write(9,'("** 下端の拘束節点")')
do i=0,m-2,2 !周方向に
 i2=i/2
 if(i/=m-2) then !周の開始から一つ前まで
  !下側中三角形の左上頂点、右上頂点、下(下から見て左回り)
  write(8,*) 3,m*3+i2*3+1,m*3+i2*3+2,m*4+m/2+i2*3+1,' ','0xaaaccc'
  !下側右三角形の左下頂点、上、右下頂点(下から見て左回り)
  write(8,*) 3,m*4+m/2+i2*3+1,m*3+i2*3+2,m*4+m/2+i2*3+2,' ','0xaaaccc'
  !上側左三角形の左下頂点、右下頂点、上
  write(8,*) 3,m*3+i2*3+2,m*3+i2*3+3,m*4+m/2+i2*3+2,' ','0xcccccc'
  !上側中三角形の左上頂点、下、右上頂点
  write(8,*) 3,m*4+m/2+i2*3+2,m*3+i2*3+3,m*4+m/2+i2*3+3,' ','0xcccccc'
  !上側右三角形の左下頂点、右下頂点、上頂点
  write(8,*) 3,m*3+i2*3+3,m*3+i2*3+4,m*4+m/2+i2*3+3,' ','0xcccccc'
  !下側左三角形の左下頂点、上頂点、右下頂点(下から見て左回り)
  write(8,*) 3,m*4+m/2+i2*3+3,m*3+i2*3+4,m*4+m/2+i2*3+4,' ','0xaaaccc'
 else !周の一番最後は、開始点の頂点とくっつくように
  !下側中三角形の左上頂点、右上頂点、下(下から見て左回り)
  write(8,*) 3,m*3+i2*3+1,m*3+i2*3+2,m*4+m/2+i2*3+1,' ','0xaaaccc'
  !下側右三角形の左下頂点、上、右下頂点(下から見て左回り)
  write(8,*) 3,m*4+m/2+i2*3+1,m*3+i2*3+2,m*4+m/2+i2*3+2,' ','0xaaaccc'
  !上側左三角形の左下頂点、右下頂点、上
  write(8,*) 3,m*3+i2*3+2,m*3,m*4+m/2+i2*3+2,' ','0xcccccc'
  !上側中三角形の左上頂点、下、右上頂点
  write(8,*) 3,m*4+m/2+i2*3+2,m*3,m*4+m/2,' ','0xcccccc'
  !上側右三角形の左下頂点、右下頂点、上頂点
  write(8,*) 3,m*3,m*3+1,m*4+m/2,' ','0xcccccc'
  !下側左三角形の左下頂点、上頂点、右下頂点(下から見て左回り)
  write(8,*) 3,m*4+m/2,m*3+1,m*4+m/2+1,' ','0xaaaccc'
 end if
end do
!
!
!
!4段め
do i=0,m-2,2 !周方向に
 i2=i/2
 if(i/=m-2) then !周の開始から一つ前まで
  !上側右三角形の左下頂点、右下頂点、上頂点
  write(8,*) 3,m*4+m/2+i2*3+1,m*4+m/2+i2*3+2,m*5+m+i2*3+2,' ','0xcccccc'
  !下側左三角形の左下頂点、上頂点、右下頂点(下から見て左回り)
  write(8,*) 3,m*5+m+i2*3+2,m*4+m/2+i2*3+2,m*5+m+i2*3+3,' ','0xaaaccc'
  !下側中三角形の左上頂点、右上頂点、下(下から見て左回り)
  write(8,*) 3,m*4+m/2+i2*3+2,m*4+m/2+i2*3+3,m*5+m+i2*3+3,' ','0xaaaccc'
  !下側右三角形の左下頂点、上、右下頂点(下から見て左回り)
  write(8,*) 3,m*5+m+i2*3+3,m*4+m/2+i2*3+3,m*5+m+i2*3+4,' ','0xaaaccc'
  !上側左三角形の左下頂点、右下頂点、上
  write(8,*) 3,m*4+m/2+i2*3+3,m*4+m/2+i2*3+4,m*5+m+i2*3+4,' ','0xcccccc'
  !上側中三角形の左上頂点、下、右上頂点
  write(8,*) 3,m*5+m+i2*3+4,m*4+m/2+i2*3+4,m*5+m+i2*3+5,' ','0xcccccc'
 else !周の一番最後は、開始点の頂点とくっつくように
  !上側右三角形の左下頂点、右下頂点、上頂点
  write(8,*) 3,m*4+m/2+i2*3+1,m*4+m/2+i2*3+2,m*5+m+i2*3+2,' ','0xcccccc'
  !下側左三角形の左下頂点、上頂点、右下頂点(下から見て左回り)
  write(8,*) 3,m*5+m+i2*3+2,m*4+m/2+i2*3+2,m*6,' ','0xaaaccc'
  !下側中三角形の左上頂点、右上頂点、下(下から見て左回り)
  write(8,*) 3,m*4+m/2+i2*3+2,m*3+i2*3+3,m*6,' ','0xaaaccc'
  !下側右三角形の左下頂点、上、右下頂点(下から見て左回り)
  write(8,*) 3,m*6,m*3+i2*3+3,m*6+1,' ','0xaaaccc'
  !上側左三角形の左下頂点、右下頂点、上
  write(8,*) 3,m*3+i2*3+3,m*3+i2*3+4,m*6+1,' ','0xcccccc'
  !上側中三角形の左上頂点、下、右上頂点
  write(8,*) 3,m*6+1,m*3+i2*3+4,m*6+2,' ','0xcccccc'
 end if
end do
!
!
!
write(8,'("#")') !daiya.objファイルから要素数を読み取るとき用
close(8)
close(9)
!
!
end
