Changeset 8435


Ignore:
Timestamp:
2017-08-14T12:19:09+02:00 (3 years ago)
Author:
isabella
Message:

Added changes after review.

Location:
branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO/LIM_SRC_2/dom_ice_2.F90

    r7065 r8435  
    11MODULE dom_ice_2 
    2 !test 
    32   !!====================================================================== 
    43   !!                   ***  MODULE  dom_ice  *** 
  • branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO/OPA_SRC/ASM/bias.F90

    r7667 r8435  
    211211      !! * Local declarations 
    212212 
    213       CHARACTER(len=100) ::  cn_dir       ! dir for location ofline bias 
     213      CHARACTER(len=100) ::  cn_dir          ! dir for location ofline bias 
    214214      INTEGER            ::  ierror 
    215       INTEGER            ::  ios          ! Local integer output status for namelist read 
    216       REAL(wp)           ::  eft_rlx,  &  ! efolding time (bias memory) in days 
    217          &                   eft_asm,  &  !     "      " 
     215      INTEGER            ::  ios             ! Local integer output status for namelist read 
     216      REAL(wp)           ::  eft_rlx,  &     ! efolding time (bias memory) in days 
     217         &                   eft_asm,  &     !     "      " 
    218218         &                   log2,     & 
    219          &                   lenscl_bias, &  !lengthscale of the pressure bias decay between minlat and maxlat. 
    220          &                   minlat_bias, &     !used in ipc 
    221          &                   maxlat_bias    !used in ipc 
     219         &                   lenscl_bias, &  ! lengthscale pressure bias decay between minlat and maxlat. 
     220         &                   minlat_bias, &  ! used in ipc 
     221         &                   maxlat_bias     ! used in ipc 
    222222       
    223223      NAMELIST/nambias/ ln_bias, ln_bias_asm, ln_bias_rlx, ln_bias_ofl,   & 
     
    396396         tbias_i(:,:,:) = 0.0_wp 
    397397         sbias_i(:,:,:) = 0.0_wp 
    398          IF(lwp) WRITE(numout,*) 'gru_pc and grv_pc are set to zero ' 
    399398         gru_pc(:,:)    = 0.0_wp 
    400399         grv_pc(:,:)    = 0.0_wp 
     400          
    401401         IF ( ln_bias_rlx ) THEN 
    402402            tbias_rlx(:,:,:) = 0.0_wp 
    403403            sbias_rlx(:,:,:) = 0.0_wp 
    404404         ENDIF 
     405 
    405406         IF ( ln_bias_asm ) THEN   !now rlx and asm bias in same file 
    406            tbias_asm(:,:,:) = 0.0_wp 
    407            sbias_asm(:,:,:) = 0.0_wp 
    408            tbias_asm_out(:,:,:) = 0.0_wp 
    409            sbias_asm_out(:,:,:) = 0.0_wp 
     407            tbias_asm(:,:,:) = 0.0_wp 
     408            sbias_asm(:,:,:) = 0.0_wp 
     409            tbias_asm_out(:,:,:) = 0.0_wp 
     410            sbias_asm_out(:,:,:) = 0.0_wp 
    410411         ENDIF 
    411412 
    412413         IF ( ln_incpc ) THEN   !incr pressure correction 
    413            tbias_asm_stscale(:,:,:) = 0.0_wp 
    414            sbias_asm_stscale(:,:,:) = 0.0_wp 
    415            tbias_asm_stscale_out(:,:,:) = 0.0_wp 
    416            sbias_asm_stscale_out(:,:,:) = 0.0_wp 
     414            tbias_asm_stscale(:,:,:) = 0.0_wp 
     415            sbias_asm_stscale(:,:,:) = 0.0_wp 
     416            tbias_asm_stscale_out(:,:,:) = 0.0_wp 
     417            sbias_asm_stscale_out(:,:,:) = 0.0_wp 
    417418         ENDIF 
    418419 
     
    458459                  ! Get the T and S bias data 
    459460                  CALL iom_get( numbias_tot, jpdom_autoglo, 'tbias_asm_stscale', tbias_asm_stscale ) 
    460                      CALL iom_get( numbias_tot, jpdom_autoglo, 'sbias_asm_stscale', sbias_asm_stscale ) 
     461                  CALL iom_get( numbias_tot, jpdom_autoglo, 'sbias_asm_stscale', sbias_asm_stscale ) 
    461462               ELSE 
    462                      CALL ctl_stop( 'Short time scale bias assim variables not found in ',cn_bias_tot ) 
     463                  CALL ctl_stop( 'Short time scale bias assim variables not found in ',cn_bias_tot ) 
    463464               ENDIF 
    464465            ENDIF   
     
    550551 
    551552      IF ( ln_incpc) THEN 
    552        ! not sure if this should be a special case of nn_lat_ramp == 1 
    553553         minlat_bias = 3.0_wp 
    554554         maxlat_bias = 8.0_wp   
     
    616616            &                              fbcoef(:,:)   * fb_t_asm_max ) 
    617617         zcof2(:,:) = ( 1.0_wp - fbcoef(:,:) ) * fb_p_asm 
    618       
    619618 
    620619         IF ( ln_itdecay ) THEN   !decay the pressure correction at each time step 
    621620     
    622        ztsday  = rday / real(rdt) 
    623  
    624             zdecay = (1-ztscale)**(1/real(ztsday)) ! used in ipc 
    625             zfrac1 = zdecay**real(kt) ! used in ipc 
     621       ztsday  = rday / REAL(rdt,wp) 
     622 
     623            zdecay = (1-ztscale)**(1/REAL(ztsday,wp)) ! used in ipc 
     624            zfrac1 = zdecay**REAL(kt,wp) ! used in ipc 
    626625            IF ( zfrac1 <= 0.0 ) zfrac1 = 0.0_wp 
    627626 
    628  
    629627            IF( ln_asmiau .and. ln_trainc ) THEN  ! nudge in increments and decay historical contribution 
    630  
    631628                
    632629               IF ( kt <= nitiaufin ) THEN  ! During IAU calculate the fraction of increments to apply at each time step 
    633630 
    634                   ztfrac = real(kt) / real(nitiaufin)  ! nudging factor during the IAU 
    635                    
     631                  ztfrac = REAL(kt,wp) / REAL(nitiaufin,wp)  ! nudging factor during the IAU 
    636632           
    637633                  IF (lwp) THEN 
     
    639635                     WRITE(numout,*) '~~~~~~~~~~~~' 
    640636                     WRITE(numout,* ) "proportion of  increment applied in pcbias ",ztfrac 
    641                      WRITE(numout,* ) "proportion of  historical pcbias applied ",t_asm_mem**(real(kt)/ztsday) 
     637                     WRITE(numout,* ) "proportion of  historical pcbias applied ",t_asm_mem**(REAL(kt,wp)/ztsday) 
    642638                  ENDIF 
    643639 
    644  
    645                   DO jk = 1, jpkm1 
    646                       
     640                  DO jk = 1, jpkm1   
    647641                     tbias(:,:,jk) = tbias(:,:,jk) +                            & 
    648                      &                ( t_asm_mem**(real(kt)/ztsday) * tbias_asm(:,:,jk)  +                    & 
     642                     &                ( t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk)  +                    & 
    649643                     &                ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof1(:,:) 
    650644                     sbias(:,:,jk) = sbias(:,:,jk) +                            & 
    651                      &               ( t_asm_mem**(real(kt)/ztsday) * sbias_asm(:,:,jk)  +                     & 
     645                     &               ( t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk)  +                     & 
    652646                     &               ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof1(:,:) 
    653647 
    654                     tbias_p(:,:,jk) = tbias_p(:,:,jk) +                        & 
    655                      &               ( t_asm_mem**(real(kt)/ztsday) * tbias_asm(:,:,jk)  +                     & 
     648                     tbias_p(:,:,jk) = tbias_p(:,:,jk) +                        & 
     649                     &               ( t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk)  +                     & 
    656650                     &               ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof2(:,:) 
    657651                     sbias_p(:,:,jk) = sbias_p(:,:,jk) +                        &   
    658                      &               ( t_asm_mem**(real(kt)/ztsday) * sbias_asm(:,:,jk)  +                     & 
     652                     &               ( t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk)  +                     & 
    659653                     &               ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk) ) * zcof2(:,:) 
    660654                  ENDDO 
    661655 
    662  
    663656                  IF (ln_incpc) THEN 
    664  
    665657 
    666658                     IF (lwp) THEN 
     
    687679 
    688680                        DO jk = 1, jpkm1 
    689                            tbias_asm_out(:,:,jk) =  t_asm_mem**(real(kt)/ztsday) * tbias_asm(:,:,jk)  +       & 
     681                           tbias_asm_out(:,:,jk) =  t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk)  +       & 
    690682                           &                     ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk) 
    691                            sbias_asm_out(:,:,jk) =  t_asm_mem**(real(kt)/ztsday) * sbias_asm(:,:,jk)  +       & 
     683                           sbias_asm_out(:,:,jk) =  t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk)  +       & 
    692684                           &                     ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk) 
    693685                         END DO 
     
    707699                  IF ( kt == nitiaufin ) THEN 
    708700                     DO jk = 1, jpkm1 
    709                         tbias_asm(:,:,jk) =  t_asm_mem**(real(kt)/ztsday) * tbias_asm(:,:,jk)  +       & 
     701                        tbias_asm(:,:,jk) =  t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk)  +       & 
    710702                        &                     ztfrac * t_asm_upd * t_bkginc(:,:,jk) * tmask(:,:,jk) 
    711                         sbias_asm(:,:,jk) =  t_asm_mem**(real(kt)/ztsday) * sbias_asm(:,:,jk)  +       & 
     703                        sbias_asm(:,:,jk) =  t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk)  +       & 
    712704                        &                     ztfrac * t_asm_upd * s_bkginc(:,:,jk) * tmask(:,:,jk) 
    713705                     END DO 
     
    721713      
    722714                  ENDIF 
    723  
    724715                
    725716               ELSE ! decay pressure correction from combined historical component and increments after IAU 
    726717 
    727                   ztfrac=t_asm_mem**(real(kt-nitiaufin)/real(nitiaufin))  ! decay from end of IAU 
    728                    
     718                  ztfrac=t_asm_mem**(REAL(kt-nitiaufin,wp)/REAL(nitiaufin,wp))  ! decay from end of IAU 
    729719                   
    730720                  DO jk = 1, jpkm1 
     
    737727                     sbias_p(:,:,jk) = sbias_p(:,:,jk) +                        &   
    738728                     &               ( ztfrac * sbias_asm(:,:,jk) ) * zcof2(:,:) 
    739  
    740729                  ENDDO 
    741730 
    742731                 IF (ln_incpc) THEN 
    743732 
    744                    zfrac  = zdecay**real(kt-nitiaufin)  
     733                   zfrac  = zdecay**REAL(kt-nitiaufin,wp)  
    745734                   IF ( zfrac <= 0.0 ) zfrac = 0.0_wp 
    746735    
    747  
    748736                   DO jk = 1, jpkm1 
    749737                      tbias_i(:,:,jk) =  tbias_asm_stscale(:,:,jk) * zfrac * (1.0 - fbcoef_stscale(:,:)) 
     
    752740 
    753741                   IF (lwp) THEN 
    754                      WRITE(numout,*) 'tra_bias : bias weights' 
    755                      WRITE(numout,*) '~~~~~~~~~~~~' 
    756                      WRITE(numout,* ) "IPC - proportion of increments + historical pcbias applied ",zfrac 
     742                      WRITE(numout,*) 'tra_bias : bias weights' 
     743                      WRITE(numout,*) '~~~~~~~~~~~~' 
     744                      WRITE(numout,* ) "IPC - proportion of increments + historical pcbias applied ",zfrac 
    757745                   ENDIF 
    758746 
    759747                 ENDIF 
    760748 
    761                   IF (lwp) THEN 
    762                      WRITE(numout,*) 'tra_bias : bias weights' 
    763                      WRITE(numout,*) '~~~~~~~~~~~~' 
    764                      WRITE(numout,* ) "proportion of increments + historical pcbias applied ",ztfrac 
    765                   ENDIF 
    766  
    767                   IF ( ln_trainc .and. .not.ln_bsyncro ) THEN  
    768                      IF ( kt == nn_bias_itwrt ) THEN 
    769                         DO jk = 1, jpkm1 
    770                            tbias_asm_out(:,:,jk) =  ztfrac * tbias_asm(:,:,jk)  
    771                            sbias_asm_out(:,:,jk) =  ztfrac * sbias_asm(:,:,jk)  
    772                         END DO 
    773  
    774                         IF ( ln_incpc ) THEN 
    775                          IF (lwp) WRITE(numout,*) 'after end of IAU - IPC - saving s/tbias_asm_stscale' 
    776                            DO jk = 1, jpkm1 
    777                               tbias_asm_stscale_out(:,:,jk) =  tbias_asm_stscale(:,:,jk) * zfrac 
    778                               sbias_asm_stscale_out(:,:,jk) =  sbias_asm_stscale(:,:,jk) * zfrac 
    779                            ENDDO 
    780                         ENDIF 
    781  
    782                      ENDIF 
    783  
    784                   ENDIF 
    785  
    786                ENDIF 
    787  
     749                 IF (lwp) THEN 
     750                    WRITE(numout,*) 'tra_bias : bias weights' 
     751                    WRITE(numout,*) '~~~~~~~~~~~~' 
     752                    WRITE(numout,* ) "proportion of increments + historical pcbias applied ",ztfrac 
     753                 ENDIF 
     754 
     755                 IF ( ln_trainc .and. .not.ln_bsyncro ) THEN  
     756                    IF ( kt == nn_bias_itwrt ) THEN 
     757                       DO jk = 1, jpkm1 
     758                          tbias_asm_out(:,:,jk) =  ztfrac * tbias_asm(:,:,jk)  
     759                          sbias_asm_out(:,:,jk) =  ztfrac * sbias_asm(:,:,jk)  
     760                       END DO 
     761 
     762                       IF ( ln_incpc ) THEN 
     763                          IF (lwp) WRITE(numout,*) 'after end of IAU - IPC - saving s/tbias_asm_stscale' 
     764                             DO jk = 1, jpkm1 
     765                                tbias_asm_stscale_out(:,:,jk) =  tbias_asm_stscale(:,:,jk) * zfrac 
     766                                sbias_asm_stscale_out(:,:,jk) =  sbias_asm_stscale(:,:,jk) * zfrac 
     767                             ENDDO 
     768                          ENDIF 
     769                       ENDIF 
     770                    ENDIF 
     771                 ENDIF 
    788772 
    789773            ELSE ! no assimilation increments, simply decay pressure correction (e.g for forecasts) 
     
    792776               DO jk = 1, jpkm1 
    793777                  tbias(:,:,jk) = tbias(:,:,jk) +                                                         & 
    794                   &               ( t_asm_mem**(real(kt)/ztsday) * tbias_asm(:,:,jk) ) * zcof1(:,:) 
     778                  &               ( t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk) ) * zcof1(:,:) 
    795779                  sbias(:,:,jk) = sbias(:,:,jk) +                                                         & 
    796                   &               ( t_asm_mem**(real(kt)/ztsday) * sbias_asm(:,:,jk) ) * zcof1(:,:) 
     780                  &               ( t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk) ) * zcof1(:,:) 
    797781                  tbias_p(:,:,jk) = tbias_p(:,:,jk) +                                                     & 
    798                   &               ( t_asm_mem**(real(kt)/ztsday) * tbias_asm(:,:,jk) ) * zcof2(:,:) 
     782                  &               ( t_asm_mem**(REAL(kt,wp)/ztsday) * tbias_asm(:,:,jk) ) * zcof2(:,:) 
    799783                  sbias_p(:,:,jk) = sbias_p(:,:,jk) +                                                     &   
    800                   &               ( t_asm_mem**(real(kt)/ztsday) * sbias_asm(:,:,jk) ) * zcof2(:,:) 
     784                  &               ( t_asm_mem**(REAL(kt,wp)/ztsday) * sbias_asm(:,:,jk) ) * zcof2(:,:) 
    801785               ENDDO 
    802786 
     
    804788                  WRITE(numout,*) 'tra_bias : bias weights' 
    805789                  WRITE(numout,*) '~~~~~~~~~~~~' 
    806                   WRITE(numout,* ) "proportion of historical pcbias applied ",t_asm_mem**(real(kt)/ztsday) 
     790                  WRITE(numout,* ) "proportion of historical pcbias applied ",t_asm_mem**(REAL(kt,wp)/ztsday) 
    807791               ENDIF 
    808  
    809  
    810792 
    811793               IF (ln_incpc) THEN 
     
    957939            tsw(:,:,:,jp_sal) = tsa(:,:,:,jp_sal) - sbias_p(:,:,:) - sbias_i(:,:,:) 
    958940            IF ( ln_incpc_only ) THEN 
    959                IF(lwp) WRITE(numout,*) 'SUM(tbias_i) is', SUM(tbias_i), 'before to be added to tsw' 
    960941               IF(lwp) WRITE(numout,*) 'ln_incpc_only =', ln_incpc_only, 'tsw updated with IPC only and ln_dynhpg_imp = ',ln_dynhpg_imp 
    961942               tsw(:,:,:,jp_tem) = tsa(:,:,:,jp_tem) - tbias_i(:,:,:) 
     
    970951            tsw(:,:,:,jp_sal) = tsb(:,:,:,jp_sal) - sbias_p(:,:,:) - sbias_i(:,:,:) 
    971952            IF ( ln_incpc_only ) THEN 
    972                IF(lwp) WRITE(numout,*) 'SUM(tbias_i) is', SUM(tbias_i), 'before to be added to tsw' 
    973953               IF(lwp) WRITE(numout,*) 'ln_incpc_only =', ln_incpc_only, 'tsw updated with IPC only and ln_dynhpg_imp = ',ln_dynhpg_imp 
    974954               tsw(:,:,:,jp_tem) = tsb(:,:,:,jp_tem) - tbias_i(:,:,:) 
     
    978958      ENDIF 
    979959 
    980       
    981  
    982960      IF(lwp) WRITE(numout,*) 'dyn_bias( kt ) calculating rhd_pc, kt =', kt 
    983       !! GO5: CALL eos( tsw, rhd_pc, rhop ) 
    984               CALL eos( tsw, rhd_pc, rhop, fsdept_n(:,:,:) ) 
     961      CALL eos( tsw, rhd_pc, rhop, fsdept_n(:,:,:) ) 
    985962       
    986       ! is this needed? 
    987       IF(lwp) WRITE(numout,*) 'CALL lbc_lnk( rhd_pc, T, 1.0_wp )' 
    988963      CALL lbc_lnk( rhd_pc, 'T', 1.0_wp ) 
    989964 
    990965      ! Partial steps: now horizontal gradient of t,s,rd  
    991966      ! at the bottom ocean level 
    992       IF(lwp) WRITE(numout,*) 'horizontal gradient of t,s,rd ' 
     967 
    993968      IF( ln_zps    )  THEN 
    994969         CALL zps_hde( kt, jpts, tsw, gtsu, gtsv,  & 
  • branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO/OPA_SRC/DYN/dynhpg.F90

    r7658 r8435  
    9999      ! 
    100100      IF ( ln_bias .AND. ln_bias_pc_app ) THEN 
    101          IF(lwp) THEN  
    102          WRITE(numout,*) " ! store orig density, use pressure corrected density using gru_pc and grv_pc" 
    103          ENDIF 
    104101         z_rhd_st(:,:,:) = rhd(:,:,:)     ! store orig density  
    105102         rhd(:,:,:)      = rhd_pc(:,:,:)  ! use pressure corrected density 
     
    110107      ENDIF 
    111108 
    112          IF(lwp) WRITE(numout,*) " computing pressure gradient with nhpg = ",nhpg 
    113109      SELECT CASE ( nhpg )      ! Hydrostatic pressure gradient computation 
    114110      CASE (  0 )   ;   CALL hpg_zco    ( kt )      ! z-coordinate 
  • branches/UKMO/dev_r5518_pcbias_ipc/NEMOGCM/NEMO/OPA_SRC/step.F90

    r7665 r8435  
    273273                             CALL tra_nxt( kstp )                ! tracer fields at next time step 
    274274            IF( ln_sto_eos ) CALL sto_pts( tsn )                 ! Random T/S fluctuations 
    275                              CALL eos    ( tsa, rhd, rhop, fsdept_n(:,:,:) )  ! Time-filtered in situ density for hpg computation 
    276                              CALL lbc_lnk( rhd, 'T', 1.0_wp ) 
     275            CALL eos    ( tsa, rhd, rhop, fsdept_n(:,:,:) )  ! Time-filtered in situ density for hpg computation 
     276            CALL lbc_lnk( rhd, 'T', 1.0_wp ) 
    277277            IF( ln_zps .AND. .NOT. ln_isfcav)                                & 
    278278               &             CALL zps_hde    ( kstp, jpts, tsa, gtsu, gtsv,  &    ! Partial steps: before horizontal gradient 
     
    283283               &                                    gtui, gtvi, grui, grvi, arui, arvi, gzui, gzvi, ge3rui, ge3rvi    ) ! of t, s, rd at the last ocean level 
    284284            IF( ln_bias )    CALL dyn_bias( kstp ) 
    285       ELSE                                                  ! centered hpg  (eos then time stepping) 
    286          IF ( .NOT. lk_dynspg_ts ) THEN                     ! eos already called in time-split case 
    287             IF( ln_sto_eos ) CALL sto_pts( tsn )    ! Random T/S fluctuations 
    288                              CALL eos    ( tsn, rhd, rhop, fsdept_n(:,:,:) )  ! now in situ density for hpg computation 
    289                              CALL lbc_lnk( rhd, 'T', 1.0_wp ) 
     285      ELSE                                                   ! centered hpg  (eos then time stepping) 
     286         IF ( .NOT. lk_dynspg_ts ) THEN                      ! eos already called in time-split case 
     287            IF( ln_sto_eos ) CALL sto_pts( tsn )             ! Random T/S fluctuations 
     288            CALL eos    ( tsn, rhd, rhop, fsdept_n(:,:,:) )  ! now in situ density for hpg computation 
     289            CALL lbc_lnk( rhd, 'T', 1.0_wp ) 
    290290         IF( ln_zps .AND. .NOT. ln_isfcav)                                   & 
    291291               &             CALL zps_hde    ( kstp, jpts, tsn, gtsu, gtsv,  &    ! Partial steps: before horizontal gradient 
Note: See TracChangeset for help on using the changeset viewer.