C FACED SOURCE PV 18/06/15 21:15:01 9859 C TRACE D'UN OBJET PAR FACE EN COMMENCANT PAR CELLES DE DERRIERE C C imod=0 trace en couleur d'effacement C imod=1 trace en couleir normale C SG 2016/07/18 Programmation comme envvo2 pour gerer les faces TRI7/QUA9 C SUBROUTINE FACED(MELEME,XPROJ,ICPR,IVU,MCOUP,KON,LNDEGR,imod) IMPLICIT INTEGER(I-N) -INC PPARAM -INC CCOPTIO -INC SMELEME -INC CCGEOME -INC CCREEL C+PP DEGRADE OU NON LOGICAL lndegr C+PP SEGMENT XPROJ(3,1) SEGMENT ICPR(1) SEGMENT IVU(ITE) SEGMENT MCOUP(0) SEGMENT KON(3,NBCON,NMAX) * * Type de faces prises en compte: T3, Q4, T6, Q8, POLY, T7, Q9 * Numero dans KDFAC 1 2 3 4 6 7 8 * Pour ne pas se prendre la tête, on numerote pareil que KDFAC * Pour les 5 (non utilisé), 6 (polygone) et >8, ca restera à 0 * NTYFAC=20, codé en dur dans CCGEOME pour KDFAC PARAMETER (NTYFAC=20) * Nb de faces de chaque type, sert également de compteur SEGMENT NBFAC(NTYFAC) * Un segment pointant sur les IFACI SEGMENT IPOFAC(NTYFAC) * Segment IFACI contenant les noeuds, la couleur et si la face d'un * type donné est vue ou non SEGMENT IFACI(NNODE+2,NFACI) * Nombre de noeuds max pour les polygones PARAMETER (NNOMAX=14) * SEGMENT NSOMP(NFACP) SEGMENT IFACOL(NFAC) SEGMENT TFAC(NFAC) SEGMENT KFAK(NFAC) SEGMENT NAUX(max(2,NFAC)) DIMENSION XTR(NNOMAX),YTR(NNOMAX),ZTR(NNOMAX),NPTR(NNOMAX) DIMENSION XTR2(3),YTR2(3),ZTR2(3) * *dbg write(ioimp,*) 'coucou faced lndegr=',LNDEGR,' imod=',imod SEGACT MELEME MELSAU=MELEME IPT1=MELEME * IF (MCOUP.NE.0) THEN * NBNN=0 * NBELEM=0 * NBREF=0 * NBSOUS=LISOUS(/1)-1 * SEGINI IPT1 * DO 5 I=1,NBSOUS * IPT1.LISOUS(I)=LISOUS(I) * 5 CONTINUE * ISCOUP=LISOUS(NBSOUS+1) * ENDIF * call ecmail(ipt1,0) CALL ECROBJ('MAILLAGE',IPT1) CALL ENVELO * IF (MCOUP.NE.0) SEGSUP IPT1 CALL LIROBJ('MAILLAGE',MELEME,1,IRETOU) IF (IERR.NE.0) RETURN SEGACT MELEME * IF (MCOUP.NE.0) THEN * IF (LISOUS(/1).NE.0) THEN * NBSOUS=LISOUS(/1)+1 * SEGADJ MELEME * LISOUS(NBSOUS)=ISCOUP * ELSE * NBSOUS=2 * SEGINI IPT1 * IPT1.LISOUS(1)=MELEME * IPT1.LISOUS(2)=ISCOUP * MELEME=IPT1 * ENDIF * ENDIF SEGACT XPROJ,ICPR SEGACT KON*MOD NBCON=KON(/2) NBCONR=NBCON-1 NBPOIN=XPROJ(/2) TMIN=xsgran TMAX=-xsgran DO 1 I=1,NBPOIN TMIN=MIN(TMIN,XPROJ(3,I)) TMAX=MAX(TMAX,XPROJ(3,I)) 1 CONTINUE TDIST=TMAX-TMIN c c on compte le nombre d elements dont les faces sont de type 1 2 3 4 c 6 7 8 dans NBFAC, attention à 6 : gestion des polygones SEGINI NBFAC IPT1=MELEME SEGACT MELEME DO 10 IOB=1,MAX(1,LISOUS(/1)) IF (LISOUS(/1).NE.0) THEN IPT1=LISOUS(IOB) SEGACT IPT1 ENDIF NBELEM=IPT1.NUM(/2) ILTEL=LTEL(1,IPT1.ITYPEL) IF (ILTEL.EQ.0) GOTO 12 ILTAD=LTEL(2,IPT1.ITYPEL) DO 13 IF=1,ILTEL IFT=LDEL(1,ILTAD+IF-1) NBFAC(IFT)=NBFAC(IFT)+NBELEM 13 CONTINUE 12 CONTINUE 10 CONTINUE * NFACP=NBFAC(6) SEGINI NSOMP c==== CREATION DES FACES ============================================== * Initialisation des IFACI SEGINI IPOFAC DO ITYFAC=1,NTYFAC NNODE=KDFAC(1,ITYFAC) * Polygone IF (ITYFAC.EQ.6) NNODE=NNOMAX IF (NNODE.GT.0) THEN NFACI=NBFAC(ITYFAC) SEGINI IFACI IPOFAC(ITYFAC)=IFACI ENDIF ENDDO c NBFAC sert maintenant de compteur DO ITYFAC=1,NTYFAC NBFAC(ITYFAC)=0 ENDDO ICOUPE=0 DO 50 IOB=1,MAX(1,LISOUS(/1)) IF (LISOUS(/1).NE.0) THEN IPT1=LISOUS(IOB) ICOUPE=0 IF (IOB.EQ.LISOUS(/1).AND.MCOUP.NE.0) ICOUPE=1 ENDIF NBELEM=IPT1.NUM(/2) ILTEL=LTEL(1,IPT1.ITYPEL) IF (ILTEL.EQ.0) GOTO 52 ILTAD=LTEL(2,IPT1.ITYPEL) DO 60 IF=1,ILTEL ITYFAC=LDEL(1,ILTAD+IF-1) IAD=LDEL(2,ILTAD+IF-1) NNODE=KDFAC(1,ITYFAC) NNODF=NNODE * Polygone IF (ITYFAC.EQ.6) THEN NBNN=IPT1.NUM(/1) * 23 1 *Erreur dans le module de trace IF (NBNN.GT.NNOMAX) THEN CALL ERREUR(23) RETURN ENDIF NNODE=NBNN NNODF=NNOMAX ENDIF IF (NNODE.GT.0) THEN IFACI=IPOFAC(ITYFAC) DO 70 IEL=1,NBELEM NBFAC(ITYFAC)=NBFAC(ITYFAC)+1 j=NBFAC(ITYFAC) * Polygone IF (ITYFAC.EQ.6) NSOMP(j)=NNODE IFACI(NNODF+2,j)=1 DO i=1,NNODE IFACI(i,j)=ICPR(IPT1.NUM(LFAC(IAD+i-1),IEL)) IF (IVU(IFACI(i,j)).NE.1) IFACI(NNODF+2,j)=0 ENDDO IFACI(NNODF+1,j)=IPT1.ICOLOR(IEL) * TRI3 cas des coupes IF (ITYFAC.EQ.1) THEN IF (ICOUPE.EQ.1) THEN IF (MCOUP(IEL)/8.EQ.1) IFACI(NNODF+2,j)=2 IF (MCOUP(IEL)/16.EQ.1) IFACI(NNODF+2,j)=3 ENDIF ENDIF 70 CONTINUE ENDIF 60 CONTINUE 52 CONTINUE IF (LISOUS(/1).NE.0) SEGDES IPT1 50 CONTINUE SEGDES MELEME C IF FAUT MAINTENANT RETASSER ET CLASSER LES TABLEAUX DES FACES C PROBLEME ELLES NE SONT PAS FORCEMENT DECRITES DE LA MEME FACON NFAC=0 DO ITYFAC=1,NTYFAC NFAC=NFAC+NBFAC(ITYFAC) ENDDO SEGINI TFAC,KFAK,IFACOL IFAC=0 DO ITYFAC=1,NTYFAC NNODE=KDFAC(1,ITYFAC) IF (ITYFAC.EQ.6) THEN NNODF=NNOMAX ELSE NNODF=NNODE ENDIF IF (NNODE.GT.0.OR.ITYFAC.EQ.6) THEN IFACI=IPOFAC(ITYFAC) DO I=1,NBFAC(ITYFAC) IFAC=IFAC+1 * Polygone IF (ITYFAC.EQ.6) NNODE=NSOMP(I) XXXX = 0. DO J=1,NNODE * ligne suivante : erreur de faced pour poly ? * XXXX = XXXX + XPROJ(3,ICPR(IFACI(J,I))) XXXX = XXXX + XPROJ(3,IFACI(J,I)) ENDDO XXXX=XXXX/NNODE TFAC(IFAC)=XXXX IF (IFACI(NNODF+2,I).EQ.1) TFAC(IFAC)=TFAC(IFAC)-TDIST IFACOL(IFAC)=IFACI(NNODF+1,I) KFAK(IFAC)=I+((ITYFAC-1)*NFAC) * TRI3/coupe IF (ITYFAC.EQ.1) THEN IF (IFACI(NNODF+2,I).EQ.2) TFAC(I)=TFAC(IFAC)-2*TDIST IF (IFACI(NNODF+2,I).EQ.3) KFAK(IFAC)=0 ENDIF ENDDO ENDIF ENDDO C IL N'Y A PLUS QU'A TRIER ET RETASSER KFAK SUIVANT TFAC SEGINI NAUX IF (IDIM.EQ.2) GOTO 209 NAUX(1)=1 NAUX(2)=NFAC IZ=2 208 IZ=IZ-1 IF (IZ.LE.0) GOTO 209 IPB=NAUX(IZ*2-1) IPH=NAUX(IZ*2) IF(IPB.GE.IPH) GOTO 208 JPB=IPB-1 JPH=IPH+1 C CALCUL DU PIVOT PV=0. * DO 207 J=IPB,IPH * PV=PV+TFAC(J) *207 CONTINUE * PV=PV/(IPH-IPB+1) PV=(TFAC(IPB)+TFAC(IPH))/2. 242 JPB=JPB+1 IF (JPH.EQ.JPB) GOTO 245 IF (TFAC(JPB).LT.PV) GOTO 243 GOTO 242 243 JPH=JPH-1 IF (JPH.EQ.JPB) GOTO 245 IF (TFAC(JPH).GT.PV) GOTO 244 GOTO 243 244 IAUX=KFAK(JPB) IAUXX=IFACOL(JPB) KFAK(JPB)=KFAK(JPH) IFACOL(JPB)=IFACOL(JPH) KFAK(JPH)=IAUX IFACOL(JPH)=IAUXX TAUX=TFAC(JPB) TFAC(JPB)=TFAC(JPH) TFAC(JPH)=TAUX GOTO 242 245 IF (JPB.GE.IPB) GOTO 247 JPB=JPB+1 JPH=JPH+2 GOTO 248 247 IF (JPH.LE.IPH) GOTO 249 JPB=JPB-2 JPH=JPH-1 GOTO 248 249 IF (TFAC(JPB).LE.PV) GOTO 250 IF (JPH.EQ.IPH) GOTO 251 252 JPH=JPH+1 GOTO 248 250 IF (JPB.EQ.IPB) GOTO 252 251 JPB=JPB-1 248 IF (JPB.EQ.IPB) GOTO 253 NAUX(2*IZ)=JPB IZ=IZ+1 253 IF (JPH.EQ.IPH) GOTO 208 NAUX(2*IZ)=IPH NAUX(2*IZ-1)=JPH IZ=IZ+1 GOTO 208 209 CONTINUE C LES FACES SONT CLASSEES DANS KFAK LES FACES ON ETE ELIMINEES PAR C ENVELO . IL NE RESTE PLUS QU'A TRACER DO 300 I=1,NFAC IFF=KFAK(I) IF (IFF.EQ.0) GOTO 300 IT=(IFF-1)/NFAC IF=IFF-IT*NFAC IT=IT+1 IOK=0 * ITYFAC=IT NNODE=KDFAC(1,ITYFAC) IF (ITYFAC.EQ.6) NNODE=NSOMP(IF) IF (NNODE.GT.0) THEN IFACI=IPOFAC(ITYFAC) DO IP=1,NNODE XTR(IP)=XPROJ(1,IFACI(IP,IF)) YTR(IP)=XPROJ(2,IFACI(IP,IF)) ZTR(IP)=XPROJ(3,IFACI(IP,IF)) NPTR(IP)=IFACI(IP,IF) IF (IVU(IFACI(IP,IF)).EQ.1) IOK=1 *!! write(ioimp,*) 'ip=',ip,' x y z n iok',xtr(ip),ytr(ip), *!! $ ztr(ip), iok ENDDO ENDIF IF (IOK.EQ.0) GOTO 300 NP=NNODE * SG Pour les TRI7,QUA9 seuls les points du contour sont pris en compte * pour le calcul de la normale et la suppression des points du contour. IF (ITYFAC.EQ.7.OR.ITYFAC.EQ.8) NP=NP-1 *!! write(ioimp,*) 'np=',np C DETERMINATION DE LA COULEUR D'APRES L'ORIENTATION DE LA NORMALE XN=0. YN=0. ZN=0. XDE=XTR(NP) YDE=YTR(NP) ZDE=ZTR(NP) XW1=XTR(1)-XDE YW1=YTR(1)-YDE ZW1=ZTR(1)-ZDE DO 900 I2=2,NP-1 XW2=XTR(I2)-XDE YW2=YTR(I2)-YDE ZW2=ZTR(I2)-ZDE XN=XN+(YW1*ZW2-ZW1*YW2) YN=YN+(ZW1*XW2-XW1*ZW2) ZN=ZN+(XW1*YW2-YW1*XW2) XW1=XW2 YW1=YW2 ZW1=ZW2 900 CONTINUE DN=SQRT(XN**2+YN**2+ZN**2) IF (DN.EQ.0.) DN=1. ZN=ASIN(ABS(ZN/DN)) C+PP DEGRADE OU NON IF (lndegr) ZN=REAL(XPI/2.D0) C+PP IFACO=IFACOL(I) if (imod.eq.0) ifaco=8 * SG Pour les TRI7 et QUA9, on choisit de decouper en TRI3 : * idealement, on aurait mis le calcul de la normale dedans IF (ITYFAC.EQ.7.OR.ITYFAC.EQ.8) THEN XTR2(1)=XTR(NP+1) YTR2(1)=YTR(NP+1) ZTR2(1)=ZTR(NP+1) DO IPO=1,NP IPOP=IPO+1 IF (IPOP.GT.NP) IPOP=1 XTR2(2)=XTR(IPO) YTR2(2)=YTR(IPO) ZTR2(2)=ZTR(IPO) XTR2(3)=XTR(IPOP) YTR2(3)=YTR(IPOP) ZTR2(3)=ZTR(IPOP) CALL TRFACE(3,XTR2,YTR2,ZTR2,ZN,IFACO,IEFF) if (imod.eq.0) ieff=0 ENDDO ELSE CALL TRFACE(NP,XTR,YTR,ZTR,ZN,IFACO,IEFF) ENDIF if (imod.eq.0) ieff=0 C SI IEFF <> 0 IL FAUT METTRE LES TRAITS EN EFFACEMENT N1=NPTR(NP) DO 450 IIP=1,NP N2=NPTR(IIP) NI=N1 NJ=N2 * 459 CONTINUE 457 DO 454 K=1,NBCONR IF (KON(1,K,NI).NE.NJ) GOTO 454 C 8 = EFFACEMENT KON(2,K,NI)=8 IF (IEFF.EQ.0) KON(2,K,NI)=IFACOL(I) GOTO 456 454 CONTINUE IF (KON(1,NBCON,NI).EQ.0) GOTO 456 NI=KON(1,NBCON,NI) GOTO 457 456 CONTINUE NI=N2 NJ=N1 * 469 CONTINUE 467 DO 464 K=1,NBCONR IF (KON(1,K,NI).NE.NJ) GOTO 464 C 8 = EFFACEMENT KON(2,K,NI)=8 IF (IEFF.EQ.0) KON(2,K,NI)=IFACOL(I) GOTO 466 464 CONTINUE IF (KON(1,NBCON,NI).EQ.0) GOTO 466 NI=KON(1,NBCON,NI) GOTO 467 466 CONTINUE N1=N2 450 CONTINUE 300 CONTINUE C'EST FINI SEGACT MELEME IF (LISOUS(/1).NE.0) THEN NBSOUS=LISOUS(/1) * IF (MCOUP.NE.0) NBSOUS=NBSOUS-1 DO 490 IO=1,NBSOUS IPT2=LISOUS(IO) segact ipt2 if (ipt2.itypel.gt.3.AND.ipt2.itypel.NE.32) SEGSUP IPT2 490 CONTINUE ENDIF if ((itypel.eq.0) .AND. (MELEME .NE. MELSAU)) SEGSUP MELEME MELEME=MELSAU SEGDES KON DO ITYFAC=1,NTYFAC IFACI=IPOFAC(ITYFAC) IF (IFACI.NE.0) THEN SEGSUP IFACI ENDIF ENDDO SEGSUP IPOFAC SEGSUP NBFAC SEGSUP NAUX,TFAC,KFAK,IFACOL,NSOMP RETURN END