/[MITgcm]/MITgcm/lsopt/lsline.F
ViewVC logotype

Contents of /MITgcm/lsopt/lsline.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph


Revision 1.3 - (show annotations) (download)
Fri Dec 6 01:42:25 2002 UTC (21 years, 3 months ago) by heimbach
Branch: MAIN
CVS Tags: checkpoint51k_post, checkpoint47e_post, checkpoint47c_post, checkpoint50c_post, checkpoint48e_post, checkpoint50c_pre, checkpoint51o_pre, checkpoint51l_post, checkpoint48i_post, checkpoint50d_pre, checkpoint51, checkpoint50, checkpoint52, checkpoint50d_post, checkpoint50b_pre, checkpoint51f_post, checkpoint48b_post, checkpoint51d_post, checkpoint48c_pre, checkpoint47d_pre, checkpoint51t_post, checkpoint51n_post, checkpoint51s_post, checkpoint48d_pre, checkpoint51j_post, checkpoint47i_post, checkpoint51n_pre, checkpoint47d_post, checkpoint48d_post, checkpoint48f_post, checkpoint51l_pre, checkpoint48h_post, checkpoint51q_post, checkpoint51b_pre, checkpoint47g_post, checkpoint51h_pre, checkpoint48a_post, checkpoint50f_post, checkpoint50a_post, checkpoint50f_pre, checkpoint47j_post, branch-exfmods-tag, branchpoint-genmake2, checkpoint51r_post, checkpoint48c_post, checkpoint51i_post, checkpoint51b_post, checkpoint51c_post, checkpoint50g_post, checkpoint52a_pre, checkpoint50h_post, checkpoint50e_pre, checkpoint50i_post, checkpoint51i_pre, checkpoint47f_post, checkpoint50e_post, checkpoint51e_post, checkpoint48, checkpoint49, checkpoint51o_post, checkpoint51f_pre, checkpoint48g_post, checkpoint47h_post, checkpoint52a_post, checkpoint51g_post, ecco_c52_e35, checkpoint50b_post, checkpoint51m_post, checkpoint51a_post, checkpoint51p_post, checkpoint51u_post
Branch point for: branch-exfmods-curt, branch-genmake2, branch-nonh, tg2-branch, checkpoint51n_branch
Changes since 1.2: +3 -3 lines
o lsopt:
  changed BLAS calls from single prec. (SDOT, SNRM2,SSCAL)
  to double prec. (DDOT, DNRM2, DSCAL)
  for compatibility with IBM SP3/SP4
o optim:
  bringing optim_readdata/optim_writedata formats up-to-date
  with latest ctrl_pack/ctrl_unpack formats.
NB: need to be merged in release1 and ecco-branch

