1 subroutine findslpw(x,ilwork,iper,iBC) 2 use pointer_data 3 use slpw 4 include "common.h" 5 include "mpif.h" 6 real*8 x(numnp,nsd) 7 integer iper(nshg) 8 integer ilwork(nlwork) 9 integer ifirstvisit(nshg) 10 integer, allocatable :: ienb(:) 11 real*8 elnrm(nsd),wn1(nsd),wn2(nsd),wn3(nsd),aedg(nsd),bedg(nsd) 12 real*8 mag 13 integer iBC(nshg) 14 integer, allocatable :: IDslpw(:) 15 real*8 AA(nsd),BB(nsd),tmp 16 integer iparall 17 18 allocate(idxslpw(nshg)) 19 allocate(mpslpw(nshg,2)) 20 allocate(wlnorm(nshg,nsd,9)) 21 22 23 nslpwnd=0 24 idxslpw=0 25 mpslpw=0 26 wlnorm=0 27 ifirstvisit(:)=1 28 do iblk=1,nelblb 29 npro=lcblkb(1,iblk+1)-lcblkb(1,iblk) 30 nenbl=lcblkb(6,iblk) 31 nshl=lcblkb(9,iblk) 32 allocate(ienb(nshl)) 33 do i=1,npro 34 isfID=miBCB(iblk)%p(i,2) 35 if(isfID.gt.709.or.isfID.lt.701)cycle 36 islpw=mod(isfID,700) 37 ienb(1:nshl)=mienb(iblk)%p(i,1:nshl) 38 do j=1,nenbl 39 nn=ienb(j) 40 if(ifirstvisit(nn).eq.1)then 41 ifirstvisit(nn)=0 42 nslpwnd=nslpwnd+1 43 idxslpw(nslpwnd)=nn 44 endif 45 if(.not.btest(mpslpw(nn,2),islpw))then 46 mpslpw(nn,1)=mpslpw(nn,1)+1 47 mpslpw(nn,2)=mpslpw(nn,2)+2**islpw 48 endif 49 if(j.ne.1.and.j.ne.nenbl)then 50 aedg(:)=x(ienb(j+1),:)-x(nn,:) 51 bedg(:)=x(ienb(j-1),:)-x(nn,:) 52 elseif(j.eq.1)then 53 aedg(:)=x(ienb(j+1),:)-x(nn,:) 54 bedg(:)=x(ienb(nenbl),:)-x(nn,:) 55 elseif(j.eq.nenbl)then 56 aedg(:)=x(ienb(1),:)-x(nn,:) 57 bedg(:)=x(ienb(j-1),:)-x(nn,:) 58 endif 59 elnrm(1)=aedg(2)*bedg(3)-aedg(3)*bedg(2) 60 elnrm(2)=aedg(3)*bedg(1)-aedg(1)*bedg(3) 61 elnrm(3)=aedg(1)*bedg(2)-aedg(2)*bedg(1) 62 wlnorm(nn,:,islpw)=wlnorm(nn,:,islpw)+elnrm(:) 63 enddo 64 enddo 65 deallocate(ienb) 66 enddo 67 68 do k=1,9 69 if(numpe.gt.1) call commu(wlnorm(:,:,k),ilwork,3,'in ') 70 call bc3per(iBC,wlnorm(:,:,k),iper,ilwork,3) 71 if(numpe.gt.1) call commu(wlnorm(:,:,k),ilwork,3,'out') 72 enddo 73cccccccccccccccccccccccccccccccccccccccccccccccccc 74 75 76 do id=1,nslpwnd 77 nn=idxslpw(id) 78 nslpw=mpslpw(nn,1) 79 ihis=mpslpw(nn,2) 80 nslpwg=nslpw 81 82 allocate(IDslpw(nslpw)) 83 jslpw=0 84 do ni=1,nslpw 85 do islpw=1+jslpw,9 86 if(btest(ihis,islpw))exit 87 enddo 88 IDslpw(ni)=islpw 89 jslpw=islpw 90 enddo 91 92 ired=0 93 do i=1,nslpwg 94 if(btest(ired,i))cycle 95 do j=i+1,nslpwg 96 AA(:)=wlnorm(nn,:,IDslpw(i)) 97 BB(:)=wlnorm(nn,:,IDslpw(j)) 98 if(iparall(AA,BB).eq.1)then 99 nslpw=nslpw-1 100 ihis=ihis-2**IDslpw(j) 101 wlnorm(nn,:,IDslpw(i))=AA(:)+BB(:) 102 ired=ired+2**j 103 endif 104 enddo 105 enddo 106 107 do k=1,9 108 wn1(:)=wlnorm(nn,:,k) 109 mag=(wn1(1)**2+wn1(2)**2+wn1(3)**2)**0.5 110 if(mag.ne.0)wlnorm(nn,:,k)=wn1(:)/mag 111 enddo 112 113 mpslpw(nn,1)=nslpw 114 mpslpw(nn,2)=ihis 115 116 deallocate(IDslpw) 117 enddo 118 119 120c call write_field(myrank,'a'//char(0), 121c & 'wnrm'//char(0),4,wlnorm(:,:,1),'d'//char(0),nshg,3,701) 122cc 7 is the number of charaters of the name of array wnrm701 123c call write_field(myrank,'a'//char(0), 124c & 'wnrm'//char(0),4,wlnorm(:,:,2),'d'//char(0),nshg,3,702) 125c stop 126 127 return 128 end 129 130cccccccccccccccccccccc 131 132 integer function iparall(A,B) 133 real*8 A(3),B(3) 134 real*8 angtol,tmp 135 136 tmp=(A(1)**2+A(2)**2+A(3)**2)**0.5 137 A(:)=A(:)/tmp 138 tmp=(B(1)**2+B(2)**2+B(3)**2)**0.5 139 B(:)=B(:)/tmp 140 pi=3.1415926535 141 angtol=10 142 angtol=sin(angtol*(pi/180)) 143 C1=A(2)*B(3)-A(3)*B(2) 144 C2=A(3)*B(1)-A(1)*B(3) 145 C3=A(1)*B(2)-A(2)*B(1) 146 Cmag=(C1**2+C2**2+C3**2)**0.5 147 dotp=A(1)*B(1)+A(2)*B(2)+A(3)*B(3) 148 if(dotp.gt.0.and.Cmag.lt.angtol)then 149 iparall=1 150 else 151 iparall=0 152 endif 153 154 return 155 end 156