!onsi21220.f90謾ｹ濶ｯ迚�
!縺溘□縺ｮ邂ｱ譁ｭ髱｢
implicit real*8(a-h,o-z)
dimension x(1000000),y(1315500),z(1315500),e66(6,6),b(3,3)&
&,meshmoku(300000,8),meshkou(300000,8)
character filename*128,filename2*48
ell=10.0d0
!縺薙％繧貞､峨∴縺溘ｉ荳九�!job縺ｮ險ｭ螳� 縺ｮi3繧ょ酔縺倥￥螟峨∴縺ｦ
write(filename, '("hako15-ell", i3.3, ".inp")') i3!
open  (7,file=filename, status='replace')
b1=8.4d-1
b2=9.0d-3
h1=0.0
h2=5.0d-1
hall=2*h1+h2
p=1.000d-3!闕ｷ驥�
do i3=100,330,10!1縺九ｉ100縺ｾ縺ｧ1蛻ｻ縺ｿ
h1=h1+1.0d-2
h2=h2-2.0d-2
!xyz譁ｹ蜷代�隕∫ｴ�蛻�牡謨ｰ
nx1=10!譛ｨ譚�
nx2=8!驪ｼ譚ｿ
ny1=10!荳企Κ蛻�→荳矩Κ蛻��鬮倥＆譁ｹ蜷代�蛻�牡謨ｰ
ny2=8!髢薙�驪ｼ譚ｿ縺ｮ鬮倥＆譁ｹ蜷代�蛻�牡謨ｰ(3莉･荳�)
nz=500
!繝､繝ｳ繧ｰ邇��蜈･蜉�
ezz=7.694d3
exx=ezz/25.0d0!x霆ｸ譁ｹ蜷代�繝､繝ｳ繧ｰ邇�
eyy=ezz/25.0d0!y霆ｸ譁ｹ蜷代�繝､繝ｳ繧ｰ邇�
ezz2=206.0d3!z霆ｸ譁ｹ蜷代�繝､繝ｳ繧ｰ邇�
!繝昴い繧ｽ繝ｳ豈斐�蜈･蜉�
poixy=0.016d0
poixz=poixy
poiyx=poixy
poiyz=poixy
poizx=0.4d0
poizy=poizx
!縺帙ｓ譁ｭ蠑ｾ諤ｧ菫よ焚縺ｮ蜈･蜉�
gxy=ezz/15.0d0
gxz=gxy
gyz=gxy
!1縺､縺ゅ◆繧翫�隕∫ｴ�縺ｮ髟ｷ縺�?
xn1=b1/nx1
xn2=b2/nx2
yn1=h1/ny1
yn2=h2/ny2
zn=ell/nz
! 遽轤ｹ逡ｪ蜿ｷ縺斐→縺ｮ遽轤ｹ蠎ｧ讓�!i~m縺ｧ蟋九∪繧九�縺慧o縺ｧ菴ｿ縺医ｋ
!荳九→荳企Κ蛻��遽轤ｹ縲驪ｼ譚ｿ!
do m=0,1
do k=0,nz
do j=0,ny1
do l=0,1
do i=0,nx2
nset=1+i+l*(nx2+nx1)+j*(nx2*2+nx1+1)&
&+k*(2*(nx1+nx2*2+1)*(ny1+1)+2*(nx2+1)*(ny2-1))&
&+m*((nx2*2+nx1+1)*(ny1+1)+(nx2+1)*2*(ny2-1))
!荳翫→荳九≠繧九�縺ｧm繧貞虚縺九＠縺ｦ陦ｨ縺�
!k縺ｯ髱｢縺ｮ荳髱｢縺吶∋縺ｦ縺ｮ遽轤ｹ繧呈寺縺代※z譁ｹ蜷代↓蜍輔°縺�

