New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 4578 for branches/2012 – NEMO

Changeset 4578 for branches/2012


Ignore:
Timestamp:
2014-03-26T10:02:56+01:00 (10 years ago)
Author:
pabouttier
Message:

Add missing initalization in dynspg_flt_tam, and correct adjoint test routines, see Ticket #1273

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2012/dev_v3_4_STABLE_2012/NEMOGCM/NEMO/OPATAM_SRC/DYN/dynspg_flt_tam.F90

    r3611 r4578  
    4646   USE dom_oce 
    4747   USE solver 
     48   USE dynspg_flt 
    4849   USE sol_oce 
    4950   USE oce_tam 
     
    6364   USE lib_fortran 
    6465   USE timing 
     66   USE iom 
    6567 
    6668 
     
    149151         spgu_tl(:,:) = 0.0_wp                    ! surface pressure gradient (i-direction) 
    150152         spgv_tl(:,:) = 0.0_wp                    ! surface pressure gradient (j-direction) 
    151          !CALL solver_init( nit000 )               ! Elliptic solver initialisation 
     153         ! Reinitialize the solver arrays 
     154         gcxb_tl(:,:) = 0.e0 
     155         gcx_tl (:,:) = 0.e0 
     156         CALL sol_mat( nit000 ) 
    152157      ENDIF 
    153158      ! Local constant initialization 
     
    168173            END DO 
    169174         END DO 
     175 
    170176         DO jk = 1, jpkm1              ! unweighted time stepping 
    171177            DO jj = 2, jpjm1 
     
    236242         END DO 
    237243      END DO 
     244 
    238245      ! apply the lateral boundary conditions 
    239246      IF( nn_solv == 2 .AND. MAX( jpr2di, jpr2dj ) > 0 ) CALL lbc_lnk_e( gcb_tl, c_solver_pt, 1.0_wp ) 
     
    291298         END DO 
    292299      END DO 
     300 
    293301      ! 
    294302      IF( nn_timing == 1 )  CALL timing_stop('dyn_spg_flt_tan') 
     
    370378         CALL sol_mat(kt)       ! initialize matrix 
    371379 
    372       ELSEIF (   neuler == 0 .AND. kt == nit000 + 1 ) THEN 
     380      ELSEIF (                   kt == nitend    ) THEN 
    373381 
    374382         z2dt = 2.0_wp * rdt    ! time step: leap-frog 
    375383         CALL sol_mat(kt)       ! reinitialize matrix 
    376384 
    377       ELSEIF (                   kt == nitend     ) THEN 
     385      ELSE 
    378386 
    379387         z2dt = 2.0_wp * rdt    ! time step: leap-frog 
    380          CALL sol_mat(kt)       ! reinitialize matrix 
    381  
    382       ELSEIF ( neuler /= 0 .AND. kt == nit000     ) THEN 
    383  
    384          z2dt = 2.0_wp * rdt    ! time step: leap-frog 
    385          CALL sol_mat(kt)       ! initialize matrix 
    386  
    387       ELSE 
    388  
    389          z2dt = 2.0_wp * rdt    ! time step: leap-frog 
    390388 
    391389      ENDIF 
    392390 
    393391      z2dtg  = grav * z2dt 
     392 
     393      ! set to zero free surface specific arrays (they are actually local variables) 
     394      spgu_ad(:,:) = 0.0_wp    ;      spgv_ad(:,:) = 0.0_wp 
    394395 
    395396      ! Add the trends multiplied by z2dt to the after velocity 
     
    546547                  ua_ad(  ji,jj,jk) = ua_ad(ji,jj,jk) * umask(ji,jj,jk) 
    547548                  ub_ad(  ji,jj,jk) = ub_ad(ji,jj,jk) + ua_ad(ji,jj,jk) 
    548                   spgu_ad(ji,jj   ) = spgu_ad(ji,jj)  + ua_ad(ji,jj,jk) * z2dt 
    549549                  ua_ad(  ji,jj,jk) = ua_ad(ji,jj,jk) * z2dt 
     550                  spgu_ad(ji,jj   ) = spgu_ad(ji,jj)  + ua_ad(ji,jj,jk) 
    550551                  va_ad(  ji,jj,jk) = va_ad(ji,jj,jk) * vmask(ji,jj,jk) 
    551552                  vb_ad(  ji,jj,jk) = vb_ad(ji,jj,jk) + va_ad(ji,jj,jk) 
    552                   spgv_ad(ji,jj   ) = spgv_ad(ji,jj)  + va_ad(ji,jj,jk) * z2dt 
    553553                  va_ad(  ji,jj,jk) = va_ad(ji,jj,jk) * z2dt 
     554                  spgv_ad(ji,jj   ) = spgv_ad(ji,jj)  + va_ad(ji,jj,jk) 
    554555               END DO 
    555556            END DO 
     
    573574         spgu_ad(:,:) = 0.0_wp                    ! surface pressure gradient (i-direction) 
    574575         spgv_ad(:,:) = 0.0_wp                    ! surface pressure gradient (j-direction) 
     576         ! Reinitialize the solver arrays 
     577         gcxb_ad(:,:) = 0.e0 
     578         gcx_ad (:,:) = 0.e0 
    575579      ENDIF 
    576580      ! 
     
    615619         & zua_tlout,   & ! Tangent output: ua_tl 
    616620         & zva_tlout,   & ! Tangent output: va_tl 
    617          & zub_tlout,   & ! Tangent output: ua_tl 
    618          & zvb_tlout,   & ! Tangent output: va_tl 
    619          & zub_adin,    & ! Tangent output: ua_ad 
    620          & zvb_adin,    & ! Tangent output: va_ad 
    621621         & zua_adin,    & ! Adjoint input: ua_ad 
    622622         & zva_adin,    & ! Adjoint input: va_ad 
     
    628628 
    629629      REAL(wp), DIMENSION(:,:), ALLOCATABLE :: & 
    630          & zgcx_tlin, zgcxb_tlin, zgcb_tlin, zgcx_tlout, zgcxb_tlout, zgcb_tlout,  & 
    631          & zgcx_adin, zgcxb_adin, zgcb_adin, zgcx_adout, zgcxb_adout, zgcb_adout,  & 
    632          & zspgu_tlout, zspgv_tlout, zspgu_adin, zspgv_adin 
     630         & zgcx_tlin, zgcxb_tlin, zgcx_tlout, zgcxb_tlout,  & 
     631         & zgcx_adin, zgcxb_adin, zgcx_adout, zgcxb_adout 
    633632 
    634633      REAL(KIND=wp), DIMENSION(:,:), ALLOCATABLE :: & 
    635634         & zsshn_tlin, & ! Tangent input: sshn_tl 
    636635         & zsshn_adout,& ! Adjoint output: sshn_ad 
    637          & zemp_tlin,  & ! Tangent input: emp_tl 
    638          & zemp_adout, & ! Adjoint output: emp_ad 
    639636         & znssh         ! 2D random field for SSH 
    640637      REAL(wp) :: & 
     
    654651      INTEGER ::             & 
    655652         & jpert 
    656       INTEGER, PARAMETER :: jpertmax = 7 
     653      INTEGER, PARAMETER :: jpertmax = 6 
    657654 
    658655      ! Allocate memory 
     
    676673         & zsshn_tlin(jpi,jpj), & 
    677674         & zsshn_adout(jpi,jpj),& 
    678          & zemp_tlin(jpi,jpj),  & 
    679          & zemp_adout(jpi,jpj), & 
    680675         & znssh(jpi,jpj)       & 
    681676         & ) 
    682677 
    683678      ALLOCATE( zgcx_tlin (jpi,jpj), zgcx_tlout (jpi,jpj), zgcx_adin (jpi,jpj), zgcx_adout (jpi,jpj),  & 
    684                 zgcxb_tlin(jpi,jpj), zgcxb_tlout(jpi,jpj), zgcxb_adin(jpi,jpj), zgcxb_adout(jpi,jpj),  & 
    685                 zgcb_tlin (jpi,jpj), zgcb_tlout (jpi,jpj), zgcb_adin (jpi,jpj), zgcb_adout (jpi,jpj)   & 
    686                & ) 
    687  
    688       ALLOCATE ( zub_tlout(jpi,jpj,jpk), zvb_tlout(jpi,jpj,jpk),  & 
    689                  zub_adin (jpi,jpj,jpk), zvb_adin (jpi,jpj,jpk)  ) 
    690  
    691       ALLOCATE( zspgu_tlout (jpi,jpj), zspgv_tlout (jpi,jpj), zspgu_adin (jpi,jpj), zspgv_adin (jpi,jpj)) 
     679                zgcxb_tlin(jpi,jpj), zgcxb_tlout(jpi,jpj), zgcxb_adin(jpi,jpj), zgcxb_adout(jpi,jpj)   ) 
    692680 
    693681      !========================================================================= 
     
    698686      ! Test for time steps nit000 and nit000 + 1 (the matrix changes) 
    699687 
    700       DO jstp = nit000, nit000 + 1 
    701          DO jpert = 1, jpertmax 
     688      DO jstp = nit000, nitend, nitend-nit000 
     689         DO jpert = jpertmax, jpertmax 
    702690            istp = jstp 
    703691 
     
    720708 
    721709            zsshn_tlin (:,:) = 0.0_wp 
    722             zemp_tlin  (:,:) = 0.0_wp 
    723710            zsshn_adout(:,:) = 0.0_wp 
    724             zemp_adout (:,:) = 0.0_wp 
    725             zspgu_adin (:,:) = 0.0_wp 
    726             zspgv_adin (:,:) = 0.0_wp 
    727             zspgu_tlout(:,:) = 0.0_wp 
    728             zspgv_tlout(:,:) = 0.0_wp 
    729711 
    730712            zgcx_tlout (:,:) = 0.0_wp ; zgcx_adin (:,:) = 0.0_wp ; zgcx_adout (:,:) = 0.0_wp 
    731713            zgcxb_tlout(:,:) = 0.0_wp ; zgcxb_adin(:,:) = 0.0_wp ; zgcxb_adout(:,:) = 0.0_wp 
    732             zgcb_tlout (:,:) = 0.0_wp ; zgcb_adin (:,:) = 0.0_wp ; zgcb_adout (:,:) = 0.0_wp 
    733714 
    734715            ub_tl(:,:,:) = 0.0_wp 
     
    737718            va_tl(:,:,:) = 0.0_wp 
    738719            sshn_tl(:,:) = 0.0_wp 
    739             emp_tl(:,:)  = 0.0_wp 
    740             gcb_tl(:,:)  = 0.0_wp 
    741720            gcx_tl(:,:)  = 0.0_wp 
    742721            gcxb_tl(:,:) = 0.0_wp 
     
    748727            va_ad(:,:,:) = 0.0_wp 
    749728            sshn_ad(:,:) = 0.0_wp 
    750             emp_ad(:,:)  = 0.0_wp 
    751729            gcb_ad(:,:)  = 0.0_wp 
    752730            gcx_ad(:,:)  = 0.0_wp 
     
    807785            ENDIF 
    808786            IF ( (jpert == 5) .OR. (jpert == jpertmax) ) THEN 
    809                CALL grid_random(  znssh, 'T', 0.0_wp, stdemp ) 
    810  
    811                DO jj = nldj, nlej 
    812                   DO ji = nldi, nlei 
    813                      zemp_tlin(ji,jj) = znssh(ji,jj) 
    814                   END DO 
    815                END DO 
    816  
    817             ENDIF 
    818             IF ( (jpert == 6) .OR. (jpert == jpertmax) ) THEN 
     787 
    819788               CALL grid_random(  znssh, 'T', 0.0_wp, stdssh ) 
    820789               DO jj = nldj, nlej 
     
    834803            ub_tl(:,:,:) = zub_tlin(:,:,:) 
    835804            vb_tl(:,:,:) = zvb_tlin(:,:,:) 
    836             emp_tl (:,:) = zemp_tlin (:,:) 
    837805            sshn_tl(:,:) = zsshn_tlin(:,:) 
    838806 
    839             gcx_tl (:,:) = 0.e0              ;   gcxb_tl(:,:) = 0.e0 
    840807            gcb_tl (:,:) = 0.e0 
    841808            gcx_tl (:,:) = zgcx_tlin (:,:)   ;   gcxb_tl(:,:) = zgcxb_tlin(:,:) 
    842809 
     810            CALL sol_mat( istp ) ! for nitend, it is not called in _tan so it is still set to the nit000 case    
    843811            CALL dyn_spg_flt_tan( istp, indic ) 
    844812 
    845813            zua_tlout(:,:,:) = ua_tl(:,:,:)   ;   zva_tlout(:,:,:) = va_tl(:,:,:) 
    846             zspgu_tlout(:,:) = spgu_tl(:,:)   ;   zspgv_tlout(:,:) = spgv_tl(:,:) 
    847             zgcb_tlout (:,:) = gcb_tl (:,:) 
     814            zgcxb_tlout(:,:) = gcxb_tl(:,:)   ;   zgcx_tlout (:,:) = gcx_tl (:,:) 
    848815 
    849816            !-------------------------------------------------------------------- 
     
    865832            DO jj = nldj, nlej 
    866833               DO ji = nldi, nlei 
    867                   zgcb_adin (ji,jj) = zgcb_tlout (ji,jj) & 
    868                      &              * e1t(ji,jj) * e2t(ji,jj) * fse3u(ji,jj,1) * tmask(ji,jj,1) 
    869                   zspgu_adin (ji,jj) = zspgu_tlout (ji,jj) & 
    870                      &              * e1u(ji,jj) * e2u(ji,jj) * fse3u(ji,jj,1) * umask(ji,jj,1) 
    871                   zspgv_adin(ji,jj) = zspgv_tlout(ji,jj) & 
    872                      &              * e1v(ji,jj) * e2v(ji,jj) * fse3v(ji,jj,1) * vmask(ji,jj,1) 
     834                  zgcx_adin (ji,jj) = zgcx_tlout (ji,jj) & 
     835                     &              * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1) * tmask(ji,jj,1) 
     836                  zgcxb_adin(ji,jj) = zgcxb_tlout(ji,jj) & 
     837                     &              * e1t(ji,jj) * e2t(ji,jj) * fse3t(ji,jj,1) * tmask(ji,jj,1) 
    873838               END DO 
    874839            END DO 
     
    878843            !-------------------------------------------------------------------- 
    879844 
    880             zsp1 =   DOT_PRODUCT( zua_tlout  , zua_adin   ) & 
    881                &   + DOT_PRODUCT( zgcb_tlout , zgcb_adin  ) & 
    882                &   + DOT_PRODUCT( zspgu_tlout , zspgu_adin  ) & 
    883                &   + DOT_PRODUCT( zspgv_tlout , zspgv_adin  ) & 
    884                &   + DOT_PRODUCT( zva_tlout  , zva_adin   ) 
     845            zsp1 =   DOT_PRODUCT( zua_tlout   , zua_adin   ) & 
     846               &   + DOT_PRODUCT( zgcx_tlout  , zgcx_adin  ) & 
     847               &   + DOT_PRODUCT( zgcxb_tlout , zgcxb_adin ) & 
     848               &   + DOT_PRODUCT( zva_tlout   , zva_adin   ) 
    885849 
    886850 
     
    892856            va_ad(:,:,:) = zva_adin(:,:,:) 
    893857 
    894             gcx_ad (:,:)   = 0.0_wp            ;   gcxb_ad(:,:) = 0.0_wp 
    895             gcb_ad (:,:)   = zgcb_adin (:,:) 
    896             spgu_ad(:,:)   = zspgu_adin(:,:) 
    897             spgv_ad(:,:)   = zspgv_adin(:,:) 
    898             ub_ad  (:,:,:) = zub_adin  (:,:,:) ;  vb_ad  (:,:,:) = zvb_adin  (:,:,:) 
     858            gcx_ad (:,:)   = zgcx_adin (:,:)   ;   gcxb_ad(:,:) = zgcxb_adin (:,:) 
     859            ub_ad  (:,:,:) = 0.0_wp ;  vb_ad  (:,:,:) = 0.0_wp 
    899860 
    900861            CALL dyn_spg_flt_adj( istp, indic ) 
     
    919880               &   + DOT_PRODUCT( zgcx_tlin , zgcx_adout  ) & 
    920881               &   + DOT_PRODUCT( zgcxb_tlin, zgcxb_adout ) & 
    921                &   + DOT_PRODUCT( zsshn_tlin, zsshn_adout ) & 
    922                &   + DOT_PRODUCT( zemp_tlin , zemp_adout  ) 
     882               &   + DOT_PRODUCT( zsshn_tlin, zsshn_adout ) 
    923883 
    924884            ! Compare the scalar products 
     
    936896                  cl_name = 'spg_flt  Vb T1' 
    937897               CASE(5) 
    938                   cl_name = 'spg_flt emp T1' 
    939                CASE(6) 
    940898                  cl_name = 'spg_flt ssh T1' 
    941899               CASE(jpertmax) 
    942900                  cl_name = 'dyn_spg_flt T1' 
    943901               END SELECT 
    944             ELSEIF ( istp == nit000 + 1 ) THEN 
     902            ELSEIF ( istp == nitend ) THEN 
    945903               SELECT CASE (jpert) 
    946904               CASE(1) 
     
    953911                  cl_name = 'spg_flt  Vb T2' 
    954912               CASE(5) 
    955                   cl_name = 'spg_flt emp T2' 
    956                CASE(6) 
    957913                  cl_name = 'spg_flt ssh T2' 
    958914               CASE(jpertmax) 
     
    965921      END DO 
    966922 
    967 !nn_nmod = kmod  ! restore initial frequency of test for the SOR solver 
     923      nitsor(:) = jp_it0adj  ! restore nitsor to avoid non reproducible results with or without the tests 
    968924 
    969925      ! Deallocate memory 
     
    986942      DEALLOCATE( & 
    987943         & zsshn_tlin, & 
    988          & zemp_tlin,  & 
    989944         & zsshn_adout,& 
    990          & zemp_adout, & 
    991945         & znssh       & 
    992946         & ) 
    993947      DEALLOCATE( zgcx_tlin , zgcx_tlout , zgcx_adin , zgcx_adout,  & 
    994          & zgcxb_tlin, zgcxb_tlout, zgcxb_adin, zgcxb_adout,  & 
    995          & zgcb_tlin , zgcb_tlout , zgcb_adin , zgcb_adout    & 
    996          & ) 
    997       DEALLOCATE ( zub_tlout, zvb_tlout, zub_adin , zvb_adin ) 
     948         & zgcxb_tlin, zgcxb_tlout, zgcxb_adin, zgcxb_adout         )  
    998949   END SUBROUTINE dyn_spg_flt_adj_tst 
    999950 
Note: See TracChangeset for help on using the changeset viewer.