/[MITgcm]/MITgcm/pkg/grdchk/grdchk_main.F
ViewVC logotype

Diff of /MITgcm/pkg/grdchk/grdchk_main.F

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

revision 1.3.6.6 by heimbach, Mon Jul 7 16:18:18 2003 UTC revision 1.16 by heimbach, Thu Jul 28 13:54:36 2005 UTC
# Line 1  Line 1 
1    C
2  C $Header$  C $Header$
3    C $Name$
4    
5    #include "AD_CONFIG.h"
6  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
7    
8  CBOI  CBOI
# Line 78  c     == global variables == Line 81  c     == global variables ==
81  #include "PARAMS.h"  #include "PARAMS.h"
82  #include "grdchk.h"  #include "grdchk.h"
83  #include "cost.h"  #include "cost.h"
84    #ifdef ALLOW_TANGENTLINEAR_RUN
85    #include "g_cost.h"
86    #endif
87    
88  C     !INPUT/OUTPUT PARAMETERS:  C     !INPUT/OUTPUT PARAMETERS:
89  c     == routine arguments ==  c     == routine arguments ==
90        integer mythid        integer mythid
91    
92  #ifdef ALLOW_GRADIENT_CHECK  #ifdef ALLOW_GRDCHK
93  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
94  c     == local variables ==  c     == local variables ==
95        integer myiter        integer myiter
# Line 94  c     == local variables == Line 100  c     == local variables ==
100        integer j,  jmin, jmax        integer j,  jmin, jmax
101        integer k        integer k
102    
       integer iobcs  
103        integer icomp        integer icomp
104        integer ichknum        integer ichknum
105        integer icvrec        integer icvrec
# Line 123  c     == local variables == Line 128  c     == local variables ==
128        _RL tmpplot2(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)        _RL tmpplot2(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)
129        _RL tmpplot3(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)        _RL tmpplot3(1-olx:snx+olx,1-oly:sny+oly,nr,nsx,nsy)
130    
131  #ifdef ALLOW_TANGENTLINEAR_RUN        CHARACTER*(MAX_LEN_MBUF) msgBuf
       _RL     g_fc  
       common /g_cost_r/ g_fc  
 #endif  
132    
133  c     == end of interface ==  c     == end of interface ==
134  CEOP  CEOP
# Line 141  c--   Set the loop ranges. Line 143  c--   Set the loop ranges.
143        imin = 1        imin = 1
144        imax = snx        imax = snx
145    
146          print *, 'ph-check entering grdchk_main '
147    
148  c--   initialise variables  c--   initialise variables
149        call grdchk_init( mythid )        call grdchk_init( mythid )
150    
# Line 150  cph   Gradient via adjoint has already b Line 154  cph   Gradient via adjoint has already b
154  cph   and so has unperturbed cost function,  cph   and so has unperturbed cost function,
155  cph   assuming all xx_ fields are initialised to zero.  cph   assuming all xx_ fields are initialised to zero.
156    
       fcref = fc  