x(nset)=i*xn2+l*(b1+b2) !i縺�0~xn2縺ｾ縺ｧ縺�▲縺ｦi=0,l=1縺ｫ縺ｪ繧�
y(nset)=j*yn1+m*(h1+h2)
z(nset)=k*zn
end do
end do
end do
end do
end do
!荳九→荳企Κ蛻��遽轤ｹ 譛ｨ譚�
do m=0,1
do k=0,nz
do j=0,ny1
do i=0,nx1-2
nset=nx2+1+1+i+j*(nx2*2+nx1+1)&
&+k*(2*(nx1+nx2*2+1)*(ny1+1)+2*(nx2+1)*(ny2-1))&
&+m*((nx2*2+nx1+1)*(ny1+1)+(nx2+1)*2*(ny2-1))
!荳翫→荳�2縺､縺ゅｋ縺ｮ縺ｧm繧貞虚縺九＠縺ｦ陦ｨ縺�
x(nset)=i*xn1+b2+xn1 !i縺�0~xn2縺ｾ縺ｧ縺�▲縺ｦi=0,l=1縺ｫ縺ｪ繧�
y(nset)=j*yn1+m*(h1+h2)
z(nset)=k*zn
end do
end do
end do
end do

!逵溘ｓ荳ｭ驛ｨ蛻��驪ｼ譚舌�遽轤ｹ
do k=0,nz
do j=0,ny2-2
do l=0,1
do i=0,nx2
nset=(1+nx1+2*nx2)*(ny1+1)+1+i+l*(1+nx2)+j*2*(nx2+1)&
&+k*(2*(nx1+nx2*2+1)*(ny1+1)+2*(nx2+1)*(ny2-1))
x(nset)=i*xn2+l*(b1+b2)
y(nset)=j*yn2+h1+yn2
z(nset)=k*zn
end do
end do
end do
end do

!i5~6縺梧紛謨ｰ,a縺梧枚蟄�,f9,6縺悟ｰ乗焚.1pe縺ｨ縺九�E縺ｮ�贋ｹ�
!謗･轤ｹ逡ｪ蜿ｷ縺ｨ蠎ｧ讓吶�蜃ｺ蜉�
write(7,'("*NODE, NSET=Nall")')
do i=1,(nz+1)*(2*(nx1+nx2*2+1)*(ny1+1)+2*(nx2+1)*(ny2-1))
write(7,'(i6,a,f9.6,a,f9.6,a,f9.6)'),&
&i,',',x(i),',',y(i),',',z(i)
nsettensuu=i
end do

