program f90_FEM_DOMAIN implicit none integer,parameter::ntmax=10 integer,parameter::mdmax=100 integer::NODT,NELT,MATEL,nod integer,allocatable::kakom(:,:) integer,allocatable::matno(:) real(8),allocatable::x(:),y(:) real(8),allocatable::xg(:),yg(:) integer::nt(1:ntmax) real(8):: x_dom(1:ntmax,1:mdmax),y_dom(1:ntmax,1:mdmax) integer::i,j,ne real(8)::sumx,sumy character(len=50)::fnameR1 character(len=50)::fnameR2 character(len=50)::fnameW character(len=50)::dummy call getarg(1,dummy) call getarg(2,fnameR1) call getarg(3,fnameR2) call getarg(4,fnameW) read(dummy,*) nod ! Input of Mesh data open(11,file=fnameR1,status='old') read(11,*) NODT,NELT allocate(kakom(1:NELT,1:nod),matno(1:NELT)) allocate(x(1:NODT),y(1:NODT)) allocate(xg(1:NELT),yg(1:NELT)) do ne=1,NELT read(11,*) (kakom(ne,j),j=1,nod) end do do i=1,NODT read(11,*) x(i),y(i) end do close(11) do ne=1,NELT sumx=0.0D0 sumy=0.0D0 do i=1,nod sumx=sumx+x(kakom(ne,i)) sumy=sumy+y(kakom(ne,i)) end do xg(ne)=sumx/dble(nod) yg(ne)=sumy/dble(nod) end do ! Input of Domain data open(11,file=fnameR2,status='old') read(11,*) MATEL do i=1,MATEL read(11,*) nt(i) do j=1,nt(i) read(11,*) x_dom(i,j),y_dom(i,j) end do end do close(11) ! Definition of Domain (material No.) call DEFDOM(NELT,MATEL,nt,x_dom,y_dom,xg,yg,matno,ntmax,mdmax) open(12, file=fnameW, status='replace') write(12,'(4(i6))') nod,NODT,NELT,MATEL if(nod==3)then do ne=1,NELT write(12,'(4(i6))') (kakom(ne,j),j=1,nod),matno(ne) end do end if if(nod==4)then do ne=1,NELT write(12,'(5(i6))') (kakom(ne,j),j=1,nod),matno(ne) end do end if do i=1,NODT write(12,'(2(e15.7))') x(i),y(i) end do close(12) do i=1,MATEL dummy(1:50)='' write (dummy,'("_mat_",i2.2,".txt")') i open(12,file=dummy,status='replace') do ne=1,NELT if(matno(ne)==i)then write(12,'(a)') '>' do j=1,nod write(12,'(2(e15.7))') x(kakom(ne,j)),y(kakom(ne,j)) end do end if end do close(12) end do stop contains subroutine DEFDOM(NELT,MATEL,nt,x_dom,y_dom,xg,yg,matno,ntmax,mdmax) integer,intent(in)::NELT,MATEL,ntmax,mdmax integer,intent(in)::nt(1:ntmax) real(8),intent(in):: x_dom(1:ntmax,1:mdmax),y_dom(1:ntmax,1:mdmax) real(8),intent(in)::xg(1:NELT),yg(1:NELT) integer,intent(out)::matno(1:NELT) real(8)::theta,dx1,dy1,dx2,dy2,r1,r2,cs,sn integer::i,j,ne,j1,j2 real(8),parameter::pi=3.14159265358979323846D0 real(8),parameter::eps=1D-10 do ne=1,NELT matno(ne)=0 do i=1,MATEL theta=0.0D0 do j=1,nt(i) j1=j j2=j+1 if(j1==nt(i))j2=1 dx1=x_dom(i,j1)-xg(ne) dy1=y_dom(i,j1)-yg(ne) dx2=x_dom(i,j2)-xg(ne) dy2=y_dom(i,j2)-yg(ne) r1=sqrt(dx1*dx1+dy1*dy1) r2=sqrt(dx2*dx2+dy2*dy2) if(eps