/[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.4 by heimbach, Fri Mar 7 04:01:59 2003 UTC revision 1.20 by heimbach, Sat Mar 24 02:25:29 2007 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 100  c     == local variables == Line 106  c     == local variables ==
106        integer jtile        integer jtile
107        integer itile        integer itile
108        integer layer        integer layer
109          integer obcspos
110        integer itilepos        integer itilepos
111        integer jtilepos        integer jtilepos
112          integer icglo
113        integer itest        integer itest
114        integer ierr        integer ierr
115        integer ierr_grdchk        integer ierr_grdchk
# Line 121  c     == local variables == Line 129  c     == local variables ==
129        _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)
130        _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)
131    
132  #ifdef ALLOW_TANGENTLINEAR_RUN        CHARACTER*(MAX_LEN_MBUF) msgBuf
       _RL     g_fc  
       common /g_cost_r/ g_fc  
 #endif  
133    
134  c     == end of interface ==  c     == end of interface ==
135  CEOP  CEOP
# Line 139  c--   Set the loop ranges. Line 144  c--   Set the loop ranges.
144        imin = 1        imin = 1
145        imax = snx        imax = snx
146    
147          print *, 'ph-check entering grdchk_main '
148    
149  c--   initialise variables  c--   initialise variables
150        call grdchk_init( mythid )        call grdchk_init( mythid )
151    
# Line 148  cph   Gradient via adjoint has already b Line 155  cph   Gradient via adjoint has already b
155  cph   and so has unperturbed cost function,  cph   and so has unperturbed cost function,
156  cph   assuming all xx_ fields are initialised to zero.  cph   assuming all xx_ fields are initialised to zero.
157    
       fcref = fc  
