MODULE trcrst !!====================================================================== !! *** MODULE trcrst *** !! TOP : create, write, read the restart files for passive tracers !!====================================================================== !! History : 1.0 ! 2007-02 (C. Ethe) adaptation from the ocean !!---------------------------------------------------------------------- #if defined key_top !!---------------------------------------------------------------------- !! 'key_top' TOP models !!---------------------------------------------------------------------- !! trc_rst_opn : open restart file !! trc_rst_read : read restart file !! trc_rst_wri : write restart file !!---------------------------------------------------------------------- USE oce_trc USE trc USE sms_lobster ! LOBSTER variables USE sms_pisces ! PISCES variables USE trcsms_cfc ! CFC variables USE trcsms_c14b ! C14 variables USE trcsms_my_trc ! MY_TRC variables USE trctrp_lec USE lib_mpp USE iom IMPLICIT NONE PRIVATE PUBLIC trc_rst_opn ! called by ??? PUBLIC trc_rst_read ! called by ??? PUBLIC trc_rst_wri ! called by ??? INTEGER, PUBLIC :: numrtr, numrtw !: logical unit for trc restart (read and write) !! * Substitutions # include "top_substitute.h90" !!---------------------------------------------------------------------- !! NEMO/TOP 1.0 , LOCEAN-IPSL (2005) !! $Id$ !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) !!---------------------------------------------------------------------- CONTAINS SUBROUTINE trc_rst_opn( kt ) !!---------------------------------------------------------------------- !! *** trc_rst_opn *** !! !! ** purpose : output of sea-trc variable in a netcdf file !!---------------------------------------------------------------------- INTEGER, INTENT(in) :: kt ! number of iteration ! CHARACTER(LEN=20) :: clkt ! ocean time-step define as a character CHARACTER(LEN=50) :: clname ! trc output restart file name !!---------------------------------------------------------------------- ! IF( kt == nit000 ) lrst_trc = .FALSE. # if defined key_off_tra IF( MOD( kt - 1, nstock ) == 0 ) THEN ! we use kt - 1 and not kt - nit000 to keep the same periodicity from the beginning of the experiment nitrst = kt + nstock - 1 ! define the next value of nitrst for restart writing IF( nitrst > nitend ) nitrst = nitend ! make sure we write a restart at the end of the run ENDIF #endif IF( kt == nitrst - ndttrc .OR. nitend - nit000 + 1 < 2 * ndttrc ) THEN ! beware if model runs less than 2*ndttrc time step ! beware of the format used to write kt (default is i8.8, that should be large enough) IF( nitrst > 1.0e9 ) THEN ; WRITE(clkt,* ) nitrst ELSE ; WRITE(clkt,'(i8.8)') nitrst ENDIF ! create the file IF(lwp) WRITE(numout,*) clname = TRIM(cexper)//"_"//TRIM(ADJUSTL(clkt))//"_"//TRIM(cn_trcrst_out) IF(lwp) WRITE(numout,*) ' open trc restart.output NetCDF file: '//clname CALL iom_open( clname, numrtw, ldwrt = .TRUE., kiolib = jprstlib ) lrst_trc = .TRUE. ENDIF ! END SUBROUTINE trc_rst_opn SUBROUTINE trc_rst_read !!---------------------------------------------------------------------- !! *** trc_rst_opn *** !! !! ** purpose : read passive tracer fields in restart files !!---------------------------------------------------------------------- INTEGER :: jn INTEGER :: iarak0 REAL(wp) :: zkt, zarak0 !!---------------------------------------------------------------------- IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) 'trc_rst_read : read the TOP restart file' IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd ) THEN ; iarak0 = 1 ELSE ; iarak0 = 0 ENDIF IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) ' the present run starts at the time step nit000 = ', nit000 IF(lwp .AND. iarak0 == 1 ) WRITE(numout,*) ' and needs previous fields for Arakawa sheme ??? ' ! Time domain : restart ! ------------------------- IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) ' *** passive tracer restart option' SELECT CASE ( nrsttr ) CASE ( 0 ) IF(lwp) WRITE(numout,*) ' nrsttr = 0 no control of nit000' CASE ( 1 ) IF(lwp) WRITE(numout,*) ' nrsttr = 1 we control the date of nit000' CASE ( 2 ) IF(lwp) WRITE(numout,*) ' nrsttr = 2 the date adatrj is read in restart file' CASE DEFAULT IF(lwp) WRITE(numout,*) ' ===>>>> nrsttr not equal 0, 1 or 2 : no control of the date' IF(lwp) WRITE(numout,*) ' ======= =========' END SELECT CALL iom_open( cn_trcrst_in, numrtr, kiolib = jprstlib ) CALL iom_get( numrtr, 'kt' , zkt ) CALL iom_get( numrtr, 'arak0', zarak0 ) IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) ' Info on the restart file read : ' IF(lwp) WRITE(numout,*) ' time-step : ', NINT( zkt ) IF(lwp) WRITE(numout,*) ' arakawa option : ', NINT( zarak0 ) IF( nittrc000 - NINT( zkt ) /= 1 .AND. nrsttr /= 0 ) & ! control of date & CALL ctl_stop( ' ===>>>> : problem with nit000 for the restart', & & ' verify the restart file or rerun with nrstdt = 0 (namelist)' ) IF( iarak0 /= NINT( zarak0 ) ) & ! Control of the scheme & CALL ctl_stop( ' ===>>>> : problem with advection scheme', & & ' it must be the same type for both restart and previous run', & & ' centered or euler ' ) ! READ prognostic variables and computes diagnostic variable DO jn = 1, jptra CALL iom_get( numrtr, jpdom_autoglo, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) ) END DO DO jn = 1, jptra CALL iom_get( numrtr, jpdom_autoglo, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) ) END DO #if defined key_lobster CALL iom_get( numrtr, jpdom_autoglo, 'SEDB'//ctrcnm(jpdet), sedpocb(:,:) ) CALL iom_get( numrtr, jpdom_autoglo, 'SEDN'//ctrcnm(jpdet), sedpocn(:,:) ) #endif #if defined key_pisces CALL iom_get( numrtr, jpdom_autoglo, 'Silicalim', xksi(:,:) ) CALL iom_get( numrtr, jpdom_autoglo, 'Silicamax', xksimax(:,:) ) CALL trc_rst_ini ! Initialisation of some variables #endif #if defined key_cfc DO jn = jp_cfc0, jp_cfc1 CALL iom_get( numrtr, jpdom_autoglo, 'qint_'//ctrcnm(jn), qint_cfc(:,:,jn) ) END DO #endif #if defined key_c14b CALL iom_get( numrtr, jpdom_autoglo, 'qint_'//ctrcnm(jn) , qint_c14(:,:) ) #endif #if defined key_my_trc #endif CALL iom_close( numrtr ) ! END SUBROUTINE trc_rst_read SUBROUTINE trc_rst_wri( kt ) !!---------------------------------------------------------------------- !! *** trc_rst_wri *** !! !! ** purpose : write passive tracer fields in restart files !!---------------------------------------------------------------------- INTEGER, INTENT( in ) :: kt ! ocean time-step index !! INTEGER :: ji, jj, jk, jn !!---------------------------------------------------------------------- IF( kt == nitrst ) THEN IF(lwp) WRITE(numout,*) IF(lwp) WRITE(numout,*) 'trc_wri : write the TOP restart file (NetCDF) at it= ', kt, ' date= ', ndastp IF(lwp) WRITE(numout,*) '~~~~~~~' ENDIF CALL iom_rstput( kt, nitrst, numrtw, 'kt' , REAL( kt, wp ) ) IF( ln_trcadv_cen2 .OR. ln_trcadv_tvd ) THEN CALL iom_rstput( kt, nitrst, numrtw, 'arak0', 1. ) ELSE CALL iom_rstput( kt, nitrst, numrtw, 'arak0', 0. ) ENDIF ! prognostic variables ! -------------------- DO jn = 1, jptra CALL iom_rstput( kt, nitrst, numrtw, 'TRN'//ctrcnm(jn), trn(:,:,:,jn) ) END DO DO jn = 1, jptra CALL iom_rstput( kt, nitrst, numrtw, 'TRB'//ctrcnm(jn), trb(:,:,:,jn) ) END DO #if defined key_lobster CALL iom_rstput( kt, nitrst, numrtw, 'SEDB'//ctrcnm(jpdet), sedpocb(:,:) ) CALL iom_rstput( kt, nitrst, numrtw, 'SEDN'//ctrcnm(jpdet), sedpocn(:,:) ) #endif #if defined key_pisces CALL iom_rstput( kt, nitrst, numrtw, 'Silicalim', xksi(:,:) ) CALL iom_rstput( kt, nitrst, numrtw, 'Silicamax', xksimax(:,:) ) #endif #if defined key_cfc DO jn = jp_cfc0, jp_cfc1 CALL iom_rstput( kt, nitrst, numrtw, 'qint_'//ctrcnm(jn), qint_cfc(:,:,jn) ) END DO #endif #if defined key_c14b CALL iom_rstput( kt, nitrst, numrtw, 'qint_'//ctrcnm(jn), qint_c14(:,:) ) #endif #if defined key_my_trc #endif IF( kt == nitrst ) THEN CALL trc_rst_stat ! statistics CALL iom_close( numrtw ) ! close the restart file (only at last time step) #if ! defined key_trdmld_trc lrst_trc = .FALSE. #endif ENDIF ! END SUBROUTINE trc_rst_wri # if defined key_pisces SUBROUTINE trc_rst_ini !!---------------------------------------------------------------------- !! *** trc_rst_ini *** !! !! ** purpose : Initialisation of some variables ( hi !!---------------------------------------------------------------------- INTEGER :: ji, jj, jk, jn REAL(wp) :: & alkmean = 2426. , & ! mean value of alkalinity ( Glodap ; for Goyet 2391. ) po4mean = 2.165 , & ! mean value of phosphates no3mean = 30.90 , & ! mean value of nitrate siomean = 91.51 ! mean value of silicate REAL(wp) :: zvol, ztrasum REAL(wp) :: caralk, bicarb, co3 IF(lwp) WRITE(numout,*) IF( cp_cfg == "orca" .AND. .NOT. lk_trc_c1d ) THEN ! ORCA condiguration (not 1D) ! ! ! --------------------------- ! ! set total alkalinity, phosphate, NO3 & silicate ! total alkalinity ztrasum = 0.e0 DO jk = 1, jpk DO jj = 1, jpj DO ji = 1, jpi zvol = cvol(ji,jj,jk) # if defined key_off_degrad zvol = zvol * facvol(ji,jj,jk) # endif ztrasum = ztrasum + trn(ji,jj,jk,jptal) * zvol END DO END DO END DO IF( lk_mpp ) CALL mpp_sum( ztrasum ) ! sum over the global domain ztrasum = ztrasum / areatot * 1.e6 IF(lwp) WRITE(numout,*) ' TALK mean : ', ztrasum trn(:,:,:,jptal) = trn(:,:,:,jptal) * alkmean / ztrasum ! phosphate ztrasum = 0.e0 DO jk = 1, jpk DO jj = 1, jpj DO ji = 1, jpi zvol = cvol(ji,jj,jk) # if defined key_off_degrad zvol = zvol * facvol(ji,jj,jk) # endif ztrasum = ztrasum + trn(ji,jj,jk,jppo4) * zvol END DO END DO END DO IF( lk_mpp ) CALL mpp_sum( ztrasum ) ! sum over the global domain ztrasum = ztrasum / areatot * 1.e6 / 122. IF(lwp) WRITE(numout,*) ' PO4 mean : ', ztrasum trn(:,:,:,jppo4) = trn(:,:,:,jppo4) * po4mean / ztrasum ! Nitrates ztrasum = 0.e0 DO jk = 1, jpk DO jj = 1, jpj DO ji = 1, jpi zvol = cvol(ji,jj,jk) # if defined key_off_degrad zvol = zvol * facvol(ji,jj,jk) # endif ztrasum = ztrasum + trn(ji,jj,jk,jpno3) * zvol END DO END DO END DO IF( lk_mpp ) CALL mpp_sum( ztrasum ) ! sum over the global domain ztrasum = ztrasum / areatot * 1.e6 / 7.6 IF(lwp) WRITE(numout,*) ' NO3 mean : ', ztrasum trn(:,:,:,jpno3) = trn(:,:,:,jpno3) * no3mean / ztrasum ! Silicate ztrasum = 0.e0 DO jk = 1, jpk DO jj = 1, jpj DO ji = 1, jpi zvol = cvol(ji,jj,jk) # if defined key_off_degrad zvol = zvol * facvol(ji,jj,jk) # endif ztrasum = ztrasum + trn(ji,jj,jk,jpsil) * zvol END DO END DO END DO IF( lk_mpp ) CALL mpp_sum( ztrasum ) ! sum over the global domain ztrasum = ztrasum / areatot * 1.e6 IF(lwp) WRITE(numout,*) ' SiO3 mean : ', ztrasum trn(:,:,:,jpsil) = MIN( 400.e-6,trn(:,:,:,jpsil) * siomean / ztrasum ) ! ENDIF !#if defined key_kriest ! !! Initialize number of particles from a standart restart file ! !! The name of big organic particles jpgoc has been only change ! !! and replace by jpnum but the values here are concentration ! trn(:,:,:,jppoc) = trn(:,:,:,jppoc) + trn(:,:,:,jpnum) ! trn(:,:,:,jpnum) = trn(:,:,:,jppoc) / ( 6. * xkr_massp ) !#endif !! Set hi (???) from total alcalinity, borat (???), akb3 (???) and ak23 (???) !! --------------------------------------------------------------------- DO jk = 1, jpk DO jj = 1, jpj DO ji = 1,jpi caralk = trn(ji,jj,jk,jptal) - borat(ji,jj,jk) / ( 1. + 1.e-8 / ( rtrn + akb3(ji,jj,jk) ) ) co3 = ( caralk - trn(ji,jj,jk,jpdic) ) * tmask(ji,jj,jk) & & + 0.5e-3 * ( 1.- tmask(ji,jj,jk) ) bicarb = 2.* trn(ji,jj,jk,jpdic) - caralk hi(ji,jj,jk) = ( ak23(ji,jj,jk) * bicarb / co3 ) * tmask(ji,jj,jk) & & + 1.0e-9 * ( 1.- tmask(ji,jj,jk) ) END DO END DO END DO END SUBROUTINE trc_rst_ini #endif !!---------------------------------------------------------------------- SUBROUTINE trc_rst_stat !!---------------------------------------------------------------------- !! *** trc_rst_stat *** !! !! ** purpose : Compute tracers statistics !!---------------------------------------------------------------------- INTEGER :: ji, jj, jk, jn REAL(wp) :: zdiag_var, zdiag_varmin, zdiag_varmax, zdiag_tot REAL(wp) :: zder, zvol !!---------------------------------------------------------------------- IF( lwp ) THEN WRITE(numout,*) WRITE(numout,*) ' ----TRACER STAT---- ' WRITE(numout,*) ENDIF zdiag_tot = 0.e0 DO jn = 1, jptra zdiag_var = 0.e0 zdiag_varmin = 0.e0 zdiag_varmax = 0.e0 DO ji = 1, jpi DO jj = 1, jpj DO jk = 1,jpk zvol = cvol(ji,jj,jk) # if defined key_off_degrad zvol = zvol * facvol(ji,jj,jk) # endif zdiag_var = zdiag_var + trn(ji,jj,jk,jn) * zvol END DO END DO END DO zdiag_varmin = MINVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) ) zdiag_varmax = MAXVAL( trn(:,:,:,jn), mask= ((tmask*SPREAD(tmask_i,DIM=3,NCOPIES=jpk).NE.0.)) ) IF( lk_mpp ) THEN CALL mpp_min( zdiag_varmin ) ! min over the global domain CALL mpp_max( zdiag_varmax ) ! max over the global domain CALL mpp_sum( zdiag_var ) ! sum over the global domain END IF zdiag_tot = zdiag_tot + zdiag_var zdiag_var = zdiag_var / areatot IF(lwp) WRITE(numout,*) ' MEAN NO ', jn, ctrcnm(jn), ' = ', zdiag_var, & & ' MIN = ', zdiag_varmin, ' MAX = ', zdiag_varmax END DO zder = ( ( zdiag_tot - trai ) / ( trai + 1.e-12 ) ) * 100._wp IF(lwp) WRITE(numout,*) ' Integral of all tracers over the full domain = ', zdiag_tot IF(lwp) WRITE(numout,*) ' Drift of the sum of all tracers =', zder, ' %' END SUBROUTINE trc_rst_stat #else !!---------------------------------------------------------------------- !! Dummy module : No passive tracer !!---------------------------------------------------------------------- CONTAINS SUBROUTINE trc_rst_read ! Empty routines END SUBROUTINE trc_rst_read SUBROUTINE trc_rst_wri( kt ) INTEGER, INTENT ( in ) :: kt WRITE(*,*) 'trc_rst_wri: You should not have seen this print! error?', kt END SUBROUTINE trc_rst_wri #endif !!====================================================================== END MODULE trcrst