Changeset 2887 for branches/2011/dev_r2802_NOCL_Smagorinsky/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilapg.F90
- Timestamp:
- 2011-10-05T11:16:29+02:00 (13 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2011/dev_r2802_NOCL_Smagorinsky/NEMOGCM/NEMO/OPA_SRC/DYN/dynldf_bilapg.F90
r2715 r2887 4 4 !! Ocean dynamics: lateral viscosity trend 5 5 !!====================================================================== 6 !! History : OPA ! 1997-07 (G. Madec) Original code7 !! NEMO 1.0 ! 2002-08 (G. Madec) F90: Free form and module8 !! 2.0 ! 2004-08 (C. Talandier) New trends organization9 !!----------------------------------------------------------------------10 6 #if defined key_ldfslp || defined key_esopa 11 7 !!---------------------------------------------------------------------- … … 16 12 !! ldfguv : 17 13 !!---------------------------------------------------------------------- 14 !! * Modules used 18 15 USE oce ! ocean dynamics and tracers 19 16 USE dom_oce ! ocean space and time domain 20 17 USE ldfdyn_oce ! ocean dynamics lateral physics 21 18 USE zdf_oce ! ocean vertical physics 19 USE in_out_manager ! I/O manager 22 20 USE trdmod ! ocean dynamics trends 23 21 USE trdmod_oce ! ocean variables trends 24 22 USE ldfslp ! iso-neutral slopes available 25 USE in_out_manager ! I/O manager26 USE lib_mpp ! MPP library27 23 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 28 24 USE prtctl ! Print control … … 31 27 PRIVATE 32 28 33 PUBLIC dyn_ldf_bilapg ! called by step.F90 34 35 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zfuw, zfvw , zdiu, zdiv ! 2D workspace (ldfguv) 36 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: zdju, zdj1u, zdjv, zdj1v ! 2D workspace (ldfguv) 29 !! * Routine accessibility 30 PUBLIC dyn_ldf_bilapg ! called by step.F90 37 31 38 32 !! * Substitutions … … 42 36 !! NEMO/OPA 3.3 , NEMO Consortium (2010) 43 37 !! $Id$ 44 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 45 !!---------------------------------------------------------------------- 38 !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 39 !!---------------------------------------------------------------------- 40 46 41 CONTAINS 47 48 INTEGER FUNCTION dyn_ldf_bilapg_alloc()49 !!----------------------------------------------------------------------50 !! *** ROUTINE dyn_ldf_bilapg_alloc ***51 !!----------------------------------------------------------------------52 ALLOCATE( zfuw(jpi,jpk) , zfvw (jpi,jpk) , zdiu(jpi,jpk) , zdiv (jpi,jpk) , &53 & zdju(jpi,jpk) , zdj1u(jpi,jpk) , zdjv(jpi,jpk) , zdj1v(jpi,jpk) , STAT=dyn_ldf_bilapg_alloc )54 !55 IF( dyn_ldf_bilapg_alloc /= 0 ) CALL ctl_warn('dyn_ldf_bilapg_alloc: failed to allocate arrays')56 END FUNCTION dyn_ldf_bilapg_alloc57 58 42 59 43 SUBROUTINE dyn_ldf_bilapg( kt ) … … 83 67 !! biharmonic mixing trend. 84 68 !! - save the trend in (zwk3,zwk4) ('key_trddyn') 69 !! 70 !! History : 71 !! 8.0 ! 97-07 (G. Madec) Original code 72 !! 8.5 ! 02-08 (G. Madec) F90: Free form and module 73 !! 9.0 ! 04-08 (C. Talandier) New trends organization 85 74 !!---------------------------------------------------------------------- 86 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 87 USE wrk_nemo, ONLY: zwk1 => wrk_3d_3 , zwk2 => wrk_3d_4 ! 3D workspace 88 USE oce , ONLY: zwk3 => ta , zwk4 => sa ! ta, sa used as 3D workspace 89 ! 75 !! * Modules used 76 USE oce, ONLY : zwk3 => ta, & ! use ta as 3D workspace 77 zwk4 => sa ! use sa as 3D workspace 78 79 !! * Arguments 90 80 INTEGER, INTENT( in ) :: kt ! ocean time-step index 91 ! 81 82 !! * Local declarations 92 83 INTEGER :: ji, jj, jk ! dummy loop indices 84 REAL(wp), DIMENSION(jpi,jpj,jpk) :: & 85 zwk1, zwk2 ! work array used for rotated biharmonic 86 ! ! operator on tracers and/or momentum 93 87 !!---------------------------------------------------------------------- 94 95 IF( wrk_in_use(3, 3,4) ) THEN96 CALL ctl_stop('dyn_ldf_bilapg: requested workspace arrays unavailable') ; RETURN97 ENDIF98 88 99 89 IF( kt == nit000 ) THEN … … 103 93 zwk1(:,:,:) = 0.e0 ; zwk3(:,:,:) = 0.e0 104 94 zwk2(:,:,:) = 0.e0 ; zwk4(:,:,:) = 0.e0 105 ! ! allocate dyn_ldf_bilapg arrays106 IF( dyn_ldf_bilapg_alloc() /= 0 ) CALL ctl_stop('STOP', 'dyn_ldf_bilapg: failed to allocate arrays')107 95 ENDIF 108 96 109 97 ! Laplacian of (ub,vb) multiplied by ahm 110 98 ! -------------------------------------- 111 CALL ldfguv( ub, vb, zwk1, zwk2, 1 ) ! rotated harmonic operator applied to (ub,vb) 112 ! ! and multiply by ahmu, ahmv (output in (zwk1,zwk2) ) 113 CALL lbc_lnk( zwk1, 'U', -1. ) ; CALL lbc_lnk( zwk2, 'V', -1. ) ! Lateral boundary conditions 99 ! rotated harmonic operator applied to (ub,vb) 100 ! and multiply by ahmu, ahmv (output in (zwk1,zwk2) ) 101 102 CALL ldfguv ( ub, vb, zwk1, zwk2, 1 ) 103 104 105 ! Lateral boundary conditions on (zwk1,zwk2) 106 CALL lbc_lnk( zwk1, 'U', -1. ) 107 CALL lbc_lnk( zwk2, 'V', -1. ) 108 114 109 115 110 ! Bilaplacian of (ub,vb) 116 111 ! ---------------------- 117 CALL ldfguv( zwk1, zwk2, zwk3, zwk4, 2 ) ! rotated harmonic operator applied to (zwk1,zwk2) 118 ! ! (output in (zwk3,zwk4) ) 119 120 ! Update the momentum trends 112 ! rotated harmonic operator applied to (zwk1,zwk2) (output in (zwk3,zwk4) ) 113 114 CALL ldfguv ( zwk1, zwk2, zwk3, zwk4, 2 ) 115 116 117 ! Update the momentum trends (j-slab : 2, jpj-1) 121 118 ! -------------------------- 122 DO jj = 2, jpjm1 ! add the diffusive trend to the general momentum trends 119 ! ! =============== 120 DO jj = 2, jpjm1 ! Vertical slab 121 ! ! =============== 123 122 DO jk = 1, jpkm1 124 123 DO ji = 2, jpim1 124 ! add the diffusive trend to the general momentum trends 125 125 ua(ji,jj,jk) = ua(ji,jj,jk) + zwk3(ji,jj,jk) 126 126 va(ji,jj,jk) = va(ji,jj,jk) + zwk4(ji,jj,jk) 127 127 END DO 128 128 END DO 129 END DO130 !131 IF( wrk_not_released(3, 3,4) ) CALL ctl_stop('dyn_ldf_bilapg: failed to release workspace arrays')132 ! 129 ! ! =============== 130 END DO ! End of slab 131 ! ! =============== 132 133 133 END SUBROUTINE dyn_ldf_bilapg 134 134 … … 174 174 !! second order vertical derivative term) 175 175 !! 'key_trddyn' defined: the trend is saved for diagnostics. 176 !! 177 !! History : 178 !! 8.0 ! 97-07 (G. Madec) Original code 179 !! 8.5 ! 02-08 (G. Madec) F90: Free form and module 176 180 !!---------------------------------------------------------------------- 177 USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 178 USE wrk_nemo, ONLY: ziut => wrk_2d_1 , zjuf => wrk_2d_2 , zjvt => wrk_2d_3 179 USE wrk_nemo, ONLY: zivf => wrk_2d_4 , zdku => wrk_2d_5 , zdk1u => wrk_2d_6 180 USE wrk_nemo, ONLY: zdkv => wrk_2d_7 , zdk1v => wrk_2d_8 181 !! 182 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT(in ) :: pu , pv ! 1st call: before horizontal velocity 183 ! ! 2nd call: ahm x these fields 184 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( out) :: plu, plv ! partial harmonic operator applied to 185 ! ! pu and pv (all the components except 186 ! ! second order vertical derivative term) 187 INTEGER , INTENT(in ) :: kahm ! =1 1st call ; =2 2nd call 188 ! 189 INTEGER :: ji, jj, jk ! dummy loop indices 190 REAL(wp) :: zabe1 , zabe2 , zcof1 , zcof2 ! local scalar 191 REAL(wp) :: zcoef0, zcoef3, zcoef4 ! - - 192 REAL(wp) :: zbur, zbvr, zmkt, zmkf, zuav, zvav ! - - 193 REAL(wp) :: zuwslpi, zuwslpj, zvwslpi, zvwslpj ! - - 181 !! * Arguments 182 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( in ) :: & 183 pu, pv ! momentum fields (before u and v for the 1st call, and 184 ! ! laplacian of these fields multiplied by ahm for the 2nd 185 REAL(wp), DIMENSION(jpi,jpj,jpk), INTENT( out ) :: & 186 plu, plv ! partial harmonic operator applied to 187 ! ! pu and pv (all the components except 188 ! ! second order vertical derivative term) 189 INTEGER, INTENT( in ) :: & 190 kahm ! =1 the laplacian is multiplied by the eddy diffusivity coef. 191 ! ! =2 no multiplication 192 193 !! * Local declarations 194 INTEGER :: ji, jj, jk ! dummy loop indices 195 REAL(wp) :: & 196 zabe1, zabe2, zcof1, zcof2, & ! temporary scalars 197 zcoef0, zcoef3, zcoef4 198 REAL(wp) :: & 199 zbur, zbvr, zmkt, zmkf, zuav, zvav, & 200 zuwslpi, zuwslpj, zvwslpi, zvwslpj 201 REAL(wp), DIMENSION(jpi,jpj) :: & 202 ziut, zjuf , zjvt, zivf, & ! workspace 203 zdku, zdk1u, zdkv, zdk1v 204 REAL(wp), DIMENSION(jpi,jpk) :: & 205 zfuw, zfvw, zdiu, zdiv, & ! workspace 206 zdju, zdj1u, zdjv, zdj1v 194 207 !!---------------------------------------------------------------------- 195 208 196 IF( wrk_in_use(2, 1,2,3,4,5,6,7,8) ) THEN197 CALL ctl_stop('dyn:ldfguv: requested workspace arrays unavailable') ; RETURN198 END IF199 209 ! ! ********** ! ! =============== 200 210 DO jk = 1, jpkm1 ! First step ! ! Horizontal slab … … 412 422 ! II.3 Divergence of vertical fluxes added to the horizontal divergence 413 423 ! --------------------------------------------------------------------- 414 424 #if defined key_dynldf_smag 415 425 IF( kahm == 1 ) THEN 416 426 ! multiply the laplacian by the eddy viscosity coefficient … … 448 458 STOP 'ldfguv' 449 459 ENDIF 460 #else 461 IF( kahm == 2 ) THEN 462 ! multiply the laplacian by the eddy viscosity coefficient 463 DO jk = 1, jpkm1 464 DO ji = 2, jpim1 465 ! eddy coef. divided by the volume element 466 zbur = fsahmu(ji,jj,jk) / ( e1u(ji,jj)*e2u(ji,jj)*fse3u(ji,jj,jk) ) 467 zbvr = fsahmv(ji,jj,jk) / ( e1v(ji,jj)*e2v(ji,jj)*fse3v(ji,jj,jk) ) 468 ! vertical divergence 469 zuav = zfuw(ji,jk) - zfuw(ji,jk+1) 470 zvav = zfvw(ji,jk) - zfvw(ji,jk+1) 471 ! harmonic operator applied to (pu,pv) and multiply by ahm 472 plu(ji,jj,jk) = ( plu(ji,jj,jk) + zuav ) * zbur 473 plv(ji,jj,jk) = ( plv(ji,jj,jk) + zvav ) * zbvr 474 END DO 475 END DO 476 ELSEIF( kahm == 1 ) THEN 477 ! second call, no multiplication 478 DO jk = 1, jpkm1 479 DO ji = 2, jpim1 480 ! inverse of the volume element 481 zbur = 1. / ( e1u(ji,jj)*e2u(ji,jj)*fse3u(ji,jj,jk) ) 482 zbvr = 1. / ( e1v(ji,jj)*e2v(ji,jj)*fse3v(ji,jj,jk) ) 483 ! vertical divergence 484 zuav = zfuw(ji,jk) - zfuw(ji,jk+1) 485 zvav = zfvw(ji,jk) - zfvw(ji,jk+1) 486 ! harmonic operator applied to (pu,pv) 487 plu(ji,jj,jk) = ( plu(ji,jj,jk) + zuav ) * zbur 488 plv(ji,jj,jk) = ( plv(ji,jj,jk) + zvav ) * zbvr 489 END DO 490 END DO 491 ELSE 492 IF(lwp)WRITE(numout,*) ' ldfguv: kahm= 1 or 2, here =', kahm 493 IF(lwp)WRITE(numout,*) ' We stop' 494 STOP 'ldfguv' 495 ENDIF 496 497 #endif 450 498 ! ! =============== 451 499 END DO ! End of slab 452 500 ! ! =============== 453 454 IF( wrk_not_released(2, 1,2,3,4,5,6,7,8) ) CALL ctl_stop('dyn:ldfguv: failed to release workspace arrays')455 !456 501 END SUBROUTINE ldfguv 457 502 … … 462 507 CONTAINS 463 508 SUBROUTINE dyn_ldf_bilapg( kt ) ! Dummy routine 464 INTEGER, INTENT(in) :: kt465 509 WRITE(*,*) 'dyn_ldf_bilapg: You should not have seen this print! error?', kt 466 510 END SUBROUTINE dyn_ldf_bilapg
Note: See TracChangeset
for help on using the changeset viewer.
