PROBLEME SUBROUTINE SVD NUMERICAL RECIPES

**BAmal** · 19/09/2013, 19h30

Bonjour ,

J'aimerais vous demander si quelqu'un peut m'expliquer veut dire quoi exactement le "NMAX" qui est un paramètre dans la subroutine donnée par Numerical recipes qui permet de calculer la SVD d une matrice, car cette subroutine arrive bien à me donner de bons resultats pour une matrice de (6*6) mais pas pour une matrice de (n*n ) où n >121 ! j'attends votre reponse (j ai pas trouvé la signification de ce paramètre sur numerical recipes) ,merci à vous.
je joins la subroutine à mon message

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
SUBROUTINE SVDCMP(A,M,N,MP,NP,W,V)
      PARAMETER (NMAX=800)
      DIMENSION A(MP,NP),W(NP),V(NP,NP),RV1(NMAX)
      G=0.0
      SCALE=0.0
      ANORM=0.0
      DO 25 I=1,N
        L=I+1
        RV1(I)=SCALE*G
        G=0.0
        S=0.0
        SCALE=0.0
        IF (I.LE.M) THEN
          DO 11 K=I,M
            SCALE=SCALE+ABS(A(K,I))
11        CONTINUE
          IF (SCALE.NE.0.0) THEN
            DO 12 K=I,M
              A(K,I)=A(K,I)/SCALE
              S=S+A(K,I)*A(K,I)
12          CONTINUE
            F=A(I,I)
            G=-SIGN(SQRT(S),F)
            H=F*G-S
            A(I,I)=F-G
            IF (I.NE.N) THEN
              DO 15 J=L,N
                S=0.0
                DO 13 K=I,M
                  S=S+A(K,I)*A(K,J)
13              CONTINUE
                F=S/H
                DO 14 K=I,M
                  A(K,J)=A(K,J)+F*A(K,I)
14              CONTINUE
15            CONTINUE
            ENDIF
            DO 16 K= I,M
              A(K,I)=SCALE*A(K,I)
16          CONTINUE
          ENDIF
        ENDIF
        W(I)=SCALE *G
        G=0.0
        S=0.0
        SCALE=0.0
        IF ((I.LE.M).AND.(I.NE.N)) THEN
          DO 17 K=L,N
            SCALE=SCALE+ABS(A(I,K))
17        CONTINUE
          IF (SCALE.NE.0.0) THEN
            DO 18 K=L,N
              A(I,K)=A(I,K)/SCALE
              S=S+A(I,K)*A(I,K)
18          CONTINUE
            F=A(I,L)
            G=-SIGN(SQRT(S),F)
            H=F*G-S
            A(I,L)=F-G
            DO 19 K=L,N
              RV1(K)=A(I,K)/H
19          CONTINUE
            IF (I.NE.M) THEN
              DO 23 J=L,M
                S=0.0
                DO 21 K=L,N
                  S=S+A(J,K)*A(I,K)
21              CONTINUE
                DO 22 K=L,N
                  A(J,K)=A(J,K)+S*RV1(K)
22              CONTINUE
23            CONTINUE
            ENDIF
            DO 24 K=L,N
              A(I,K)=SCALE*A(I,K)
24          CONTINUE
          ENDIF
        ENDIF
        ANORM=MAX(ANORM,(ABS(W(I))+ABS(RV1(I))))
25    CONTINUE
      DO 32 I=N,1,-1
        IF (I.LT.N) THEN
          IF (G.NE.0.0) THEN
            DO 26 J=L,N
              V(J,I)=(A(I,J)/A(I,L))/G
26          CONTINUE
            DO 29 J=L,N
              S=0.0
              DO 27 K=L,N
                S=S+A(I,K)*V(K,J)
27            CONTINUE
              DO 28 K=L,N
                V(K,J)=V(K,J)+S*V(K,I)
28            CONTINUE
29          CONTINUE
          ENDIF
          DO 31 J=L,N
            V(I,J)=0.0
            V(J,I)=0.0
31        CONTINUE
        ENDIF
        V(I,I)=1.0
        G=RV1(I)
        L=I
32    CONTINUE
      DO 39 I=N,1,-1
        L=I+1
        G=W(I)
        IF (I.LT.N) THEN
          DO 33 J=L,N
            A(I,J)=0.0
33        CONTINUE
        ENDIF
        IF (G.NE.0.0) THEN
          G=1.0/G
          IF (I.NE.N) THEN
            DO 36 J=L,N
              S=0.0
              DO 34 K=L,M
                S=S+A(K,I)*A(K,J)
34            CONTINUE
              F=(S/A(I,I))*G
              DO 35 K=I,M
                A(K,J)=A(K,J)+F*A(K,I)
35            CONTINUE
36          CONTINUE
          ENDIF
          DO 37 J=I,M
            A(J,I)=A(J,I)*G
37        CONTINUE
        ELSE
          DO 38 J= I,M
            A(J,I)=0.0
38        CONTINUE
        ENDIF
        A(I,I)=A(I,I)+1.0
39    CONTINUE
      DO 49 K=N,1,-1
        DO 48 ITS=1,30
          DO 41 L=K,1,-1
            NM=L-1
            IF ((ABS(RV1(L))+ANORM).EQ.ANORM)  GO TO 2
            IF ((ABS(W(NM))+ANORM).EQ.ANORM)  GO TO 1
