Changeset 10822
- Timestamp:
- 2019-04-01T17:50:07+02:00 (5 years ago)
- Location:
- NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE
- Files:
-
- 4 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE/BDY/bdytides.F90
r10811 r10822 115 115 WRITE(numout,*) ' Tidal components: ' 116 116 DO itide = 1, nb_harmo 117 WRITE(numout,*) ' ', tide_ components(itide)%cname_tide117 WRITE(numout,*) ' ', tide_harmonics(itide)%cname_tide 118 118 END DO 119 119 ENDIF … … 156 156 igrd = 1 ! Everything is at T-points here 157 157 DO itide = 1, nb_harmo 158 CALL iom_get( inum, jpdom_autoglo, TRIM(tide_ components(itide)%cname_tide)//'_z1', ztr(:,:) )159 CALL iom_get( inum, jpdom_autoglo, TRIM(tide_ components(itide)%cname_tide)//'_z2', zti(:,:) )158 CALL iom_get( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_z1', ztr(:,:) ) 159 CALL iom_get( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_z2', zti(:,:) ) 160 160 DO ib = 1, ilen0(igrd) 161 161 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) … … 172 172 igrd = 2 ! Everything is at U-points here 173 173 DO itide = 1, nb_harmo 174 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_ components(itide)%cname_tide)//'_u1', ztr(:,:) )175 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_ components(itide)%cname_tide)//'_u2', zti(:,:) )174 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_u1', ztr(:,:) ) 175 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_u2', zti(:,:) ) 176 176 DO ib = 1, ilen0(igrd) 177 177 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) … … 188 188 igrd = 3 ! Everything is at V-points here 189 189 DO itide = 1, nb_harmo 190 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_ components(itide)%cname_tide)//'_v1', ztr(:,:) )191 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_ components(itide)%cname_tide)//'_v2', zti(:,:) )190 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_v1', ztr(:,:) ) 191 CALL iom_get ( inum, jpdom_autoglo, TRIM(tide_harmonics(itide)%cname_tide)//'_v2', zti(:,:) ) 192 192 DO ib = 1, ilen0(igrd) 193 193 ii = idx_bdy(ib_bdy)%nbi(ib,igrd) … … 217 217 DO itide = 1, nb_harmo 218 218 ! ! SSH fields 219 clfile = TRIM(filtide)//TRIM(tide_ components(itide)%cname_tide)//'_grid_T.nc'219 clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_T.nc' 220 220 CALL iom_open( clfile, inum ) 221 221 CALL fld_map( inum, 'z1' , dta_read(1:ilen0(1),1:1,1:1) , 1, ibmap_ptr(1) ) … … 225 225 CALL iom_close( inum ) 226 226 ! ! U fields 227 clfile = TRIM(filtide)//TRIM(tide_ components(itide)%cname_tide)//'_grid_U.nc'227 clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_U.nc' 228 228 CALL iom_open( clfile, inum ) 229 229 CALL fld_map( inum, 'u1' , dta_read(1:ilen0(2),1:1,1:1) , 1, ibmap_ptr(2) ) … … 233 233 CALL iom_close( inum ) 234 234 ! ! V fields 235 clfile = TRIM(filtide)//TRIM(tide_ components(itide)%cname_tide)//'_grid_V.nc'235 clfile = TRIM(filtide)//TRIM(tide_harmonics(itide)%cname_tide)//'_grid_V.nc' 236 236 CALL iom_open( clfile, inum ) 237 237 CALL fld_map( inum, 'v1' , dta_read(1:ilen0(3),1:1,1:1) , 1, ibmap_ptr(3) ) … … 335 335 336 336 DO itide = 1, nb_harmo 337 z_sarg = z_arg * omega_tide(itide)337 z_sarg = z_arg * tide_harmonics(itide)%omega 338 338 z_cost(itide) = COS( z_sarg ) 339 339 z_sist(itide) = SIN( z_sarg ) … … 440 440 DO itide = 1, nb_harmo 441 441 ! 442 z_sarg = (z_arg + zoff) * omega_tide(itide)442 z_sarg = (z_arg + zoff) * tide_harmonics(itide)%omega 443 443 z_cost = zramp * COS( z_sarg ) 444 444 z_sist = zramp * SIN( z_sarg ) … … 501 501 END DO 502 502 DO ib = 1 , ilen0(igrd) 503 mod_tide(ib)=mod_tide(ib)* ftide(itide)504 phi_tide(ib)=phi_tide(ib)+ v0tide(itide)+utide(itide)503 mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f 504 phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0+tide_harmonics(itide)%u 505 505 ENDDO 506 506 DO ib = 1 , ilen0(igrd) … … 540 540 END DO 541 541 DO ib = 1, ilen0(igrd) 542 mod_tide(ib)=mod_tide(ib)* ftide(itide)543 phi_tide(ib)=phi_tide(ib)+ v0tide(itide)+utide(itide)542 mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f 543 phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0 + tide_harmonics(itide)%u 544 544 ENDDO 545 545 DO ib = 1, ilen0(igrd) … … 561 561 END DO 562 562 DO ib = 1, ilen0(igrd) 563 mod_tide(ib)=mod_tide(ib)* ftide(itide)564 phi_tide(ib)=phi_tide(ib)+ v0tide(itide)+utide(itide)563 mod_tide(ib)=mod_tide(ib)*tide_harmonics(itide)%f 564 phi_tide(ib)=phi_tide(ib)+tide_harmonics(itide)%v0 + tide_harmonics(itide)%u 565 565 ENDDO 566 566 DO ib = 1, ilen0(igrd) -
NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE/DIA/diaharm.F90
r10811 r10822 39 39 40 40 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) :: ana_temp 41 REAL(wp), ALLOCATABLE, DIMENSION(:) :: ana_freq, ut , vt , ft42 41 REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) :: out_eta , out_u, out_v 43 42 … … 73 72 INTEGER :: ios ! Local integer output status for namelist read 74 73 TYPE(tide), PUBLIC, DIMENSION(:), POINTER :: tide_components ! Selected tidal components 74 TYPE(tide), DIMENSION(:), POINTER :: tide_harmonics ! Oscillation parameters of selected tidal components 75 75 76 76 NAMELIST/nam_diaharm/ nit000_han, nitend_han, nstep_han, tname … … 127 127 ENDIF 128 128 129 ! Initialize frequency array: 130 ! --------------------------- 131 ALLOCATE( ana_freq(nb_ana), ut(nb_ana), vt(nb_ana), ft(nb_ana) ) 132 133 CALL tide_harmo( ana_freq, vt, ut, ft, nb_ana ) 129 ! Initialize oscillation parameters of selected tidal components 130 ! -------------------------------------------------------------- 131 CALL tide_init_harmonics( tide_components, tide_harmonics) 134 132 135 133 IF(lwp) WRITE(numout,*) 'Analysed frequency : ',nb_ana ,'Frequency ' 136 134 137 135 DO jh = 1, nb_ana 138 IF(lwp) WRITE(numout,*) ' : ',tname(jh),' ', ana_freq(jh)136 IF(lwp) WRITE(numout,*) ' : ',tname(jh),' ',tide_harmonics(jh)%omega 139 137 END DO 140 138 … … 173 171 DO jc = 1, 2 174 172 nhc = nhc+1 175 ztemp =( MOD(jc,2) * ft(jh) *COS(ana_freq(jh)*ztime + vt(jh) + ut(jh)) & 176 & +(1.-MOD(jc,2))* ft(jh) *SIN(ana_freq(jh)*ztime + vt(jh) + ut(jh))) 173 ztemp =( MOD(jc,2) * tide_harmonics(jh)%f *COS(tide_harmonics(jh)%omega*ztime + & 174 & tide_harmonics(jh)%v0 + tide_harmonics(jh)%u) & 175 & +(1.-MOD(jc,2))* tide_harmonics(jh)%f *SIN(tide_harmonics(jh)%omega*ztime + & 176 & tide_harmonics(jh)%v0 + tide_harmonics(jh)%u)) 177 177 ! 178 178 DO jj = 1,jpj … … 235 235 nisparse(ksp) = keq 236 236 njsparse(ksp) = kun 237 valuesparse(ksp) = ( MOD(jc,2) * ft(jh) * COS(ana_freq(jh)*ztime + vt(jh) + ut(jh)) & 238 & + (1.-MOD(jc,2))* ft(jh) * SIN(ana_freq(jh)*ztime + vt(jh) + ut(jh)) ) 237 valuesparse(ksp) = ( MOD(jc,2) * tide_harmonics(jh)%f * COS(tide_harmonics(jh)%omega*ztime & 238 & + tide_harmonics(jh)%v0 + tide_harmonics(jh)%u) & 239 & + (1.-MOD(jc,2))* tide_harmonics(jh)%f * SIN(tide_harmonics(jh)%omega*ztime & 240 & + tide_harmonics(jh)%v0 + tide_harmonics(jh)%u) ) 239 241 END DO 240 242 END DO -
NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE/SBC/sbctide.F90
r10811 r10822 61 61 kt_tide = kt 62 62 ENDIF 63 CALL tide_harmo( omega_tide, v0tide, utide, ftide, nb_harmo )63 CALL tide_harmo(tide_components, tide_harmonics) ! Update oscillation parameters of tidal components 64 64 ! 65 65 ! … … 69 69 WRITE(numout,*) '~~~~~~~~ ' 70 70 DO jk = 1, nb_harmo 71 WRITE(numout,*) tide_components(jk)%cname_tide, utide(jk), ftide(jk), v0tide(jk), omega_tide(jk) 71 WRITE(numout,*) tide_harmonics(jk)%cname_tide, tide_harmonics(jk)%u, & 72 & tide_harmonics(jk)%f,tide_harmonics(jk)%v0, tide_harmonics(jk)%omega 72 73 END DO 73 74 ENDIF -
NEMO/branches/2019/dev_r10742_ENHANCE-12_SimonM-Tides/src/OCE/TDE/tide_mod.F90
r10811 r10822 22 22 PUBLIC tide_harmo ! called by tideini and diaharm modules 23 23 PUBLIC tide_init_components ! called internally and by module diaharm 24 PUBLIC tide_init_harmonics ! called internally and by module diaharm 24 25 PUBLIC tide_init_load 25 26 PUBLIC tide_init_potential … … 39 40 TYPE(tide), PUBLIC, DIMENSION(:), POINTER :: tide_components !: Array of selected tidal component parameters 40 41 41 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:) :: omega_tide !: 42 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:) :: v0tide !: 43 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:) :: utide !: 44 REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:) :: ftide !: 42 TYPE, PUBLIC :: tide_harmonic !: Oscillation parameters of harmonic tidal components 43 CHARACTER(LEN=4) :: cname_tide ! Name of component 44 REAL(wp) :: equitide ! Amplitude of equilibrium tide 45 REAL(wp) :: f ! Node factor 46 REAL(wp) :: omega ! Angular velocity 47 REAL(wp) :: v0 ! Initial phase at prime meridian 48 REAL(wp) :: u ! Phase correction 49 END type tide_harmonic 50 51 TYPE(tide_harmonic), PUBLIC, DIMENSION(:), POINTER :: tide_harmonics !: Oscillation parameters of selected tidal components 45 52 46 53 LOGICAL , PUBLIC :: ln_tide !: … … 140 147 & CALL ctl_stop('rn_tide_ramp_dt must be positive') 141 148 ! 142 ALLOCATE( omega_tide(nb_harmo), v0tide (nb_harmo), & 143 & utide (nb_harmo), ftide (nb_harmo) ) 149 ! Initialise array used to store tidal oscillation parameters (frequency, 150 ! amplitude, phase) 151 CALL tide_init_harmonics(tide_components, tide_harmonics) 152 ! 153 ! Reference time step for time-dependent tidal parameters 144 154 kt_tide = nit000 145 155 ! … … 209 219 210 220 221 SUBROUTINE tide_init_harmonics(ptide_comp, ptide_harmo) 222 !!---------------------------------------------------------------------- 223 !! *** ROUTINE tide_init_harmonics *** 224 !! 225 !! Returns pointer to array of variables of type 'tide_harmonics' that 226 !! contain oscillation parameters of the selected harmonic tidal 227 !! components 228 !! ---------------------------------------------------------------------- 229 TYPE(tide), POINTER, DIMENSION(:) :: ptide_comp ! Selected components 230 TYPE(tide_harmonic), POINTER, DIMENSION(:) :: ptide_harmo ! Oscillation parameters of selected components 231 232 ! Allocate and populate array of oscillation parameters 233 ALLOCATE(ptide_harmo(size(ptide_comp))) 234 ptide_harmo(:)%cname_tide = ptide_comp(:)%cname_tide 235 ptide_harmo(:)%equitide = ptide_comp(:)%equitide 236 CALL tide_harmo(ptide_comp, ptide_harmo) 237 238 END SUBROUTINE tide_init_harmonics 239 240 211 241 SUBROUTINE tide_init_potential 212 242 !!---------------------------------------------------------------------- … … 218 248 219 249 DO jk = 1, nb_harmo 220 zcons = rn_tide_gamma * tide_components(jk)%equitide * ftide(jk)250 zcons = rn_tide_gamma * tide_components(jk)%equitide * tide_harmonics(jk)%f 221 251 DO ji = 1, jpi 222 252 DO jj = 1, jpj 223 ztmp1 = ftide(jk) * amp_pot(ji,jj,jk) * COS( phi_pot(ji,jj,jk) + v0tide(jk) + utide(jk) ) 224 ztmp2 = -ftide(jk) * amp_pot(ji,jj,jk) * SIN( phi_pot(ji,jj,jk) + v0tide(jk) + utide(jk) ) 253 ztmp1 = tide_harmonics(jk)%f * amp_pot(ji,jj,jk) * COS( phi_pot(ji,jj,jk) & 254 & + tide_harmonics(jk)%v0 + tide_harmonics(jk)%u ) 255 ztmp2 = -tide_harmonics(jk)%f * amp_pot(ji,jj,jk) * SIN( phi_pot(ji,jj,jk) & 256 & + tide_harmonics(jk)%v0 + tide_harmonics(jk)%u ) 225 257 zlat = gphit(ji,jj)*rad !! latitude en radian 226 258 zlon = glamt(ji,jj)*rad !! longitude en radian 227 ztmp = v0tide(jk) + utide(jk)+ tide_components(jk)%nutide * zlon259 ztmp = tide_harmonics(jk)%v0 + tide_harmonics(jk)%u + tide_components(jk)%nutide * zlon 228 260 ! le potentiel est composé des effets des astres: 229 261 IF ( tide_components(jk)%nutide == 1 ) THEN ; zcs = zcons * SIN( 2._wp*zlat ) … … 277 309 278 310 279 SUBROUTINE tide_harmo( pomega, pvt, put , pcor, kc) 280 !!---------------------------------------------------------------------- 281 !!---------------------------------------------------------------------- 282 INTEGER , INTENT(in ) :: kc ! Total number of tidal constituents 283 REAL(wp), DIMENSION(kc), INTENT(out) :: pomega ! pulsation in radians/s 284 REAL(wp), DIMENSION(kc), INTENT(out) :: pvt, put, pcor ! 285 !!---------------------------------------------------------------------- 311 SUBROUTINE tide_harmo( ptide_comp, ptide_harmo ) 312 ! 313 TYPE(tide), DIMENSION(:), POINTER :: ptide_comp ! Array of selected tidal component parameters 314 TYPE(tide_harmonic), DIMENSION(:), POINTER :: ptide_harmo ! Oscillation parameters of selected tidal components 286 315 ! 287 316 CALL astronomic_angle 288 CALL tide_pulse( p omega, kc)289 CALL tide_vuf ( pvt, put, pcor, kc)317 CALL tide_pulse( ptide_comp, ptide_harmo ) 318 CALL tide_vuf( ptide_comp, ptide_harmo ) 290 319 ! 291 320 END SUBROUTINE tide_harmo … … 383 412 384 413 385 SUBROUTINE tide_pulse( p omega, kc)414 SUBROUTINE tide_pulse( ptide_comp, ptide_harmo ) 386 415 !!---------------------------------------------------------------------- 387 416 !! *** ROUTINE tide_pulse *** … … 389 418 !! ** Purpose : Compute tidal frequencies 390 419 !!---------------------------------------------------------------------- 391 INTEGER , INTENT(in ) :: kc ! Total number of tidal constituents392 REAL(wp), DIMENSION(kc), INTENT(out) :: pomega ! pulsation in radians/s420 TYPE(tide), DIMENSION(:), POINTER :: ptide_comp ! Array of selected tidal component parameters 421 TYPE(tide_harmonic), DIMENSION(:), POINTER :: ptide_harmo ! Oscillation parameters of selected tidal components 393 422 ! 394 423 INTEGER :: jh … … 404 433 zscale = rad / ( 36525._wp * 86400._wp ) 405 434 ! 406 DO jh = 1, kc407 p omega(jh) = ( zomega_T * tide_components( jh )%nT &408 & + zomega_s * tide_components( jh )%ns &409 & + zomega_h * tide_components( jh )%nh &410 & + zomega_p * tide_components( jh )%np &411 & + zomega_p1* tide_components( jh )%np1 ) * zscale435 DO jh = 1, nb_harmo 436 ptide_harmo(jh)%omega = ( zomega_T * ptide_comp( jh )%nT & 437 & + zomega_s * ptide_comp( jh )%ns & 438 & + zomega_h * ptide_comp( jh )%nh & 439 & + zomega_p * ptide_comp( jh )%np & 440 & + zomega_p1* ptide_comp( jh )%np1 ) * zscale 412 441 END DO 413 442 ! … … 415 444 416 445 417 SUBROUTINE tide_vuf( p vt, put, pcor, kc)446 SUBROUTINE tide_vuf( ptide_comp, ptide_harmo ) 418 447 !!---------------------------------------------------------------------- 419 448 !! *** ROUTINE tide_vuf *** … … 425 454 !! ft: Nodal correction factor 426 455 !!---------------------------------------------------------------------- 427 INTEGER , INTENT(in ) :: kc ! Total number of tidal constituents428 REAL(wp), DIMENSION(kc), INTENT(out) :: pvt, put, pcor !456 TYPE(tide), DIMENSION(:), POINTER :: ptide_comp ! Array of selected tidal component parameters 457 TYPE(tide_harmonic), DIMENSION(:), POINTER :: ptide_harmo ! Oscillation parameters of selected tidal components 429 458 ! 430 459 INTEGER :: jh ! dummy loop index 431 460 !!---------------------------------------------------------------------- 432 461 ! 433 DO jh = 1, kc462 DO jh = 1, nb_harmo 434 463 ! Phase of the tidal potential relative to the Greenwhich 435 464 ! meridian (e.g. the position of the fictuous celestial body). Units are radian: 436 p vt(jh) = sh_T * tide_components( jh )%nT &437 & + sh_s * tide_components( jh )%ns &438 & + sh_h * tide_components( jh )%nh &439 & + sh_p * tide_components( jh )%np &440 & + sh_p1* tide_components( jh )%np1 &441 & + tide_components( jh )%shift * rad465 ptide_harmo(jh)%v0 = sh_T * ptide_comp( jh )%nT & 466 & + sh_s * ptide_comp( jh )%ns & 467 & + sh_h * ptide_comp( jh )%nh & 468 & + sh_p * ptide_comp( jh )%np & 469 & + sh_p1* ptide_comp( jh )%np1 & 470 & + ptide_comp( jh )%shift * rad 442 471 ! 443 472 ! Phase correction u due to nodal motion. Units are radian: 444 p ut(jh) = sh_xi * tide_components( jh )%nksi &445 & + sh_nu * tide_components( jh )%nnu0 &446 & + sh_nuprim * tide_components( jh )%nnu1 &447 & + sh_nusec * tide_components( jh )%nnu2 &448 & + sh_R * tide_components( jh )%R473 ptide_harmo(jh)%u = sh_xi * ptide_comp( jh )%nksi & 474 & + sh_nu * ptide_comp( jh )%nnu0 & 475 & + sh_nuprim * ptide_comp( jh )%nnu1 & 476 & + sh_nusec * ptide_comp( jh )%nnu2 & 477 & + sh_R * ptide_comp( jh )%R 449 478 450 479 ! Nodal correction factor: 451 p cor(jh) = nodal_factort( tide_components( jh )%nformula )480 ptide_harmo(jh)%f = nodal_factort( ptide_comp( jh )%nformula ) 452 481 END DO 453 482 ! … … 672 701 ENDIF 673 702 ! 674 zwt(:) = omega_tide(:)* zt703 zwt(:) = tide_harmonics(:)%omega * zt 675 704 676 705 pot_astro(:,:) = 0._wp ! update tidal potential (sum of all harmonics)
Note: See TracChangeset
for help on using the changeset viewer.