C PSIPHI SOURCE CB215821 20/11/25 13:37:58 10792 SUBROUTINE PSIPHI *********************************************************************** c OPERATEUR : PSIP c c APPEL : c CHP1 (CHP2) (CHP3) = PSIP MAIL1 MAIL2 (CRIT1) (|'DEUX' | | P1 (P2) | ); c |'TROI' | | MAIL3 | c c FONCTION : calcule les fonctions distances signées (level set) c relatives aux maillages MAIL2 (surface de fissure) c et MAIL3 en 3D (front de fissure) c ou P1 en 2D (pointe de fissure) c aux noeuds de MAIL1 (pas trop éloignés). c c CREATION : chat (16/07/2007) c MODIFS : bp (2009 -> 2012) : diverses corrections/ameliorations c repertoriees par les fiches d'anomalies/developpement c bp 14/03/2012 : on split psiphi en psiphi psip2d psip3d et c zonag2 (evolué seulement le 18/12/2013) c c TO DO : cas 3d avec plusieurs fronts (ou 1 front discontinu) c *********************************************************************** IMPLICIT REAL*8(A-H,O-Z) IMPLICIT INTEGER(I-N) -INC PPARAM -INC CCOPTIO -INC SMCOORD -INC SMELEME -INC SMCHPOI -INC CCREEL PARAMETER(NMOCLE=2) CHARACTER*4 MOCLE(NMOCLE) DATA MOCLE/'DEUX','TROI'/ if(iimpi.ge.1) write(ioimp,*)'==== ENTREE DANS PSIPHI =====' segact mcoord *********************************************************************** * INITIALISATIONs et LECTUREs *********************************************************************** * * calcule t'on phi? phi ET psi? phi, psi et tau? ideux=0 call lirmot (mocle,nmocle,ideux,0) if (ideux.eq.2.and.idim.lt.3) then write(*,*) ' ON NE PEUT CALCULER TROIS LEVEL SET QU EN 3D !' write(*,*) ' ON CONTINUE AVEC LE CALCUL DE 2 LEVEL SET ...' ideux=1 endif * lecture des points pour lesquels on veut calculer phi et psi call LIROBJ('MAILLAGE',ipt1,1,iretou) melmai=ipt1 if(ierr.ne.0) return call change(ipt1,1) c segact ipt1 c rem: inutile car change laisse ipt1 actif * lecture du maillage de la fissure call LIROBJ('MAILLAGE',meleme,1,iretou) melfis=meleme if(ierr.ne.0) return * lecture de la pointe de fissure (objet de type point) ou du front ip1=0 ip2=0 melfro=0 1 continue if (ideux.ge.1) then if (idim.eq.2) then call LIROBJ('POINT',ipt,0,iretou) if (iretou.ne.0) then if (ip1.eq.0) then ip1 = ipt go to 1 else ip2 = ipt endif endif else call LIROBJ('MAILLAGE', melfro ,1,iretou) if(iretou.eq.0) write(ioimp,*) 'Il manque le maillage du front' if(IERR.NE.0) return endif endif * lecture facultative d'une longueur max (=plus grande taille des elements * concernés par la fissure) xcrit=0.D0 call lirree(xcrit,0,ircrit) *********************************************************************** * CRÉATION DES MCHPOI DE SORTIE (1 POUR PSI , 2 POUR PHI, 3 POUR TAU) *********************************************************************** mpova1=0 mpova2=0 mpova3=0 *-----PHI----- c if(idebug.eq.1) write(6,*) '----Creation de PHI----' nat=1 nsoupo=1 nc=1 n=ipt1.num(/2) segini,mchpo2 segini,msoup2 mchpo2.jattri(1)=1 mchpo2.ipchp(1)=msoup2 mchpo2.ifopoi=ifour segdes mchpo2 msoup2.igeoc=ipt1 msoup2.nocomp(1)='PHI' msoup2.noharm(1)=nifour segini,mpova2 msoup2.ipoval=mpova2 c segdes,msoup2 *-----PSI----- if (ideux.ge.1) then c if(idebug.eq.1) write(6,*) '----Creation de PSI----' segini,mchpo1 segini,msoup1 mchpo1.jattri(1)=1 mchpo1.ipchp(1)=msoup1 mchpo1.ifopoi=ifour segdes,mchpo1 msoup1.igeoc=ipt1 msoup1.nocomp(1)='PSI' msoup1.noharm(1)=nifour segini,mpova1 msoup1.ipoval=mpova1 c segdes,msoup1 else msoup1 = 0 endif *-----TAU----- if (ideux.ge.2) then c if(idebug.eq.1) write(6,*) '----Creation de TAU----' segini,mchpo3 segini,msoup3 mchpo3.jattri(1)=1 mchpo3.ipchp(1)=msoup3 mchpo3.ifopoi=ifour segdes,mchpo3 msoup3.igeoc=ipt1 msoup3.nocomp(1)='TAU' msoup3.noharm(1)=nifour segini,mpova3 msoup3.ipoval=mpova3 c segdes,msoup3 else msoup3 = 0 endif *********************************************************************** * PSIPHI 3D/2D *********************************************************************** if (idim.eq.3) then c write(ioimp,*)'appel PSIP3D (',ideux,ipt1,melfis,melfro,xcrit call PSIP3D(ideux,ipt1,melfis,melfro,xcrit, & msoup1,msoup2,msoup3) else c write(ioimp,*)'appel PSIP3D (',ideux,ipt1,melfis,ip1,ip2,xcrit call PSIP2D(ideux,ipt1,melfis,ip1,ip2,xcrit, & msoup1,msoup2) endif *********************************************************************** * ECRITURE ET FIN DU PROGRAMME *********************************************************************** c-----ecriture de (PSI) PHI ((TAU)) -------------- if(ideux.ge.2) then call ACTOBJ('CHPOINT ',mchpo3,1) call ECROBJ('CHPOINT ',mchpo3) endif call ACTOBJ('CHPOINT ',mchpo2,1) call ECROBJ('CHPOINT ',mchpo2) if(ideux.ge.1) then call ACTOBJ('CHPOINT',mchpo1,1) call ECROBJ('CHPOINT',mchpo1) endif end