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
|
subroutine degl2(tg,ng,neq,imina,iminb,iminc,ta,tb,tc,iback,
& ntback,tback,mxsol,w,iw,ierr)
C!but
C Cette procedure a pour objectif de determiner le plus grand
C nombre de minimums de degre "neq".
C!liste d'appel
C subroutine degre (neq,imina,iminb,iminc,ta,tb,tc,
C & iback,ntback,tback)
C
C Entree :
C -neq. est le degre des polynomes parmi lesquels ont
C recherche les minimums.
C -imina. est le nombre de minimums de degre "neq-1"
C contenus dans ta.
C -iminb. est le nombre de minimums de degre "neq-2"
C contenus dans tb.
C -iminc. est le nombre de minimums de degre "neq-3"
C contenus dans tc.
C -ta. est le tableau contenant donc les minimums de degre
C "neq-1"
C -tb. est le tableau contenant donc les minimums de degre
C "neq-2"
C -tc. est le tableau contenant donc les minimums de degre
C "neq-3"
C -iback. est le nombre de minimums obtenus apres une
C intersection avec la frontiere
C -ntback est un tableau d'entier qui contient les degre
C de ces minimums
C -tback. est le tableau qui contient leurs coefficients,
C ou ils sont ordonnes degre par degre.
C
C Sortie :
C -imina. est le nombre de minimums de degre neq que l'on
C vient de determiner
C -iminb. est le nombre de minimums de degre "neq-1"
C -iminc. est le nombre de minimums de degre "neq-2"
C -ta. contient les mins de degre neq, -tb. ceux de degre
C neq-1 et tc ceux de degre neq-2
C -iback,ntback,tback ont pu etre augmente des mins obtenus
C apres intersection eventuelle avec la frontiere.
C
C tableaux de travail
C w : 33+33*neq+7*ng+neq*ng+neq**2*(ng+2)
C iw :29+neq**2+4*neq
c
c Copyright INRIA
C!
implicit double precision (a-h,o-y)
dimension ta(mxsol,*), tb(mxsol,*), tc(mxsol,*),tg(ng+1),
& ntback(mxsol), tback(mxsol,*)
dimension w(*), iw(*), xx(1)
C
dimension tps(0:1), tms(0:1)
double precision x,phi0,phi,gnrm
C
external feq, jacl2
common /comall/ nall1
common /sortie/ io,info,ll
common /no2f/ gnrm
C
tps(0) = 1.0d+0
tps(1) = 1.0d+0
tms(0) = -1.0d+0
tms(1) = 1.0d+0
C
C
C -------- Reinitialisation des tableaux --------
C
if (neq .eq. 1) goto 111
C
do 110 j = 1,iminb
call dcopy(neq,tb(j,1),mxsol,tc(j,1),mxsol)
110 continue
iminc = iminb
C
111 do 120 j = 1,imina
call dcopy(neq,ta(j,1),mxsol,tb(j,1),mxsol)
120 continue
iminb = imina
imina = 0
neq = neq + 1
neqbac = neq
c
lrw = neq**2 + 9*neq + 22
liw = 20+neq
C decoupage du tableau de travail w
ltq = 1
lwopt = ltq+6+6*neq+6*ng+neq*ng+neq**2*(ng+1)
ltr = lwopt +25+26*neq+ng+neq**2
lfree=ltr+neq+1
c
c les lrw elements de w suivant w(lwopt) ne doivent pas etre modifies
c d'un appel de optml2 a l'autre
lw=lwopt+lrw
ltg=ltq+neq+1
call dcopy(ng+1,tg,1,w(ltg),1)
C decoupage du tableau de travail iw
lneq=1
liwopt=lneq+3+(neq+1)*(neq+2)
lifree =liwopt + 20+neq
c
iw(lneq)=neq
iw(lneq+1)=ng
iw(lneq+2)=neq
if (info .gt. 0) call outl2(51,neq,neq,xx,xx,x,x)
C
C -------- Boucles de calculs --------
C
do 190 k = 1,iminb
C
call dcopy(neq-1,tb(k,1),mxsol,w(ltr),1)
w(ltr+neq-1) = 1.0d+0
C
do 180 imult = 1,2
C
if (imult .eq. 1) then
call dpmul1(w(ltr),neq-1,tps,1,w(ltq))
elseif (imult .eq. 2) then
call dpmul1(w(ltr),neq-1,tms,1,w(ltq))
endif
C
140 continue
C
nch = 1
call optml2(feq,jacl2,iw(lneq),w(ltq),nch,w(lwopt),iw(liwopt))
neq=iw(lneq)
if(info.gt.1) call outl2(nch,iw(lneq),neqbac,w(ltq),xx,x,x)
if (info .gt. 0) then
call lq(neq,w(ltq),w(lw),w(ltg),ng)
x=sqrt(gnrm)
call dscal(neq,x,w(lw),1)
call outl2(nch,neq,neqbac,w(ltq),w(lw),x,x)
phi0= abs(phi(w(ltq),neq,w(ltg),ng,w(lw)))
lqdot=lw
call feq(iw(lneq),t,w(ltq),w(lqdot))
call outl2(17,neq,neq,w(ltq),w(lqdot),phi0,x)
endif
if (nch.eq.15 .and. nall1.eq.0) then
ierr = 4
return
endif
C
if (nch .eq. -1) goto 140
if (nch .eq. -2) goto 140
C
nch = 2
call optml2(feq,jacl2,iw(lneq),w(ltq),nch,w(lwopt),iw(liwopt))
neq=iw(lneq)
if (info .gt. 1) then
call lq(neq,w(ltq),w(lw),w(ltg),ng)
x=sqrt(gnrm)
call dscal(neq,x,w(lw),1)
call outl2(nch,neq,neqbac,w(ltq),w(lw),x,x)
phi0= abs(phi(w(ltq),neq,w(ltg),ng,w(lw)))
lqdot=lw
call feq(iw(lneq),t,w(ltq),w(lqdot))
call outl2(17,neq,neq,w(ltq),w(lqdot),phi0,x)
endif
if (nch.eq.15 .and. nall1.eq.0) then
ierr = 4
return
endif
C
C
if (nch .eq. -1) goto 140
if (nch .eq. -2) goto 140
C
if (nch .eq. 15) then
if (info .gt. 0) call outl2(50,neq,neq,xx,xx,x,x)
goto 170
endif
C
nch = neq - neqbac
if (nch .eq. -2) then
call storl2(neq,w(ltq),w(ltg),ng,iminc,tc,iback,ntback,
& tback,nch,mxsol,w(lwopt),ierr)
elseif (nch .eq. -1) then
call storl2(neq,w(ltq),w(ltg),ng,iminb,tb,iback,ntback,
& tback,nch,mxsol,w(lwopt),ierr)
else
call storl2(neq,w(ltq),w(ltg),ng,imina,ta,iback,ntback,
& tback,nch,mxsol,w(lwopt),ierr)
endif
C
170 neq = neqbac
iw(lneq)=neq
C
180 continue
190 continue
if (info .gt. 0) then
x = real(mxsol)
call outl2(53,neq,imina,ta,xx,x,x)
endif
return
end
|
