C KCOT SOURCE GOUNAND 25/11/12 21:15:16 12399 SUBROUTINE KCOT IMPLICIT INTEGER(I-N) IMPLICIT REAL*8 (A-H,O-Z) C************************************************************************* C C OBJET : Cree un CHAMPOINT CENTRE contenant les informations C décrites ci dessous, sur les éléments du domaine. C C SYNTAXE : CHPC = KCOT OBJDOM ; C C OBJDOM : TABLE de SOUSTYPE DOMAINE C C AVIS IMPORTANT ; C C Ces informations sont destinées aux opérateurs de discrétisation C et sont rangées dans un ordre particulier. C C Les dimensions du tableau sont : (DIME NBEL) DIME=11 en 2D C DIME=13 en 3D C et il contient : C C pour un SEG2 la longueur de l'element (XML) et la matrice P(2,2) en 2D C P(3,3) en 3D C pour un TRI3 XML,XMH,AJ1/XML,AJ2/XML, et la matrice P(2,2) en 2D C P(3,3) en 3D C pour un QUA4 IDEM C pour un CUB8 XML XMH XME et la matrice P C pour un PRI6 IDEM C C*********************************************************************** C MATRICE P C LA MATRICE DE ROTATION DU REPERE GLOBALE VERS LE REPERE LOCAL C DEFINI PAR DEUX OU TROIS POINTS PRIS DANS XYZ SUIVANT QU'ON EST C EN 2D OU EN 3D C ON PREND P1 P2 ET PNP C C U TEL QUE T SOIT DIRIGE SUIVANT P1P2 V TOURNE VERS C . .V P1PNP ET U = T VECTORIEL V C . . C . . __ __ C (P1). . . . .T (P2) | tx ty tz | C | | C ON A ALORS WL= P WG P = | vx vy vz | C | | C | ux uy uz | C |__ __| C************************************************************************* -INC PPARAM -INC CCOPTIO -INC CCGEOME -INC SMTABLE -INC SMELEME POINTEUR MELEMC.MELEME -INC SMCOORD -INC SMCHPOI -INC SIZFFB PARAMETER (NLTS=8) CHARACTER*8 LISTS(NLTS) C*** REAL*8 AAJ(3,3,9),U,XC(3),T(13),XQ(81) PARAMETER (NTB=1) CHARACTER*8 LTAB(NTB),TYPE,NOM0,TYPC DIMENSION KTAB(NTB) DATA LTAB/'DOMAINE '/ DATA LISTS/'QUA4 ','CUB8 ','TRI3 ','PRI6 ','SEG2 ', &'SEG3 ','QUA9 ','TRI7 '/ C*** C write(6,*)' DBUT KCOT ' NTO=NTB CALL LITABS(LTAB,KTAB,NTB,NTO,IRET) IF (IERR.NE.0) RETURN MTABD=KTAB(1) CALL INITD(AAJ,81,0.D0) TYPE=' ' CALL ACMO(MTABD,'CENTRE',TYPE,MELEMC) IF(TYPE.NE.'MAILLAGE')GO TO 90 IF(IDIM.EQ.2)NC=11 IF(IDIM.EQ.3)NC=13 TYPE='CENTRE' CALL CRCHPT(TYPE,MELEMC,NC,1,MCHPOI) CALL LICHT(MCHPOI,MPOVAL,TYPC,IGEOM) TYPE=' ' CALL ACMO(MTABD,'MAILLAGE',TYPE,MELEME) IF(TYPE.NE.'MAILLAGE')GO TO 90 SEGACT MELEME IAXI=0 IF(IFOMOD.EQ.0)IAXI=2 NPGI=1 NBSOUS=LISOUS(/1) IF(NBSOUS.EQ.0)NBSOUS=1 K0=0 DO 1 L=1,NBSOUS IPT1=MELEME IF(NBSOUS.NE.1)IPT1=LISOUS(L) SEGACT IPT1 NP=IPT1.NUM(/1) NEL=IPT1.NUM(/2) NOM0=NOMS(IPT1.ITYPEL)//' ' CALL OPTLI(IP,LISTS,NOM0,NLTS) IF(IP.EQ.0)THEN WRITE(6,*)' Il y a des types d''éléments pour lesquels on ne sait' &,'pas faire' GO TO 1 ENDIF CALL KALPBG(NOM0,'FONFORM ',IZFFM) SEGACT IZFFM*MOD IZHR=KZHR(1) SEGACT IZHR*MOD NES=GR(/1) NPG=GR(/3) DO 10 K=1,NEL KK=K0+K C C REMPLISSAGE DE XYZ C ------------------ C DO 12 I=1,NP J=IPT1.NUM(I,K) DO 121 N=1,IDIM XYZ(N,I)=XCOOR((J-1)*(IDIM+1)+N) 121 CONTINUE 12 CONTINUE IF(IP.EQ.6)THEN IN=0 DO 13 I=1,NP,2 J=IPT1.NUM(I,K) DO 131 N=1,IDIM IN=IN+1 XQ(IN)=XCOOR((J-1)*(IDIM+1)+N) 131 CONTINUE 13 CONTINUE ELSEIF(IP.GT.6)THEN IN=0 DO 14 I=1,NP-1,2 J=IPT1.NUM(I,K) DO 141 N=1,IDIM IN=IN+1 XQ(IN)=XCOOR((J-1)*(IDIM+1)+N) 141 CONTINUE 14 CONTINUE ENDIF IF(IP.LE.5)CALL CALJAJ(GR,XYZ,NES,IDIM,NP,NPGI,AAJ) C CALL CALJBC(FN,GR,PG,XYZ,HR,PGSQ,RPG,NES,IDIM,NP,NPG,IAXI,AIRE) IF(AIRE.LE.0)AIRE=ABS(AIRE) C CALCUL DES LONGUEURS DES COTES DE L'ELEMENT C UNIQUEMENT POUR LES SEG2 QUA4 ET CUB8 XML=0.D0 DO 231 N=1,IDIM XML=XML+(XYZ(N,1)-XYZ(N,2))*(XYZ(N,1)-XYZ(N,2)) 231 CONTINUE XML=SQRT(XML) C C IF(XML.EQ.0.D0) &WRITE(6,*)' ELEMENT DEGENERE ( AYANT UN COTE DE LONGUEUR NULLE) ' C C IF(NOM0.EQ.'SEG2 ')THEN T(1)=XML CALL CALJQB(XYZ,T(5),IDIM,IDIM) GOTO 244 ELSEIF(NOM0.EQ.'SEG3 ')THEN T(1)=XML CALL CALJQB(XQ,T(5),IDIM,IDIM) GOTO 244 ENDIF C C XMH=0.D0 DO 232 N=1,IDIM XMH=XMH+(XYZ(N,2)-XYZ(N,3))*(XYZ(N,2)-XYZ(N,3)) 232 CONTINUE XMH=SQRT(XMH) C WRITE(6,*)'........XML XMH NP NES ',XML,XMH,NP,NES C C IF(XMH.EQ.0.D0) &WRITE(6,*)' ELEMENT DEGENERE ( AYANT UN COTE DE LONGUEUR NULLE) ' C C IF(NOM0.EQ.'QUA4 '.OR.NOM0.EQ.'TRI3 ')THEN T(1)=XML T(2)=XMH T(3)=AAJ(1,1,1)/XML T(4)=AAJ(2,1,1)/XML CALL CALJQB(XYZ,T(5),IDIM,IDIM) GOTO 244 ELSEIF(NOM0.EQ.'QUA9 '.OR.NOM0.EQ.'TRI7 ')THEN T(1)=XML T(2)=XMH T(3)=AAJ(1,1,1)/XML T(4)=AAJ(2,1,1)/XML CALL CALJQB(XQ,T(5),IDIM,IDIM) GOTO 244 ENDIF IF(NOM0.EQ.'CUB8 '.OR.NOM0.EQ.'PRI6 ')THEN XME=0.D0 IF(NOM0.EQ.'CUB8 ')THEN DO 242 N=1,IDIM XME=XME+(XYZ(N,1)-XYZ(N,5))*(XYZ(N,1)-XYZ(N,5)) 242 CONTINUE ENDIF IF(NOM0.EQ.'PRI6 ')THEN DO 243 N=1,IDIM XME=XME+(XYZ(N,1)-XYZ(N,4))*(XYZ(N,1)-XYZ(N,4)) 243 CONTINUE ENDIF XME=SQRT(XME) C WRITE(6,*)'........XML XMH XME NP NES ',XML,XMH,XME,NP,NES C IF(XME.EQ.0.D0) &WRITE(6,*)' ELEMENT DEGENERE ( AYANT UN COTE DE LONGUEUR NULLE) ' C C CALL CALJQB(XYZ,T(5),IDIM,4) T(1)=XML T(2)=XMH T(3)=XME GO TO 244 ENDIF 244 CONTINUE CALL RSETIV(VPOCHA,NC,KK,T,NC) C write(6,*)' KK=',KK C write(6,1002)T 10 CONTINUE K0=K0+NEL SEGDES IPT1 1 CONTINUE SEGDES MELEME,MPOVAL,MCHPOI C CALL ECROBJ('CHPOINT ',MCHPOI) C write(6,*)' FIN KCOT ' RETURN 90 CONTINUE WRITE(6,*)' Interruption anormale de KCOT' C Option %m1:8 incompatible avec les données CALL ERREUR(803) RETURN 1001 FORMAT(20(1X,I5)) 1002 FORMAT(10(1X,1PE11.4)) END