!隕∫ｴ�逡ｪ蜿ｷ縺ｨ縺昴�隕∫ｴ�縺ｮ8謗･轤ｹ縺ｮ謗･轤ｹ逡ｪ蜿ｷ
!荳九→荳企Κ蛻�惠譚舌�隕∫ｴ�
nyou=0
do m=0,1
do k=0,nz-1
do j=0,ny1-1
do i=1,nx1
nyou=nyou+1
meshmoku(nyou,1)=nx2+i+(2*nx2+nx1+1)*j&
&+k*((2*nx2+nx1+1)*(ny1+1)*2+2*(nx2+1)*(ny2-1))&
&+m*((2*nx2+nx1+1)*(ny1+1)+2*(nx2+1)*(ny2-1))
meshmoku(nyou,2)=meshmoku(nyou,1)+(2*nx2+nx1+1)
meshmoku(nyou,3)=meshmoku(nyou,2)+(2*nx2+nx1+1)*(ny1+1)*2+2*(nx2+1)*(ny2-1)
meshmoku(nyou,4)=meshmoku(nyou,1)+(2*nx2+nx1+1)*(ny1+1)*2+2*(nx2+1)*(ny2-1)
meshmoku(nyou,5)=meshmoku(nyou,1)+1
meshmoku(nyou,6)=meshmoku(nyou,5)+(2*nx2+nx1+1)
meshmoku(nyou,7)=meshmoku(nyou,6)+(2*nx2+nx1+1)*(ny1+1)*2+2*(nx2+1)*(ny2-1)
meshmoku(nyou,8)=meshmoku(nyou,5)+(2*nx2+nx1+1)*(ny1+1)*2+2*(nx2+1)*(ny2-1)
end do
end do
end do
end do
!荳企Κ蛻�→荳矩Κ蛻��驪ｼ譚ｿ
nyou=0
do l=0,1
do m=0,1
do k=0,nz-1
do j=0,ny1-1
do i=1,nx2
nyou=nyou+1
meshkou(nyou,1)=i+(2*nx2+nx1+1)*j&
&+k*((2*nx2+nx1+1)*(ny1+1)*2+2*(nx2+1)*(ny2-1))&
&+m*((2*nx2+nx1+1)*(ny1+1)+2*(nx2+1)*(ny2-1))&
&+l*(nx2+nx1)
meshkou(nyou,2)=meshkou(nyou,1)+(2*nx2+nx1+1)
meshkou(nyou,3)=meshkou(nyou,2)+(2*nx2+nx1+1)*(ny1+1)*2+2*(nx2+1)*(ny2-1)
meshkou(nyou,4)=meshkou(nyou,1)+(2*nx2+nx1+1)*(ny1+1)*2+2*(nx2+1)*(ny2-1)
meshkou(nyou,5)=meshkou(nyou,1)+1
meshkou(nyou,6)=meshkou(nyou,5)+(2*nx2+nx1+1)
meshkou(nyou,7)=meshkou(nyou,6)+(2*nx2+nx1+1)*(ny1+1)*2+2*(nx2+1)*(ny2-1)
meshkou(nyou,8)=meshkou(nyou,5)+(2*nx2+nx1+1)*(ny1+1)*2+2*(nx2+1)*(ny2-1)
end do
end do
end do
end do
end do
!逵溘ｓ荳ｭ縺ｮ驪ｼ譚ｿ(譛蛻昴�荳谿ｵ逶ｮ)
do m=0,1
do k=0,nz-1
do i=1,nx2
nyou=nyou+1
meshkou(nyou,1)=i+(nx2*2+nx1+1)*ny1&
&+k*((nx2*2+1+nx1)*(ny1+1)*2+2*(nx2+1)*(ny2-1))&
&+m*(nx2+nx1)
meshkou(nyou,2)=(2*nx2+nx1+1)*(ny1+1)+i&
&+k*((nx2*2+1+nx1)*(ny1+1)*2+2*(nx2+1)*(ny2-1))&
&+m*(nx2+1)
meshkou(nyou,3)=meshkou(nyou,2)+(nx2*2+1+nx1)*(ny1+1)*2+2*(nx2+1)*(ny2-1)
meshkou(nyou,4)=meshkou(nyou,1)+(nx2*2+1+nx1)*(ny1+1)*2+2*(nx2+1)*(ny2-1)
meshkou(nyou,5)=meshkou(nyou,1)+1
meshkou(nyou,6)=meshkou(nyou,2)+1
meshkou(nyou,7)=meshkou(nyou,3)+1
meshkou(nyou,8)=meshkou(nyou,4)+1
end do
end do
end do
!逵溘ｓ荳ｭ縺ｮ驪ｼ譚ｿ(2谿ｵ逶ｮ縺九ｉ譛蠕後�蜑阪∪縺ｧ)
do m=0,1
do k=0,nz-1
do j=0,ny2-3!ny2>3
do i=1,nx2
nyou=nyou+1
meshkou(nyou,1)=(2*nx2+nx1+1)*(ny1+1)+i&
&+k*((nx2*2+1+nx1)*(ny1+1)*2+2*(nx2+1)*(ny2-1))&
&+m*(nx2+1)+j*(nx2+1+nx2+1)
meshkou(nyou,2)=meshkou(nyou,1)+nx2+1+nx2+1
meshkou(nyou,3)=meshkou(nyou,2)+(nx2+1)*(ny2-1)*2+(nx2*2+1+nx1)*(ny1+1)*2
meshkou(nyou,4)=meshkou(nyou,1)+(nx2+1)*(ny2-1)*2+(nx2*2+1+nx1)*(ny1+1)*2
meshkou(nyou,5)=meshkou(nyou,1)+1
meshkou(nyou,6)=meshkou(nyou,2)+1
meshkou(nyou,7)=meshkou(nyou,3)+1
meshkou(nyou,8)=meshkou(nyou,4)+1
end do
end do
end do
end do
!
!逵溘ｓ荳ｭ縺ｮ驪ｼ譚ｿ(譛蠕後�谿ｵ)
do m=0,1
do k=0,nz-1
do i=1,nx2
nyou=nyou+1
meshkou(nyou,1)=(nx2*2+nx1+1)*(ny1+1)+(nx2+1)*(ny2-2)*2+i&
&+k*((nx2*2+1+nx1)*(ny1+1)*2+2*(nx2+1)*(ny2-1))&
&+m*(nx2+1)
meshkou(nyou,2)=(nx2*2+nx1+1)*(ny1+1)+(nx2+1)*(ny2-1)*2+i&
&+k*((nx2*2+1+nx1)*(ny1+1)*2+2*(nx2+1)*(ny2-1))&
&+m*(nx2+nx1)
meshkou(nyou,3)=meshkou(nyou,2)+(nx2*2+nx1+1)*2*(ny1+1)+(nx2+1)*(ny2-1)*2
meshkou(nyou,4)=meshkou(nyou,1)+(nx2*2+nx1+1)*2*(ny1+1)+(nx2+1)*(ny2-1)*2
meshkou(nyou,5)=meshkou(nyou,1)+1
meshkou(nyou,6)=meshkou(nyou,2)+1
meshkou(nyou,7)=meshkou(nyou,3)+1
meshkou(nyou,8)=meshkou(nyou,4)+1
end do
end do
end do
write(7,'("*ELEMENT, TYPE=C3D8, ELSET=Ekou")')
do i=1,2*nx2*(ny1*2+ny2)*nz
write(7,'(i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6)'),&
&i,',',&
&meshkou(i,1),',',&
&meshkou(i,2),',',&
&meshkou(i,3),',',&
&meshkou(i,4),',',&
&meshkou(i,5),',',&
&meshkou(i,6),',',&
&meshkou(i,7),',',&
&meshkou(i,8)
end do
nyousosuukou=2*nx2*(ny1*2+ny2)*nz
write(7,'("*ELEMENT, TYPE=C3D8, ELSET=Emoku")')
do i=1,2*nx1*ny1*nz
write(7,'(i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6,a,i6)'),&
&i+nyousosuukou,',',&
&meshmoku(i,1),',',&
&meshmoku(i,2),',',&
&meshmoku(i,3),',',&
&meshmoku(i,4),',',&
&meshmoku(i,5),',',&
&meshmoku(i,6),',',&
&meshmoku(i,7),',',&
&meshmoku(i,8)
end do

