/[MITgcm]/MITgcm_contrib/mlosch/optim_m1qn3/optim_sub.F
ViewVC logotype

Contents of /MITgcm_contrib/mlosch/optim_m1qn3/optim_sub.F

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


Revision 1.10 - (show annotations) (download)
Tue Feb 12 15:50:49 2019 UTC (5 years, 1 month ago) by mlosch
Branch: MAIN
CVS Tags: HEAD
Changes since 1.9: +2 -2 lines
fix a bug reported by Andrew McRae: reduce memory allocation
from ndz = 3*nn+nupdate*(3*nn+1)to ndz = 4*nn+nupdate*(2*nn+1)

1 C $Header: /u/gcmpack/MITgcm_contrib/mlosch/optim_m1qn3/optim_sub.F,v 1.9 2018/05/02 10:06:42 mlosch Exp $
2 C $Name: $
3
4 C ECCO_CPPOPTIONS used to affect maxcvars and defined ALLOW_OBCS?_CONTROL
5 C#include "ECCO_CPPOPTIONS.h"
6 C now:
7 C CTRL_OPTIONS affects maxcvars and may define ALLOW_OBCS?_CONTROL
8 #include "CTRL_OPTIONS.h"
9
10 subroutine optim_sub(
11 I nn, ff
12 & )
13
14 c ==================================================================
15 c subroutine optim_sub
16 c ==================================================================
17 c
18 c o This is the main driver routine for the offline version of
19 c m1qn3 (m1qn3_offline). It
20 c - sets all m1qn3 relevant parameters
21 c - reads the model state (control vector, cost function,
22 c and gradient)
23 c - reads the state of m1qn3_offline
24 c - calls m1qn3_offline
25 c - saves model state (control vector) and state of m1qn3_offline
26 c The routine is somewhat lengthy and could be split into separate
27 c subroutines, but I felt that it was easier to write and test it
28 c in this form
29 c Martin Losch (Martin.Losch@awi.de), Apr, 2012
30 c
31 c ==================================================================
32
33 implicit none
34
35 c == global variables ==
36
37 #include "EEPARAMS.h"
38 #include "SIZE.h"
39
40 #if (defined (ALLOW_GENARR2D_CONTROL) || defined (ALLOW_GENARR3D_CONTROL) || defined (ALLOW_GENTIM2D_CONTROL))
41 # include "CTRL_SIZE.h"
42 #endif
43 #include "ctrl.h"
44 #include "optim.h"
45 #include "m1qn3_common.h"
46
47 c == routine arguments ==
48
49 integer nn
50 _RL ff
51
52 c == local variables ==
53
54 _RL objf
55
56 #ifdef DYNAMIC
57 _RL, dimension(:), allocatable :: xx, adxx
58 #else
59 integer nmax
60 parameter( nmax = MAX_INDEPEND )
61 _RL xx(nmax)
62 _RL adxx(nmax)
63 #endif
64 _RL xxmean
65
66 CML logical coldStart
67 c formal parameters of m1qn3
68 integer reverse
69 integer impres,imode(3),omode,niter,nsim,iz(5),indic
70 _RL dxmin,df1
71 character*3 normtype
72 c work arrays
73 integer ndz
74 CML _RL dz(ndz)
75 double precision, dimension(:), allocatable :: dz
76 c extra dummy variables
77 integer izs(1)
78 _RS rzs(1)
79 _RL dzs(1)
80 integer, parameter :: io = 60
81 character*(*), parameter :: fname_m1qn3 = 'm1qn3_output.txt'
82 c end of m1qn3 parameters
83
84 integer i
85
86 c == external ==
87
88 external simul_rc,euclid,ctonbe,ctcabe
89
90 c == end of interface ==
91
92 c-- Allocate memory for the control variables and the gradient vector.
93 #if defined(DYNAMIC)
94 allocate( xx(nn) )
95 allocate( adxx(nn) )
96 #endif
97
98 #ifndef DYNAMIC
99 if (nn .gt. nmax) then
100 print*,' OPTIMUM: Not enough space.'
101 print*,' nmax = ',nmax
102 print*,' nn = ',nn
103 print*
104 print*,' Set MAX_INDEPEND in Makefile .ge. ',nn
105 print*
106 stop ' ... stopped in OPTIMUM.'
107 endif
108 #endif
109
110 print*, ' OPTIM_SUB: Calling m1qn3_optim for iteration: ',
111 & optimcycle
112 print*, ' OPTIM_SUB: with nn, REAL_BYTE = ', nn, REAL_BYTE
113
114 c can be 'two','sup','dfn', see m1qn3 documentation for details
115 normtype='two'
116 c after reading data.optim some of these parameter values can be guessed
117 c impres=6, impres determines the amount of m1qn3-output see documentation
118 impres=iprint
119 c these should strictly be different (nsim>niter), but in practice
120 c it does not matter
121 niter = numiter
122 nsim = nfunc*niter
123 c epsg=1.d-8
124 dxmin=epsx
125 c will be set later
126 df1=-UNSET_RL
127 c
128 imode=(/0,1,0/)
129 omode=-1
130 c initialise work array
131 ndz = 4*nn+nupdate*(2*nn+1)
132 do i=1,5
133 iz(i)=0
134 enddo
135 allocate(dz(ndz))
136 do i=1,ndz
137 dz(i) = 0.
138 enddo
139 c these alway have to be set like this
140 reverse=1
141 indic=4
142 c initialise the dummy arguments that are not used
143 izs(1)=UNSET_I
144 rzs(1)=UNSET_RS
145 dzs(1)=UNSET_RL
146
147 c-- first read the model output into xx, adxx, and cost function
148 c value into objf
149 do i = 1,nn
150 xx(i) = 0.
151 adxx(i) = 0.
152 enddo
153 c
154 print *, ' OPTIM_SUB: read model state'
155 call optim_readdata( nn, ctrlname, .false., objf, xx )
156 call optim_readdata( nn, costname, .false., objf, adxx )
157 print *, ' OPTIM_SUB after reading ',
158 & ctrlname, ' and ', costname, ':'
159 print *, ' OPTIM_SUB nn = ', nn
160 print *, ' OPTIM_SUB objf = ', objf
161 print *, ' OPTIM_SUB xx(1) = ', xx(1)
162 print *, ' OPTIM_SUB adxx(1) = ', adxx(1)
163
164 if ( coldStart ) then
165 c-- cold start
166 print *, ' OPTIM_SUB: cold start, optimcycle =', optimcycle
167 imode(2) = 0
168 c this variable is the only one of the m1qn3-related common blocks
169 c that needs to be initialized here to make sure that we have a
170 c clean start
171 reentry = 0
172 c compute expected decrease of cost function from objf and fmin;
173 c this value is only used for a cold start of m1qn3_offline, for a
174 c warm start df1 is overwritten with data from a restart file
175 df1=objf-fmin
176 if ( df1 .le. 0. ) then
177 print *, ' OPTIM_SUB: df1 = objf-fmin = ', df1,
178 & ' but should be > 0.'
179 stop 'ABNORMAL in S/R OPTIM_SUB'
180 endif
181
182 c open output file for m1qn3
183 open(io,file=fname_m1qn3,status='unknown')
184 else
185 c-- warm restart
186 c requires restoring the state of m1qn3 with pickup file
187 print *, ' OPTIM_SUB: warm start, optimcycle =', optimcycle
188 imode(2) = 1
189 call optim_store_m1qn3(ndz,iz,dz,niter,nsim,epsg,df1,
190 I optimcycle,
191 I .false.)
192 c re-open output file for m1qn3
193 open(io,file=fname_m1qn3,status='old',position='append')
194 endif
195
196 c-- call the minimizer, a slightly modified version of m1qn3 v3.3
197 c (Gilbert & Lemarechal, 1989), downloaded in April 2012.
198 c simul_rc,euclid,ctonbe,ctcabe are external subroutines that
199 c are provided within m1qn3. euclid, ctonbe, ctcabe can be replaced
200 c by something more efficient, simul_rc is a dummy routine for
201 c the reverse communication mode and should not be changed.
202 print *, ' OPTIM_SUB: call m1qn3_offline ........'
203 call m1qn3_offline (simul_rc,euclid,ctonbe,ctcabe,
204 & nn,xx,objf,adxx,dxmin,df1,
205 & epsg,normtype,impres,io,imode,omode,niter,nsim,
206 & iz,dz,ndz,reverse,indic,izs,rzs,dzs)
207 close(io)
208 print *, ' OPTIM_SUB: ...........................'
209 print *, ' OPTIM_SUB: returned from m1qn3_offline'
210 print *, ' OPTIM_SUB: nn = ', nn
211 print *, ' OPTIM_SUB: xx(1) = ', xx(1), xx(2)
212 print *, ' OPTIM_SUB: adxx(1) = ', adxx(1), adxx(2)
213 print *, ' OPTIM_SUB: omode = ', omode
214 print *, ' OPTIM_SUB: niter = ', niter
215 print *, ' OPTIM_SUB: nsim = ', nsim
216 print *, ' OPTIM_SUB: reverse = ', reverse
217
218 c compute min/max/mean/std of output vector see if it is within
219 c reasonable bounds (prior to scaling)
220 xxmean = sum(xx)/dble(nn)
221 print *
222 print *,' OPTIM_SUB: mean(xx) =', xxmean
223 print *,' OPTIM_SUB: max(xx) =', maxval(xx)
224 print *,' OPTIM_SUB: min(xx) =', minval(xx)
225 print *,' OPTIM_SUB: std(xx) =',sum((xx-xxmean)**2)/dble(nn)
226 print *
227 c write state of m1qn3 into pickup file for warm restart
228 call optim_store_m1qn3(ndz,iz,dz,niter,nsim,epsg,df1,
229 I optimcycle,
230 I .true.)
231 c write model control vector
232 print *,' OPTIM_SUB: writing ', nn,' sized control to file ',
233 & ctrlname
234 c give the cost function a funny value to make sure that nobody
235 c mistakes it for the real one
236 call optim_writedata( nn, ctrlname, .false., -9999. _d 0, xx )
237
238 c clean up
239 #ifdef DYNAMIC
240 deallocate(xx, adxx)
241 #endif /* DYNAMIC */
242 deallocate(dz)
243
244 c stopping criterion
245 if ( reverse .lt. 0 ) then
246 print *, ' OPTIM_SUB: reverse = ', reverse
247 print *, ' OPTIM_SUB: optimization terminated with omode = ',
248 & omode
249 stop 'ABNORMAL in S/R OPTIM_SUB'
250 endif
251
252 return
253 end

  ViewVC Help
Powered by ViewVC 1.1.22