Changeset 508 for trunk/NEMO/OPA_SRC/istate.F90
- Timestamp:
- 2006-10-03T17:58:55+02:00 (18 years ago)
- File:
-
- 1 edited
-
trunk/NEMO/OPA_SRC/istate.F90 (modified) (23 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/istate.F90
r479 r508 4 4 !! Ocean state : initial state setting 5 5 !!===================================================================== 6 !! History : 4.0 ! 89-12 (P. Andrich) Original code 7 !! 5.0 ! 91-11 (G. Madec) rewritting 8 !! 6.0 ! 96-01 (G. Madec) terrain following coordinates 9 !! 8.0 ! 01-09 (M. Levy, M. Ben Jelloul) istate_eel 10 !! 8.0 ! 01-09 (M. Levy, M. Ben Jelloul) istate_uvg 11 !! 9.0 ! 03-08 (G. Madec) F90: Free form, modules 12 !! 9.0 ! 03-09 (G. Madec, C. Talandier) add EEL R5 13 !! 9.0 ! 04-05 (A. Koch-Larrouy) istate_gyre 14 !! 9.0 ! 06-07 (S. Masson) distributed restart using iom 15 !!---------------------------------------------------------------------- 6 16 7 17 !!---------------------------------------------------------------------- … … 13 23 !! istate_uvg : initial velocity in geostropic balance 14 24 !!---------------------------------------------------------------------- 15 !! * Modules used16 25 USE oce ! ocean dynamics and active tracers 17 26 USE dom_oce ! ocean space and time domain … … 19 28 USE ldftra_oce ! ocean active tracers: lateral physics 20 29 USE zdf_oce ! ocean vertical physics 21 USE in_out_manager ! I/O manager22 30 USE phycst ! physical constants 23 31 USE wzvmod ! verctical velocity (wzv routine) … … 26 34 USE restart ! ocean restart (rst_read routine) 27 35 USE solisl ! ??? 28 36 USE in_out_manager ! I/O manager 37 USE iom 38 29 39 IMPLICIT NONE 30 40 PRIVATE 31 41 32 !! * Routine accessibility 33 PUBLIC istate_init ! routine called by step.F90 42 PUBLIC istate_init ! routine called by step.F90 34 43 35 44 !! * Substitutions … … 37 46 # include "vectopt_loop_substitute.h90" 38 47 !!---------------------------------------------------------------------- 39 !! OPA 9.0 , LOCEAN-IPSL (200 5)48 !! OPA 9.0 , LOCEAN-IPSL (2006) 40 49 !! $Header$ 41 !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt50 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 42 51 !!---------------------------------------------------------------------- 43 52 … … 48 57 !! *** ROUTINE istate_init *** 49 58 !! 50 !! ** Purpose : Initialization of the dynamics and tracers. 51 !! 52 !! ** Method : 53 !! 54 !! History : 55 !! 4.0 ! 91-03 () Original code 56 !! ! 91-11 (G. Madec) 57 !! 9.0 ! 03-09 (G. Madec) F90: Free form, modules, orthogonality 58 !!---------------------------------------------------------------------- 59 USE iom 60 !! * Local declarations 61 !CT INTEGER :: inum 62 !!---------------------------------------------------------------------- 63 64 65 ! Initialization to zero 66 ! ---------------------- 67 68 ! before fields ! now fields ! after fields ! 69 ; ub (:,:,:) = 0.e0 ; un (:,:,:) = 0.e0 ; ua (:,:,:) = 0.e0 70 ; vb (:,:,:) = 0.e0 ; vn (:,:,:) = 0.e0 ; va (:,:,:) = 0.e0 71 ; ; wn (:,:,:) = 0.e0 ; 72 ; rotb (:,:,:) = 0.e0 ; rotn (:,:,:) = 0.e0 ; 73 ; hdivb(:,:,:) = 0.e0 ; hdivn(:,:,:) = 0.e0 ; 74 75 ; tb (:,:,:) = 0.e0 ; tn (:,:,:) = 0.e0 ; ta (:,:,:) = 0.e0 76 ; sb (:,:,:) = 0.e0 ; sn (:,:,:) = 0.e0 ; sa (:,:,:) = 0.e0 59 !! ** Purpose : Initialization of the dynamics and tracer fields. 60 !!---------------------------------------------------------------------- 61 62 IF(lwp) WRITE(numout,*) 63 IF(lwp) WRITE(numout,*) 'istate_ini : Initialization of the dynamics and tracers' 64 IF(lwp) WRITE(numout,*) '~~~~~~~~~~' 77 65 78 66 rhd (:,:,:) = 0.e0 79 67 rhop (:,:,:) = 0.e0 80 68 rn2 (:,:,:) = 0.e0 81 82 #if defined key_dynspg_rl83 ! rigid-lid formulation84 bsfb(:,:) = 0.e0 ! before barotropic stream-function85 bsfn(:,:) = 0.e0 ! now barotropic stream-function86 bsfd(:,:) = 0.e0 ! barotropic stream-function trend87 #endif88 ! free surface formulation89 sshb(:,:) = 0.e0 ! before sea-surface height90 sshn(:,:) = 0.e0 ! now sea-surface height91 92 69 93 70 IF( ln_rstart ) THEN ! Restart from a file … … 100 77 neuler = 0 ! Set time-step indicator at nit000 (euler forward) 101 78 adatrj = 0._wp 79 ! ! Initialization of ocean to zero 80 ! before fields ! now fields 81 ; ub (:,:,:) = 0.e0 ; un (:,:,:) = 0.e0 82 ; vb (:,:,:) = 0.e0 ; vn (:,:,:) = 0.e0 83 ; rotb (:,:,:) = 0.e0 ; rotn (:,:,:) = 0.e0 84 ; hdivb(:,:,:) = 0.e0 ; hdivn(:,:,:) = 0.e0 85 ! 102 86 IF( cp_cfg == 'eel' ) THEN 103 87 CALL istate_eel ! EEL configuration : start from pre-defined … … 107 91 ! ! and salinity fields 108 92 ELSE 109 ! ! Other configurations: Initial temperature and salinity fields 110 111 !CT CALL iom_open ('ssh', inum) 112 !CT CALL iom_get( inum, jpdom_local, 'sshb', sshb ) ! free surface formulation (ssh) 113 !CT sshn(:,:) = sshb(:,:) 114 !CT CALL iom_close (inum) 115 93 ! ! Other configurations: Initial temperature and salinity fields 116 94 #if defined key_dtatem 117 95 CALL dta_tem( nit000 ) ! read 3D temperature data … … 120 98 #else 121 99 IF(lwp) WRITE(numout,*) ! analytical temperature profile 122 IF(lwp) WRITE(numout,*)' Temperature initialization using an analytic profile'100 IF(lwp) WRITE(numout,*)' Temperature initialization using an analytic profile' 123 101 CALL istate_tem 124 102 #endif … … 130 108 ! No salinity data 131 109 IF(lwp)WRITE(numout,*) ! analytical salinity profile 132 IF(lwp)WRITE(numout,*)' Salinity initialisation using a constant value'110 IF(lwp)WRITE(numout,*)' Salinity initialisation using a constant value' 133 111 CALL istate_sal 134 112 #endif … … 139 117 ! ! ----------------- 140 118 CALL wzv( nit000 ) ! from horizontal divergence 141 119 ! 142 120 END SUBROUTINE istate_init 143 121 … … 153 131 !! 154 132 !! References : Philander ??? 155 !! 156 !! History : 157 !! 4.0 ! 89-12 (P. Andrich) Original code 158 !! 6.0 ! 96-01 (G. Madec) terrain following coordinates 159 !! 9.0 ! 02-09 (G. Madec) F90: Free form 160 !!---------------------------------------------------------------------- 161 !! * Local declarations 133 !!---------------------------------------------------------------------- 162 134 INTEGER :: ji, jj, jk 163 135 !!---------------------------------------------------------------------- 164 136 ! 165 137 IF(lwp) WRITE(numout,*) 166 138 IF(lwp) WRITE(numout,*) 'istate_tem : initial temperature profile' … … 181 153 & 1 , jpi , 5 , 1 , jpk , & 182 154 & 1 , 1. , numout ) 183 155 ! 184 156 END SUBROUTINE istate_tem 185 157 … … 194 166 !! 195 167 !! ** Action : Initialize sn and sb 196 !! 197 !! History : 198 !! 4.0 ! 89-12 (P. Andrich) Original code 199 !! 8.5 ! 02-09 (G. Madec) F90: Free form 200 !!---------------------------------------------------------------------- 201 !! * Local declarations 168 !!---------------------------------------------------------------------- 202 169 REAL(wp) :: zsal = 35.50_wp 203 170 !!---------------------------------------------------------------------- … … 224 191 !! - set velocity field including horizontal divergence 225 192 !! and relative vorticity fields 226 !! 227 !! History : 228 !! 8.0 ! 01-09 (M. Levy, M. Ben Jelloul) read file for EEL 2 & 6 229 !! 9.0 ! 03-09 (G. Madec, C. Talandier) EEL 5 230 !! 9.0 ! 05-11 (V. Garnier) Surface pressure gradient organization 231 !!---------------------------------------------------------------------- 232 !! * Modules used 193 !!---------------------------------------------------------------------- 233 194 USE eosbn2 ! eq. of state, Brunt Vaisala frequency (eos routine) 234 195 USE divcur ! hor. divergence & rel. vorticity (div_cur routine) 235 196 USE iom 236 197 237 !! * Local declarations238 198 INTEGER :: inum ! temporary logical unit 239 199 INTEGER :: ji, jj, jk ! dummy loop indices 240 200 INTEGER :: ijloc 241 REAL(wp) :: & 242 zh1, zh2, zslope, zcst, & ! temporary scalars 243 zfcor 244 REAL(wp) :: & 245 zt1 = 12._wp, & ! surface temperature value (EEL R5) 246 zt2 = 2._wp, & ! bottom temperature value (EEL R5) 247 zsal = 35.5_wp, & ! constant salinity (EEL R2, R5 and R6) 248 zueel = 0.1_wp ! constant uniform zonal velocity (EEL R5) 201 REAL(wp) :: zh1, zh2, zslope, zcst, zfcor ! temporary scalars 202 REAL(wp) :: zt1 = 12._wp, & ! surface temperature value (EEL R5) 203 & zt2 = 2._wp, & ! bottom temperature value (EEL R5) 204 & zsal = 35.5_wp, & ! constant salinity (EEL R2, R5 and R6) 205 & zueel = 0.1_wp ! constant uniform zonal velocity (EEL R5) 249 206 # if ! defined key_dynspg_rl 250 REAL(wp), DIMENSION(jpiglo,jpjglo) :: & 251 zssh ! initial ssh over the global domain 207 REAL(wp), DIMENSION(jpiglo,jpjglo) :: zssh ! initial ssh over the global domain 252 208 # endif 253 209 !!---------------------------------------------------------------------- … … 389 345 !! ** Method : - set temprature field 390 346 !! - set salinity field 391 !! 392 !! ** History : 393 !! 9.0 ! 04-05 (A. Koch-Larrouy) Original code 394 !!---------------------------------------------------------------------- 395 !! * Modules used 396 USE iom 397 398 !! * Local variables 399 INTEGER :: inum ! temporary logical unit 400 INTEGER, PARAMETER :: & 401 ntsinit = 0 ! (0/1) (analytical/input data files) T&S initialization 402 347 !!---------------------------------------------------------------------- 403 348 INTEGER :: ji, jj, jk ! dummy loop indices 349 INTEGER :: inum ! temporary logical unit 350 INTEGER, PARAMETER :: ntsinit = 0 ! (0/1) (analytical/input data files) T&S initialization 404 351 !!---------------------------------------------------------------------- 405 352 … … 469 416 ENDIF 470 417 471 472 418 END SUBROUTINE istate_gyre 473 474 419 475 420 … … 484 429 !! surface to the bottom. 485 430 !! p=integral [ rau*g dz ] 486 !! 487 !! History : 488 !! 8.1 ! 01-09 (M. Levy, M. Ben Jelloul) Original code 489 !! 8.5 ! 02-09 (G. Madec) F90: Free form 490 !! 9.0 ! 05-11 (V. Garnier) Surface pressure gradient organization 491 !!---------------------------------------------------------------------- 492 !! * Modules used 431 !!---------------------------------------------------------------------- 493 432 USE eosbn2 ! eq. of state, Brunt Vaisala frequency (eos routine) 494 433 USE dynspg ! surface pressure gradient (dyn_spg routine) … … 496 435 USE lbclnk ! ocean lateral boundary condition (or mpp link) 497 436 498 !! * Local declarations499 437 INTEGER :: ji, jj, jk ! dummy loop indices 500 438 INTEGER :: indic ! ??? 501 REAL(wp) :: & 502 zmsv, zphv, zmsu, zphu, & ! temporary scalars 503 zalfg 504 REAL(wp), DIMENSION (jpi,jpj,jpk) :: & 505 zprn ! workspace 439 REAL(wp) :: zmsv, zphv, zmsu, zphu, zalfg ! temporary scalars 440 REAL(wp), DIMENSION (jpi,jpj,jpk) :: zprn ! workspace 506 441 !!---------------------------------------------------------------------- 507 442 … … 514 449 515 450 zalfg = 0.5 * grav * rau0 516 ! Surface value 517 zprn(:,:,1) = zalfg * fse3w(:,:,1) * ( 1 + rhd(:,:,1) ) 518 519 ! Vertical integration from the surface 520 DO jk = 2, jpkm1 451 452 zprn(:,:,1) = zalfg * fse3w(:,:,1) * ( 1 + rhd(:,:,1) ) ! Surface value 453 454 DO jk = 2, jpkm1 ! Vertical integration from the surface 521 455 zprn(:,:,jk) = zprn(:,:,jk-1) & 522 456 & + zalfg * fse3w(:,:,jk) * ( 2. + rhd(:,:,jk) + rhd(:,:,jk-1) ) … … 525 459 ! Compute geostrophic balance 526 460 ! --------------------------- 527 528 461 DO jk = 1, jpkm1 529 462 DO jj = 2, jpjm1 … … 571 504 ! after initializing u and v, we need to calculate the initial streamfunction bsf. 572 505 ! Otherwise, only the trend will be computed and the model will blow up (inconsistency). 573 574 506 ! to do that, we call dyn_spg with a special trick: 575 ! we fill ua and va with the velocities divided by dt, 576 ! and the streamfunction will be brought to the right 577 ! value assuming the velocities have been set up in 578 ! one time step. 579 ! we then set bsfd to zero (first guess for next step 580 ! is d(psi)/dt = 0.) 581 582 ! sets up s false trend to calculate the barotropic 583 ! streamfunction. 507 ! we fill ua and va with the velocities divided by dt, and the streamfunction will be brought to the 508 ! right value assuming the velocities have been set up in one time step. 509 ! we then set bsfd to zero (first guess for next step is d(psi)/dt = 0.) 510 ! sets up s false trend to calculate the barotropic streamfunction. 584 511 585 512 ua(:,:,:) = ub(:,:,:) / rdt … … 612 539 hdivb(:,:,:) = hdivn(:,:,:) ! set the before to the now value 613 540 rotb (:,:,:) = rotn (:,:,:) ! set the before to the now value 614 541 ! 615 542 END SUBROUTINE istate_uvg 616 543
Note: See TracChangeset
for help on using the changeset viewer.