!霈芽差謗･轤ｹ縲∵鋸譚滓磁轤ｹ縺ｮ險ｭ螳�
write(7,'("*NSET,NSET=kousoku")')
do i=1,(nx2*2+nx1+1)*2*(ny1+1)+(nx2+1)*(ny2-1)*2
nkousoku=i
write(7,'(i6,a)'),nkousoku,','
end do
write(7,'("*NSET,NSET=saika")')
nsaikasettensuu=0
do i=1,nsettensuu    !遽轤ｹ謨ｰ
if (ell-0.00002<z(i).and.z(i)<ell+0.00002)then
write(7,'(i6,a)'),i,','
nsaikasettensuu=nsaikasettensuu+1
end if
end do
write(7,'("*NSET,NSET=saikasita")')
do i=1,nsettensuu    !遽轤ｹ謨ｰ
if (ell-0.00002<z(i).and.z(i)<ell+0.00002&
    &.and.y(i)>-0.00001d0.and.y(i)<0.00001d0&
        &.and.b2+0.0001d0<x(i).and.x(i)<b1+b2-0.0001d0)then
write(7,'(i6,a)'),i,','
end if
end do
! 蠅�阜譚｡莉ｶ
write(7,'("*BOUNDARY")')
write(7,'("kousoku,1,3")')
! 譛ｨ譚舌�險ｭ螳�
write(7,'("*MATERIAL,NAME=ELmoku")')
write(7,'("*ELASTIC,TYPE=ORTHOROPIC")')
call ortho(e66,exx,eyy,ezz,poixy,poixz,poiyz,gyz,&
& gxz,gxy,poiyx,poizx,poizy)
write(7,'(f12.4,a,f12.4,a,f12.4,a,&
&f12.4,a,f12.4,a,f12.4,a,f12.4,&
& a ,f12.4,a)'),&
&e66(1,1),',',e66(1,2),',',e66(2,2),',',e66(1,3),','&
&,e66(2,3),',',e66(3,3),',',e66(6,6),',',e66(5,5),','
write(7,'(f12.4,a,a)'),e66(4,4),',','293'
write(7,'("*MATERIAL,NAME=ELkou")')
write(7,'("*ELASTIC")')
write(7,'("206.d3, 0.3")')
write(7,'("*SOLID SECTION,ELSET=Emoku,MATERIAL=ELmoku")')
write(7,'("*SOLID SECTION,ELSET=Ekou,MATERIAL=ELkou")')
! 霈芽差譚｡莉ｶ
write(7,'("*STEP")')
write(7,'("*STATIC,SOLVER=SPOOLES")')
write(7,'("*CLOAD")')
write(7,'(A,a,i6,a,1pe13.5)'),'saika',',',2,',',p/nsaikasettensuu
! dat縺ｮ繧ｪ繝励す繝ｧ繝ｳ
write(7,'("*NODE print,NSET=saika")')
write(7,'("U")')
write(7,'("*NODE FILE")')
write(7,'("U")')
write(7,'("*END STEP")')
!
close (7)
!
end do
!

