1 : C GAUSS SOURCE CHAT 05/01/13 00:16:08 5004
2 : c Sous-Programme du module MISTRAL0
3 : C --------------------------------------------------------------------------
4 : SUBROUTINE GAUSS (A,B,N,NSY,X,ISOL)
5 : C --------------------------------------------------------------------------
6 : C Resolution de NSY systemes de N equations lineaires a N inconnues :
7 : C A.X = B
8 : C de meme matrice carree A et de matrices colonnes B differentes.
9 : C ..........................................................................
10 : C Pour tout ISY compris entre 1 et NSY inclus,
11 : C tout J compris entre 1 et N inclus :
12 : C N
13 : C SOMME A(I,J) * X(I,ISY) = B(J,ISY)
14 : C I=1
15 : C ..........................................................................
16 : C Methode de Gauss :
17 : C (version 2.2 FORTRAN du 17/03/1998)
18 : C transformation en matrice triangulaire inferieure (virtuellement :
19 : C permutation interne des indices uniquement).
20 : C
21 : C Pour tout M compris entre 1 et N inclus,
22 : C il existe (au debut de l'etape M) :
23 : C - une permutation IM de {1,N},
24 : C - une matrice carree AM (stockee dans A),
25 : C - des matrices colonnes BM (stockees dans B),
26 : C telles que, pour tout ISY compris entre 1 et NSY inclus :
27 : C - pour tout J compris entre 1 et J0=N-M+1 inclus :
28 : C J0
29 : C SOMME AM(IM(K),J) * X(IM(K),ISY) = BM(J,ISY)
30 : C K=1
31 : C - pour tout J compris entre J0+1 et N inclus :
32 : C J-1
33 : C X(IM(J),ISY) = BM(J,ISY) - SOMME AM(IM(K),J) * X(IM(K),ISY)
34 : C K=1
35 : C
36 : C Les matrices A et B sont donc modifiees par GAUSS.
37 : C --------------------------------------------------------------------------
38 : IMPLICIT INTEGER(I-N)
39 : IMPLICIT REAL*8 (A-H, O-Z)
40 : c std :
41 : PARAMETER ( NMAX = 7 )
42 : DIMENSION A(1:N,1:*),B(1:N,1:*),IM(1:NMAX),X(1:N,1:*),ISOL(1:*)
43 : c std.
44 : C
45 : C Initialisation de la permutation
46 : C
47 : DO K = 1,N
48 : IM(K) = K
49 : END DO
50 : C ------------------------------------------------------------------
51 : C Debut de la boucle sur M de reduction du nombre d'equations
52 : C -----------------------------------------------------------
53 : DO M = 1,N-1
54 : J0M1 = N-M
55 : J0 = J0M1+1
56 : C ----------------------------------
57 : C Recherche du pivot sur la ligne J0
58 : C
59 : K0 = J0
60 : AI0J0 = A(IM(K0),J0)
61 : DO K = J0M1,1,-1
62 : AIJ0 = A(IM(K),J0)
63 : IF (ABS(AIJ0).GT.ABS(AI0J0)) THEN
64 : K0 = K
65 : AI0J0 = AIJ0
66 : END IF
67 : END DO
68 : I0 = IM(K0)
69 : C
70 : C Systemes singuliers
71 : C
72 : IF (AI0J0.EQ.0.) THEN
73 : GOTO 1
74 : END IF
75 : C ------------------------------------
76 : C Permutation eventuelle des inconnues
77 : C
78 : IF (K0.LT.J0) THEN
79 : DO K = K0,J0M1
80 : IM(K) = IM(K+1)
81 : END DO
82 : IM(J0) = I0
83 : END IF
84 : C -----------------------------
85 : C Coefficients de X(IM(J0),ISY)
86 : C
87 : DO K = 1,J0M1
88 : I = IM(K)
89 : A(I,J0) = A(I,J0)/AI0J0
90 : END DO
91 : DO ISY = 1,NSY
92 : B(J0,ISY) = B(J0,ISY)/AI0J0
93 : END DO
94 : C -----------------------------------------
95 : C Etablissement des N-M equations restantes
96 : C
97 : DO J = 1,J0M1
98 : AI0J = A(I0,J)
99 : DO K = 1,J0M1
100 : I = IM(K)
101 : A(I,J) = A(I,J)-AI0J*A(I,J0)
102 : END DO
103 : DO ISY = 1,NSY
104 : B(J,ISY) = B(J,ISY)-AI0J*B(J0,ISY)
105 : END DO
106 : END DO
107 : C ------
108 : END DO
109 : C ---------------------------------------------------------
110 : C Fin de la boucle sur M de reduction du nombre d'equations
111 : C ------------------------------------------------------------------
112 : IM1 = IM(1)
113 : AIM11 = A(IM1,1)
114 : C
115 : C Systemes singuliers
116 : C
117 : IF (AIM11.EQ.0.) THEN
118 : J0 = 1
119 : GOTO 1
120 : END IF
121 : C ------------------------------------------------------------------
122 : C Calculs des solutions X(IM(J),ISY)
123 : C
124 : DO ISY = 1,NSY
125 : X(IM1,ISY) = B(1,ISY)/AIM11
126 : DO J = 2,N
127 : XINJ = B(J,ISY)
128 : DO K = 1,J-1
129 : I = IM(K)
130 : XINJ = XINJ-A(I,J)*X(I,ISY)
131 : END DO
132 : X(IM(J),ISY) = XINJ
133 : END DO
134 : ISOL(ISY) = 1
135 : END DO
136 : RETURN
137 : C ------------------------------------------------------------------
138 : C Systemes singuliers
139 : C
140 : 1 DO ISY = 1,NSY
141 : IF (B(J0,ISY).EQ.0.) THEN
142 : C Infinite de solutions
143 : ISOL(ISY) = 2
144 : ELSE
145 : C Aucune solution
146 : ISOL(ISY) = 0
147 : END IF
148 : END DO
149 : RETURN
150 : END
151 :
152 :
© Cast3M 2003 - All rights reserved.
Disclaimer