157        ierr_grdchk = 0        ierr_grdchk = 0
158    cphadmtlm(
159          fcref = fc
160    cphadmtlm      fcref = objf_state_final(45,4,1,1,1)
161    cphadmtlm)
162    
163          print *, 'ph-check fcref = ', fcref
164    
165        do bj = jtlo, jthi        do bj = jtlo, jthi
166           do bi = itlo, ithi           do bi = itlo, ithi
# Line 173  cph   assuming all xx_ fields are initia Line 182  cph   assuming all xx_ fields are initia
182           grdchk_epsfac = 1. _d 0           grdchk_epsfac = 1. _d 0
183        end if        end if
184    
185        print *, 'ph-grd 3 -------------------------------'        print *, 'grad-res -------------------------------'
186        print ('(2a)'),        print ('(2a)'),
187       &     ' ph-grd 3  proc    #    i    j    k            fc ref',       &     ' grad-res  proc    #    i    j    k            fc ref',
188       &     '        fc + eps        fc - eps'       &     '        fc + eps        fc - eps'
189  #ifdef ALLOW_TANGENTLINEAR_RUN  #ifdef ALLOW_TANGENTLINEAR_RUN
190        print ('(2a)'),        print ('(2a)'),
191       &     ' ph-grd 3  proc    #    i    j    k          tlm grad',       &     ' grad-res  proc    #    i    j    k          tlm grad',
192       &     '         fd grad      1 - fd/tlm'       &     '         fd grad      1 - fd/tlm'
193  #else  #else
194        print ('(2a)'),        print ('(2a)'),
195       &     ' ph-grd 3  proc    #    i    j    k          adj grad',       &     ' grad-res  proc    #    i    j    k          adj grad',
196       &     '         fd grad      1 - fd/adj'       &     '         fd grad      1 - fd/adj'
197  #endif  #endif
198    
# Line 191  c--   Compute the finite difference appr Line 200  c--   Compute the finite difference appr
200  c--   Cycle through all processes doing NINT(nend-nbeg+1)/nstep  c--   Cycle through all processes doing NINT(nend-nbeg+1)/nstep
201  c--   gradient checks.  c--   gradient checks.
202    
203          if ( nbeg .EQ. 0 ) call grdchk_get_position( mythid )
204    
205        do icomp = nbeg, nend, nstep        do icomp = nbeg, nend, nstep
206    
207           ichknum = (icomp - nbeg)/nstep + 1           ichknum = (icomp - nbeg)/nstep + 1
# Line 242  c-- Line 253  c--
253  c--   2. perform tangent linear run  c--   2. perform tangent linear run
254              mytime = starttime              mytime = starttime
255              myiter = niter0              myiter = niter0
256    cphadmtlm(
257              g_fc = 0.              g_fc = 0.
258    cphadmtlm            do j=1,sny
259    cphadmtlm               do i=1,snx
260    cphadmtlm                  g_objf_state_final(i,j,1,1,1) = 0.
261    cphadmtlm                  g_objf_state_final(i,j,1,1,2) = 0.
262    cphadmtlm               enddo
263    cphadmtlm            enddo
264    cphadmtlm)
265              call g_the_main_loop( mytime, myiter, mythid )              call g_the_main_loop( mytime, myiter, mythid )
266    cphadmtlm(
267              ftlxxmemo = g_fc              ftlxxmemo = g_fc
268    cphadmtlm            ftlxxmemo = g_objf_state_final(45,4,1,1,1)
269    cphadmtlm)
270              _BARRIER              _BARRIER
271  c--  c--
272  c--   3. reset control vector  c--   3. reset control vector
# Line 282  c--   forward run with perturbed control Line 304  c--   forward run with perturbed control
304              mytime = starttime              mytime = starttime
305              myiter = niter0              myiter = niter0
306              call the_main_loop( mytime, myiter, mythid )              call the_main_loop( mytime, myiter, mythid )
307    cphadmtlm(
308              fcpertplus = fc              fcpertplus = fc
309    cphadmtlm            fcpertplus = objf_state_final(45,4,1,1,1)
310    cphadmtlm)
311                print *, 'ph-check fcpertplus = ', fcpertplus
312              _BARRIER              _BARRIER
313                                        
314  c--   Reset control vector.  c--   Reset control vector.
# Line 296  c--   Reset control vector. Line 322  c--   Reset control vector.
322              _BARRIER              _BARRIER
323    
324              fcpertminus = fcref              fcpertminus = fcref
325                print *, 'ph-check fcpertminus = ', fcpertminus
326    
327              if ( useCentralDiff ) then              if ( useCentralDiff ) then
328    
# Line 389  c*************************************** Line 416  c***************************************
416                 itestmem  ( ichknum ) = itest                 itestmem  ( ichknum ) = itest
417                 ierrmem   ( ichknum ) = ierr                 ierrmem   ( ichknum ) = ierr
418                                        
419                 print *, 'ph-grd 3 -------------------------------'                 print *, 'grad-res -------------------------------'
420                 print '(a,5I5,2x,3(1x,E15.9))', ' ph-grd 3 ',                 print '(a,5I5,2x,3(1x,E15.9))', ' grad-res ',
421       &              myprocid,ichknum,itilepos,jtilepos,layer,       &              myprocid,ichknum,itilepos,jtilepos,layer,
422       &              fcref, fcpertplus, fcpertminus       &              fcref, fcpertplus, fcpertminus
423  #ifdef ALLOW_TANGENTLINEAR_RUN  #ifdef ALLOW_TANGENTLINEAR_RUN
424                 print '(a,5I5,2x,3(1x,E15.9))', ' ph-grd 3 ',                 print '(a,5I5,2x,3(1x,E15.9))', ' grad-res ',
425       &              myprocid,ichknum,ichkmem(ichknum),       &              myprocid,ichknum,ichkmem(ichknum),
426       &              icompmem(ichknum),itestmem(ichknum),       &              icompmem(ichknum),itestmem(ichknum),
427       &              ftlxxmemo, gfd, ratio_ftl       &              ftlxxmemo, gfd, ratio_ftl
428                   WRITE(msgBuf,'(A34,2(1PE24.14,X))')
429         &              'precision_grdchk_result TLM ', fcref, ftlxxmemo
430                   CALL PRINT_MESSAGE
431         &              (msgBuf,standardMessageUnit,SQUEEZE_RIGHT,1)
432  #else  #else
433                 print '(a,5I5,2x,3(1x,E15.9))', ' ph-grd 3 ',                 print '(a,5I5,2x,3(1x,E15.9))', ' grad-res ',
434       &              myprocid,ichknum,ichkmem(ichknum),       &              myprocid,ichknum,ichkmem(ichknum),
435       &              icompmem(ichknum),itestmem(ichknum),       &              icompmem(ichknum),itestmem(ichknum),
436       &              adxxmemo, gfd, ratio_ad       &              adxxmemo, gfd, ratio_ad
437                   WRITE(msgBuf,'(A34,2(1PE24.14,X))')
438         &              'precision_grdchk_result ADM ', fcref, adxxmemo
439                   CALL PRINT_MESSAGE
440         &              (msgBuf,standardMessageUnit,SQUEEZE_RIGHT,1)
441  #endif  #endif
442    
443              endif              endif
# Line 438  c--   Everyone has to wait for the compo Line 473  c--   Everyone has to wait for the compo
473  c--   Print the results of the gradient check.  c--   Print the results of the gradient check.
474        call grdchk_print( ichknum, ierr_grdchk, mythid )        call grdchk_print( ichknum, ierr_grdchk, mythid )
475    
476  #endif /* ALLOW_GRADIENT_CHECK */  #endif /* ALLOW_GRDCHK */
477    
478        end        end
Legend:
Removed from v.1.3.6.6  
changed lines
  Added in v.1.16
  ViewVC Help
Powered by ViewVC 1.1.22