!job縺ｮ險ｭ螳�
do i=100,330,10
write(filename2, '("ccx_2.5 hako15-ell", i3.3,">hako15-ell", i3.3,".kekka" )') i,i
open  (8,file='job')
write(8,*) filename2
end do
close(8)
end
!######## 逶ｴ莠､逡ｰ譁ｹ諤ｧ縺ｮ螳壽焚繧貞�縺吶�繝ｭ繧ｰ繝ｩ繝� ########
subroutine ortho(e66,exx,eyy,ezz,poixy,poixz,poiyz,gyz,&
&gxz,gxy,poiyx,poizx,poizy)
implicit real*8(a-h,o-z)
dimension e66(6,6),b(3,3)
!蠢懷鴨�阪�縺壹∩陦悟��啼66
do i=1,6
do j=1,6
e66(i,j)=0.d0
end do
end do
!
e66(1,1)=1.d0/exx
e66(1,2)=-poixy/exx
e66(1,3)=-poixz/exx
e66(2,1)=-poiyx/eyy
e66(2,2)=1.d0/eyy
e66(2,3)=-poiyz/eyy
e66(3,1)=-poizx/ezz
e66(3,2)=-poizy/ezz
e66(3,3)=1.d0/ezz
!
!****** 騾�｡悟�繧貞叙縺｣縺ｦ縲∝ｿ懷鴨�阪�縺壹∩陦悟�縺ｫ螟画鋤縺吶ｋ******
det=e66(1,1)*e66(2,2)*e66(3,3)+e66(2,1)*e66(3,2)*e66(1,3)&
&+e66(3,1)*e66(1,2)*e66(2,3)-e66(2,1)*e66(1,2)*e66(3,3)&
&-e66(3,1)*e66(2,2)*e66(1,3)-e66(1,1)*e66(3,2)*e66(2,3)
b(1,1)=(-1.d0)**(1+1)*(e66(2,2)*e66(3,3)-e66(3,2)*e66(2,3))
b(1,2)=(-1.d0)**(1+2)*(e66(2,1)*e66(3,3)-e66(3,1)*e66(2,3))
b(1,3)=(-1.d0)**(1+3)*(e66(2,1)*e66(3,2)-e66(3,1)*e66(2,2))
b(2,1)=(-1.d0)**(2+1)*(e66(1,2)*e66(3,3)-e66(3,2)*e66(1,3))
b(2,2)=(-1.d0)**(2+2)*(e66(1,1)*e66(3,3)-e66(3,1)*e66(1,3))
b(2,3)=(-1.d0)**(2+3)*(e66(1,1)*e66(3,2)-e66(3,1)*e66(1,2))
b(3,1)=(-1.d0)**(3+1)*(e66(1,2)*e66(2,3)-e66(2,2)*e66(1,3))
b(3,2)=(-1.d0)**(3+2)*(e66(1,1)*e66(2,3)-e66(2,1)*e66(1,3))
b(3,3)=(-1.d0)**(3+3)*(e66(1,1)*e66(2,2)-e66(2,1)*e66(1,2))
do i=1,3
do j=1,3
e66(i,j)=b(j,i)/det
end do
end do
!
!****** 縺帙ｓ譁ｭ陬懈ｭ｣菫よ焚繧貞�蜉� ******
!print*,'xy髱｢縺ｮ縺帙ｓ譁ｭ蠑ｾ諤ｧ菫よ焚繧貞�繧後※荳九＆縺�'
!read*, gxy
e66(4,4)=gyz
!print*,'xz髱｢縺ｮ縺帙ｓ譁ｭ蠑ｾ諤ｧ菫よ焚繧貞�繧後※荳九＆縺�'
!read*, gxz
e66(5,5)=gxz
!print*,'yz髱｢縺ｮ縺帙ｓ譁ｭ蠑ｾ諤ｧ菫よ焚繧貞�繧後※荳九＆縺�'
!read*, gyz
e66(6,6)=gxy
!
!
return
end