[4715] | 1 | MODULE diacfl |
---|
[9019] | 2 | !!====================================================================== |
---|
[4715] | 3 | !! *** MODULE diacfl *** |
---|
| 4 | !! Output CFL diagnostics to ascii file |
---|
[9019] | 5 | !!====================================================================== |
---|
| 6 | !! History : 3.4 ! 2010-03 (E. Blockley) Original code |
---|
| 7 | !! 3.6 ! 2014-06 (T. Graham) Removed CPP key & Updated to vn3.6 |
---|
| 8 | !! 4.0 ! 2017-09 (G. Madec) style + comments |
---|
[4715] | 9 | !!---------------------------------------------------------------------- |
---|
| 10 | !! dia_cfl : Compute and output Courant numbers at each timestep |
---|
| 11 | !!---------------------------------------------------------------------- |
---|
| 12 | USE oce ! ocean dynamics and active tracers |
---|
| 13 | USE dom_oce ! ocean space and time domain |
---|
[9019] | 14 | USE domvvl ! |
---|
| 15 | ! |
---|
[4715] | 16 | USE lib_mpp ! distribued memory computing |
---|
| 17 | USE lbclnk ! ocean lateral boundary condition (or mpp link) |
---|
| 18 | USE in_out_manager ! I/O manager |
---|
[10824] | 19 | USE iom ! |
---|
[4715] | 20 | USE timing ! Performance output |
---|
| 21 | |
---|
| 22 | IMPLICIT NONE |
---|
| 23 | PRIVATE |
---|
| 24 | |
---|
[9019] | 25 | CHARACTER(LEN=50) :: clname="cfl_diagnostics.ascii" ! ascii filename |
---|
| 26 | INTEGER :: numcfl ! outfile unit |
---|
| 27 | ! |
---|
| 28 | INTEGER, DIMENSION(3) :: nCu_loc, nCv_loc, nCw_loc ! U, V, and W run max locations in the global domain |
---|
| 29 | REAL(wp) :: rCu_max, rCv_max, rCw_max ! associated run max Courant number |
---|
[4715] | 30 | |
---|
| 31 | PUBLIC dia_cfl ! routine called by step.F90 |
---|
| 32 | PUBLIC dia_cfl_init ! routine called by nemogcm |
---|
| 33 | |
---|
| 34 | !! * Substitutions |
---|
| 35 | # include "vectopt_loop_substitute.h90" |
---|
| 36 | !!---------------------------------------------------------------------- |
---|
[10068] | 37 | !! NEMO/OCE 4.0 , NEMO Consortium (2018) |
---|
[4715] | 38 | !! $Id$ |
---|
[10068] | 39 | !! Software governed by the CeCILL license (see ./LICENSE) |
---|
[4715] | 40 | !!---------------------------------------------------------------------- |
---|
| 41 | CONTAINS |
---|
| 42 | |
---|
| 43 | SUBROUTINE dia_cfl ( kt ) |
---|
| 44 | !!---------------------------------------------------------------------- |
---|
| 45 | !! *** ROUTINE dia_cfl *** |
---|
| 46 | !! |
---|
| 47 | !! ** Purpose : Compute the Courant numbers Cu=u*dt/dx and Cv=v*dt/dy |
---|
| 48 | !! and output to ascii file 'cfl_diagnostics.ascii' |
---|
| 49 | !!---------------------------------------------------------------------- |
---|
[9019] | 50 | INTEGER, INTENT(in) :: kt ! ocean time-step index |
---|
| 51 | ! |
---|
[11532] | 52 | INTEGER :: ji, jj, jk ! dummy loop indices |
---|
| 53 | REAL(wp) :: z2dt, zCu_max, zCv_max, zCw_max ! local scalars |
---|
| 54 | INTEGER , DIMENSION(3) :: iloc_u , iloc_v , iloc_w , iloc ! workspace |
---|
| 55 | REAL(wp), DIMENSION(jpi,jpj,jpk) :: zCu_cfl, zCv_cfl, zCw_cfl ! workspace |
---|
[9019] | 56 | !!---------------------------------------------------------------------- |
---|
| 57 | ! |
---|
[9124] | 58 | IF( ln_timing ) CALL timing_start('dia_cfl') |
---|
[9019] | 59 | ! |
---|
| 60 | ! ! setup timestep multiplier to account for initial Eulerian timestep |
---|
| 61 | IF( neuler == 0 .AND. kt == nit000 ) THEN ; z2dt = rdt |
---|
| 62 | ELSE ; z2dt = rdt * 2._wp |
---|
| 63 | ENDIF |
---|
| 64 | ! |
---|
| 65 | ! |
---|
| 66 | DO jk = 1, jpk ! calculate Courant numbers |
---|
| 67 | DO jj = 1, jpj |
---|
[11532] | 68 | DO ji = 1, jpi |
---|
[9019] | 69 | zCu_cfl(ji,jj,jk) = ABS( un(ji,jj,jk) ) * z2dt / e1u (ji,jj) ! for i-direction |
---|
| 70 | zCv_cfl(ji,jj,jk) = ABS( vn(ji,jj,jk) ) * z2dt / e2v (ji,jj) ! for j-direction |
---|
| 71 | zCw_cfl(ji,jj,jk) = ABS( wn(ji,jj,jk) ) * z2dt / e3w_n(ji,jj,jk) ! for k-direction |
---|
| 72 | END DO |
---|
| 73 | END DO |
---|
| 74 | END DO |
---|
| 75 | ! |
---|
[10824] | 76 | ! write outputs |
---|
| 77 | IF( iom_use('cfl_cu') ) CALL iom_put( 'cfl_cu', MAXVAL( zCu_cfl, dim=3 ) ) |
---|
| 78 | IF( iom_use('cfl_cv') ) CALL iom_put( 'cfl_cv', MAXVAL( zCv_cfl, dim=3 ) ) |
---|
| 79 | IF( iom_use('cfl_cw') ) CALL iom_put( 'cfl_cw', MAXVAL( zCw_cfl, dim=3 ) ) |
---|
| 80 | |
---|
[9019] | 81 | ! ! calculate maximum values and locations |
---|
| 82 | IF( lk_mpp ) THEN |
---|
[10425] | 83 | CALL mpp_maxloc( 'diacfl', zCu_cfl, umask, zCu_max, iloc_u ) |
---|
| 84 | CALL mpp_maxloc( 'diacfl', zCv_cfl, vmask, zCv_max, iloc_v ) |
---|
| 85 | CALL mpp_maxloc( 'diacfl', zCw_cfl, wmask, zCw_max, iloc_w ) |
---|
[9019] | 86 | ELSE |
---|
| 87 | iloc = MAXLOC( ABS( zcu_cfl(:,:,:) ) ) |
---|
| 88 | iloc_u(1) = iloc(1) + nimpp - 1 |
---|
| 89 | iloc_u(2) = iloc(2) + njmpp - 1 |
---|
| 90 | iloc_u(3) = iloc(3) |
---|
| 91 | zCu_max = zCu_cfl(iloc(1),iloc(2),iloc(3)) |
---|
| 92 | ! |
---|
| 93 | iloc = MAXLOC( ABS( zcv_cfl(:,:,:) ) ) |
---|
| 94 | iloc_v(1) = iloc(1) + nimpp - 1 |
---|
| 95 | iloc_v(2) = iloc(2) + njmpp - 1 |
---|
| 96 | iloc_v(3) = iloc(3) |
---|
| 97 | zCv_max = zCv_cfl(iloc(1),iloc(2),iloc(3)) |
---|
| 98 | ! |
---|
| 99 | iloc = MAXLOC( ABS( zcw_cfl(:,:,:) ) ) |
---|
| 100 | iloc_w(1) = iloc(1) + nimpp - 1 |
---|
| 101 | iloc_w(2) = iloc(2) + njmpp - 1 |
---|
| 102 | iloc_w(3) = iloc(3) |
---|
| 103 | zCw_max = zCw_cfl(iloc(1),iloc(2),iloc(3)) |
---|
| 104 | ENDIF |
---|
| 105 | ! |
---|
| 106 | ! ! write out to file |
---|
| 107 | IF( lwp ) THEN |
---|
[11532] | 108 | WRITE(numcfl,FMT='(2x,i6,3x,a6,4x,f7.4,1x,i4,1x,i4,1x,i4)') kt, 'Max Cu', zCu_max, iloc_u(1), iloc_u(2), iloc_u(3) |
---|
[9019] | 109 | WRITE(numcfl,FMT='(11x, a6,4x,f7.4,1x,i4,1x,i4,1x,i4)') 'Max Cv', zCv_max, iloc_v(1), iloc_v(2), iloc_v(3) |
---|
| 110 | WRITE(numcfl,FMT='(11x, a6,4x,f7.4,1x,i4,1x,i4,1x,i4)') 'Max Cw', zCw_max, iloc_w(1), iloc_w(2), iloc_w(3) |
---|
| 111 | ENDIF |
---|
| 112 | ! |
---|
| 113 | ! ! update maximum Courant numbers from whole run if applicable |
---|
| 114 | IF( zCu_max > rCu_max ) THEN ; rCu_max = zCu_max ; nCu_loc(:) = iloc_u(:) ; ENDIF |
---|
| 115 | IF( zCv_max > rCv_max ) THEN ; rCv_max = zCv_max ; nCv_loc(:) = iloc_v(:) ; ENDIF |
---|
| 116 | IF( zCw_max > rCw_max ) THEN ; rCw_max = zCw_max ; nCw_loc(:) = iloc_w(:) ; ENDIF |
---|
[4715] | 117 | |
---|
[9019] | 118 | ! ! at end of run output max Cu and Cv and close ascii file |
---|
| 119 | IF( kt == nitend .AND. lwp ) THEN |
---|
| 120 | ! to ascii file |
---|
| 121 | WRITE(numcfl,*) '******************************************' |
---|
| 122 | WRITE(numcfl,FMT='(3x,a12,6x,f7.4,1x,i4,1x,i4,1x,i4)') 'Run Max Cu', rCu_max, nCu_loc(1), nCu_loc(2), nCu_loc(3) |
---|
| 123 | WRITE(numcfl,FMT='(3x,a8,11x,f15.1)') ' => dt/C', z2dt/rCu_max |
---|
| 124 | WRITE(numcfl,*) '******************************************' |
---|
| 125 | WRITE(numcfl,FMT='(3x,a12,6x,f7.4,1x,i4,1x,i4,1x,i4)') 'Run Max Cv', rCv_max, nCv_loc(1), nCv_loc(2), nCv_loc(3) |
---|
| 126 | WRITE(numcfl,FMT='(3x,a8,11x,f15.1)') ' => dt/C', z2dt/rCv_max |
---|
| 127 | WRITE(numcfl,*) '******************************************' |
---|
| 128 | WRITE(numcfl,FMT='(3x,a12,6x,f7.4,1x,i4,1x,i4,1x,i4)') 'Run Max Cw', rCw_max, nCw_loc(1), nCw_loc(2), nCw_loc(3) |
---|
| 129 | WRITE(numcfl,FMT='(3x,a8,11x,f15.1)') ' => dt/C', z2dt/rCw_max |
---|
| 130 | CLOSE( numcfl ) |
---|
| 131 | ! |
---|
| 132 | ! to ocean output |
---|
| 133 | WRITE(numout,*) |
---|
| 134 | WRITE(numout,*) 'dia_cfl : Maximum Courant number information for the run ' |
---|
| 135 | WRITE(numout,*) '~~~~~~~' |
---|
| 136 | WRITE(numout,*) ' Max Cu = ', rCu_max, ' at (i,j,k) = (',nCu_loc(1),nCu_loc(2),nCu_loc(3),') => dt/C = ', z2dt/rCu_max |
---|
| 137 | WRITE(numout,*) ' Max Cv = ', rCv_max, ' at (i,j,k) = (',nCv_loc(1),nCv_loc(2),nCv_loc(3),') => dt/C = ', z2dt/rCv_max |
---|
| 138 | WRITE(numout,*) ' Max Cw = ', rCw_max, ' at (i,j,k) = (',nCw_loc(1),nCw_loc(2),nCw_loc(3),') => dt/C = ', z2dt/rCw_max |
---|
[4715] | 139 | ENDIF |
---|
[9019] | 140 | ! |
---|
[9124] | 141 | IF( ln_timing ) CALL timing_stop('dia_cfl') |
---|
[9019] | 142 | ! |
---|
[4715] | 143 | END SUBROUTINE dia_cfl |
---|
| 144 | |
---|
[9019] | 145 | |
---|
[4715] | 146 | SUBROUTINE dia_cfl_init |
---|
| 147 | !!---------------------------------------------------------------------- |
---|
| 148 | !! *** ROUTINE dia_cfl_init *** |
---|
| 149 | !! |
---|
| 150 | !! ** Purpose : create output file, initialise arrays |
---|
| 151 | !!---------------------------------------------------------------------- |
---|
[9019] | 152 | ! |
---|
| 153 | IF(lwp) THEN |
---|
| 154 | WRITE(numout,*) |
---|
| 155 | WRITE(numout,*) 'dia_cfl : Outputting CFL diagnostics to ',TRIM(clname), ' file' |
---|
| 156 | WRITE(numout,*) '~~~~~~~' |
---|
| 157 | WRITE(numout,*) |
---|
| 158 | ! |
---|
| 159 | ! create output ascii file |
---|
| 160 | CALL ctl_opn( numcfl, clname, 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL', 1, numout, lwp, 1 ) |
---|
| 161 | WRITE(numcfl,*) 'Timestep Direction Max C i j k' |
---|
| 162 | WRITE(numcfl,*) '******************************************' |
---|
[4715] | 163 | ENDIF |
---|
[9019] | 164 | ! |
---|
| 165 | rCu_max = 0._wp |
---|
| 166 | rCv_max = 0._wp |
---|
| 167 | rCw_max = 0._wp |
---|
| 168 | ! |
---|
[4715] | 169 | END SUBROUTINE dia_cfl_init |
---|
| 170 | |
---|
[9019] | 171 | !!====================================================================== |
---|
[4715] | 172 | END MODULE diacfl |
---|