!ipk last update sep 20 2013 add more output of progress and flushing of messages SUBROUTINE ADDTOMESH(IADDFIL,ISWT) ! iswt = 0 ADD TO MESH ! ISWT = 1 MERGE MESHES USE WINTERACTER USE BLK1MOD USE BLK2MOD INCLUDE 'D.INC' ! INCLUDE 'BLK1.COM' INCLUDE 'BFILES.I90' IADD=IADDFIL+50 CALL RDTOCLIP(IADD) IF(ISWT .EQ. 1) THEN CALL OUTLINES(1) ISWT1=0 ! IF(NOUTLST(2) .EQ. 0) THEN ISWT2=1 ! ELSE ! ISWT2=0 ! ENDIF CALL MERGEMESH1(ISWT1,ISWT2) write(90,*) 'finished mergemesh1' IF(ISWT2 .EQ. 0) CALL MERGEMESH ! CALL MERGEMESH write(90,*) 'finished mergemesh' flush(90) ENDIF CALL ADDMESH(0) write(90,*) 'finished addmesh' IF(ISWT .EQ. 1 ) THEN CALL WMessageBox(YesNo,QuestionIcon,CommonOK,'Do you wish to delete unused nodes?'//& CHAR(13)//' ','Delete unused nodes?') ! ! If answer 'No', return ! IF (WInfoDialog(4).EQ.2) return ! ! Delete all unused nodes ! CALL DELETM(2) ENDIF RETURN END SUBROUTINE RDTOCLIP(IUNIT) USE BLK1MOD ! INCLUDE 'BLK1.COM' CHARACTER*80 ALINE REWIND IUNIT READ(IUNIT) TITLE,NPSTO(1),NESTO(1) WRITE(90,*) 'IN RDTOCLIP',IUNIT WRITE(90,*) TITLE,NPSTO(1),NESTO(1) READ(IUNIT) ISLP,IPRT,IPNN,IPEN,IPO,IRO,IPP,IRFN & & ,IGEN,NXZL,NITST,ISCTXT,IFILL,IALTGM,NLAYD,xadded,yadded,ntempinc WRITE(90,*) ISLP,IPRT,IPNN,IPEN,IPO,IRO,IPP,IRFN & & ,IGEN,NXZL,NITST,ISCTXT,IFILL,IALTGM,NLAYD,xadded,yadded,ntempinc READ(IUNIT) HORIZ,VERT,XSALE,YSALE,XFACT,YFACT,AR,ANG WRITE(90,*) HORIZ,VERT,XSALE,YSALE,XFACT,YFACT,AR,ANG IF(IPP .GT. 0) READ(IIN) ALINE READ(IUNIT) ((NOPSTO(J,K,1),K=1,8),IMATSTO(J,1),THTASTO(J,1),J=1,NESTO(1)) READ(IUNIT) & & (XUSRSTO(J,1),YUSRSTO(J,1),WDSTO(J,1),WIDTHSTO(J,1),SS1STO(J,1),SS2STO(J,1),WIDSSTO(J,1), & & WIDBSSTO(J,1),SSOSTO(J,1),BS1STO(J,1),J=1,NPSTO(1)) READ(IUNIT) NLSTSTO(1) IF(NLSTSTO(1) .GT. 0) THEN READ(IUNIT) (LLISTSTO(J,1),J=1,NLSTSTO(1)), & ((ILISTSTO(J,I,1),I=1,LLISTSTO(J,1)),J=1,NLSTSTO(1)) ENDIF READ(IUNIT) NENTRYC,NLAYDC,NCLMSTO(1) IF(NENTRYC .GT. 0) THEN READ(IUNIT) ((NEFC,J=1,3),I=1,NENTRYC) ENDIF IF(NLAYDC .GT. 0) THEN READ(IUNIT) (LAYC,I=1,NPSTO(1)) ENDIF IF(NCLMSTO(1) .GT. 0) THEN READ(IUNIT) ((ICCLNSTO(I,J,1),J=1,350),I=1,NCLMSTO(1)) ENDIF REWIND IUNIT RETURN END SUBROUTINE ADDMESH(ISWT) USE BLK1MOD ! INCLUDE 'BLK1.COM' ALLOCATABLE NODETRAN(:) DATA VDX9/-9.E9/ ! Loop through nodes assigning new number and adding to list IF(.NOT. ALLOCATED(NODETRAN)) ALLOCATE (NODETRAN(maxp)) IF(ISWT .EQ. 0) THEN DO N=1,NPSTO(1) IF(XUSRSTO(N,1) .GT. VDX9) THEN CALL GETNOD(J) NODETRAN(N)=J XUSR(J)=XUSRSTO(N,1) YUSR(J)=YUSRSTO(N,1) WD(J)=WDSTO(N,1) WIDTH(J)=WIDTHSTO(N,1) SS1(J)=SS1STO(N,1) SS2(J)=SS2STO(N,1) WIDS(J)=WIDSSTO(N,1) WIDBS(J)=WIDBSSTO(N,1) SSO(J)=SSOSTO(N,1) BS1(J)=BS1STO(N,1) INSKP(J) = 0 INEW(J) = 1 ENDIF ENDDO ELSE DO N=1,NPSTO(1) NODETRAN(N)=N ENDDO ENDIF ! Loop through elements assigning new number and adding to list DO N=1,NESTO(1) IF(IMATSTO(N,1) .GT. 0) THEN CALL GETELM(M) DO K=1,8 IF(NOPSTO(N,K,1) .GT. 0) THEN J=NODETRAN(NOPSTO(N,K,1)) NOP(M,K)=J ELSE NOP(M,K)=0 ENDIF ENDDO IMAT(M)=IMATSTO(N,1) THTA(M)=THTASTO(N,1) IESKP(M)=0 NCN = 2 IF (NOP(M,3) .NE. 0) NCN = 3 IF (NOP(M,4) .NE. 0) NCN = 4 IF (NOP(M,5) .NE. 0 .AND. NOP(M,4) .NE. 0) NCN = 5 IF (NOP(M,5) .NE. 0 .AND. NOP(M,4) .EQ. 0) NCN = 6 IF (NOP(M,6) .NE. 0) NCN = 6 IF (NOP(M,7) .NE. 0) NCN = 8 NCORN(M) = NCN ENDIF ENDDO if(iswt .eq. 0) CALL RESCAL CALL HEDR RETURN END SUBROUTINE MERGEMESH1(ISWT1,ISWT2) USE BLK1MOD USE BLK2MOD USE WINTERACTER ! INCLUDE 'BLK1.COM' REAL*8 ELXMIN,ELXMAX,ELYMIN,ELYMAX,XLC,YLC,XXX,YYY LOGICAL LSTAT ALLOCATABLE ELXMIN(:),ELXMAX(:),ELYMIN(:),ELYMAX(:),KEY(:),NKEY(:) DIMENSION XOUT1(1000),YOUT1(1000) IF(.NOT. ALLOCATED(ELXMIN)) & ALLOCATE (ELXMIN(MAXE),ELXMAX(MAXE),ELYMIN(MAXE),ELYMAX(MAXE),KEY(MAXE),NKEY(MAXP)) IF(ISWT2 .EQ. 0) GO TO 110 ! first eliminate any elements inside outline CALL KCONST(0) NKEP=0 DO K=1,10 IF(NOUTLST(K) .LE. 0) THEN DO J=1,NPSTO(1) XXXX=XUSRSTO(J,1) YYYY=YUSRSTO(J,1) LSTAT=IGrInsidePolygon(XOUT(1,K),YOUT(1,K),-NOUTLST(K),XXXX,YYYY) IF(LSTAT) THEN NKEP(J)=1 ENDIF ENDDO ENDIF ENDDO DO K=1,10 IF(NOUTLST(K) .GT. 0) THEN DO J=1,NPSTO(1) IF(NKEP(J) .EQ. 1) CYCLE XXXX=XUSRSTO(J,1) YYYY=YUSRSTO(J,1) ! WRITE(155,*) J,XXXX,YYYY LSTAT=IGrInsidePolygon(XOUT(1,K),YOUT(1,K),NOUTLST(K),XXXX,YYYY) ! WRITE(155,*) J,LSTAT IF(LSTAT) THEN DO L=1,NDELM(J) NCAN=NECON(J,L) CALL DELEM(NCAN) ENDDO ENDIF ENDDO ENDIF 100 CONTINUE ENDDO IF(ISWT2 .EQ. 1) RETURN ! First sort coordinates for min of element connection ! List all limiting values 110 CONTINUE DO N=1,NE IF(IMAT(N) .GT. 0) THEN ELXMIN(N)=XUSR(NOP(N,1)) ELXMAX(N)=XUSR(NOP(N,1)) ELYMIN(N)=YUSR(NOP(N,1)) ELYMAX(N)=YUSR(NOP(N,1)) DO M=2,8 IF(NOP(N,M) .NE. 0) THEN ELXMIN(N)=MIN(ELXMIN(N),XUSR(NOP(N,M))) ELXMAX(N)=MAX(ELXMAX(N),XUSR(NOP(N,M))) ELYMIN(N)=MIN(ELYMIN(N),YUSR(NOP(N,M))) ELYMAX(N)=MAX(ELYMAX(N),YUSR(NOP(N,M))) ENDIF ENDDO ELSE ELXMIN(N)=VOID ELXMAX(N)=VOID ELYMIN(N)=VOID ELYMAX(N)=VOID ENDIF ENDDO CALL SORTDB(XUSRSTO,NKEY,NPSTO(1)) CALL SORTDB(ELXMIN,KEY,NE) ! Loop on elements to check for overlap DO KK=1,NESTO(1) IF (NOPSTO(KK,6,1) .EQ. 0) CYCLE IF(IMATSTO(KK,1) .GT. 0) THEN if(mod(kk,1000) .eq. 0) write(90,*) 'merged',kk flush(90) KL=1 200 CONTINUE IF(ISWT1 .EQ. 0) THEN DO K=KL,8 J=NOPSTO(KK,K,1) IF(J .GT. 0) THEN KLL=KL XXX=XUSRSTO(J,1) YYY=YUSRSTO(J,1) GO TO 220 ENDIF ENDDO KLL=8 GO TO 400 220 CONTINUE ELSE XXX=0. YYY=0. DO K=1,7,2 JJ=NOPSTO(KK,K,1) IF(JJ .GT. 0) THEN XXX=XXX+XUSRSTO(JJ,1) YYY=YYY+YUSRSTO(JJ,1) ENDIF ENDDO IF(JJ .EQ. 0) THEN XXX=XXX/3. YYY=YYY/3. ELSE XXX=XXX/4. YYY=YYY/4. ENDIF ENDIF ! Search on elements to find a startin point DO NN=1,NE N=KEY(NN) IF(IMAT(N) .GT. 0) THEN !- !...... DETERMINE ELEMENT TYPE !- NCN=8 IT=1 IF(NOP(N,7) .EQ. 0) THEN NCN=6 IT=2 ENDIF IF(NOP(N,6) .EQ. 0) THEN GOTO 350 ENDIF ! Test for point inside an element ! Test for max and min within IF(XXX .GT. ELXMIN(N)) THEN IF(XXX .GT. ELXMAX(N)) GO TO 350 IF(YYY .GT. ELYMIN(N)) THEN IF(YYY .GT. ELYMAX(N)) GO TO 350 ! Now get local coordinate as final test CALL GPTEV(N,XXX,YYY,XLC,YLC,IT,NCN) IF(IT .EQ. 2) THEN IF(XLC .LT. 0. .OR. YLC .LT. 0. .OR. XLC+YLC .GT. 1.) THEN GO TO 350 ELSE CALL DELEM(KK) GO TO 400 ENDIF ELSE IF(XLC .LT. -1. .OR. YLC .LT. -1. .OR. & XLC .GT. 1. .OR. YLC .GT. 1.) THEN GO TO 350 ELSE CALL DELEM(KK) GO TO 400 ENDIF ENDIF ENDIF ENDIF ENDIF 350 CONTINUE ENDDO KL=KLL+1 IF(KL .LT. 8 .AND. ISWT1 .EQ. 0) GO TO 200 ENDIF ! Finished test 400 CONTINUE ENDDO RETURN END SUBROUTINE GPTEV(N,XSW,YSW,XG,YG,IT,NCN) !- !......EVALUATE FUNCTION AT GRID POINTS !- !- N = ELEMENT NUMBER !_ XSW = X COORDINATE OF DESIRED POINT !_ YSW = Y COORDINATE OF DESIRED POINT ! XG = X LOCAL COORDINATE ! YG = Y LOCAL COORDINATE ! IT = SWITCH FOR CHOICE BETWEEN LINEAR AND QUADRATIC WEIGHTING ! = 1 FOR LINEAR ! = 2 FOR QUADRATIC ! FROM COMMON ! NOP = LIST OF NODAL CONNECTIONS AROUND AN ELEMET ! XUSR = REAL*8 ARRAY OF NODAL COORDINATES ! USE BLK1MOD ! INCLUDE 'BLK1.COM' REAL*8 XN,DNX,DNY,XSW,YSW DOUBLE PRECISION XG,YG,XK,YK,XP,YP !- DIMENSION X(9),Y(9),WGT(8) !- DATA TOL/0.01/ !- !- !......ESTABLISH LOCAL COORDINATES FOR EACH NODE POINT OF ELEMENT !- K1=NOP(N,1) X(1)=0. Y(1)=0. DO 300 K=3,NCN,2 K2=NOP(N,K) X(K)=XUSR(K2)-XUSR(K1) Y(K)=YUSR(K2)-YUSR(K1) 300 END DO X(2)=X(3)/2. Y(2)=Y(3)/2. X(4)=(X(3)+X(5))/2. Y(4)=(Y(3)+Y(5))/2. IF(IT .EQ. 2) THEN X(6)=X(5)/2. Y(6)=Y(5)/2. xminl=min(x(1),x(3),x(5)) yminl=min(y(1),y(3),y(5)) xmaxl=max(x(1),x(3),x(5)) ymaxl=max(y(1),y(3),y(5)) ELSE X(6)=(X(5)+X(7))/2. Y(6)=(Y(5)+Y(7))/2. X(8)=X(7)/2. Y(8)=Y(7)/2. xminl=min(x(1),x(3),x(5),x(7)) yminl=min(y(1),y(3),y(5),y(7)) xmaxl=max(x(1),x(3),x(5),x(7)) ymaxl=max(y(1),y(3),y(5),y(7)) ENDIF !- !......ESTABLISH LOCAL COORDINATES OF DESIRED POINT !- XP=XSW-XUSR(K1) YP=YSW-YUSR(K1) XG=0. YG=0. !- !......ITERATE TO FIND LOCAL COORDINATE !- DO ITER=1,10 DXKDX=0. DXKDY=0. DYKDX=0. DYKDY=0. XK=-XP YK=-YP DO K=2,NCN XK=XK+XN(IT,K,XG,YG)*X(K) YK=YK+XN(IT,K,XG,YG)*Y(K) DXKDX=DXKDX+DNX(IT,K,XG,YG)*X(K) DYKDX=DYKDX+DNX(IT,K,XG,YG)*Y(K) DXKDY=DXKDY+DNY(IT,K,XG,YG)*X(K) DYKDY=DYKDY+DNY(IT,K,XG,YG)*Y(K) END DO DET=DXKDX*DYKDY-DXKDY*DYKDX DX=(-DYKDY*XK+DXKDY*YK)/DET DY=( DYKDX*XK-DXKDX*YK)/DET XG=XG+DX YG=YG+DY IF(ABS(DX).LT.TOL .AND. ABS(DY).LT.TOL) GO TO 420 END DO !- !......NOW GET WEIGHTING FUNCTIONS FOR QUAD FUNCTION !- 420 CONTINUE RETURN END SUBROUTINE DELEM(J) ! USE BLK1MOD ! INCLUDE 'BLK1.COM' ! !- !......DELETE ELEMENT ! ! Search for elements that attach to node J and remove them ! IMATSTO(J,1)=0 DO KK=1,8 NOPSTO(J,KK,1)=0 ENDDO ! RETURN END SUBROUTINE MERGEMESH USE BLK1MOD LOGICAL LSTAT ! INCLUDE 'BLK1.COM' ! Loop on element to be added DO N=1,NESTO(1) IF(IMATSTO(N,1) .NE. 0) THEN if(mod(n,1000) .eq. 0) write(90,*) 'adding',n,nesto(1) flush(90) IF(IMATSTO(N,1) .GT. 900 .AND. IMATSTO(N,1) .LT. 904) THEN X1=XUSRSTO(NOPSTO(N,1,1),1) Y1=YUSRSTO(NOPSTO(N,1,1),1) CALL CHECKIN(X1,Y1,LSTAT) IF(ISTATUS .EQ. 5) THEN CALL DELEM(N) GO TO 400 ENDIF GO TO 400 ENDIF ! loop on sides DO M=1,7,2 N1=NOPSTO(N,M,1) IF(M .EQ. 3 .AND. NOPSTO(N,5,1) .EQ. 0) GO TO 400 IF(N1 .GT. 0) THEN IF((M .EQ. 5 .AND. NOPSTO(N,7,1) .EQ. 0) .OR. (M .EQ. 7)) THEN N2=NOPSTO(N,1,1) ELSE N2=NOPSTO(N,M+2,1) ENDIF IF(NKEP(N1) .EQ. 1 .AND. NKEP(N2) .EQ. 1) GO TO 380 ! Now loop trough existing elements DO I=1,NE IF(IMAT(I) .NE. 0) THEN DO J=1,7,2 M1=NOP(I,J) IF(J .EQ. 3 .AND. NOP(I,5) .EQ. 0) GO TO 360 IF(M1 .GT. 0) THEN IF((J .EQ. 5 .AND. NOP(I,7) .EQ. 0) .OR. (J .EQ. 7)) THEN M2=NOP(I,1) ELSE M2=NOP(I,J+2) ENDIF if(m2 .eq. 0) cycle X1=XUSRSTO(N1,1) X2=XUSRSTO(N2,1) Y1=YUSRSTO(N1,1) Y2=YUSRSTO(N2,1) X3=XUSR(M1) X4=XUSR(M2) Y3=YUSR(M1) Y4=YUSR(M2) CALL IGrIntersectLine(X1,Y1,X2,Y2,X3,Y3,X4,Y4,XINTER,YINTER,ISTATUS) IF(ISTATUS .EQ. 5) THEN CALL DELEM(N) GO TO 400 ENDIF ENDIF ENDDO ENDIF 360 CONTINUE ENDDO ENDIF 380 CONTINUE ENDDO ENDIF 400 CONTINUE ENDDO RETURN END SUBROUTINE CHECKIN(X1,Y1,LSTAT) USE BLK1MOD LOGICAL LSTAT DIMENSION XP(4),YP(4) ! Now loop trough existing elements DO I=1,NE IF(IMAT(I) .NE. 0) THEN JJ=0 DO J=1,7,2 INODE=NOP(I,J) IF(INODE .GT. 0) THEN JJ=JJ+1 XP(JJ)=XUSR(INODE) YP(JJ)=YUSR(INODE) ENDIF ENDDO LSTAT=IGrInsidePolygon(XP,YP,JJ,X1,Y1) IF(LSTAT) RETURN ENDIF ENDDO RETURN END SUBROUTINE KCONST(isw1) ! ! ESTABLISH ELEMENT CONNECTED TO ELEMENT TABLE ! USE BLK1MOD USE BLK2MOD ! INCLUDE 'BLK1.COM' ! INCLUDE 'BLK2.COM' ! ! INITIALIZE ! NCM=11 DO 200 J=1,NCM DO 200 N=1,NPSTO(1) 200 NECON(N,J)=0 DO 230 N=1,NPSTO(1) 230 NDELM(N)=0 ! ! FORM TABLE OF ELEMENTS CONNECTED TO EACH NODE ! DO 300 M=1,NESTO(1) IF(IMATSTO(M,1) .EQ. 0) GO TO 300 if(isw1 .eq. 1) then if(imat(m) .eq. 999) go to 300 endif DO 280 K=1,8 N=NOPSTO(M,K,1) IF (N .GT. 0) THEN NDELM(N)=NDELM(N)+1 J=NDELM(N) NECON(N,J)=M !ipkoct93 ELSE !ipkoct93 GO TO 300 ENDIF 280 CONTINUE 300 END DO RETURN END