1 : C ETOILE SOURCE MAGN 09/07/30 21:16:51 6458 2 : SUBROUTINE ETOILE(M1,M2,M3,IRET) 3 : IMPLICIT INTEGER(I-N) 4 : IMPLICIT REAL*8 (A-H,O-Z) 5 : C ************************************************ 6 : C * Routine effectuant la multiplication de 2 * 7 : C * MATRIK sous stockage morse * 8 : C * M1 : entree * 9 : C * M2 : entree * 10 : C * M3 : Sortie sous stockage morse * 11 : C * M3 = M1 * M2 * 12 : C ************************************************ 13 : 14 : -INC SMELEME 15 : POINTEUR MELEMD.MELEME 16 : 17 : -INC SMMATRIK 18 : POINTEUR M1.MATRIK,M2.MATRIK,M3.MATRIK 19 : POINTEUR PMS3.PMORS, IZA3.IZA 20 : POINTEUR IMAT1.IMATRI,IMAT2.IMATRI,IMAT3.IMATRI 21 : 22 : INTEGER CI 23 : REAL*8 TEMP 24 : 25 : IRET=0 26 : 27 : SEGACT M1 28 : 29 : C Si les matrices ne sont pas morse, on les crees 30 : C en morse 31 : IF (M1.IRIGEL(7,1).NE.6) THEN 31 : CALL ELMORS(M1,1) 32 : END IF 33 : 34 : SEGACT M2 35 : 36 : IF (M2.IRIGEL(7,1).NE.6) THEN 36 : CALL ELMORS(M2,1) 37 : END IF 38 : 39 : SEGACT M1,M2 40 : 41 : C *********************** 42 : MELEMD=M1.IRIGEL(2,1) 43 : MELEME=M2.IRIGEL(1,1) 44 : 45 : SEGACT MELEMD,MELEME 46 : NBSOUP=MELEME.LISOUS(/1) 47 : NBSOUD=MELEMD.LISOUS(/1) 48 : SEGDES MELEMD,MELEME 49 : NMATR1=M1.IRIGEL(/2) 50 : 51 : IF (NMATR1.GT.1) THEN 52 : WRITE(6,*) 'ETOILE: Impossible de mutiplier la' 53 : WRITE(6,*) 'matrice morse 1: NMATRI =',NMATR1,'>1.' 54 : IRET=1 55 : RETURN 56 : END IF 57 : 58 : NMATR2=M2.IRIGEL(/2) 59 : 60 : IF (NMATR2.GT.1) THEN 61 : WRITE(6,*) 'ETOILE: Impossible de mutiplier la' 62 : WRITE(6,*) 'matrice morse 2: NMATRI =',NMATR2,'>1.' 63 : IRET=1 64 : RETURN 65 : END IF 66 : 67 : C *********************** 68 : C On recupere les segments necessaire pour le produit 69 : 70 : PMS1=M1.IRIGEL(5,1) 71 : PMS2=M2.IRIGEL(5,1) 72 : IZA1=M1.IRIGEL(6,1) 73 : IZA2=M2.IRIGEL(6,1) 74 : IMAT1=M1.IRIGEL(4,1) 75 : IMAT2=M2.IRIGEL(4,1) 76 : 77 : SEGACT IMAT1,IMAT2 78 : SEGACT PMS1,PMS2 79 : SEGACT IZA1,IZA2 80 : 81 : C **************************************************** 82 : 83 : NBME1=IMAT1.LIZAFM(/2) 84 : NBSOU1=IMAT1.LIZAFM(/1) 85 : 86 : NBME2=IMAT2.LIZAFM(/2) 87 : NBSOU2=IMAT2.LIZAFM(/1) 88 : 89 : IF (NBME1.NE.NBME2) THEN 90 : WRITE(6,*) 'ETOILE : Produit entre ces deux MATRIK' 91 : WRITE(6,*) 'incompatible' 92 : IRET=1 93 : RETURN 94 : ELSE 95 : IFLAG=0 96 : DO I=1,NBME1 97 : IF (IMAT1.LISPRI(I).NE.IMAT2.LISDUA(I)) THEN 98 : IFLAG=1 99 : END IF 100 : END DO 101 : IF (IFLAG.EQ.1) THEN 102 : WRITE(6,*) 'ETOILE : Produit entre ces deux MATRIK' 103 : WRITE(6,*) 'incompatible' 104 : IRET=1 105 : RETURN 106 : END IF 107 : END IF 108 : 109 : NBME=NBME1 110 : NBSOUS=MAX(NBSOUP,NBSOUD) 111 : IF (NBSOUS.EQ.0) NBSOUS=1 112 : 113 : SEGINI IMAT3 114 : 115 : C Le primal de M3 est le primal de M2 et 116 : C le dual de M3 est le dual de M1 117 : DO I=1,NBME 118 : IMAT3.LISPRI(I)=IMAT2.LISPRI(I) 119 : IMAT3.LISDUA(I)=IMAT1.LISDUA(I) 120 : END DO 121 : 122 : C Si Apres le produit on a plus qu'une composante en primal 123 : C ET en dual avec NBME>1 alors on ajuste a NBME=1. 124 : IF (NBME.NE.1) THEN 125 : IF ((IMAT3.LISPRI(1).EQ.IMAT3.LISPRI(2)).AND. 126 : & (IMAT3.LISDUA(1).EQ.IMAT3.LISDUA(2))) THEN 127 : NBME=1 128 : SEGADJ IMAT3 129 : END IF 130 : END IF 131 : 132 : C On initialise M3 133 : NMATRI=1 134 : NRIGE=7 135 : NKID=9 136 : NKMT=7 137 : SEGINI M3 138 : 139 : NTT=PMS1.IA(/1)-1 140 : NTTP=M2.KNTTP 141 : 142 : C On initialise les segments morses de M3 143 : NJA=1 144 : NBVA=1 145 : SEGINI PMS3,IZA3 146 : M=0 147 : 148 : C On calcule le produit: 149 : DO I=1,NTT 150 : LI1=PMS1.IA(I) 151 : NB1=PMS1.IA(I+1) - PMS1.IA(I) 152 : 153 : DO K=1,NTTP 154 : 155 : TEMP=0.0D0 156 : DO L=1,NB1 157 : CI=PMS1.JA(LI1+L-1) 158 : 159 : LI2=PMS2.IA(CI) 160 : NB2=PMS2.IA(CI+1)-PMS2.IA(CI) 161 : 162 : DO J=1,NB2 163 : IF (PMS2.JA(LI2+J-1).EQ.K) THEN 164 : TEMP=TEMP + (IZA1.A(LI1+L-1)* 165 : & IZA2.A(LI2+J-1)) 166 : ELSEIF (PMS2.JA(LI2+J-1).GT.K) THEN 167 : GOTO 10 168 : END IF 169 : END DO 170 : 10 CONTINUE 171 : 172 : 173 : END DO 174 : 175 : C On optimise en considerant que ce qui est < E-12 est = 0. 176 : IF (ABS(TEMP).LT.1.E-15) TEMP=0.0 177 : IF (TEMP.NE.0.0D0) THEN 178 : M=M+1 179 : 20 IF (NBVA.LT.M) THEN 180 : NBVA=NBVA+100 181 : NJA=NBVA 182 : SEGADJ IZA3,PMS3 183 : GOTO 20 184 : END IF 185 : 186 : C remplissage ligne 187 : IF (PMS3.IA(I).EQ.0) THEN 188 : PMS3.IA(I)=M 189 : END IF 190 : C remplissage colonne 191 : PMS3.JA(M)=K 192 : 193 : IZA3.A(M)=TEMP 194 : END IF 195 : C On fait attention aux lignes vide. 196 : IF ((PMS3.IA(I-1).EQ.0).AND.(I.NE.1)) PMS3.IA(I-1) = M+1 197 : END DO 198 : END DO 199 : 200 : C On oublie pas de remplir la derniere ligne. 201 : PMS3.IA(NTT+1)=M+1 202 : 203 : C on fait attention si la derniere ligne est vide. 204 : IF (PMS3.IA(NTT).EQ.0) PMS3.IA(NTT) = M+1 205 : 206 : C On ajuste les segments morse de M3 si necessaire 207 : IF (NBVA.GT.M) THEN 208 : NBVA=M 209 : NJA=NBVA 210 : SEGADJ IZA3,PMS3 211 : END IF 212 : 213 : C On remplit les pointeurs de M3 214 : M3.IRIGEL(5,1)=PMS3 215 : M3.IRIGEL(6,1)=IZA3 216 : 217 : IMAT3.KSPGP=IMAT2.KSPGP 218 : IMAT3.KSPGD=IMAT1.KSPGD 219 : 220 : M3.IRIGEL(1,1)=M2.IRIGEL(1,1) 221 : M3.IRIGEL(2,1)=M1.IRIGEL(2,1) 222 : M3.IRIGEL(4,1)=IMAT3 223 : M3.IRIGEL(7,1)=6 224 : 225 : M3.KSYM=2 226 : M3.KNTTP=M2.KNTTP 227 : M3.KNTTD=M1.KNTTD 228 : M3.KMINCP=M2.KMINCP 229 : M3.KMINCD=M1.KMINCD 230 : 231 : SEGDES PMS3,IZA3 232 : SEGDES PMS1,PMS2 233 : SEGDES IZA1,IZA2 234 : SEGDES M1,M2,M3 235 : SEGDES IMAT1,IMAT2,IMAT3 236 : RETURN 237 : END 238 : 239 : 240 :
© Cast3M 2003 - All rights reserved.
Disclaimer