41        CONTINUE
1         C=0.0
          S=1.0
          DO 43 I=L,K
            F=S*RV1(I)
            IF ((ABS(F)+ANORM).NE.ANORM) THEN
              G=W(I)
              H=SQRT(F*F+G*G)
              W(I)=H
              H=1.0/H
              C= (G*H)
              S=-(F*H)
              DO 42 J=1,M
                Y=A(J,NM)
                Z=A(J,I)
                A(J,NM)=(Y*C)+(Z*S)
                A(J,I)=-(Y*S)+(Z*C)
42            CONTINUE
            ENDIF
43        CONTINUE
2         Z=W(K)
          IF (L.EQ.K) THEN
            IF (Z.LT.0.0) THEN
              W(K)=-Z
              DO 44 J=1,N
                V(J,K)=-V(J,K)
44            CONTINUE
            ENDIF
            GO TO 3
          ENDIF
          IF (ITS.EQ.30)  stop 
		  !PAUSE 'No convergence in 30 iterations'
		  !pause
          X=W(L)
          NM=K-1
          Y=W(NM)
          G=RV1(NM)
          H=RV1(K)
          F=((Y-Z)*(Y+Z)+(G-H)*(G+H))/(2.0*H*Y)
          G=SQRT(F*F+1.0)
          F=((X-Z)*(X+Z)+H*((Y/(F+SIGN(G,F)))-H))/X
          C=1.0
          S=1.0
          DO 47 J=L,NM
            I=J+1
            G=RV1(I)
            Y=W(I)
            H=S*G
            G=C*G
            Z=SQRT(F*F+H*H)
            RV1(J)=Z
            C=F/Z
            S=H/Z
            F= (X*C)+(G*S)
            G=-(X*S)+(G*C)
            H=Y*S
            Y=Y*C
            DO 45 NM=1,N
              X=V(NM,J)
              Z=V(NM,I)
              V(NM,J)= (X*C)+(Z*S)
              V(NM,I)=-(X*S)+(Z*C)
45          CONTINUE
            Z=SQRT(F*F+H*H)
            W(J)=Z
            IF (Z.NE.0.0) THEN
              Z=1.0/Z
              C=F*Z
              S=H*Z
            ENDIF
            F= (C*G)+(S*Y)
            X=-(S*G)+(C*Y)
            DO 46 NM=1,M
              Y=A(NM,J)
              Z=A(NM,I)
              A(NM,J)= (Y*C)+(Z*S)
              A(NM,I)=-(Y*S)+(Z*C)
46          CONTINUE
47        CONTINUE
          RV1(L)=0.0
          RV1(K)=F
          W(K)=X
48      CONTINUE
3       CONTINUE
49    continue
      RETURN
      END

**LittleOwl** · 20/09/2013, 15h27

Bonjour,

NMAX est une constante (parameter) entière locale, dont la valeur fixe la taille du tableau à une dimension RV.

Afin de pouvoir mieux répondre à votre question, pouvez-vous préciser ce que vous entendez par "ne pas donner de bons résultats" ? Cela signifie t'il que, dans le cas d'une matrice de grande taille, le programme ne s'éxécute pas correctement et produit une erreur d'exécution, auquel cas, il faudrait connaitre non seulement le corps de la subroutine mais également le contexte de l'appel à la subroutine et le message donné par le compilateur pour répondre à votre question, ou voulez-vous dire que le programme s'exécute et se termine mais produit un résultat erroné?

**BAmal** · 20/09/2013, 16h46

Je vous remercie pour votre intérêt.
Alors ,non le programme s exécute normalement (c à dire y a pas des erreurs) mais le problème c est qu il m affiche pas des résultats, quand je lui demande de m afficher les matrices u,w et v ,il m affiche rien ,donc je me suis dit peut etre il arrive pas à faire une svd pour une matrice de grande taille parce que sur le meme programme je fais appel à la la même subroutine (pour une matrice carrée de 6*6) et j obtiens de bons résultats.

**FR119492** · 23/09/2013, 23h59

Salut!
Si tu ne trouves pas la solution autrement, essaie avec la routine ssvdc, de la bibliothèque linpack, que tu trouveras sur le site www.netlib.org
Jean-Marc Blanc

**Ehouarn** · 25/09/2013, 07h43

Bonjour,

Envoyé par BAmal

Alors ,non le programme s exécute normalement (c à dire y a pas des erreurs)

Même si le programme est compilé avec des options usuelles de débogage?

**BAmal** · 25/09/2013, 15h10

Vous voulez dire quoi par '' le programme est compilé avec des options usuelles de débogage?'' ?

**Ehouarn** · 26/09/2013, 07h45

Bonjour,

Et bien je pense aux options de compilations qui effectuent des vérifications complémentaires lors de la compilation et l'exécution du code, du genre: vérification de dépassement d'indices de tableaux, utilisation de variables dont les valeurs ne sont pas initialisées, etc.

Par exemple avec le compilateur gfortran (chaque compilateur a sa propre nomenclature, voir les docs respectives):

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
 
 gfortran -ffpe-trap=invalid,zero,overflow -fbounds-check -g3 -O0 -fstack-protector-all monprogramme.f

ces options sont coûteuses (en terme de temps de calcul); il est donc habituel de les utiliser lors de la mise au point d'un programme, mais pas lorsque l'on passe à l'utilisation du code pour produire des résultats.

Bonne continuation.

**BAmal** · 28/10/2013, 16h47

Si quelq un rencontre un souci lors de l utilisation de la subroutine SVD merci de me contacter

PROBLEME SUBROUTINE SVD NUMERICAL RECIPES

Fortran

Discussions similaires

Partager

Partager