implicit real*8(a-h,o-z)
dimension x(20000),y(20000),z(20000),&
&mesh6kaku(1000,6),mesh3kaku(20000,6),mesh2(20000,6),&
&narabi(20000,20000)
!高さの入力(y方向の長さ)
h=10.0d-3
!六角形1辺の長さ(板の中心)
b=10.0d-3
!板の厚さ
atusa=1.0d-3
!分割数（ハニカムの数）を代入する
m=5!x方向のハニカムの数(多いほう)
n=5!z方向のハニカムの数(絶対奇数 3以上)
nsettensuu=(m*2+(m+1)*2)*(n+1)/2
ny=2
nbunkatu=2
ell=2.0d0*b*(n+1)/2+b*(n-1)/2+b
!板の厚さの中心の線の節点の入力
!端の出っ張り
do k=1,ny+1
do j=1,2
do i=1,m
nset=i+(j-1)*(nsettensuu+m)+(k-1)*(nsettensuu+m*2)
x(nset)=(i-1)*(b*sqrt(3.0d0))+(b*sqrt(3.0d0))/2.0d0
y(nset)=(k-1)*(h/ny)
z(nset)=(j-1)*ell
end do
end do
end do
!1列目の節点
do k=1,ny+1
do j=1,(n+1)/2
do i=1,m
nset=i+(m*2+(m+1)*2)*(j-1)+(k-1)*(nsettensuu+m*2)+m
x(nset)=(i-1)*(b*sqrt(3.0d0))+(b*sqrt(3.0d0))/2.0d0
y(nset)=(k-1)*(h/ny)
z(nset)=(j-1)*3.0d0*b+b/2.0d0
end do
end do
end do
!2列目の節点
do k=1,ny+1
do j=1,(n+1)/2
do i=1,m+1
nset=i+(m*2+(m+1)*2)*(j-1)+m+(k-1)*(nsettensuu+m*2)+m
x(nset)=(i-1)*(b*sqrt(3.0d0))
y(nset)=(k-1)*(h/ny)
z(nset)=(j-1)*3.0d0*b+b
end do
end do
end do
!3列目の節点
do k=1,ny+1
do j=1,(n+1)/2
do i=1,m+1
nset=i+(m*2+(m+1)*2)*(j-1)+m*2+1+(k-1)*(nsettensuu+m*2)+m
x(nset)=(i-1)*(b*sqrt(3.0d0))
y(nset)=(k-1)*(h/ny)
z(nset)=(j-1)*3.0d0*b+b*2.0d0
end do
end do
end do
!4列目の節点
do k=1,ny+1
do j=1,(n+1)/2
do i=1,m
nset=i+(m*2+(m+1)*2)*(j-1)+m+(m+1)*2+(k-1)*(nsettensuu+m*2)+m
x(nset)=(i-1)*(b*sqrt(3.0d0))+(b*sqrt(3.0d0))/2.0d0
y(nset)=(k-1)*(h/ny)
z(nset)=(j-1)*3.0d0*b+b*2.0d0+b/2.0d0
end do
end do
end do
!print'("*NODE, NSET=Nall")'
!do i=1,(nsettensuu+m*2)*(ny+1)
!print'(i6,a,f12.8,a,f12.8,a,f12.8)',&
!&i,',',x(i),',',y(i),',',z(i)
!end do
!六角形のグループを作る
!六角形がm個並んでる列
do j=1,(n+1)/2
do i=1,m
nyou=i+(m+m-1)*(j-1)
mesh6kaku(nyou,1)=i+(m*2+(m+1)*2)*(j-1)+m
mesh6kaku(nyou,2)=mesh6kaku(nyou,1)+m
mesh6kaku(nyou,3)=mesh6kaku(nyou,1)+m+1
mesh6kaku(nyou,4)=mesh6kaku(nyou,1)+m*2+1
mesh6kaku(nyou,5)=mesh6kaku(nyou,1)+m*2+2
mesh6kaku(nyou,6)=mesh6kaku(nyou,1)+m+(m+1)*2
end do
end do
!六角形がm-1個並んでる列
do j=1,(n-1)/2
do i=1,m-1
nyou=i+(m+m-1)*(j-1)+m
mesh6kaku(nyou,1)=i+(4*m+2)*(j-1)+m+m+(m+1)+1
mesh6kaku(nyou,2)=mesh6kaku(nyou,1)+m
mesh6kaku(nyou,3)=mesh6kaku(nyou,1)+m+1
mesh6kaku(nyou,4)=mesh6kaku(nyou,1)+m*2
mesh6kaku(nyou,5)=mesh6kaku(nyou,1)+m*2+1
mesh6kaku(nyou,6)=mesh6kaku(nyou,1)+m*3+1
end do
end do
!3角形の要素を作る
!最初の出っ張り
do k=1,ny
do i=1,m
mesh3kaku((i-1)*2+1+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&i+m+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+1+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&i+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+1+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&i+m+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+2+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&i+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+2+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&i+m+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+2+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&i+(k-1)*(nsettensuu+m*2)
end do 
end do
!6角形の数がm個の列
do k=1,ny
do j=1,(n+1)/2
do i=1,m
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+1+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+(j-1)*(m+m-1),2)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+1+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+(j-1)*(m+m-1),1)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+1+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+(j-1)*(m+m-1),2)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+2+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+(j-1)*(m+m-1),1)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+2+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+(j-1)*(m+m-1),2)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+2+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+(j-1)*(m+m-1),1)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+3+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+(j-1)*(m+m-1),4)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+3+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+(j-1)*(m+m-1),2)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+3+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+(j-1)*(m+m-1),4)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+4+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+(j-1)*(m+m-1),2)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+4+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+(j-1)*(m+m-1),4)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+4+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+(j-1)*(m+m-1),2)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+5+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+(j-1)*(m+m-1),6)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+5+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+(j-1)*(m+m-1),4)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+5+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+(j-1)*(m+m-1),6)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+6+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+(j-1)*(m+m-1),4)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+6+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+(j-1)*(m+m-1),6)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+6+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+(j-1)*(m+m-1),4)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+7+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+(j-1)*(m+m-1),6)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+7+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+(j-1)*(m+m-1),5)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+7+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+(j-1)*(m+m-1),6)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+8+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+(j-1)*(m+m-1),5)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+8+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+(j-1)*(m+m-1),6)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+8+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+(j-1)*(m+m-1),5)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+9+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+(j-1)*(m+m-1),3)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+9+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+(j-1)*(m+m-1),1)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+9+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+(j-1)*(m+m-1),3)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+10+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+(j-1)*(m+m-1),1)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+10+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+(j-1)*(m+m-1),3)+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*10+m*2+(j-1)*(m*10+2+m*2)&
&+10+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+(j-1)*(m+m-1),1)+(k-1)*(nsettensuu+m*2)
end do
end do
end do
!6角形の数がm個の列 端
do k=1,ny
do j=1,(n+1)/2
mesh3kaku(1+(j-1)*(m*10+2+m*2)+m*2+m*10&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(m+(j-1)*(m+m-1),5)+(k-1)*(nsettensuu+m*2)
mesh3kaku(1+(j-1)*(m*10+2+m*2)+m*2+m*10&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(m+(j-1)*(m+m-1),3)+(k-1)*(nsettensuu+m*2)
mesh3kaku(1+(j-1)*(m*10+2+m*2)+m*2+m*10&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(m+(j-1)*(m+m-1),5)+(k-1)*(nsettensuu+m*2)+nsettensuu+m*2
mesh3kaku(2+(j-1)*(m*10+2+m*2)+m*2+m*10&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(m+(j-1)*(m+m-1),3)+(k-1)*(nsettensuu+m*2)+nsettensuu+m*2
mesh3kaku(2+(j-1)*(m*10+2+m*2)+m*2+m*10&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(m+(j-1)*(m+m-1),5)+(k-1)*(nsettensuu+m*2)+nsettensuu+m*2
mesh3kaku(2+(j-1)*(m*10+2+m*2)+m*2+m*10&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(m+(j-1)*(m+m-1),3)+(k-1)*(nsettensuu+m*2)
end do
end do
!6角形の数がm-1個の列
do k=1,ny
do j=1,(n-1)/2
do i=1,m-1
mesh3kaku((i-1)*2+1+m*2+m*10+2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+m+(j-1)*(m+m-1),4)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+1+m*2+m*10+2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+m+(j-1)*(m+m-1),2)+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+1+m*2+m*10+2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+m+(j-1)*(m+m-1),4)+(k-1)*(nsettensuu+m*2)+nsettensuu+m*2
mesh3kaku((i-1)*2+2+m*2+m*10+2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(i+m+(j-1)*(m+m-1),2)+(k-1)*(nsettensuu+m*2)+nsettensuu+m*2
mesh3kaku((i-1)*2+2+m*2+m*10+2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(i+m+(j-1)*(m+m-1),4)+(k-1)*(nsettensuu+m*2)+nsettensuu+m*2
mesh3kaku((i-1)*2+2+m*2+m*10+2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(i+m+(j-1)*(m+m-1),2)+(k-1)*(nsettensuu+m*2)
end do
end do
end do
!print'("*ELEMENT, TYPE=S6, ELSET=Eall")'
!do i=1,ny*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2)
!print'(i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6)',&
!&i&
!&,',',mesh3kaku(i,1)&
!&,',',mesh3kaku(i,2)&
!&,',',mesh3kaku(i,3)&
!&,',',mesh3kaku(i,4)&
!&,',',mesh3kaku(i,5)&
!&,',',mesh3kaku(i,6)
!end do
!6角形の数がm-1個の列 端
do k=1,ny
do j=1,(n-1)/2
mesh3kaku(1+m*2+m*10+2+(m-1)*2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(m+m-1+(j-1)*(m+m-1),5)+(k-1)*(nsettensuu+m*2)
mesh3kaku(1+m*2+m*10+2+(m-1)*2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(m+m-1+(j-1)*(m+m-1),3)+(k-1)*(nsettensuu+m*2)
mesh3kaku(1+m*2+m*10+2+(m-1)*2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(m+m-1+(j-1)*(m+m-1),5)+(k-1)*(nsettensuu+m*2)+nsettensuu+m*2
mesh3kaku(2+m*2+m*10+2+(m-1)*2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&mesh6kaku(m+m-1+(j-1)*(m+m-1),3)+(k-1)*(nsettensuu+m*2)+nsettensuu+m*2
mesh3kaku(2+m*2+m*10+2+(m-1)*2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&mesh6kaku(m+m-1+(j-1)*(m+m-1),5)+(k-1)*(nsettensuu+m*2)+nsettensuu+m*2
mesh3kaku(2+m*2+m*10+2+(m-1)*2+(j-1)*(m*10+2+m*2)&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&mesh6kaku(m+m-1+(j-1)*(m+m-1),3)+(k-1)*(nsettensuu+m*2)
end do
end do
!最後の出っ張り
do k=1,ny
do i=1,m
mesh3kaku((i-1)*2+1+(m*10+2)*(n+1)/2+m*2*(n-1)/2+m*2&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&i+m+nsettensuu+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+1+(m*10+2)*(n+1)/2+m*2*(n-1)/2+m*2&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&i+nsettensuu+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+1+(m*10+2)*(n+1)/2+m*2*(n-1)/2+m*2&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&i+m+nsettensuu+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+2+(m*10+2)*(n+1)/2+m*2*(n-1)/2+m*2&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),1)=&
&i+nsettensuu+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+2+(m*10+2)*(n+1)/2+m*2*(n-1)/2+m*2&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),2)=&
&i+m+nsettensuu+nsettensuu+m*2+(k-1)*(nsettensuu+m*2)
mesh3kaku((i-1)*2+2+(m*10+2)*(n+1)/2+m*2*(n-1)/2+m*2&
&+(k-1)*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2),3)=&
&i+nsettensuu+(k-1)*(nsettensuu+m*2)
end do 
end do
!中間節点を入れて6節点要素にする
do i=1,(nsettensuu+2*m)*(ny+1)
do j=1,(nsettensuu+2*m)*(ny+1)
narabi(i,j)=-1
end do
end do
iset=(nsettensuu+2*m)*(ny+1)
!print *,iset
do i=1,ny*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2)
!1-2
if(narabi(mesh3kaku(i,1),mesh3kaku(i,2))==-1) then
iset=iset+1
narabi(mesh3kaku(i,1),mesh3kaku(i,2))=iset !頂点1-2間と
narabi(mesh3kaku(i,2),mesh3kaku(i,1))=iset !頂点2-1間とに節点番号isetをふり
mesh3kaku(i,4)=iset
x(iset)=(x(mesh3kaku(i,1))+x(mesh3kaku(i,2)))/2.d0
y(iset)=(y(mesh3kaku(i,1))+y(mesh3kaku(i,2)))/2.d0
z(iset)=(z(mesh3kaku(i,1))+z(mesh3kaku(i,2)))/2.d0
else
mesh3kaku(i,4)=narabi(mesh3kaku(i,1),mesh3kaku(i,2))
end if
!2-3
if(narabi(mesh3kaku(i,2),mesh3kaku(i,3))==-1) then
iset=iset+1
narabi(mesh3kaku(i,2),mesh3kaku(i,3))=iset !頂点2-3間と
narabi(mesh3kaku(i,3),mesh3kaku(i,2))=iset !頂点3-2間とに節点番号isetをふり
mesh3kaku(i,5)=iset
x(iset)=(x(mesh3kaku(i,2))+x(mesh3kaku(i,3)))/2.d0
y(iset)=(y(mesh3kaku(i,2))+y(mesh3kaku(i,3)))/2.d0
z(iset)=(z(mesh3kaku(i,2))+z(mesh3kaku(i,3)))/2.d0
else
mesh3kaku(i,5)=narabi(mesh3kaku(i,2),mesh3kaku(i,3))
end if
!3-1
if(narabi(mesh3kaku(i,3),mesh3kaku(i,1))==-1) then
iset=iset+1
narabi(mesh3kaku(i,3),mesh3kaku(i,1))=iset !頂点3-1間と
narabi(mesh3kaku(i,1),mesh3kaku(i,3))=iset !頂点1-3間とに節点番号isetをふり
mesh3kaku(i,6)=iset
x(iset)=(x(mesh3kaku(i,1))+x(mesh3kaku(i,3)))/2.d0
y(iset)=(y(mesh3kaku(i,1))+y(mesh3kaku(i,3)))/2.d0
z(iset)=(z(mesh3kaku(i,1))+z(mesh3kaku(i,3)))/2.d0
else
mesh3kaku(i,6)=narabi(mesh3kaku(i,1),mesh3kaku(i,3))
end if
end do
!増やす時ためのmesh2(i,j)の設定
nyou=ny*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2)
nyou1=ny*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2)
!narabiに-1を入れる
do i=1,iset
do j=1,iset
narabi(i,j)=-1
end do
end do
!mesh2にmeshを入れる
do i=1,nyou
do j=1,6
mesh2(i,j)=mesh3kaku(i,j)
end do 
end do
!
do 
if (nbunkatu==0) exit
!
do i=1,nyou
mesh2(1+(i-1)*4,1)=mesh3kaku(i,1)
mesh2(1+(i-1)*4,2)=mesh3kaku(i,4)
mesh2(1+(i-1)*4,3)=mesh3kaku(i,6)
mesh2(2+(i-1)*4,1)=mesh3kaku(i,4)
mesh2(2+(i-1)*4,2)=mesh3kaku(i,2)
mesh2(2+(i-1)*4,3)=mesh3kaku(i,5)
mesh2(3+(i-1)*4,1)=mesh3kaku(i,6)
mesh2(3+(i-1)*4,2)=mesh3kaku(i,5)
mesh2(3+(i-1)*4,3)=mesh3kaku(i,3)
mesh2(4+(i-1)*4,1)=mesh3kaku(i,5)
mesh2(4+(i-1)*4,2)=mesh3kaku(i,6)
mesh2(4+(i-1)*4,3)=mesh3kaku(i,4)
end do
!中間節点を入れて要素を増やす
do i=1,nyou
!1-4
if(narabi(mesh3kaku(i,1),mesh3kaku(i,4))==-1) then
iset=iset+1
narabi(mesh3kaku(i,1),mesh3kaku(i,4))=iset !頂点1-4間と
narabi(mesh3kaku(i,4),mesh3kaku(i,1))=iset !頂点4-1間とに節点番号isetをふり
mesh2(1+(i-1)*4,4)=iset
x(iset)=(x(mesh3kaku(i,1))+x(mesh3kaku(i,4)))/2.d0
y(iset)=(y(mesh3kaku(i,1))+y(mesh3kaku(i,4)))/2.d0
z(iset)=(z(mesh3kaku(i,1))+z(mesh3kaku(i,4)))/2.d0
else
mesh2(1+(i-1)*4,4)=narabi(mesh3kaku(i,1),mesh3kaku(i,4))
end if
!4-6
if(narabi(mesh3kaku(i,4),mesh3kaku(i,6))==-1) then
iset=iset+1
narabi(mesh3kaku(i,4),mesh3kaku(i,6))=iset !頂点1-4間と
narabi(mesh3kaku(i,6),mesh3kaku(i,4))=iset !頂点4-1間とに節点番号isetをふり
mesh2(1+(i-1)*4,5)=iset
mesh2(4+(i-1)*4,5)=iset
x(iset)=(x(mesh3kaku(i,4))+x(mesh3kaku(i,6)))/2.d0
y(iset)=(y(mesh3kaku(i,4))+y(mesh3kaku(i,6)))/2.d0
z(iset)=(z(mesh3kaku(i,4))+z(mesh3kaku(i,6)))/2.d0
else
mesh2(1+(i-1)*4,5)=narabi(mesh3kaku(i,4),mesh3kaku(i,6))
mesh2(4+(i-1)*4,5)=narabi(mesh3kaku(i,4),mesh3kaku(i,6))
end if
!1-6
if(narabi(mesh3kaku(i,1),mesh3kaku(i,6))==-1) then
iset=iset+1
narabi(mesh3kaku(i,1),mesh3kaku(i,6))=iset !頂点1-6間と
narabi(mesh3kaku(i,6),mesh3kaku(i,1))=iset !頂点6-1間とに節点番号isetをふり
mesh2(1+(i-1)*4,6)=iset
x(iset)=(x(mesh3kaku(i,1))+x(mesh3kaku(i,6)))/2.d0
y(iset)=(y(mesh3kaku(i,1))+y(mesh3kaku(i,6)))/2.d0
z(iset)=(z(mesh3kaku(i,1))+z(mesh3kaku(i,6)))/2.d0
else
mesh2(1+(i-1)*4,6)=narabi(mesh3kaku(i,1),mesh3kaku(i,6))
end if
!2-4
if(narabi(mesh3kaku(i,2),mesh3kaku(i,4))==-1) then
iset=iset+1
narabi(mesh3kaku(i,2),mesh3kaku(i,4))=iset !頂点1-4間と
narabi(mesh3kaku(i,4),mesh3kaku(i,2))=iset !頂点4-1間とに節点番号isetをふり
mesh2(2+(i-1)*4,4)=iset
x(iset)=(x(mesh3kaku(i,2))+x(mesh3kaku(i,4)))/2.d0
y(iset)=(y(mesh3kaku(i,2))+y(mesh3kaku(i,4)))/2.d0
z(iset)=(z(mesh3kaku(i,2))+z(mesh3kaku(i,4)))/2.d0
else
mesh2(2+(i-1)*4,4)=narabi(mesh3kaku(i,2),mesh3kaku(i,4))
end if
!2-5
if(narabi(mesh3kaku(i,2),mesh3kaku(i,5))==-1) then
iset=iset+1
narabi(mesh3kaku(i,2),mesh3kaku(i,5))=iset !頂点1-4間と
narabi(mesh3kaku(i,5),mesh3kaku(i,2))=iset !頂点4-1間とに節点番号isetをふり
mesh2(2+(i-1)*4,5)=iset
x(iset)=(x(mesh3kaku(i,2))+x(mesh3kaku(i,5)))/2.d0
y(iset)=(y(mesh3kaku(i,2))+y(mesh3kaku(i,5)))/2.d0
z(iset)=(z(mesh3kaku(i,2))+z(mesh3kaku(i,5)))/2.d0
else
mesh2(2+(i-1)*4,5)=narabi(mesh3kaku(i,2),mesh3kaku(i,5))
end if
!4-5
if(narabi(mesh3kaku(i,4),mesh3kaku(i,5))==-1) then
iset=iset+1
narabi(mesh3kaku(i,4),mesh3kaku(i,5))=iset !頂点1-4間と
narabi(mesh3kaku(i,5),mesh3kaku(i,4))=iset !頂点4-1間とに節点番号isetをふり
mesh2(2+(i-1)*4,6)=iset
mesh2(4+(i-1)*4,6)=iset
x(iset)=(x(mesh3kaku(i,4))+x(mesh3kaku(i,5)))/2.d0
y(iset)=(y(mesh3kaku(i,4))+y(mesh3kaku(i,5)))/2.d0
z(iset)=(z(mesh3kaku(i,4))+z(mesh3kaku(i,5)))/2.d0
else
mesh2(2+(i-1)*4,6)=narabi(mesh3kaku(i,4),mesh3kaku(i,5))
mesh2(4+(i-1)*4,6)=narabi(mesh3kaku(i,4),mesh3kaku(i,5))
end if
!5-6
if(narabi(mesh3kaku(i,5),mesh3kaku(i,6))==-1) then
iset=iset+1
narabi(mesh3kaku(i,5),mesh3kaku(i,6))=iset !頂点1-4間と
narabi(mesh3kaku(i,6),mesh3kaku(i,5))=iset !頂点4-1間とに節点番号isetをふり
mesh2(3+(i-1)*4,4)=iset
mesh2(4+(i-1)*4,4)=iset
x(iset)=(x(mesh3kaku(i,5))+x(mesh3kaku(i,6)))/2.d0
y(iset)=(y(mesh3kaku(i,5))+y(mesh3kaku(i,6)))/2.d0
z(iset)=(z(mesh3kaku(i,5))+z(mesh3kaku(i,6)))/2.d0
else
mesh2(3+(i-1)*4,4)=narabi(mesh3kaku(i,5),mesh3kaku(i,6))
mesh2(4+(i-1)*4,4)=narabi(mesh3kaku(i,5),mesh3kaku(i,6))
end if
!3-5
if(narabi(mesh3kaku(i,3),mesh3kaku(i,5))==-1) then
iset=iset+1
narabi(mesh3kaku(i,3),mesh3kaku(i,5))=iset !頂点1-4間と
narabi(mesh3kaku(i,5),mesh3kaku(i,3))=iset !頂点4-1間とに節点番号isetをふり
mesh2(3+(i-1)*4,5)=iset
x(iset)=(x(mesh3kaku(i,3))+x(mesh3kaku(i,5)))/2.d0
y(iset)=(y(mesh3kaku(i,3))+y(mesh3kaku(i,5)))/2.d0
z(iset)=(z(mesh3kaku(i,3))+z(mesh3kaku(i,5)))/2.d0
else
mesh2(3+(i-1)*4,5)=narabi(mesh3kaku(i,3),mesh3kaku(i,5))
end if
!3-6
if(narabi(mesh3kaku(i,3),mesh3kaku(i,6))==-1) then
iset=iset+1
narabi(mesh3kaku(i,3),mesh3kaku(i,6))=iset !頂点1-4間と
narabi(mesh3kaku(i,6),mesh3kaku(i,3))=iset !頂点4-1間とに節点番号isetをふり
mesh2(3+(i-1)*4,6)=iset
x(iset)=(x(mesh3kaku(i,3))+x(mesh3kaku(i,6)))/2.d0
y(iset)=(y(mesh3kaku(i,3))+y(mesh3kaku(i,6)))/2.d0
z(iset)=(z(mesh3kaku(i,3))+z(mesh3kaku(i,6)))/2.d0
else
mesh2(3+(i-1)*4,6)=narabi(mesh3kaku(i,3),mesh3kaku(i,6))
end if
end do
!mesh(i,j)をmesh2(i,j)に入れる
do i=1,iset
do j=1,iset
narabi(i,j)=-1
end do
end do
nyou=nyou*4
do i=1,nyou
do j=1,6
mesh3kaku(i,j)=mesh2(i,j)
end do
end do
if (nyou==nyou1*4**ny) exit
!
end do
!節点番号と座標
print'("*NODE, NSET=Nall")'
do i=1,iset
print'(i6,a,f12.8,a,f12.8,a,f12.8)',&
&i,',',x(i),',',y(i),',',z(i)
end do
!要素番号と節点番号
print'("*ELEMENT, TYPE=S6, ELSET=Eall")'
!do i=1,ny*(m*2*2+(m*2*2*2+(m+1)*2)*(n+1)/2+m*2*(n-1)/2)
do i=1,nyou
print'(i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6)',&
&i&
&,',',mesh2(i,1)&
&,',',mesh2(i,2)&
&,',',mesh2(i,3)&
&,',',mesh2(i,4)&
&,',',mesh2(i,5)&
&,',',mesh2(i,6)
end do
print'("*NSET, NSET=kotei")'
do i=1,iset
if (z(i)<0.00001.and.z(i)>-0.00001) then
print '(i6,a)',i,','
end if
end do
print'("*NSET, NSET=saika")'
do i=1,iset
if (z(i)<ell+0.00001.and.z(i)>ell-0.00001) then
print '(i6,a)',i,','
end if
end do
!境界条件の設定
print'("*BOUNDARY")'
print'("kotei,1,6")'
!材料の設定
print'("*MATERIAL,NAME=EL")'
print'("*ELASTIC")'
print'("3011.7,0.32907")'
!要素の設定,どの材料を使うか等
print'("*SHELL SECTION,ELSET=Eall,MATERIAL=EL")'
print'(1pe13.5)',atusa
!計算のオプション
print'("*STEP")'
print'("*STATIC,SOLVER=SPOOLES")'
!荷重
print'("*CLOAD")'
print'(A,a,i6,a,1pe13.5)','saika ',',',2,',',-1.00E-003/(m*(ny+1))
!出力のオプション
print'("*NODE PRINT,NSET=Nall")'
print'("U")'
print'("*EL PRINT,ELSET=Eall")'
print'("S")'
print'("*NODE FILE")'
print'("U")'
print'("*EL FILE,POSITION=AVERAGED AT NODES")'
print'("S")'
print'("*END STEP")'

end
!

