1 : C TRIAKS SOURCE CHAT 05/01/13 03:46:35 5004 2 : SUBROUTINE TRIAKS(IZL) 3 : C************************************************************************ 4 : C 5 : C TRIANGULATION D'UNE MATRICE SYMETRIQUE LIGNE A LIGNE 6 : C 7 : C POINTEUR : EN ENTREE IZL CONTIENT LA MATRICE TRIANGULEE 8 : C 9 : C SUR LES CONSEILS DE PV, LORSQUE L'ON CALCULE LIJ*LKJ*DJ POUR 10 : C J VARIANT DE JMIN A J-1, ON UTILISE ZKJ=LKJ*DJ QUE L'ON DIVISE 11 : C ENSUITE PAR DJ POUR RATTRAPER LKJ. ON REMPLACE L'APPEL A SDT3 12 : C PAR UN APPEL A SDOT ET ON ECONOMISE 1 OPERATION SUR 3. 13 : C 14 : C ON POURRAIT AMELIORER L'ALGORITHME EN CODANT PAR BLOCS. CETTE 15 : C PARTIE N'EST PAS CRUCIALE CAR ELLE N'EST EFFECTUEE QU'UNE SEULE 16 : C FOIS LORS D'UN CALCUL COMPLET. 17 : C 18 : C************************************************************************ 19 : IMPLICIT INTEGER(I-N) 20 : IMPLICIT REAL*8 (A-H,O-Z) 21 : REAL*8 SDOT,SDT3,AUX,SOM 22 : 23 : -INC CCREEL 24 : *- 25 : C-INC SMMATRAKANC 26 : C************************************************************************* 27 : C 28 : C REPERAGE ET STOKAGE DES MATRICES ELEMENTAIRES puis assemblees 29 : C 30 : 31 : * LGEOC SPG de la pression et/ou des multiplicateurs de Lagrange 32 : * (points CENTRE ) pour chaque operateur de contrainte 33 : * KGEOC SPG pour la totalite des points CENTRE. 34 : * KGEOS SPG pour la totalite des points SOMMET (Diagonale vitesse) 35 : * KLEMC Connectivites de l'ensemble des contraintes 36 : * LIZAFM(NBSOUS) contient les pointeurs IZAFM des sous-zones 37 : 38 : SEGMENT MATRAK 39 : INTEGER LGEOC(NBOP),IDEBS(NBOP),IFINS(NBOP) 40 : INTEGER LIZAFM(NBSOUS) 41 : INTEGER IKAM0 (NBSOUS) 42 : INTEGER IMEM (NBELC) 43 : INTEGER KLEMC,KGEOS,KGEOC,KDIAG,KCAC,KIZCL,KIZGC 44 : ENDSEGMENT 45 : 46 : SEGMENT IZAFM 47 : REAL*8 AM(NNELP,NP,IESP),RPGI(NELAX) 48 : ENDSEGMENT 49 : 50 : POINTEUR IPMJ.IZAFM,IPMK.IZAFM 51 : 52 : C******************************************************************* 53 : -INC SMMATRK1 54 : POINTEUR IBLKK.IDBLK,IBLKI.IDBLK 55 : 56 : SEGMENT/IZAK/(AK(1)*D) 57 : SEGMENT/IZAI/(AI(1)*D) 58 : SEGMENT/IZD/(D(1)*D) 59 : C 60 : C 61 : GGIS=XPETIT*1.E-4 62 : SIG=XGRAND 63 : 64 : SEGACT IZL 65 : IDMAT=KZA1 66 : SEGACT IDMAT 67 : NBLK=IDESCR(/1) 68 : C 69 : C CALCUL DE L(I,J) A LA LIGNE J 70 : C LA MATRICE EST SYMETRIQUE 71 : 72 : N=KZA(/1) 73 : IZD=IDIAG 74 : SEGACT IZD*MOD 74 : IF(N.NE.D(/1))CALL ARRET(0) 75 : C 76 : DO 100 IBLK=1,NBLK 77 : C NUMÉRO DES LIGNES DE DÉBUT ET FIN DE BLOC KJD,KJF 78 : KJD=NLDBLK(IBLK) 79 : KJF=NLDBLK(IBLK+1)-1 80 : C ACTIVATION DE LA PORTION DE MATRICE 81 : IBLKK=IDESCR(IBLK) 82 : SEGACT IBLKK 83 : IZAK=IBLKK.IMAT 84 : SEGACT IZAK*MOD 85 : DO 1 K=KJD,KJF 86 : KLL=K-1 87 : C LAK=LONGUEUR DE LA LIGNE K 88 : LAK=IBLKK.IDEBLK(K-KJD+2)-IBLKK.IDEBLK(K-KJD+1) 89 : C DÉCALAGE DANS LE TABLEAU POUR ACCÉDER À LA LIGNE K 90 : C - NUMÉRO DE LA LIGNE K DANS LE BLOC IBLKK=K-KJD+1 91 : IDECK=IBLKK.IDEBLK(K-KJD+1)-1 92 : JK=K-LAK 93 : L2=1-JK 94 : IF(LAK.EQ.0)GO TO 1 95 : 96 : DO 4 I=JK,KLL 97 : C A QUEL BLOC APPARTIENT I ? 98 : C ATTENTION AU CAS PARTICULIER OÙ I VAUT 1 99 : IF(I.EQ.1) THEN 100 : IZAI=IZD 101 : LAI=N 102 : IDECI=0 103 : ELSE 104 : IBLI=IBLK 105 : DO 9 IBL=1,NBLK 106 : IF(I.GE.NLDBLK(IBL).AND.I.LT.NLDBLK(IBL+1)) THEN 107 : IBLI=IBL 108 : GOTO 99 109 : ENDIF 110 : 9 CONTINUE 111 : 99 CONTINUE 112 : IJD=KJD 113 : IJF=KJF 114 : IZAI=IZAK 115 : IBLKI=IBLKK 116 : C SI I APPARTIENT À UN BLOC DIFFÉRENT DE CELUI DE K 117 : IF(IBLI.NE.IBLK) THEN 118 : IJD=NLDBLK(IBLI) 119 : IJF=NLDBLK(IBLI+1)-1 120 : IBLKI=IDESCR(IBLI) 121 : SEGACT IBLKI 122 : IZAI=IBLKI.IMAT 123 : SEGACT IZAI*MOD 124 : ENDIF 125 : C LAI=LONGUEUR DE LA LIGNE I 126 : LAI=IBLKI.IDEBLK(I-IJD+2)-IBLKI.IDEBLK(I-IJD+1) 127 : C DÉCALAGE DANS LE TABLEAU POUR ACCÉDER À LA LIGNE I 128 : C - NUMÉRO DE LA LIGNE I DANS LE BLOC IBLKI=I-IJD+1 129 : IDECI=IBLKI.IDEBLK(I-IJD+1)-1 130 : ENDIF 131 : JI=I-LAI 132 : L1=1-JI 133 : I1=I-1 134 : L=L2+I 135 : AUX=0.D0 136 : JMIN=MAX0(JI,JK) 137 : IF (JMIN.GT.I1) GO TO 6 138 : JMI1=I1-JMIN+1 139 : AUX=SDOT(JMI1,AI(IDECI+L1+JMIN),1,AK(IDECK+L2+JMIN),1) 140 : CLENT AUX=SDT3(JMI1,AI(IDECI+L1+JMIN),AK(IDECK+L2+JMIN),D(JMIN)) 141 : 6 SOM=AUX 142 : AK(IDECK+L)=(AK(IDECK+L)-SOM) 143 : CLENT AK(IDECK+L)=(AK(IDECK+L)-SOM)/D(I) 144 : IF(IZAI.NE.IZAK.AND.IZAI.NE.IZD) THEN 145 : SEGDES IZAI 146 : SEGDES IBLKI 147 : ENDIF 148 : 4 CONTINUE 149 : C IL FAUT DIVISER CHAQUE TERME DE LA LIGNE QUI VIENT D'ÊTRE CALCULÉE PAR 150 : C LE TERME DIAGONAL CORRESPONDANT. 151 : DO 5 I=JK,KLL 152 : L=L2+I 153 : AK(IDECK+L)=AK(IDECK+L)/D(I) 154 : 5 CONTINUE 155 : C ON TRAITE LE TERME DIAGONAL, ON NE PEUT FAIRE LA RUSE DU TERME GÉNÉRAL 156 : AUX=SDT3(LAK,AK(IDECK+1),AK(IDECK+1),D(K-LAK)) 157 : SOM=AUX 158 : D(K)=D(K)-SOM 159 : 160 : C! DEBUT 161 : ADK=ABS(D(K)) 162 : IF(ADK.LE.GGIS)THEN 163 : WRITE(6,*)'SUB TRIAWS : DIAGONALE NULLE ADK=',ADK,' LIGNE K=', 164 : &NUNA(K) 165 : D(K)=SIG 166 : ENDIF 167 : C! FIN 168 : 169 : 1 CONTINUE 170 : 171 : SEGDES IZAK 172 : SEGDES IBLKK 173 : 100 CONTINUE 174 : 175 : SEGDES IZD 176 : SEGDES IDMAT 177 : SEGDES IZL 178 : 179 : RETURN 180 : 1002 FORMAT(10(1X,1PE11.4)) 181 : END 182 : 183 : 184 : 185 : 186 : 187 :
© Cast3M 2003 - All rights reserved.
Disclaimer