158        ierr_grdchk = 0        ierr_grdchk = 0
159    cphadmtlm(
160          fcref = fc
161    cphadmtlm      fcref = objf_state_final(45,4,1,1,1)
162    cphadmtlm)
163    
164          print *, 'ph-check fcref = ', fcref
165    
166        do bj = jtlo, jthi        do bj = jtlo, jthi
167           do bi = itlo, ithi           do bi = itlo, ithi
# Line 171  cph   assuming all xx_ fields are initia Line 183  cph   assuming all xx_ fields are initia
183           grdchk_epsfac = 1. _d 0           grdchk_epsfac = 1. _d 0
184        end if        end if
185    
186        print *, 'ph-grd 3 -------------------------------'        print *, 'grad-res -------------------------------'
187        print ('(2a)'),        print '(2a)',
188       &     ' ph-grd 3  proc    #    i    j    k            fc ref',       &    ' grad-res  proc    #    i    j    k   bi   bj iobc',
189       &     '        fc + eps        fc - eps'       &    '               fc ref           fc + eps           fc - eps'
190  #ifdef ALLOW_TANGENTLINEAR_RUN  #ifdef ALLOW_TANGENTLINEAR_RUN
191        print ('(2a)'),        print '(2a)',
192       &     ' ph-grd 3  proc    #    i    j    k          tlm grad',       &    ' grad-res  proc    #    i    j    k   bi   bj iobc',
193       &     '         fd grad      1 - fd/tlm'       &    '             tlm grad            fd grad         1 - fd/tlm'
194  #else  #else
195        print ('(2a)'),        print '(2a)',
196       &     ' ph-grd 3  proc    #    i    j    k          adj grad',       &    ' grad-res  proc    #    i    j    k   bi   bj iobc',
197       &     '         fd grad      1 - fd/adj'       &    '             adj grad            fd grad         1 - fd/adj'
198  #endif  #endif
199    
200  c--   Compute the finite difference approximations.  c--   Compute the finite difference approximations.
201  c--   Cycle through all processes doing NINT(nend-nbeg+1)/nstep  c--   Cycle through all processes doing NINT(nend-nbeg+1)/nstep
202  c--   gradient checks.  c--   gradient checks.
203    
204          if ( nbeg .EQ. 0 )
205         &     call grdchk_get_position( mythid )
206    
207        do icomp = nbeg, nend, nstep        do icomp = nbeg, nend, nstep
208    
209           ichknum = (icomp - nbeg)/nstep + 1           ichknum = (icomp - nbeg)/nstep + 1
210    
211    cph(
212    cph-print         print *, 'ph-grd _main: nbeg, icomp, ichknum ',
213    cph-print     &        nbeg, icomp, ichknum
214    cph)
215           if (ichknum .le. maxgrdchecks ) then           if (ichknum .le. maxgrdchecks ) then
216    
217  c--         Determine the location of icomp on the grid.  c--         Determine the location of icomp on the grid.
218              if ( myProcId .EQ. grdchkwhichproc ) then              if ( myProcId .EQ. grdchkwhichproc ) then
219                 call grdchk_loc( icomp, ichknum,                 call grdchk_loc( icomp, ichknum,
220       &              icvrec, itile, jtile, layer,       &              icvrec, itile, jtile, layer, obcspos,
221       &              itilepos, jtilepos, itest, ierr,       &              itilepos, jtilepos, icglo, itest, ierr,
222       &              mythid )       &              mythid )
223    cph(
224    cph-print               print *, 'ph-grd ----- back from loc -----',
225    cph-print     &             icvrec, itilepos, jtilepos, layer, obcspos
226    cph)
227              endif              endif
228              _BARRIER              _BARRIER
229                                
# Line 240  c-- Line 263  c--
263  c--   2. perform tangent linear run  c--   2. perform tangent linear run
264              mytime = starttime              mytime = starttime
265              myiter = niter0              myiter = niter0
266    cphadmtlm(
267              g_fc = 0.              g_fc = 0.
268    cphadmtlm            do j=1,sny
269    cphadmtlm               do i=1,snx
270    cphadmtlm                  g_objf_state_final(i,j,1,1,1) = 0.
271    cphadmtlm                  g_objf_state_final(i,j,1,1,2) = 0.
272    cphadmtlm               enddo
273    cphadmtlm            enddo
274    cphadmtlm)
275              call g_the_main_loop( mytime, myiter, mythid )              call g_the_main_loop( mytime, myiter, mythid )
276    cphadmtlm(
277              ftlxxmemo = g_fc              ftlxxmemo = g_fc
278    cphadmtlm            ftlxxmemo = g_objf_state_final(45,4,1,1,1)
279    cphadmtlm)
280              _BARRIER              _BARRIER
281  c--  c--
282  c--   3. reset control vector  c--   3. reset control vector
# Line 280  c--   forward run with perturbed control Line 314  c--   forward run with perturbed control
314              mytime = starttime              mytime = starttime
315              myiter = niter0              myiter = niter0
316              call the_main_loop( mytime, myiter, mythid )              call the_main_loop( mytime, myiter, mythid )
317    cphadmtlm(
318              fcpertplus = fc              fcpertplus = fc
319    cphadmtlm            fcpertplus = objf_state_final(45,4,1,1,1)
320    cphadmtlm)
321                print *, 'ph-check fcpertplus = ', fcpertplus
322              _BARRIER              _BARRIER
323                                        
324  c--   Reset control vector.  c--   Reset control vector.
# Line 294  c--   Reset control vector. Line 332  c--   Reset control vector.
332              _BARRIER              _BARRIER
333    
334              fcpertminus = fcref              fcpertminus = fcref
335                print *, 'ph-check fcpertminus = ', fcpertminus
336    
337              if ( useCentralDiff ) then              if ( useCentralDiff ) then
338    
# Line 381  c*************************************** Line 420  c***************************************
420                 ilocmem   ( ichknum ) = itilepos                 ilocmem   ( ichknum ) = itilepos
421                 jlocmem   ( ichknum ) = jtilepos                 jlocmem   ( ichknum ) = jtilepos
422                 klocmem   ( ichknum ) = layer                 klocmem   ( ichknum ) = layer
423                   iobcsmem  ( ichknum ) = obcspos
424                 icompmem  ( ichknum ) = icomp                 icompmem  ( ichknum ) = icomp
425                 ichkmem   ( ichknum ) = ichknum                 ichkmem   ( ichknum ) = ichknum
426                 itestmem  ( ichknum ) = itest                 itestmem  ( ichknum ) = itest
427                 ierrmem   ( ichknum ) = ierr                 ierrmem   ( ichknum ) = ierr
428                                     icglomem  ( ichknum ) = icglo
429                 print *, 'ph-grd 3 -------------------------------'  
430                 print '(a,5I5,2x,3(1x,E15.9))', ' ph-grd 3 ',                 print *, 'grad-res -------------------------------'
431       &              myprocid,ichknum,itilepos,jtilepos,layer,                 print '(a,8I5,2x,3(1x,E18.12))', ' grad-res ',
432         &              myprocid,ichknum,itilepos,jtilepos,
433         &              layer,itile,jtile,obcspos,
434       &              fcref, fcpertplus, fcpertminus       &              fcref, fcpertplus, fcpertminus
435  #ifdef ALLOW_TANGENTLINEAR_RUN  #ifdef ALLOW_TANGENTLINEAR_RUN
436                 print '(a,5I5,2x,3(1x,E15.9))', ' ph-grd 3 ',                 print '(a,8I5,2x,3(1x,E18.12))', ' grad-res ',
437       &              myprocid,ichknum,ichkmem(ichknum),       &              myprocid,ichknum,ichkmem(ichknum),
438       &              icompmem(ichknum),itestmem(ichknum),       &              icompmem(ichknum),itestmem(ichknum),
439         &              bimem(ichknum),bjmem(ichknum),iobcsmem(ichknum),
440       &              ftlxxmemo, gfd, ratio_ftl       &              ftlxxmemo, gfd, ratio_ftl
441                   WRITE(msgBuf,'(A34,2(1PE24.14,X))')
442         &              'precision_grdchk_result TLM ', fcref, ftlxxmemo
443                   CALL PRINT_MESSAGE
444         &              (msgBuf,standardMessageUnit,SQUEEZE_RIGHT,1)
445  #else  #else
446                 print '(a,5I5,2x,3(1x,E15.9))', ' ph-grd 3 ',                 print '(a,8I5,2x,3(1x,E18.12))', ' grad-res ',
447       &              myprocid,ichknum,ichkmem(ichknum),       &              myprocid,ichknum,ichkmem(ichknum),
448       &              icompmem(ichknum),itestmem(ichknum),       &              icompmem(ichknum),itestmem(ichknum),
449         &              bimem(ichknum),bjmem(ichknum),obcspos,
450       &              adxxmemo, gfd, ratio_ad       &              adxxmemo, gfd, ratio_ad
451                   WRITE(msgBuf,'(A34,2(1PE24.14,X))')
452         &              'precision_grdchk_result ADM ', fcref, adxxmemo
453                   CALL PRINT_MESSAGE
454         &              (msgBuf,standardMessageUnit,SQUEEZE_RIGHT,1)
455  #endif  #endif
456    
457              endif              endif
# Line 435  c--   Everyone has to wait for the compo Line 487  c--   Everyone has to wait for the compo
487  c--   Print the results of the gradient check.  c--   Print the results of the gradient check.
488        call grdchk_print( ichknum, ierr_grdchk, mythid )        call grdchk_print( ichknum, ierr_grdchk, mythid )
489    
490  #endif /* ALLOW_GRADIENT_CHECK */  #endif /* ALLOW_GRDCHK */
491    
492        end        end
Legend:
Removed from v.1.3.6.4  
changed lines
  Added in v.1.20
  ViewVC Help
Powered by ViewVC 1.1.22