1
2 subroutine lsline(
3 & simul
4 & , nn, ifail, lphprint
5 & , ifunc, nfunc
6 & , ff, dotdg
7 & , tmin, tmax, tact, epsx
8 & , dd, gg, xx, xdiff
9 & )
10
11 c ==================================================================
12 c SUBROUTINE lsline
13 c ==================================================================
14 c
15 c o line search algorithm for determining control vector update;
16 c After computing updated control vector from given gradient,
17 c a forward and adjoint model run are performed (simul.F)
18 c using the updated control vector.
19 c Tests are then applied to see whether solution has improved.
20 c
21 c o Reference: J.C. Gilbert & C. Lemarechal
22 c Some numerical experiments with variable-storage
23 c quasi-Newton algorithms
24 c Mathematical Programming 45 (1989), pp. 407-435
25 c
26 c o started: ??? not reproducible
27 c
28 c o changed: Patrick Heimbach, MIT/EAPS
29 c
30 c o Version: 2.0, 24-Feb-2000: Patrick Heimbach, MIT/EAPS
31 c - severe changes in structure including various
32 c shifts of variables which are only used in this
33 c routine
34 c - number of 3 control flags for error handling
35 c (indic, moderl, ifail) reduced to 1 (ifail)
36 c and homogenized with routine lsoptv
37 c
38 c o Version: 2.1.0, 02-Mar-2000: Patrick Heimbach, MIT/EAPS
39 c - initial computation of tact and
40 c xdiff = xx + tact*dd
41 c moved to new routine lsupdxx
42 c tmin, tmax, tact needed as parameters
43 c
44 c ==================================================================
45 c SUBROUTINE lsline
46 c ==================================================================
47
48 #include <blas1.h>
49
50 implicit none
51
52 c----------------------------------
53 c declare arguments
54 c----------------------------------
55 integer nn, ifail, ifunc, nfunc
56 double precision ff, dotdg, tmin, tmax, tact, epsx
57 double precision xx(nn), dd(nn), gg(nn), xdiff(nn)
58
59 logical lphprint
60
61 external simul
62
63 c----------------------------------
64 c declare local variables
65 c----------------------------------
66
67 double precision xpara1, xpara2
68 parameter( xpara1 = 0.0001, xpara2=0.9 )
69
70 double precision factor
71 parameter( factor = 0.2 )
72
73 double precision barmax
74 parameter( barmax = 0.3 )
75 double precision barmul
76 parameter( barmul = 5.0 )
77 double precision barmin
78 parameter( barmin = 0.01 )
79
80 integer i, indic
81
82 double precision tg, fg, td, ta
83 double precision fa, fpa, fp
84 double precision fnew, fdiff
85 double precision z, test, barr
86 double precision left, right, told
87
88 external DDOT
89 double precision DDOT
90
91 c----------------------------------
92 c check parameters
93 c----------------------------------
94
95 if ( (nn.le.0)
96 & .or. (dotdg.ge.0.0)
97 & .or. (xpara1.le.0.0) .or. (xpara1.ge.0.5)
98 & .or. (xpara2.le.xpara1) .or. (xpara2.ge.1.0) ) then
99 ifail = 9
100 go to 999
101 endif
102
103 c----------------------------------
104 c initialization
105 c----------------------------------
106 indic = 0
107
108 barr = barmin
109 fg = ff
110 fa = ff
111 fpa = dotdg
112
113 td = 0.0
114 tg = 0.0
115 ta = 0.0
116
117 c=======================================================================
118 c begin of simulation iter.
119 c=======================================================================
120
121 do ifunc = 1, nfunc
122
123 if (lphprint)
124 & print *, 'pathei-lsopt: ', ifunc, ' simul.'
125
126 c------------------------------------
127 c compute cost function and gradient
128 c------------------------------------
129 call simul( indic, nn, xdiff, fnew, gg )
130
131 fp = DDOT( nn, dd, 1, gg, 1 )
132 fdiff = fnew - ff
133
134 c-----------------------------------------
135 c apply 1st Wolfe test
136 c-----------------------------------------
137 if (fdiff .gt. tact*xpara1*dotdg) then
138 td = tact
139 ifail = 0
140 go to 500
141 endif
142
143 c-----------------------------------------
144 c 1st Wolfe test 1 ok, apply 2nd Wolf test
145 c-----------------------------------------
146 if (fp .gt. xpara2*dotdg) then
147 ifail = 0
148 go to 320
149 endif
150
151 if (ifail.eq.0) go to 350
152
153 c-----------------------------------------
154 c 2nd Wolfe test 2 ok, donc pas serieux, on sort
155 c-----------------------------------------
156
157 320 continue
158
159 ff = fnew
160 do i = 1, nn
161 xx(i) = xdiff(i)
162 end do
163 cph(
164 if (lphprint)
165 & print *, 'pathei-lsopt: no inter-/extrapolation in lsline'
166 cph)
167 go to 999
168
169
170 c-----------------------------------------
171 c extrapolation
172 c-----------------------------------------
173
174 350 continue
175
176 tg = tact
177 fg = fnew
178 if (td .ne. 0.0) go to 500
179
180 told = tact
181 left = (1.0+barmin)*tact
182 right = 10.0*tact
183 call cubic( tact, fnew, fp, ta, fa, fpa, left, right )
184 ta = told
185 if (tact.ge.tmax) then
186 ifail = 7
187 tact = tmax
188 endif
189
190 if (lphprint)
191 & print *, 'pathei-lsopt: extrapolation: ',
192 & 'td, tg, tact, ifail = ', td, tg, tact, ifail
193
194 go to 900
195
196 c-----------------------------------------
197 c interpolation
198 c-----------------------------------------
199
200 500 continue
201
202 test = barr*(td-tg)
203 left = tg+test
204 right = td-test
205 told = tact
206 call cubic( tact, fnew, fp, ta, fa, fpa, left, right )
207 ta = told
208 if (tact.gt.left .and. tact.lt.right) then
209 barr = dmax1( barmin, barr/barmul )
210 else
211 barr = dmin1( barmul*barr, barmax )
212 endif
213
214 if (lphprint)
215 & print *, 'pathei-lsopt: interpolation: ',
216 & 'td, tg, tact, ifail = ', td, tg, tact, ifail
217
218 c-----------------------------------------
219 c end of iteration loop
220 c - tact peut etre bloque sur tmax
221 c (venant de lextrapolation avec ifail=7)
222 c-----------------------------------------
223
224 900 continue
225
226 fa = fnew
227 fpa = fp
228
229 c-----------------------------------------
230 c --- faut-il continuer ?
231 c-----------------------------------------
232 if (td .eq. 0.0) go to 950
233 if (td-tg .lt. tmin) go to 920
234
235 c-----------------------------------------
236 c limit due to machine precision
237 c-----------------------------------------
238 do i = 1, nn
239 z = xx(i) + tact*dd(i)
240 if ((z.ne.xx(i)) .and. (z.ne.xdiff(i))) go to 950
241 end do
242
243 c-----------------------------------------
244 c arret sur dxmin ou de secours
245 c-----------------------------------------
246 920 continue
247 ifail = 8
248
249 c-----------------------------------------
250 c si tg=0, xx = xx_depart,
251 c sinon on prend xx=x_left qui fait decroitre f
252 c-----------------------------------------
253 if (tg .ne. 0.0) then
254 ff = fg
255 do i = 1, nn
256 xx(i) = xx(i) + tg*dd(i)
257 end do
258 endif
259
260 go to 999
261
262 c-----------------------------------------
263 c update vector for new simulation
264 c-----------------------------------------
265 950 continue
266
267 do i = 1, nn
268 xdiff(i) = xx(i) + tact*dd(i)
269 end do
270
271 c=======================================================================
272 c end of simulation iter.
273 c=======================================================================
274 end do
275
276 c-----------------------------------------
277 c too many function evaluations
278 c-----------------------------------------
279 ifail = 6
280 ff = fg
281 do i = 1, nn
282 xx(i) = xx(i) + tg*dd(i)
283 end do
284
285 c-----------------------------------------
286 c end of routine
287 c-----------------------------------------
288 999 continue
289
290 return
291
292 end

  ViewVC Help
Powered by ViewVC 1.1.22