Changeset 3294 for trunk/NEMOGCM/NEMO/OPA_SRC/DOM
- Timestamp:
- 2012-01-28T17:44:18+01:00 (12 years ago)
- Location:
- trunk/NEMOGCM/NEMO/OPA_SRC/DOM
- Files:
-
- 12 edited
- 1 copied
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/daymod.F90
r2715 r3294 34 34 USE restart ! 35 35 USE trc_oce, ONLY : lk_offline ! offline flag 36 USE timing ! Timing 36 37 37 38 IMPLICIT NONE … … 70 71 REAL(wp) :: zjul 71 72 !!---------------------------------------------------------------------- 72 73 ! 73 74 ! all calendar staff is based on the fact that MOD( rday, rdttra(1) ) == 0 74 75 IF( MOD( rday , rdttra(1) ) /= 0. ) CALL ctl_stop( 'the time step must devide the number of second of in a day' ) … … 127 128 ! call day to set the calendar parameters at the begining of the current simulaton. needed by iom_init 128 129 CALL day( nit000 ) 129 130 ! 130 131 END SUBROUTINE day_init 131 132 … … 204 205 REAL(wp) :: zprec ! fraction of day corresponding to 0.1 second 205 206 !!---------------------------------------------------------------------- 207 ! 208 IF( nn_timing == 1 ) CALL timing_start('day') 209 ! 206 210 zprec = 0.1 / rday 207 211 ! ! New time-step … … 255 259 IF( .NOT. lk_offline ) CALL rst_opn( kt ) ! Open the restart file if needed and control lrst_oce 256 260 IF( lrst_oce ) CALL day_rst( kt, 'WRITE' ) ! write day restart information 261 ! 262 IF( nn_timing == 1 ) CALL timing_stop('day') 257 263 ! 258 264 END SUBROUTINE day -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/dom_oce.F90
r2715 r3294 150 150 LOGICAL, PUBLIC, PARAMETER :: lk_vvl = .FALSE. !: fixed grid flag 151 151 #endif 152 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hur , hvr !: inverse of u and v-points ocean depth (1/m)153 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu , hv !: depth at u- and v-points (meters)154 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu_0 , hv_0 !: refernce depth at u- and v-points (meters)152 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:), TARGET :: hur , hvr !: inverse of u and v-points ocean depth (1/m) 153 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu , hv !: depth at u- and v-points (meters) 154 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: hu_0 , hv_0 !: refernce depth at u- and v-points (meters) 155 155 156 156 INTEGER, PUBLIC :: nla10 !: deepest W level Above ~10m (nlb10 - 1) -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domain.F90
r2528 r3294 35 35 USE c1d ! 1D vertical configuration 36 36 USE dyncor_c1d ! Coriolis term (c1d case) (cor_c1d routine) 37 USE timing ! Timing 37 38 38 39 IMPLICIT NONE … … 70 71 !!---------------------------------------------------------------------- 71 72 ! 73 IF( nn_timing == 1 ) CALL timing_start('dom_init') 74 ! 72 75 IF(lwp) THEN 73 76 WRITE(numout,*) … … 102 105 IF( nmsh /= 0 ) CALL dom_wri ! Create a domain file 103 106 IF( .NOT.ln_rstart ) CALL dom_ctl ! Domain control 107 ! 108 IF( nn_timing == 1 ) CALL timing_stop('dom_init') 104 109 ! 105 110 END SUBROUTINE dom_init -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domcfg.F90
r2715 r3294 15 15 USE in_out_manager ! I/O manager 16 16 USE lib_mpp ! distributed memory computing library 17 USE timing ! Timing 17 18 18 19 IMPLICIT NONE … … 35 36 !! 36 37 !!---------------------------------------------------------------------- 37 38 ! 39 IF( nn_timing == 1 ) CALL timing_start('dom_cfg') 40 ! 38 41 IF(lwp) THEN ! Control print 39 42 WRITE(numout,*) … … 56 59 ! 57 60 CALL dom_glo ! global domain versus zoom and/or local domain 61 ! 62 IF( nn_timing == 1 ) CALL timing_stop('dom_cfg') 58 63 ! 59 64 END SUBROUTINE dom_cfg -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domhgr.F90
r2715 r3294 25 25 USE in_out_manager ! I/O manager 26 26 USE lib_mpp ! MPP library 27 USE timing ! Timing 27 28 28 29 IMPLICIT NONE … … 105 106 REAL(wp) :: zphi1, zsin_alpha, zim05, zjm05 106 107 !!---------------------------------------------------------------------- 107 108 ! 109 IF( nn_timing == 1 ) CALL timing_start('dom_hgr') 110 ! 108 111 IF(lwp) THEN 109 112 WRITE(numout,*) … … 568 571 IF( znorme > 1.e-13 ) CALL ctl_stop( ' ===>>>> : symmetrical condition: rerun with good equator line' ) 569 572 ENDIF 570 573 ! 574 IF( nn_timing == 1 ) CALL timing_stop('dom_hgr') 575 ! 571 576 END SUBROUTINE dom_hgr 572 577 -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/dommsk.F90
r2715 r3294 25 25 USE oce ! ocean dynamics and tracers 26 26 USE dom_oce ! ocean space and time domain 27 USE obc_oce ! ocean open boundary conditions28 27 USE in_out_manager ! I/O manager 29 28 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 30 29 USE lib_mpp 31 30 USE dynspg_oce ! choice/control of key cpp for surface pressure gradient 31 USE wrk_nemo ! Memory allocation 32 USE timing ! Timing 32 33 33 34 IMPLICIT NONE … … 38 39 39 40 ! !!* Namelist namlbc : lateral boundary condition * 40 REAL(wp) :: rn_shlat = 2. ! type of lateral boundary condition on velocity 41 REAL(wp) :: rn_shlat = 2. ! type of lateral boundary condition on velocity 42 LOGICAL, PUBLIC :: ln_vorlat = .false. ! consistency of vorticity boundary condition 43 ! with analytical eqs. 44 41 45 42 46 INTEGER, ALLOCATABLE, SAVE, DIMENSION(:,:) :: icoord ! Workspace for dom_msk_nsa() … … 125 129 !! tmask_i : interior ocean mask 126 130 !!---------------------------------------------------------------------- 127 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released, iwrk_in_use, iwrk_not_released128 USE wrk_nemo, ONLY: zwf => wrk_2d_1 ! 2D real workspace129 USE wrk_nemo, ONLY: imsk => iwrk_2d_1 ! 2D integer workspace130 131 ! 131 132 INTEGER :: ji, jj, jk ! dummy loop indices 132 133 INTEGER :: iif, iil, ii0, ii1, ii ! local integers 133 134 INTEGER :: ijf, ijl, ij0, ij1 ! - - 134 !! 135 NAMELIST/namlbc/ rn_shlat 135 INTEGER , POINTER, DIMENSION(:,:) :: imsk 136 REAL(wp), POINTER, DIMENSION(:,:) :: zwf 137 !! 138 NAMELIST/namlbc/ rn_shlat, ln_vorlat 136 139 !!--------------------------------------------------------------------- 137 138 IF( wrk_in_use(2, 1) .OR. iwrk_in_use(2, 1) ) THEN 139 CALL ctl_stop('dom_msk: requested workspace arrays unavailable') ; RETURN 140 ENDIF 141 140 ! 141 IF( nn_timing == 1 ) CALL timing_start('dom_msk') 142 ! 143 CALL wrk_alloc( jpi, jpj, imsk ) 144 CALL wrk_alloc( jpi, jpj, zwf ) 145 ! 142 146 REWIND( numnam ) ! Namelist namlbc : lateral momentum boundary condition 143 147 READ ( numnam, namlbc ) … … 148 152 WRITE(numout,*) '~~~~~~' 149 153 WRITE(numout,*) ' Namelist namlbc' 150 WRITE(numout,*) ' lateral momentum boundary cond. rn_shlat = ',rn_shlat 154 WRITE(numout,*) ' lateral momentum boundary cond. rn_shlat = ',rn_shlat 155 WRITE(numout,*) ' consistency with analytical form ln_vorlat = ',ln_vorlat 151 156 ENDIF 152 157 … … 390 395 CALL lbc_lnk( fmask, 'F', 1._wp ) ! Lateral boundary conditions on fmask 391 396 397 ! CAUTION : The fmask may be further modified in dyn_vor_init ( dynvor.F90 ) 392 398 393 399 IF( nprint == 1 .AND. lwp ) THEN ! Control print … … 436 442 ENDIF 437 443 ! 438 IF( wrk_not_released(2, 1) .OR. & 439 iwrk_not_released(2, 1) ) CALL ctl_stop('dom_msk: failed to release workspace arrays') 444 CALL wrk_dealloc( jpi, jpj, imsk ) 445 CALL wrk_dealloc( jpi, jpj, zwf ) 446 ! 447 IF( nn_timing == 1 ) CALL timing_stop('dom_msk') 440 448 ! 441 449 END SUBROUTINE dom_msk … … 460 468 REAL(wp) :: zaa 461 469 !!--------------------------------------------------------------------- 462 470 ! 471 IF( nn_timing == 1 ) CALL timing_start('dom_msk_nsa') 472 ! 463 473 IF(lwp) WRITE(numout,*) 464 474 IF(lwp) WRITE(numout,*) 'dom_msk_nsa : noslip accurate boundary condition' … … 620 630 ENDIF 621 631 ! 632 IF( nn_timing == 1 ) CALL timing_stop('dom_msk_nsa') 633 ! 622 634 END SUBROUTINE dom_msk_nsa 623 635 -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domngb.F90
r2715 r3294 11 11 !!---------------------------------------------------------------------- 12 12 USE dom_oce ! ocean space and time domain 13 USE in_out_manager ! I/O manager 13 14 USE lib_mpp ! for mppsum 15 USE wrk_nemo ! Memory allocation 16 USE timing ! Timing 14 17 15 18 IMPLICIT NONE … … 34 37 !! -> not good if located at too high latitude... 35 38 !!---------------------------------------------------------------------- 36 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released37 USE wrk_nemo, ONLY: zglam => wrk_2d_2 , zgphi => wrk_2d_3 , zmask => wrk_2d_4 , zdist => wrk_2d_538 39 ! 39 40 REAL(wp) , INTENT(in ) :: plon, plat ! longitude,latitude of the point … … 43 44 INTEGER , DIMENSION(2) :: iloc 44 45 REAL(wp) :: zlon, zmini 46 REAL(wp), POINTER, DIMENSION(:,:) :: zglam, zgphi, zmask, zdist 45 47 !!-------------------------------------------------------------------- 46 48 ! 47 IF( wrk_in_use(2, 2,3,4,5) ) CALL ctl_stop('dom_ngb: Requested workspaces already in use') 49 IF( nn_timing == 1 ) CALL timing_start('dom_ngb') 50 ! 51 CALL wrk_alloc( jpi, jpj, zglam, zgphi, zmask, zdist ) 48 52 ! 49 53 zmask(:,:) = 0._wp … … 72 76 ENDIF 73 77 ! 74 IF( wrk_not_released(2, 2,3,4,5) ) CALL ctl_stop('dom_ngb: error releasing workspaces') 78 CALL wrk_dealloc( jpi, jpj, zglam, zgphi, zmask, zdist ) 79 ! 80 IF( nn_timing == 1 ) CALL timing_stop('dom_ngb') 75 81 ! 76 82 END SUBROUTINE dom_ngb -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domvvl.F90
r2779 r3294 20 20 USE lib_mpp ! distributed memory computing library 21 21 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 22 USE wrk_nemo ! Memory allocation 23 USE timing ! Timing 22 24 23 25 IMPLICIT NONE 24 26 PRIVATE 25 27 26 PUBLIC dom_vvl ! called by domain.F9027 PUBLIC dom_vvl_ alloc ! called by nemogcm.F9028 29 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: ee_t, ee_u, ee_v, ee_f !: ??? 30 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: mut , muu , muv , muf !: ???28 PUBLIC dom_vvl ! called by domain.F90 29 PUBLIC dom_vvl_2 ! called by domain.F90 30 PUBLIC dom_vvl_alloc ! called by nemogcm.F90 31 32 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: mut , muu , muv , muf !: 1/H_0 at t-,u-,v-,f-points 31 33 32 34 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: r2dt ! vertical profile time-step, = 2 rdttra … … 49 51 ! 50 52 ALLOCATE( mut (jpi,jpj,jpk) , muu (jpi,jpj,jpk) , muv (jpi,jpj,jpk) , muf (jpi,jpj,jpk) , & 51 & ee_t(jpi,jpj) , ee_u(jpi,jpj) , ee_v(jpi,jpj) , ee_f(jpi,jpj) , &52 53 & r2dt (jpk) , STAT=dom_vvl_alloc ) 53 54 ! … … 62 63 !! *** ROUTINE dom_vvl *** 63 64 !! 64 !! ** Purpose : compute coefficients muX at T-U-V-F points to spread65 !! ssh over the whole water column (scale factors)66 !! ----------------------------------------------------------------------67 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released68 USE wrk_nemo, ONLY: zs_t => wrk_2d_1 , zs_u_1 => wrk_2d_2 , zs_v_1 => wrk_2d_3 ! 2D workspace65 !! ** Purpose : compute mu coefficients at t-, u-, v- and f-points to 66 !! spread ssh over the whole water column (scale factors) 67 !! set the before and now ssh at u- and v-points 68 !! (also f-point in now case) 69 !!---------------------------------------------------------------------- 69 70 ! 70 71 INTEGER :: ji, jj, jk ! dummy loop indices 71 REAL(wp) :: zcoefu , zcoefv , zcoeff ! local scalars 72 REAL(wp) :: zv_t_ij, zv_t_ip1j, zv_t_ijp1, zv_t_ip1jp1 ! - - 73 !!---------------------------------------------------------------------- 74 75 IF( wrk_in_use(2, 1,2,3) ) THEN 76 CALL ctl_stop('dom_vvl: requested workspace arrays unavailable') ; RETURN 77 ENDIF 78 72 REAL(wp) :: zcoefu, zcoefv , zcoeff ! local scalars 73 REAL(wp) :: zvt , zvt_ip1, zvt_jp1, zvt_ip1jp1 ! - - 74 REAL(wp), POINTER, DIMENSION(:,:) :: zee_t, zee_u, zee_v, zee_f ! 2D workspace 75 !!---------------------------------------------------------------------- 76 ! 77 IF( nn_timing == 1 ) CALL timing_start('dom_vvl') 78 ! 79 CALL wrk_alloc( jpi, jpj, zee_t, zee_u, zee_v, zee_f ) 80 ! 79 81 IF(lwp) THEN 80 82 WRITE(numout,*) … … 97 99 98 100 ! !== mu computation ==! 99 ee_t(:,:) = fse3t_0(:,:,1) ! Lower bound : thickness of the first model level100 ee_u(:,:) = fse3u_0(:,:,1)101 ee_v(:,:) = fse3v_0(:,:,1)102 ee_f(:,:) = fse3f_0(:,:,1)101 zee_t(:,:) = fse3t_0(:,:,1) ! Lower bound : thickness of the first model level 102 zee_u(:,:) = fse3u_0(:,:,1) 103 zee_v(:,:) = fse3v_0(:,:,1) 104 zee_f(:,:) = fse3f_0(:,:,1) 103 105 DO jk = 2, jpkm1 ! Sum of the masked vertical scale factors 104 ee_t(:,:) =ee_t(:,:) + fse3t_0(:,:,jk) * tmask(:,:,jk)105 ee_u(:,:) =ee_u(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk)106 ee_v(:,:) =ee_v(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk)106 zee_t(:,:) = zee_t(:,:) + fse3t_0(:,:,jk) * tmask(:,:,jk) 107 zee_u(:,:) = zee_u(:,:) + fse3u_0(:,:,jk) * umask(:,:,jk) 108 zee_v(:,:) = zee_v(:,:) + fse3v_0(:,:,jk) * vmask(:,:,jk) 107 109 DO jj = 1, jpjm1 ! f-point : fmask=shlat at coasts, use the product of umask 108 ee_f(:,jj) =ee_f(:,jj) + fse3f_0(:,jj,jk) * umask(:,jj,jk) * umask(:,jj+1,jk)110 zee_f(:,jj) = zee_f(:,jj) + fse3f_0(:,jj,jk) * umask(:,jj,jk) * umask(:,jj+1,jk) 109 111 END DO 110 112 END DO 111 113 ! ! Compute and mask the inverse of the local depth at T, U, V and F points 112 ee_t(:,:) = 1. /ee_t(:,:) * tmask(:,:,1)113 ee_u(:,:) = 1. /ee_u(:,:) * umask(:,:,1)114 ee_v(:,:) = 1. /ee_v(:,:) * vmask(:,:,1)114 zee_t(:,:) = 1._wp / zee_t(:,:) * tmask(:,:,1) 115 zee_u(:,:) = 1._wp / zee_u(:,:) * umask(:,:,1) 116 zee_v(:,:) = 1._wp / zee_v(:,:) * vmask(:,:,1) 115 117 DO jj = 1, jpjm1 ! f-point case fmask cannot be used 116 ee_f(:,jj) = 1. /ee_f(:,jj) * umask(:,jj,1) * umask(:,jj+1,1)117 END DO 118 CALL lbc_lnk( ee_f, 'F', 1. )! lateral boundary condition on ee_f118 zee_f(:,jj) = 1._wp / zee_f(:,jj) * umask(:,jj,1) * umask(:,jj+1,1) 119 END DO 120 CALL lbc_lnk( zee_f, 'F', 1. ) ! lateral boundary condition on ee_f 119 121 ! 120 122 DO jk = 1, jpk ! mu coefficients 121 mut(:,:,jk) = ee_t(:,:) * tmask(:,:,jk) ! T-point at T levels122 muu(:,:,jk) = ee_u(:,:) * umask(:,:,jk) ! U-point at T levels123 muv(:,:,jk) = ee_v(:,:) * vmask(:,:,jk) ! V-point at T levels123 mut(:,:,jk) = zee_t(:,:) * tmask(:,:,jk) ! T-point at T levels 124 muu(:,:,jk) = zee_u(:,:) * umask(:,:,jk) ! U-point at T levels 125 muv(:,:,jk) = zee_v(:,:) * vmask(:,:,jk) ! V-point at T levels 124 126 END DO 125 127 DO jk = 1, jpk ! F-point : fmask=shlat at coasts, use the product of umask 126 128 DO jj = 1, jpjm1 127 muf(:,jj,jk) = ee_f(:,jj) * umask(:,jj,jk) * umask(:,jj+1,jk) ! at T levels128 END DO 129 muf(:,jpj,jk) = 0. e0129 muf(:,jj,jk) = zee_f(:,jj) * umask(:,jj,jk) * umask(:,jj+1,jk) ! at T levels 130 END DO 131 muf(:,jpj,jk) = 0._wp 130 132 END DO 131 133 CALL lbc_lnk( muf, 'F', 1. ) ! lateral boundary condition … … 139 141 END DO 140 142 141 ! surface at t-points and inverse surface at (u/v)-points used in surface averaging computations142 ! for ssh and scale factors143 zs_t (:,:) = e1t(:,:) * e2t(:,:)144 zs_u_1(:,:) = 0.5 / ( e1u(:,:) * e2u(:,:) )145 zs_v_1(:,:) = 0.5 / ( e1v(:,:) * e2v(:,:) )146 147 143 DO jj = 1, jpjm1 ! initialise before and now Sea Surface Height at u-, v-, f-points 148 144 DO ji = 1, jpim1 ! NO vector opt. 149 zcoefu = umask(ji,jj,1) * zs_u_1(ji,jj)150 zcoefv = vmask(ji,jj,1) * zs_v_1(ji,jj)151 zcoeff = 0. 5 * umask(ji,jj,1) * umask(ji,jj+1,1) / ( e1f(ji,jj) * e2f(ji,jj))152 ! before fields153 zv _t_ij = zs_t(ji ,jj ) * sshb(ji ,jj )154 zv _t_ip1j = zs_t(ji+1,jj ) * sshb(ji+1,jj )155 zv _t_ijp1 = zs_t(ji ,jj+1) * sshb(ji ,jj+1)156 sshu_b(ji,jj) = zcoefu * ( zv _t_ij + zv_t_ip1j)157 sshv_b(ji,jj) = zcoefv * ( zv _t_ij + zv_t_ijp1 )158 ! now fields159 zv _t_ij = zs_t(ji ,jj ) * sshn(ji ,jj )160 zv _t_ip1j = zs_t(ji+1,jj ) * sshn(ji+1,jj )161 zv _t_ijp1 = zs_t(ji ,jj+1) * sshn(ji ,jj+1)162 zv _t_ip1jp1 = zs_t(ji ,jj+1) * sshn(ji,jj+1)163 sshu_n(ji,jj) = zcoefu * ( zv _t_ij + zv_t_ip1j)164 sshv_n(ji,jj) = zcoefv * ( zv _t_ij + zv_t_ijp1 )165 sshf_n(ji,jj) = zcoeff * ( zv _t_ij + zv_t_ip1j + zv_t_ijp1 + zv_t_ip1jp1 )145 zcoefu = 0.50_wp / ( e1u(ji,jj) * e2u(ji,jj) ) * umask(ji,jj,1) 146 zcoefv = 0.50_wp / ( e1v(ji,jj) * e2v(ji,jj) ) * vmask(ji,jj,1) 147 zcoeff = 0.25_wp / ( e1f(ji,jj) * e2f(ji,jj) ) * umask(ji,jj,1) * umask(ji,jj+1,1) 148 ! 149 zvt = e1e2t(ji ,jj ) * sshb(ji ,jj ) ! before fields 150 zvt_ip1 = e1e2t(ji+1,jj ) * sshb(ji+1,jj ) 151 zvt_jp1 = e1e2t(ji ,jj+1) * sshb(ji ,jj+1) 152 sshu_b(ji,jj) = zcoefu * ( zvt + zvt_ip1 ) 153 sshv_b(ji,jj) = zcoefv * ( zvt + zvt_jp1 ) 154 ! 155 zvt = e1e2t(ji ,jj ) * sshn(ji ,jj ) ! now fields 156 zvt_ip1 = e1e2t(ji+1,jj ) * sshn(ji+1,jj ) 157 zvt_jp1 = e1e2t(ji ,jj+1) * sshn(ji ,jj+1) 158 zvt_ip1jp1 = e1e2t(ji+1,jj+1) * sshn(ji+1,jj+1) 159 sshu_n(ji,jj) = zcoefu * ( zvt + zvt_ip1 ) 160 sshv_n(ji,jj) = zcoefv * ( zvt + zvt_jp1 ) 161 sshf_n(ji,jj) = zcoeff * ( zvt + zvt_ip1 + zvt_jp1 + zvt_ip1jp1 ) 166 162 END DO 167 163 END DO … … 169 165 CALL lbc_lnk( sshv_n, 'V', 1. ) ; CALL lbc_lnk( sshv_b, 'V', 1. ) 170 166 CALL lbc_lnk( sshf_n, 'F', 1. ) 171 172 ! initialise before scale factors at (u/v)-points 173 ! Scale factor anomaly at (u/v)-points: surface averaging of scale factor at t-points 174 DO jk = 1, jpkm1 175 DO jj = 1, jpjm1 176 DO ji = 1, jpim1 177 zv_t_ij = zs_t(ji ,jj ) * fse3t_b(ji ,jj ,jk) 178 zv_t_ip1j = zs_t(ji+1,jj ) * fse3t_b(ji+1,jj ,jk) 179 zv_t_ijp1 = zs_t(ji ,jj+1) * fse3t_b(ji ,jj+1,jk) 180 fse3u_b(ji,jj,jk) = umask(ji,jj,jk) * ( zs_u_1(ji,jj) * ( zv_t_ij + zv_t_ip1j ) - fse3u_0(ji,jj,jk) ) 181 fse3v_b(ji,jj,jk) = vmask(ji,jj,jk) * ( zs_v_1(ji,jj) * ( zv_t_ij + zv_t_ijp1 ) - fse3v_0(ji,jj,jk) ) 182 END DO 183 END DO 184 END DO 185 CALL lbc_lnk( fse3u_b(:,:,:), 'U', 1. ) ! lateral boundary conditions 186 CALL lbc_lnk( fse3v_b(:,:,:), 'V', 1. ) 187 ! Add initial scale factor to scale factor anomaly 188 fse3u_b(:,:,:) = fse3u_b(:,:,:) + fse3u_0(:,:,:) 189 fse3v_b(:,:,:) = fse3v_b(:,:,:) + fse3v_0(:,:,:) 190 ! 191 IF( wrk_not_released(2, 1,2,3) ) CALL ctl_stop('dom_vvl: failed to release workspace arrays') 167 ! 168 CALL wrk_dealloc( jpi, jpj, zee_t, zee_u, zee_v, zee_f ) 169 ! 170 IF( nn_timing == 1 ) CALL timing_stop('dom_vvl') 192 171 ! 193 172 END SUBROUTINE dom_vvl 194 173 174 175 SUBROUTINE dom_vvl_2( kt, pe3u_b, pe3v_b ) 176 !!---------------------------------------------------------------------- 177 !! *** ROUTINE dom_vvl_2 *** 178 !! 179 !! ** Purpose : compute the vertical scale factors at u- and v-points 180 !! in variable volume case. 181 !! 182 !! ** Method : In variable volume case (non linear sea surface) the 183 !! the vertical scale factor at velocity points is computed 184 !! as the average of the cell surface weighted e3t. 185 !! It uses the sea surface heigth so it have to be initialized 186 !! after ssh is read/set 187 !!---------------------------------------------------------------------- 188 INTEGER , INTENT(in ) :: kt ! ocean time-step index 189 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: pe3u_b, pe3v_b ! before vertical scale factor at u- & v-pts 190 ! 191 INTEGER :: ji, jj, jk ! dummy loop indices 192 INTEGER :: iku, ikv ! local integers 193 INTEGER :: ii0, ii1, ij0, ij1 ! temporary integers 194 REAL(wp) :: zvt ! local scalars 195 !!---------------------------------------------------------------------- 196 ! 197 IF( nn_timing == 1 ) CALL timing_start('dom_vvl_2') 198 ! 199 IF( lwp .AND. kt == nit000 ) THEN 200 WRITE(numout,*) 201 WRITE(numout,*) 'dom_vvl_2 : Variable volume, fse3t_b initialization' 202 WRITE(numout,*) '~~~~~~~~~ ' 203 pe3u_b(:,:,jpk) = fse3u_0(:,:,jpk) 204 pe3v_b(:,:,jpk) = fse3u_0(:,:,jpk) 205 ENDIF 206 207 DO jk = 1, jpkm1 ! set the before scale factors at u- & v-points 208 DO jj = 2, jpjm1 209 DO ji = fs_2, fs_jpim1 210 zvt = fse3t_b(ji,jj,jk) * e1e2t(ji,jj) 211 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1e2t(ji+1,jj) ) / ( e1u(ji,jj) * e2u(ji,jj) ) 212 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e1e2t(ji,jj+1) ) / ( e1v(ji,jj) * e2v(ji,jj) ) 213 END DO 214 END DO 215 END DO 216 217 ! Correct scale factors at locations that have been individually modified in domhgr 218 ! Such modifications break the relationship between e1e2t and e1u*e2u etc. Recompute 219 ! scale factors ignoring the modified metric. 220 ! ! ===================== 221 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN ! ORCA R2 configuration 222 ! ! ===================== 223 IF( nn_cla == 0 ) THEN 224 ! 225 ii0 = 139 ; ii1 = 140 ! Gibraltar Strait (e2u was modified) 226 ij0 = 102 ; ij1 = 102 227 DO jk = 1, jpkm1 ! set the before scale factors at u-points 228 DO jj = mj0(ij0), mj1(ij1) 229 DO ji = mi0(ii0), mi1(ii1) 230 zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 231 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 232 END DO 233 END DO 234 END DO 235 ! 236 ii0 = 160 ; ii1 = 160 ! Bab el Mandeb (e2u and e1v were modified) 237 ij0 = 88 ; ij1 = 88 238 DO jk = 1, jpkm1 ! set the before scale factors at u-points 239 DO jj = mj0(ij0), mj1(ij1) 240 DO ji = mi0(ii0), mi1(ii1) 241 zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 242 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 243 END DO 244 END DO 245 END DO 246 DO jk = 1, jpkm1 ! set the before scale factors at v-points 247 DO jj = mj0(ij0), mj1(ij1) 248 DO ji = mi0(ii0), mi1(ii1) 249 zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 250 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 251 END DO 252 END DO 253 END DO 254 ENDIF 255 256 ii0 = 145 ; ii1 = 146 ! Danish Straits (e2u was modified) 257 ij0 = 116 ; ij1 = 116 258 DO jk = 1, jpkm1 ! set the before scale factors at u-points 259 DO jj = mj0(ij0), mj1(ij1) 260 DO ji = mi0(ii0), mi1(ii1) 261 zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 262 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 263 END DO 264 END DO 265 END DO 266 ! 267 ENDIF 268 ! ! ===================== 269 IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN ! ORCA R1 configuration 270 ! ! ===================== 271 272 ii0 = 281 ; ii1 = 282 ! Gibraltar Strait (e2u was modified) 273 ij0 = 200 ; ij1 = 200 274 DO jk = 1, jpkm1 ! set the before scale factors at u-points 275 DO jj = mj0(ij0), mj1(ij1) 276 DO ji = mi0(ii0), mi1(ii1) 277 zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 278 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 279 END DO 280 END DO 281 END DO 282 283 ii0 = 314 ; ii1 = 315 ! Bhosporus Strait (e2u was modified) 284 ij0 = 208 ; ij1 = 208 285 DO jk = 1, jpkm1 ! set the before scale factors at u-points 286 DO jj = mj0(ij0), mj1(ij1) 287 DO ji = mi0(ii0), mi1(ii1) 288 zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 289 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 290 END DO 291 END DO 292 END DO 293 294 ii0 = 44 ; ii1 = 44 ! Lombok Strait (e1v was modified) 295 ij0 = 124 ; ij1 = 125 296 DO jk = 1, jpkm1 ! set the before scale factors at v-points 297 DO jj = mj0(ij0), mj1(ij1) 298 DO ji = mi0(ii0), mi1(ii1) 299 zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 300 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 301 END DO 302 END DO 303 END DO 304 305 ii0 = 48 ; ii1 = 48 ! Sumba Strait (e1v was modified) [closed from bathy_11 on] 306 ij0 = 124 ; ij1 = 125 307 DO jk = 1, jpkm1 ! set the before scale factors at v-points 308 DO jj = mj0(ij0), mj1(ij1) 309 DO ji = mi0(ii0), mi1(ii1) 310 zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 311 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 312 END DO 313 END DO 314 END DO 315 316 ii0 = 53 ; ii1 = 53 ! Ombai Strait (e1v was modified) 317 ij0 = 124 ; ij1 = 125 318 DO jk = 1, jpkm1 ! set the before scale factors at v-points 319 DO jj = mj0(ij0), mj1(ij1) 320 DO ji = mi0(ii0), mi1(ii1) 321 zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 322 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 323 END DO 324 END DO 325 END DO 326 327 ii0 = 56 ; ii1 = 56 ! Timor Passage (e1v was modified) 328 ij0 = 124 ; ij1 = 125 329 DO jk = 1, jpkm1 ! set the before scale factors at v-points 330 DO jj = mj0(ij0), mj1(ij1) 331 DO ji = mi0(ii0), mi1(ii1) 332 zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 333 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 334 END DO 335 END DO 336 END DO 337 338 ii0 = 55 ; ii1 = 55 ! West Halmahera Strait (e1v was modified) 339 ij0 = 141 ; ij1 = 142 340 DO jk = 1, jpkm1 ! set the before scale factors at v-points 341 DO jj = mj0(ij0), mj1(ij1) 342 DO ji = mi0(ii0), mi1(ii1) 343 zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 344 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 345 END DO 346 END DO 347 END DO 348 349 ii0 = 58 ; ii1 = 58 ! East Halmahera Strait (e1v was modified) 350 ij0 = 141 ; ij1 = 142 351 DO jk = 1, jpkm1 ! set the before scale factors at v-points 352 DO jj = mj0(ij0), mj1(ij1) 353 DO ji = mi0(ii0), mi1(ii1) 354 zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 355 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 356 END DO 357 END DO 358 END DO 359 360 ! 361 ENDIF 362 ! ! ====================== 363 IF( cp_cfg == "orca" .AND. jp_cfg == 05 ) THEN ! ORCA R05 configuration 364 ! ! ====================== 365 ii0 = 563 ; ii1 = 564 ! Gibraltar Strait (e2u was modified) 366 ij0 = 327 ; ij1 = 327 367 DO jk = 1, jpkm1 ! set the before scale factors at u-points 368 DO jj = mj0(ij0), mj1(ij1) 369 DO ji = mi0(ii0), mi1(ii1) 370 zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 371 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 372 END DO 373 END DO 374 END DO 375 ! 376 ii0 = 627 ; ii1 = 628 ! Bosphore Strait (e2u was modified) 377 ij0 = 343 ; ij1 = 343 378 DO jk = 1, jpkm1 ! set the before scale factors at u-points 379 DO jj = mj0(ij0), mj1(ij1) 380 DO ji = mi0(ii0), mi1(ii1) 381 zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 382 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 383 END DO 384 END DO 385 END DO 386 ! 387 ii0 = 93 ; ii1 = 94 ! Sumba Strait (e2u was modified) 388 ij0 = 232 ; ij1 = 232 389 DO jk = 1, jpkm1 ! set the before scale factors at u-points 390 DO jj = mj0(ij0), mj1(ij1) 391 DO ji = mi0(ii0), mi1(ii1) 392 zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 393 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 394 END DO 395 END DO 396 END DO 397 ! 398 ii0 = 103 ; ii1 = 103 ! Ombai Strait (e2u was modified) 399 ij0 = 232 ; ij1 = 232 400 DO jk = 1, jpkm1 ! set the before scale factors at u-points 401 DO jj = mj0(ij0), mj1(ij1) 402 DO ji = mi0(ii0), mi1(ii1) 403 zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 404 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 405 END DO 406 END DO 407 END DO 408 ! 409 ii0 = 15 ; ii1 = 15 ! Palk Strait (e2u was modified) 410 ij0 = 270 ; ij1 = 270 411 DO jk = 1, jpkm1 ! set the before scale factors at u-points 412 DO jj = mj0(ij0), mj1(ij1) 413 DO ji = mi0(ii0), mi1(ii1) 414 zvt = fse3t_b(ji,jj,jk) * e1t(ji,jj) 415 pe3u_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji+1,jj,jk) * e1t(ji+1,jj) ) / ( e1u(ji,jj) ) 416 END DO 417 END DO 418 END DO 419 ! 420 ii0 = 87 ; ii1 = 87 ! Lombok Strait (e1v was modified) 421 ij0 = 232 ; ij1 = 233 422 DO jk = 1, jpkm1 ! set the before scale factors at v-points 423 DO jj = mj0(ij0), mj1(ij1) 424 DO ji = mi0(ii0), mi1(ii1) 425 zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 426 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 427 END DO 428 END DO 429 END DO 430 ! 431 ii0 = 662 ; ii1 = 662 ! Bab el Mandeb (e1v was modified) 432 ij0 = 276 ; ij1 = 276 433 DO jk = 1, jpkm1 ! set the before scale factors at v-points 434 DO jj = mj0(ij0), mj1(ij1) 435 DO ji = mi0(ii0), mi1(ii1) 436 zvt = fse3t_b(ji,jj,jk) * e2t(ji,jj) 437 pe3v_b(ji,jj,jk) = 0.5_wp * ( zvt + fse3t_b(ji,jj+1,jk) * e2t(ji,jj+1) ) / ( e2v(ji,jj) ) 438 END DO 439 END DO 440 END DO 441 ! 442 ENDIF 443 ! End of individual corrections to scale factors 444 445 IF( ln_zps ) THEN ! minimum of the e3t at partial cell level 446 DO jj = 2, jpjm1 447 DO ji = fs_2, fs_jpim1 448 iku = mbku(ji,jj) 449 ikv = mbkv(ji,jj) 450 pe3u_b(ji,jj,iku) = MIN( fse3t_b(ji,jj,iku), fse3t_b(ji+1,jj ,iku) ) 451 pe3v_b(ji,jj,ikv) = MIN( fse3t_b(ji,jj,ikv), fse3t_b(ji ,jj+1,ikv) ) 452 END DO 453 END DO 454 ENDIF 455 456 pe3u_b(:,:,:) = pe3u_b(:,:,:) - fse3u_0(:,:,:) ! anomaly to avoid zero along closed boundary/extra halos 457 pe3v_b(:,:,:) = pe3v_b(:,:,:) - fse3v_0(:,:,:) 458 CALL lbc_lnk( pe3u_b(:,:,:), 'U', 1. ) ! lateral boundary conditions 459 CALL lbc_lnk( pe3v_b(:,:,:), 'V', 1. ) 460 pe3u_b(:,:,:) = pe3u_b(:,:,:) + fse3u_0(:,:,:) ! recover the full scale factor 461 pe3v_b(:,:,:) = pe3v_b(:,:,:) + fse3v_0(:,:,:) 462 ! 463 IF( nn_timing == 1 ) CALL timing_stop('dom_vvl_2') 464 ! 465 END SUBROUTINE dom_vvl_2 466 195 467 #else 196 468 !!---------------------------------------------------------------------- … … 200 472 SUBROUTINE dom_vvl 201 473 END SUBROUTINE dom_vvl 474 SUBROUTINE dom_vvl_2(kdum, pudum, pvdum ) 475 USE par_kind 476 INTEGER , INTENT(in ) :: kdum 477 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: pudum, pvdum 478 END SUBROUTINE dom_vvl_2 202 479 #endif 203 480 -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domwri.F90
r2715 r3294 20 20 USE lbclnk ! lateral boundary conditions - mpp exchanges 21 21 USE lib_mpp ! MPP library 22 USE wrk_nemo ! Memory allocation 23 USE timing ! Timing 22 24 23 25 IMPLICIT NONE … … 63 65 !! masks, depth and vertical scale factors 64 66 !!---------------------------------------------------------------------- 65 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released66 USE wrk_nemo, ONLY: zprt => wrk_2d_1 , zprw => wrk_2d_2 ! 2D workspace67 USE wrk_nemo, ONLY: zdepu => wrk_3d_1 , zdepv => wrk_3d_2 ! 3D -68 67 !! 69 68 INTEGER :: inum0 ! temprary units for 'mesh_mask.nc' file … … 78 77 CHARACTER(len=21) :: clnam4 ! filename (vertical mesh informations) 79 78 INTEGER :: ji, jj, jk ! dummy loop indices 80 !!---------------------------------------------------------------------- 81 82 IF( wrk_in_use(2, 1,2) .OR. wrk_in_use(3, 1,2) )THEN 83 CALL ctl_stop('dom_wri: requested workspace arrays unavailable') ; RETURN 84 END IF 85 79 ! ! workspaces 80 REAL(wp), POINTER, DIMENSION(:,: ) :: zprt, zprw 81 REAL(wp), POINTER, DIMENSION(:,:,:) :: zdepu, zdepv 82 !!---------------------------------------------------------------------- 83 ! 84 IF( nn_timing == 1 ) CALL timing_start('dom_wri') 85 ! 86 CALL wrk_alloc( jpi, jpj, zprt, zprw ) 87 CALL wrk_alloc( jpi, jpj, jpk, zdepu, zdepv ) 88 ! 86 89 IF(lwp) WRITE(numout,*) 87 90 IF(lwp) WRITE(numout,*) 'dom_wri : create NetCDF mesh and mask information file(s)' … … 260 263 END SELECT 261 264 ! 262 IF( wrk_not_released(2, 1,2) .OR. & 263 wrk_not_released(3, 1,2) ) CALL ctl_stop('dom_wri: failed to release workspace arrays') 265 CALL wrk_dealloc( jpi, jpj, zprt, zprw ) 266 CALL wrk_dealloc( jpi, jpj, jpk, zdepu, zdepv ) 267 ! 268 IF( nn_timing == 1 ) CALL timing_stop('dom_wri') 264 269 ! 265 270 END SUBROUTINE dom_wri … … 275 280 !! 2) check which elements have been changed 276 281 !!---------------------------------------------------------------------- 277 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released278 USE wrk_nemo, ONLY: ztstref => wrk_2d_3 ! array with different values for each element279 282 ! 280 283 CHARACTER(len=1) , INTENT(in ) :: cdgrd ! … … 284 287 INTEGER :: ji ! dummy loop indices 285 288 LOGICAL, DIMENSION(SIZE(puniq,1),SIZE(puniq,2),1) :: lldbl ! store whether each point is unique or not 286 !!---------------------------------------------------------------------- 287 288 IF( wrk_in_use(2, 3) ) THEN 289 CALL ctl_stop('dom_uniq: requested workspace array unavailable') ; RETURN 290 ENDIF 291 289 REAL(wp), POINTER, DIMENSION(:,:) :: ztstref 290 !!---------------------------------------------------------------------- 291 ! 292 IF( nn_timing == 1 ) CALL timing_start('dom_uniq') 293 ! 294 CALL wrk_alloc( jpi, jpj, ztstref ) 295 ! 292 296 ! build an array with different values for each element 293 297 ! in mpp: make sure that these values are different even between process … … 304 308 puniq(nldi:nlei,nldj:nlej) = REAL( COUNT( lldbl(nldi:nlei,nldj:nlej,:), dim = 3 ) , wp ) 305 309 ! 306 IF( wrk_not_released(2, 3) ) CALL ctl_stop('dom_uniq: failed to release workspace array') 310 CALL wrk_dealloc( jpi, jpj, ztstref ) 311 ! 312 IF( nn_timing == 1 ) CALL timing_stop('dom_uniq') 307 313 ! 308 314 END SUBROUTINE dom_uniq -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/domzgr.F90
r2950 r3294 38 38 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 39 39 USE lib_mpp ! distributed memory computing library 40 USE wrk_nemo ! Memory allocation 41 USE timing ! Timing 40 42 41 43 IMPLICIT NONE … … 86 88 NAMELIST/namzgr/ ln_zco, ln_zps, ln_sco 87 89 !!---------------------------------------------------------------------- 88 90 ! 91 IF( nn_timing == 1 ) CALL timing_start('dom_zgr') 92 ! 89 93 REWIND( numnam ) ! Read Namelist namzgr : vertical coordinate' 90 94 READ ( numnam, namzgr ) … … 139 143 ENDIF 140 144 ! 145 IF( nn_timing == 1 ) CALL timing_stop('dom_zgr') 146 ! 141 147 END SUBROUTINE dom_zgr 142 148 … … 170 176 REAL(wp) :: za2, zkth2, zacr2 ! Values for optional double tanh function set from parameters 171 177 !!---------------------------------------------------------------------- 172 178 ! 179 IF( nn_timing == 1 ) CALL timing_start('zgr_z') 180 ! 173 181 ! Set variables from parameters 174 182 ! ------------------------------ … … 280 288 END DO 281 289 ! 290 IF( nn_timing == 1 ) CALL timing_stop('zgr_z') 291 ! 282 292 END SUBROUTINE zgr_z 283 293 … … 319 329 REAL(wp) :: r_bump , h_bump , h_oce ! bump characteristics 320 330 REAL(wp) :: zi, zj, zh, zhmin ! local scalars 321 INTEGER , DIMENSION(jpidta,jpjdta) :: idta ! global domain integer data 322 REAL(wp), DIMENSION(jpidta,jpjdta) :: zdta ! global domain scalar data 323 !!---------------------------------------------------------------------- 324 331 INTEGER , POINTER, DIMENSION(:,:) :: idta ! global domain integer data 332 REAL(wp), POINTER, DIMENSION(:,:) :: zdta ! global domain scalar data 333 !!---------------------------------------------------------------------- 334 ! 335 IF( nn_timing == 1 ) CALL timing_start('zgr_bat') 336 ! 337 CALL wrk_alloc( jpidta, jpjdta, idta ) 338 CALL wrk_alloc( jpidta, jpjdta, zdta ) 339 ! 325 340 IF(lwp) WRITE(numout,*) 326 341 IF(lwp) WRITE(numout,*) ' zgr_bat : defines level and meter bathymetry' … … 440 455 CALL iom_close( inum ) 441 456 ! ! ===================== 442 IF( cp_cfg == "orca" .AND. jp_cfg == 1 ) THEN ! ORCA R1 configuration443 ii0 = 51 ; ii1 = 53444 ij0 = 142 ; ij1 = 142 ! =====================445 DO ji = mi0(ii0), mi1(ii1) ! Close Halmera Strait446 DO jj = mj0(ij0), mj1(ij1)447 bathy(ji,jj) = 0._wp448 END DO449 END DO450 IF(lwp) WRITE(numout,*)451 IF(lwp) WRITE(numout,*) ' orca_r1: Halmera strait closed at i=',ii0,' j=',ij0,'->',ij1452 ENDIF453 ! ! =====================454 457 IF( cp_cfg == "orca" .AND. jp_cfg == 2 ) THEN ! ORCA R2 configuration 455 458 ! ! ===================== … … 512 515 ENDIF 513 516 ! 517 CALL wrk_dealloc( jpidta, jpjdta, idta ) 518 CALL wrk_dealloc( jpidta, jpjdta, zdta ) 519 ! 520 IF( nn_timing == 1 ) CALL timing_stop('zgr_bat') 521 ! 514 522 END SUBROUTINE zgr_bat 515 523 … … 589 597 !! - update bathy : meter bathymetry (in meters) 590 598 !!---------------------------------------------------------------------- 591 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released592 USE wrk_nemo, ONLY: zbathy => wrk_2d_1593 599 !! 594 600 INTEGER :: ji, jj, jl ! dummy loop indices 595 601 INTEGER :: icompt, ibtest, ikmax ! temporary integers 596 !!---------------------------------------------------------------------- 597 598 IF( wrk_in_use(2, 1) ) THEN 599 CALL ctl_stop('zgr_bat_ctl: requested workspace array unavailable') ; RETURN 600 ENDIF 601 602 REAL(wp), POINTER, DIMENSION(:,:) :: zbathy 603 !!---------------------------------------------------------------------- 604 ! 605 IF( nn_timing == 1 ) CALL timing_start('zgr_bat_ctl') 606 ! 607 CALL wrk_alloc( jpi, jpj, zbathy ) 608 ! 602 609 IF(lwp) WRITE(numout,*) 603 610 IF(lwp) WRITE(numout,*) ' zgr_bat_ctl : check the bathymetry' … … 702 709 ENDIF 703 710 ! 704 IF( wrk_not_released(2, 1) ) CALL ctl_stop('zgr_bat_ctl: failed to release workspace array') 711 CALL wrk_dealloc( jpi, jpj, zbathy ) 712 ! 713 IF( nn_timing == 1 ) CALL timing_stop('zgr_bat_ctl') 705 714 ! 706 715 END SUBROUTINE zgr_bat_ctl … … 719 728 !! (min value = 1 over land) 720 729 !!---------------------------------------------------------------------- 721 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released722 USE wrk_nemo, ONLY: zmbk => wrk_2d_1723 730 !! 724 731 INTEGER :: ji, jj ! dummy loop indices 725 !!---------------------------------------------------------------------- 726 ! 727 IF( wrk_in_use(2, 1) ) THEN 728 CALL ctl_stop('zgr_bot_level: requested 2D workspace unavailable') ; RETURN 729 ENDIF 732 REAL(wp), POINTER, DIMENSION(:,:) :: zmbk 733 !!---------------------------------------------------------------------- 734 ! 735 IF( nn_timing == 1 ) CALL timing_start('zgr_bot_level') 736 ! 737 CALL wrk_alloc( jpi, jpj, zmbk ) 730 738 ! 731 739 IF(lwp) WRITE(numout,*) … … 745 753 zmbk(:,:) = REAL( mbkv(:,:), wp ) ; CALL lbc_lnk(zmbk,'V',1.) ; mbkv (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 746 754 ! 747 IF( wrk_not_released(2, 1) ) CALL ctl_stop('zgr_bot_level: failed to release workspace array') 755 CALL wrk_dealloc( jpi, jpj, zmbk ) 756 ! 757 IF( nn_timing == 1 ) CALL timing_stop('zgr_bot_level') 748 758 ! 749 759 END SUBROUTINE zgr_bot_level … … 761 771 !!---------------------------------------------------------------------- 762 772 ! 773 IF( nn_timing == 1 ) CALL timing_start('zgr_zco') 774 ! 763 775 DO jk = 1, jpk 764 fsdept(:,:,jk) = gdept_0(jk) 765 fsdepw(:,:,jk) = gdepw_0(jk) 766 fsde3w(:,:,jk) = gdepw_0(jk) 767 fse3t (:,:,jk) = e3t_0(jk) 768 fse3u (:,:,jk) = e3t_0(jk) 769 fse3v (:,:,jk) = e3t_0(jk) 770 fse3f (:,:,jk) = e3t_0(jk) 771 fse3w (:,:,jk) = e3w_0(jk) 772 fse3uw(:,:,jk) = e3w_0(jk) 773 fse3vw(:,:,jk) = e3w_0(jk) 774 END DO 776 gdept(:,:,jk) = gdept_0(jk) 777 gdepw(:,:,jk) = gdepw_0(jk) 778 gdep3w(:,:,jk) = gdepw_0(jk) 779 e3t (:,:,jk) = e3t_0(jk) 780 e3u (:,:,jk) = e3t_0(jk) 781 e3v (:,:,jk) = e3t_0(jk) 782 e3f (:,:,jk) = e3t_0(jk) 783 e3w (:,:,jk) = e3w_0(jk) 784 e3uw(:,:,jk) = e3w_0(jk) 785 e3vw(:,:,jk) = e3w_0(jk) 786 END DO 787 ! 788 IF( nn_timing == 1 ) CALL timing_stop('zgr_zco') 775 789 ! 776 790 END SUBROUTINE zgr_zco … … 822 836 !! Reference : Pacanowsky & Gnanadesikan 1997, Mon. Wea. Rev., 126, 3248-3270. 823 837 !!---------------------------------------------------------------------- 824 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released825 USE wrk_nemo, ONLY: zprt => wrk_3d_1826 838 !! 827 839 INTEGER :: ji, jj, jk ! dummy loop indices … … 833 845 REAL(wp) :: zdiff ! temporary scalar 834 846 REAL(wp) :: zrefdep ! temporary scalar 847 REAL(wp), POINTER, DIMENSION(:,:,:) :: zprt 835 848 !!--------------------------------------------------------------------- 836 ! 837 IF( wrk_in_use(3, 1) ) THEN838 CALL ctl_stop('zgr_zps: requested workspace unavailable.') ; RETURN839 ENDIF840 849 ! 850 IF( nn_timing == 1 ) CALL timing_start('zgr_zps') 851 ! 852 CALL wrk_alloc( jpi, jpj, jpk, zprt ) 853 ! 841 854 IF(lwp) WRITE(numout,*) 842 855 IF(lwp) WRITE(numout,*) ' zgr_zps : z-coordinate with partial steps' … … 1017 1030 WRITE(numout,*) 'domzgr e3t(mbathy)' ; CALL prihre(zprt(:,:,1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1018 1031 WRITE(numout,*) 1019 WRITE(numout,*) 'domzgr e3w(mbathy)' ; CALL prihre(zprt(:,:, 1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1032 WRITE(numout,*) 'domzgr e3w(mbathy)' ; CALL prihre(zprt(:,:,2),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1020 1033 WRITE(numout,*) 1021 WRITE(numout,*) 'domzgr e3u(mbathy)' ; CALL prihre(zprt(:,:, 1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1034 WRITE(numout,*) 'domzgr e3u(mbathy)' ; CALL prihre(zprt(:,:,3),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1022 1035 WRITE(numout,*) 1023 WRITE(numout,*) 'domzgr e3v(mbathy)' ; CALL prihre(zprt(:,:, 1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1036 WRITE(numout,*) 'domzgr e3v(mbathy)' ; CALL prihre(zprt(:,:,4),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1024 1037 WRITE(numout,*) 1025 WRITE(numout,*) 'domzgr e3f(mbathy)' ; CALL prihre(zprt(:,:, 1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1038 WRITE(numout,*) 'domzgr e3f(mbathy)' ; CALL prihre(zprt(:,:,5),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1026 1039 WRITE(numout,*) 1027 WRITE(numout,*) 'domzgr gdep3w(mbathy)' ; CALL prihre(zprt(:,:, 1),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout)1040 WRITE(numout,*) 'domzgr gdep3w(mbathy)' ; CALL prihre(zprt(:,:,6),jpi,jpj,1,jpi,1,1,jpj,1,1.e-3,numout) 1028 1041 ENDIF 1029 1042 ! 1030 IF( wrk_not_released(3, 1) ) CALL ctl_stop('zgr_zps: failed to release workspace') 1043 CALL wrk_dealloc( jpi, jpj, jpk, zprt ) 1044 ! 1045 IF( nn_timing == 1 ) CALL timing_stop('zgr_zps') 1031 1046 ! 1032 1047 END SUBROUTINE zgr_zps … … 1116 1131 !! Reference : Madec, Lott, Delecluse and Crepon, 1996. JPO, 26, 1393-1408. 1117 1132 !!---------------------------------------------------------------------- 1118 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released1119 USE wrk_nemo, ONLY: zenv => wrk_2d_1 , ztmp => wrk_2d_2 , zmsk => wrk_2d_31120 USE wrk_nemo, ONLY: zri => wrk_2d_4 , zrj => wrk_2d_5 , zhbat => wrk_2d_61121 USE wrk_nemo, ONLY: gsigw3 => wrk_3d_11122 USE wrk_nemo, ONLY: gsigt3 => wrk_3d_21123 USE wrk_nemo, ONLY: gsi3w3 => wrk_3d_31124 USE wrk_nemo, ONLY: esigt3 => wrk_3d_41125 USE wrk_nemo, ONLY: esigw3 => wrk_3d_51126 USE wrk_nemo, ONLY: esigtu3 => wrk_3d_61127 USE wrk_nemo, ONLY: esigtv3 => wrk_3d_71128 USE wrk_nemo, ONLY: esigtf3 => wrk_3d_81129 USE wrk_nemo, ONLY: esigwu3 => wrk_3d_91130 USE wrk_nemo, ONLY: esigwv3 => wrk_3d_101131 1133 ! 1132 1134 INTEGER :: ji, jj, jk, jl ! dummy loop argument … … 1134 1136 REAL(wp) :: zcoeft, zcoefw, zrmax, ztaper ! temporary scalars 1135 1137 ! 1138 REAL(wp), POINTER, DIMENSION(:,: ) :: zenv, ztmp, zmsk, zri, zrj, zhbat 1139 REAL(wp), POINTER, DIMENSION(:,:,:) :: gsigw3, gsigt3, gsi3w3 1140 REAL(wp), POINTER, DIMENSION(:,:,:) :: esigt3, esigw3, esigtu3, esigtv3, esigtf3, esigwu3, esigwv3 1136 1141 1137 1142 NAMELIST/namzgr_sco/ rn_sbot_max, rn_sbot_min, rn_theta, rn_thetb, rn_rmax, ln_s_sigma, rn_bb, rn_hc 1138 1143 !!---------------------------------------------------------------------- 1139 1140 IF( wrk_in_use(2, 1,2,3,4,5,6) .OR. wrk_in_use(3, 1,2,3,4,5,6,7,8,9,10) ) THEN 1141 CALL ctl_stop('zgr_sco: ERROR - requested workspace arrays unavailable') ; RETURN 1142 ENDIF 1143 1144 ! 1145 IF( nn_timing == 1 ) CALL timing_start('zgr_sco') 1146 ! 1147 CALL wrk_alloc( jpi, jpj, zenv, ztmp, zmsk, zri, zrj, zhbat ) 1148 CALL wrk_alloc( jpi, jpj, jpk, gsigw3, gsigt3, gsi3w3 ) 1149 CALL wrk_alloc( jpi, jpj, jpk, esigt3, esigw3, esigtu3, esigtv3, esigtf3, esigwu3, esigwv3 ) 1150 ! 1144 1151 REWIND( numnam ) ! Read Namelist namzgr_sco : sigma-stretching parameters 1145 1152 READ ( numnam, namzgr_sco ) … … 1494 1501 1495 1502 ! 1496 !! H. Liu, POL. April 2009. Added for passing the scale check for the new released vvl code. 1503 where (e3t (:,:,:).eq.0.0) e3t(:,:,:) = 1.0 1504 where (e3u (:,:,:).eq.0.0) e3u(:,:,:) = 1.0 1505 where (e3v (:,:,:).eq.0.0) e3v(:,:,:) = 1.0 1506 where (e3f (:,:,:).eq.0.0) e3f(:,:,:) = 1.0 1507 where (e3w (:,:,:).eq.0.0) e3w(:,:,:) = 1.0 1508 where (e3uw (:,:,:).eq.0.0) e3uw(:,:,:) = 1.0 1509 where (e3vw (:,:,:).eq.0.0) e3vw(:,:,:) = 1.0 1510 1497 1511 1498 1512 fsdept(:,:,:) = gdept (:,:,:) … … 1590 1604 !!gm bug #endif 1591 1605 ! 1592 IF( wrk_not_released(2, 1,2,3,4,5,6) .OR. wrk_not_released(3, 1,2,3,4,5,6,7,8,9,10) ) & 1593 & CALL ctl_stop('dom:zgr_sco: failed to release workspace arrays') 1606 CALL wrk_dealloc( jpi, jpj, zenv, ztmp, zmsk, zri, zrj, zhbat ) 1607 CALL wrk_dealloc( jpi, jpj, jpk, gsigw3, gsigt3, gsi3w3 ) 1608 CALL wrk_dealloc( jpi, jpj, jpk, esigt3, esigw3, esigtu3, esigtv3, esigtf3, esigwu3, esigwv3 ) 1609 ! 1610 IF( nn_timing == 1 ) CALL timing_stop('zgr_sco') 1594 1611 ! 1595 1612 END SUBROUTINE zgr_sco -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/istate.F90
r2777 r3294 13 13 !! 2.0 ! 2006-07 (S. Masson) distributed restart using iom 14 14 !! 3.3 ! 2010-10 (C. Ethe) merge TRC-TRA 15 !! 3.4 ! 2011-04 (G. Madec) Merge of dtatem and dtasal & suppression of tb,tn/sb,sn 15 16 !!---------------------------------------------------------------------- 16 17 … … 30 31 USE zdf_oce ! ocean vertical physics 31 32 USE phycst ! physical constants 32 USE dtatem ! temperature data (dta_tem routine) 33 USE dtasal ! salinity data (dta_sal routine) 33 USE dtatsd ! data temperature and salinity (dta_tsd routine) 34 34 USE restart ! ocean restart (rst_read routine) 35 35 USE in_out_manager ! I/O manager … … 42 42 USE dynspg_exp ! pressure gradient schemes 43 43 USE dynspg_ts ! pressure gradient schemes 44 USE traswp ! Swap arrays (tra_swp routine)45 44 USE lib_mpp ! MPP library 45 USE wrk_nemo ! Memory allocation 46 USE timing ! Timing 46 47 47 48 IMPLICIT NONE … … 68 69 ! - ML - needed for initialization of e3t_b 69 70 INTEGER :: jk ! dummy loop indice 71 !!---------------------------------------------------------------------- 72 ! 73 IF( nn_timing == 1 ) CALL timing_start('istate_init') 74 ! 70 75 71 76 IF(lwp) WRITE(numout,*) … … 73 78 IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 74 79 75 rhd (:,:,:) = 0.e0 76 rhop (:,:,:) = 0.e0 77 rn2 (:,:,:) = 0.e0 78 ta (:,:,:) = 0.e0 79 sa (:,:,:) = 0.e0 80 CALL dta_tsd_init ! Initialisation of T & S input data 81 82 rhd (:,:,: ) = 0.e0 83 rhop (:,:,: ) = 0.e0 84 rn2 (:,:,: ) = 0.e0 85 tsa (:,:,:,:) = 0.e0 80 86 81 87 IF( ln_rstart ) THEN ! Restart from a file … … 83 89 neuler = 1 ! Set time-step indicator at nit000 (leap-frog) 84 90 CALL rst_read ! Read the restart file 85 CALL tra_swap ! swap 3D arrays (t,s) in a 4D array (ts) 91 ! ! define e3u_b, e3v_b from e3t_b read in restart file 92 CALL dom_vvl_2( nit000, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 86 93 CALL day_init ! model calendar (using both namelist and restart infos) 87 94 ELSE … … 92 99 CALL day_init ! model calendar (using both namelist and restart infos) 93 100 ! ! Initialization of ocean to zero 94 ! before fields ! now fields 95 sshb (:,:) = 0.e0 ; sshn (:,:) = 0.e0 96 ub (:,:,:) = 0.e0 ; un (:,:,:) = 0.e0 97 vb (:,:,:) = 0.e0 ; vn (:,:,:) = 0.e0 98 rotb (:,:,:) = 0.e0 ; rotn (:,:,:) = 0.e0 99 hdivb(:,:,:) = 0.e0 ; hdivn(:,:,:) = 0.e0 101 ! before fields ! now fields 102 sshb (:,:) = 0._wp ; sshn (:,:) = 0._wp 103 ub (:,:,:) = 0._wp ; un (:,:,:) = 0._wp 104 vb (:,:,:) = 0._wp ; vn (:,:,:) = 0._wp 105 rotb (:,:,:) = 0._wp ; rotn (:,:,:) = 0._wp 106 hdivb(:,:,:) = 0._wp ; hdivn(:,:,:) = 0._wp 107 ! 108 ! ! define e3u_b, e3v_b from e3t_b initialized in domzgr 109 CALL dom_vvl_2( nit000, fse3u_b(:,:,:), fse3v_b(:,:,:) ) 100 110 ! 101 111 IF( cp_cfg == 'eel' ) THEN … … 103 113 ELSEIF( cp_cfg == 'gyre' ) THEN 104 114 CALL istate_gyre ! GYRE configuration : start from pre-defined T-S fields 105 ELSE 106 ! ! Other configurations: Initial T-S fields 107 #if defined key_dtatem 108 CALL dta_tem( nit000 ) ! read 3D temperature data 109 tb(:,:,:) = t_dta(:,:,:) ; tn(:,:,:) = t_dta(:,:,:) 110 111 #else 112 IF(lwp) WRITE(numout,*) ! analytical temperature profile 113 IF(lwp) WRITE(numout,*)' Temperature initialization using an analytic profile' 114 CALL istate_tem 115 #endif 116 #if defined key_dtasal 117 CALL dta_sal( nit000 ) ! read 3D salinity data 118 sb(:,:,:) = s_dta(:,:,:) ; sn(:,:,:) = s_dta(:,:,:) 119 #else 120 ! No salinity data 121 IF(lwp)WRITE(numout,*) ! analytical salinity profile 122 IF(lwp)WRITE(numout,*)' Salinity initialisation using a constant value' 123 CALL istate_sal 124 #endif 115 ELSEIF( ln_tsd_init ) THEN ! Initial T-S fields read in files 116 CALL dta_tsd( nit000, tsb ) ! read 3D T and S data at nit000 117 tsn(:,:,:,:) = tsb(:,:,:,:) 118 ! 119 ELSE ! Initial T-S fields defined analytically 120 CALL istate_t_s 125 121 ENDIF 126 122 ! 127 CALL tra_swap ! swap 3D arrays (tb,sb,tn,sn) in a 4D array128 123 CALL eos( tsb, rhd, rhop ) ! before potential and in situ densities 129 124 #if ! defined key_c1d … … 148 143 ENDIF 149 144 ! 145 IF( nn_timing == 1 ) CALL timing_stop('istate_init') 146 ! 150 147 END SUBROUTINE istate_init 151 148 152 153 SUBROUTINE istate_tem 149 SUBROUTINE istate_t_s 154 150 !!--------------------------------------------------------------------- 155 !! *** ROUTINE istate_t em***151 !! *** ROUTINE istate_t_s *** 156 152 !! 157 153 !! ** Purpose : Intialization of the temperature field with an 158 154 !! analytical profile or a file (i.e. in EEL configuration) 159 155 !! 160 !! ** Method : Use Philander analytic profile of temperature 156 !! ** Method : - temperature: use Philander analytic profile 157 !! - salinity : use to a constant value 35.5 161 158 !! 162 159 !! References : Philander ??? 163 160 !!---------------------------------------------------------------------- 164 INTEGER :: ji, jj, jk 161 INTEGER :: ji, jj, jk 162 REAL(wp) :: zsal = 35.50 165 163 !!---------------------------------------------------------------------- 166 164 ! 167 165 IF(lwp) WRITE(numout,*) 168 IF(lwp) WRITE(numout,*) 'istate_t em :initial temperature profile'169 IF(lwp) WRITE(numout,*) '~~~~~~~~~~ '166 IF(lwp) WRITE(numout,*) 'istate_t_s : Philander s initial temperature profile' 167 IF(lwp) WRITE(numout,*) '~~~~~~~~~~ and constant salinity (',zsal,' psu)' 170 168 ! 171 169 DO jk = 1, jpk 172 DO jj = 1, jpj 173 DO ji = 1, jpi 174 tn(ji,jj,jk) = ( ( ( 7.5 - 0.*ABS(gphit(ji,jj))/30. ) & 175 & *( 1.-TANH((fsdept(ji,jj,jk)-80.)/30.) ) & 176 & + 10.*(5000.-fsdept(ji,jj,jk))/5000.) ) * tmask(ji,jj,jk) 177 tb(ji,jj,jk) = tn(ji,jj,jk) 178 END DO 179 END DO 170 tsn(:,:,jk,jp_tem) = ( ( ( 7.5 - 0. * ABS( gphit(:,:) )/30. ) * ( 1.-TANH((fsdept(:,:,jk)-80.)/30.) ) & 171 & + 10. * ( 5000. - fsdept(:,:,jk) ) /5000.) ) * tmask(:,:,jk) 172 tsb(:,:,jk,jp_tem) = tsn(:,:,jk,jp_tem) 180 173 END DO 181 ! 182 IF(lwp) CALL prizre( tn , jpi , jpj , jpk , jpj/2 , & 183 & 1 , jpi , 5 , 1 , jpk , & 184 & 1 , 1. , numout ) 185 ! 186 END SUBROUTINE istate_tem 187 188 189 SUBROUTINE istate_sal 190 !!--------------------------------------------------------------------- 191 !! *** ROUTINE istate_sal *** 192 !! 193 !! ** Purpose : Intialize the salinity field with an analytic profile 194 !! 195 !! ** Method : Use to a constant value 35.5 196 !! 197 !! ** Action : Initialize sn and sb 198 !!---------------------------------------------------------------------- 199 REAL(wp) :: zsal = 35.50_wp 200 !!---------------------------------------------------------------------- 201 ! 202 IF(lwp) WRITE(numout,*) 203 IF(lwp) WRITE(numout,*) 'istate_sal : initial salinity : ', zsal 204 IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 205 ! 206 sn(:,:,:) = zsal * tmask(:,:,:) 207 sb(:,:,:) = sn(:,:,:) 208 ! 209 END SUBROUTINE istate_sal 174 tsn(:,:,:,jp_sal) = zsal * tmask(:,:,:) 175 tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 176 ! 177 END SUBROUTINE istate_t_s 210 178 211 179 … … 254 222 ! 255 223 DO jk = 1, jpk 256 t n(:,:,jk) = ( zt2 + zt1 * exp( - fsdept(:,:,jk) / 1000 ) ) * tmask(:,:,jk)257 t b(:,:,jk) = tn(:,:,jk)224 tsn(:,:,jk,jp_tem) = ( zt2 + zt1 * exp( - fsdept(:,:,jk) / 1000 ) ) * tmask(:,:,jk) 225 tsb(:,:,jk,jp_tem) = tsn(:,:,jk,jp_tem) 258 226 END DO 259 227 ! 260 IF(lwp) CALL prizre( t n, jpi , jpj , jpk , jpj/2 , &261 & 1 , jpi , 5 , 1 , jpk , &262 & 1 , 1. , numout )228 IF(lwp) CALL prizre( tsn(:,:,:,jp_tem), jpi , jpj , jpk , jpj/2 , & 229 & 1 , jpi , 5 , 1 , jpk , & 230 & 1 , 1. , numout ) 263 231 ! 264 232 ! set salinity field to a constant value … … 268 236 IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 269 237 ! 270 sn(:,:,:) = zsal * tmask(:,:,:)271 sb(:,:,:) = sn(:,:,:)238 tsn(:,:,:,jp_sal) = zsal * tmask(:,:,:) 239 tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 272 240 ! 273 241 ! set the dynamics: U,V, hdiv, rot (and ssh if necessary) … … 323 291 ! 324 292 CALL iom_open ( 'eel.initemp', inum ) 325 CALL iom_get ( inum, jpdom_data, 'initemp', t b) ! read before temprature (tb)293 CALL iom_get ( inum, jpdom_data, 'initemp', tsb(:,:,:,jp_tem) ) ! read before temprature (tb) 326 294 CALL iom_close( inum ) 327 295 ! 328 t n(:,:,:) = tb(:,:,:) ! set nox temperature to tb329 ! 330 IF(lwp) CALL prizre( t n, jpi , jpj , jpk , jpj/2 , &331 & 1 , jpi , 5 , 1 , jpk , &332 & 1 , 1. , numout )296 tsn(:,:,:,jp_tem) = tsb(:,:,:,jp_tem) ! set nox temperature to tb 297 ! 298 IF(lwp) CALL prizre( tsn(:,:,:,jp_tem), jpi , jpj , jpk , jpj/2 , & 299 & 1 , jpi , 5 , 1 , jpk , & 300 & 1 , 1. , numout ) 333 301 ! 334 302 ! set salinity field to a constant value … … 338 306 IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 339 307 ! 340 sn(:,:,:) = zsal * tmask(:,:,:)341 sb(:,:,:) = sn(:,:,:)308 tsn(:,:,:,jp_sal) = zsal * tmask(:,:,:) 309 tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 342 310 ! 343 311 ! ! =========================== … … 377 345 DO jj = 1, jpj 378 346 DO ji = 1, jpi 379 t n(ji,jj,jk) = ( 16. - 12. * TANH( (fsdept(ji,jj,jk) - 400) / 700 ) ) &347 tsn(ji,jj,jk,jp_tem) = ( 16. - 12. * TANH( (fsdept(ji,jj,jk) - 400) / 700 ) ) & 380 348 & * (-TANH( (500-fsdept(ji,jj,jk)) / 150 ) + 1) / 2 & 381 349 & + ( 15. * ( 1. - TANH( (fsdept(ji,jj,jk)-50.) / 1500.) ) & … … 383 351 & + 7. * (1500. - fsdept(ji,jj,jk)) / 1500. ) & 384 352 & * (-TANH( (fsdept(ji,jj,jk) - 500) / 150) + 1) / 2 385 t n(ji,jj,jk) = tn(ji,jj,jk) * tmask(ji,jj,jk)386 t b(ji,jj,jk) = tn(ji,jj,jk)387 388 sn(ji,jj,jk) = ( 36.25 - 1.13 * TANH( (fsdept(ji,jj,jk) - 305) / 460 ) ) &353 tsn(ji,jj,jk,jp_tem) = tsn(ji,jj,jk,jp_tem) * tmask(ji,jj,jk) 354 tsb(ji,jj,jk,jp_tem) = tsn(ji,jj,jk,jp_tem) 355 356 tsn(ji,jj,jk,jp_sal) = ( 36.25 - 1.13 * TANH( (fsdept(ji,jj,jk) - 305) / 460 ) ) & 389 357 & * (-TANH((500 - fsdept(ji,jj,jk)) / 150) + 1) / 2 & 390 358 & + ( 35.55 + 1.25 * (5000. - fsdept(ji,jj,jk)) / 5000. & … … 393 361 & + 0.2 * TANH( (fsdept(ji,jj,jk) - 1000.) / 5000.) ) & 394 362 & * (-TANH((fsdept(ji,jj,jk) - 500) / 150) + 1) / 2 395 sn(ji,jj,jk) = sn(ji,jj,jk) * tmask(ji,jj,jk)396 sb(ji,jj,jk) = sn(ji,jj,jk)363 tsn(ji,jj,jk,jp_sal) = tsn(ji,jj,jk,jp_sal) * tmask(ji,jj,jk) 364 tsb(ji,jj,jk,jp_sal) = tsn(ji,jj,jk,jp_sal) 397 365 END DO 398 366 END DO … … 408 376 ! ---------------------- 409 377 CALL iom_open ( 'data_tem', inum ) 410 CALL iom_get ( inum, jpdom_data, 'votemper', t n)378 CALL iom_get ( inum, jpdom_data, 'votemper', tsn(:,:,:,jp_tem) ) 411 379 CALL iom_close( inum ) 412 380 413 t n(:,:,:) = tn(:,:,:) * tmask(:,:,:)414 t b(:,:,:) = tn(:,:,:)381 tsn(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) * tmask(:,:,:) 382 tsb(:,:,:,jp_tem) = tsn(:,:,:,jp_tem) 415 383 416 384 ! Read salinity field 417 385 ! ------------------- 418 386 CALL iom_open ( 'data_sal', inum ) 419 CALL iom_get ( inum, jpdom_data, 'vosaline', sn)387 CALL iom_get ( inum, jpdom_data, 'vosaline', tsn(:,:,:,jp_sal) ) 420 388 CALL iom_close( inum ) 421 389 422 sn(:,:,:) = sn(:,:,:) * tmask(:,:,:)423 sb(:,:,:) = sn(:,:,:)390 tsn(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) * tmask(:,:,:) 391 tsb(:,:,:,jp_sal) = tsn(:,:,:,jp_sal) 424 392 425 393 END SELECT … … 429 397 WRITE(numout,*) ' Initial temperature and salinity profiles:' 430 398 WRITE(numout, "(9x,' level gdept_0 temperature salinity ')" ) 431 WRITE(numout, "(10x, i4, 3f10.2)" ) ( jk, gdept_0(jk), t n(2,2,jk), sn(2,2,jk), jk = 1, jpk )399 WRITE(numout, "(10x, i4, 3f10.2)" ) ( jk, gdept_0(jk), tsn(2,2,jk,jp_tem), tsn(2,2,jk,jp_sal), jk = 1, jpk ) 432 400 ENDIF 433 401 … … 446 414 !! p=integral [ rau*g dz ] 447 415 !!---------------------------------------------------------------------- 448 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released449 USE wrk_nemo, ONLY: zprn => wrk_3d_1 ! 3D workspace450 451 416 USE dynspg ! surface pressure gradient (dyn_spg routine) 452 417 USE divcur ! hor. divergence & rel. vorticity (div_cur routine) … … 456 421 INTEGER :: indic ! ??? 457 422 REAL(wp) :: zmsv, zphv, zmsu, zphu, zalfg ! temporary scalars 458 !!---------------------------------------------------------------------- 459 460 IF(wrk_in_use(3, 1) ) THEN 461 CALL ctl_stop('istate_uvg: requested workspace array unavailable') ; RETURN 462 ENDIF 463 423 REAL(wp), POINTER, DIMENSION(:,:,:) :: zprn 424 !!---------------------------------------------------------------------- 425 ! 426 CALL wrk_alloc( jpi, jpj, jpk, zprn) 427 ! 464 428 IF(lwp) WRITE(numout,*) 465 429 IF(lwp) WRITE(numout,*) 'istate_uvg : Start from Geostrophy' … … 557 521 rotb (:,:,:) = rotn (:,:,:) ! set the before to the now value 558 522 ! 559 IF( wrk_not_released(3, 1) ) THEN 560 CALL ctl_stop('istate_uvg: failed to release workspace array') 561 ENDIF 523 CALL wrk_dealloc( jpi, jpj, jpk, zprn) 562 524 ! 563 525 END SUBROUTINE istate_uvg -
trunk/NEMOGCM/NEMO/OPA_SRC/DOM/phycst.F90
r2528 r3294 9 9 !! - ! 2006-08 (G. Madec) style 10 10 !! 3.2 ! 2006-08 (S. Masson, G. Madec) suppress useless variables + style 11 !! 3.4 ! 2011-11 (C. Harris) minor changes for CICE constants 11 12 !!---------------------------------------------------------------------- 12 13 … … 48 49 #endif 49 50 51 #if defined key_cice 52 REAL(wp), PUBLIC :: rau0 = 1026._wp !: reference volumic mass (density) (kg/m3) 53 #else 50 54 REAL(wp), PUBLIC :: rau0 = 1035._wp !: reference volumic mass (density) (kg/m3) 55 #endif 51 56 REAL(wp), PUBLIC :: rau0r !: reference specific volume (m3/kg) 52 57 REAL(wp), PUBLIC :: rcp = 4.e+3_wp !: ocean specific heat 53 58 REAL(wp), PUBLIC :: ro0cpr !: = 1. / ( rau0 * rcp ) 54 59 55 #if defined key_lim3 60 #if defined key_lim3 || defined key_cice 56 61 REAL(wp), PUBLIC :: rcdsn = 0.31_wp !: thermal conductivity of snow 57 62 REAL(wp), PUBLIC :: rcdic = 2.034396_wp !: thermal conductivity of fresh ice … … 100 105 rsiyea = 365.25 * rday * 2. * rpi / 6.283076 101 106 rsiday = rday / ( 1. + rday / rsiyea ) 107 #if defined key_cice 108 omega = 7.292116e-05 109 #else 102 110 omega = 2. * rpi / rsiday 111 #endif 103 112 104 113 rau0r = 1. / rau0 … … 155 164 WRITE(numout,*) ' latent heat of fusion of fresh ice / snow = ', lfus , ' J/kg' 156 165 WRITE(numout,*) ' latent heat of subl. of fresh ice / snow = ', lsub , ' J/kg' 166 #elif defined key_cice 167 WRITE(numout,*) ' latent heat of fusion of fresh ice / snow = ', lfus , ' J/kg' 157 168 #else 158 169 WRITE(numout,*) ' density times specific heat for snow = ', rcpsn , ' J/m^3/K'
Note: See TracChangeset
for help on using the changeset viewer.