Changeset 714 for trunk/NEMO
- Timestamp:
- 2007-10-11T11:58:02+02:00 (16 years ago)
- Location:
- trunk/NEMO
- Files:
-
- 7 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/C1D_SRC/diawri1d.F90
r712 r714 210 210 CALL histdef( nid_T, "sosalflx", "Surface Salt Flux" , "Kg/m2/s", & ! emps * sn 211 211 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 212 CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux" , "W/m2" , & ! q t212 CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux" , "W/m2" , & ! qsr + qns 213 213 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 214 214 CALL histdef( nid_T, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! qsr … … 373 373 zw2d(:,:) = emps(:,:) * sn(:,:,1) * tmask(:,:,1) 374 374 CALL histwrite( nid_T, "sosalflx", it, zw2d , ndim_hT, ndex_hT ) ! c/d salt flux 375 CALL histwrite( nid_T, "sohefldo", it, q t, ndim_hT, ndex_hT ) ! total heat flux375 CALL histwrite( nid_T, "sohefldo", it, qsr + qns , ndim_hT, ndex_hT ) ! total heat flux 376 376 CALL histwrite( nid_T, "soshfldo", it, qsr , ndim_hT, ndex_hT ) ! solar heat flux 377 377 CALL histwrite( nid_T, "somxl010", it, hmlp , ndim_hT, ndex_hT ) ! mixed layer depth -
trunk/NEMO/C1D_SRC/icestp1d.F90
r710 r714 214 214 IF( kt == nit000 ) THEN 215 215 qsr (:,:) = 0.e0 216 q t(:,:) = 0.e0216 qns (:,:) = 0.e0 217 217 qrp (:,:) = 0.e0 218 218 emp (:,:) = 0.e0 … … 232 232 ! ----------------- 233 233 234 q t (:,:) = fnsolar(:,:) + fsolar(:,:) ! non solar heat flux + solarflux234 qns (:,:) = fnsolar(:,:) ! non solar heat flux 235 235 qsr (:,:) = fsolar(:,:) ! solar flux 236 236 -
trunk/NEMO/C1D_SRC/step1d.F90
r710 r714 152 152 CALL prt_ctl(tab2d_1=emp , clinfo1=' emp - : ', mask1=tmask, ovlap=1) 153 153 CALL prt_ctl(tab2d_1=emps , clinfo1=' emps - : ', mask1=tmask, ovlap=1) 154 CALL prt_ctl(tab2d_1=q t , clinfo1=' qt- : ', mask1=tmask, ovlap=1)154 CALL prt_ctl(tab2d_1=qns , clinfo1=' qns - : ', mask1=tmask, ovlap=1) 155 155 CALL prt_ctl(tab2d_1=qsr , clinfo1=' qsr - : ', mask1=tmask, ovlap=1) 156 156 CALL prt_ctl(tab2d_1=runoff , clinfo1=' runoff : ', mask1=tmask, ovlap=1) -
trunk/NEMO/OPA_SRC/DIA/diawri.F90
r712 r714 260 260 CALL histdef( nid_T, "sosalflx", "Surface Salt Flux" , "Kg/m2/s", & ! emps * sn 261 261 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 262 CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux" , "W/m2" , & ! q t262 CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux" , "W/m2" , & ! qns + qsr 263 263 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 264 264 CALL histdef( nid_T, "soshfldo", "Shortwave Radiation" , "W/m2" , & ! qsr … … 430 430 zw2d(:,:) = emps(:,:) * sn(:,:,1) * tmask(:,:,1) 431 431 CALL histwrite( nid_T, "sosalflx", it, zw2d , ndim_hT, ndex_hT ) ! c/d salt flux 432 CALL histwrite( nid_T, "sohefldo", it, q t, ndim_hT, ndex_hT ) ! total heat flux432 CALL histwrite( nid_T, "sohefldo", it, qns + qsr , ndim_hT, ndex_hT ) ! total heat flux 433 433 CALL histwrite( nid_T, "soshfldo", it, qsr , ndim_hT, ndex_hT ) ! solar heat flux 434 434 CALL histwrite( nid_T, "somxl010", it, hmlp , ndim_hT, ndex_hT ) ! mixed layer depth … … 638 638 639 639 ! Write all fields on T grid 640 CALL histwrite( id_i, "votemper", 1, tn , jpi*jpj*jpk, idex ) ! now temperature641 CALL histwrite( id_i, "vosaline", 1, sn , jpi*jpj*jpk, idex ) ! now salinity642 #if defined key_dynspg_rl 643 CALL histwrite( id_i, "sobarstf", 1, bsfn , jpi*jpj , idex ) ! barotropic streamfunction640 CALL histwrite( id_i, "votemper", 1, tn , jpi*jpj*jpk, idex ) ! now temperature 641 CALL histwrite( id_i, "vosaline", 1, sn , jpi*jpj*jpk, idex ) ! now salinity 642 #if defined key_dynspg_rl 643 CALL histwrite( id_i, "sobarstf", 1, bsfn , jpi*jpj , idex ) ! barotropic streamfunction 644 644 #else 645 CALL histwrite( id_i, "sossheig", 1, sshn , jpi*jpj , idex ) ! sea surface height646 #endif 647 CALL histwrite( id_i, "vozocrtx", 1, un , jpi*jpj*jpk, idex ) ! now i-velocity648 CALL histwrite( id_i, "vomecrty", 1, vn , jpi*jpj*jpk, idex ) ! now j-velocity649 CALL histwrite( id_i, "vovecrtz", 1, wn , jpi*jpj*jpk, idex ) ! now k-velocity650 CALL histwrite( id_i, "sowaflup", 1, emp , jpi*jpj , idex ) ! freshwater budget651 CALL histwrite( id_i, "sohefldo", 1, q t, jpi*jpj , idex ) ! total heat flux652 CALL histwrite( id_i, "soshfldo", 1, qsr , jpi*jpj , idex ) ! totalheat flux653 CALL histwrite( id_i, "soicecov", 1, freeze , jpi*jpj , idex ) ! ice cover654 CALL histwrite( id_i, "sozotaux", 1, utau , jpi*jpj , idex ) ! i-wind stress655 CALL histwrite( id_i, "sometauy", 1, vtau , jpi*jpj , idex ) ! j-wind stress645 CALL histwrite( id_i, "sossheig", 1, sshn , jpi*jpj , idex ) ! sea surface height 646 #endif 647 CALL histwrite( id_i, "vozocrtx", 1, un , jpi*jpj*jpk, idex ) ! now i-velocity 648 CALL histwrite( id_i, "vomecrty", 1, vn , jpi*jpj*jpk, idex ) ! now j-velocity 649 CALL histwrite( id_i, "vovecrtz", 1, wn , jpi*jpj*jpk, idex ) ! now k-velocity 650 CALL histwrite( id_i, "sowaflup", 1, emp , jpi*jpj , idex ) ! freshwater budget 651 CALL histwrite( id_i, "sohefldo", 1, qsr + qns, jpi*jpj , idex ) ! total heat flux 652 CALL histwrite( id_i, "soshfldo", 1, qsr , jpi*jpj , idex ) ! solar heat flux 653 CALL histwrite( id_i, "soicecov", 1, freeze , jpi*jpj , idex ) ! ice cover 654 CALL histwrite( id_i, "sozotaux", 1, utau , jpi*jpj , idex ) ! i-wind stress 655 CALL histwrite( id_i, "sometauy", 1, vtau , jpi*jpj , idex ) ! j-wind stress 656 656 657 657 ! 3. Close the file -
trunk/NEMO/OPA_SRC/DIA/diawri_dimg.h90
r710 r714 43 43 !! level 1: utau(:,:) * umask(:,:,1) zonal stress in N.m-2 44 44 !! level 2: vtau(:,:) * vmask(:,:,1) meridional stress in N. m-2 45 !! level 3: qt (:,:)total heat flux (W/m2)46 !! level 4: 45 !! level 3: qsr + qns total heat flux (W/m2) 46 !! level 4: emp (:,:) E-P flux (mm/day) 47 47 !! level 5: tb (:,:,1)-sst model SST -forcing sst (degree C) 48 48 !! level 6: bsfb(:,:) streamfunction (m**3/s) … … 169 169 fsel(:,:,1 ) = fsel(:,:,1 ) + utau(:,:) * umask(:,:,1) 170 170 fsel(:,:,2 ) = fsel(:,:,2 ) + vtau(:,:) * vmask(:,:,1) 171 fsel(:,:,3 ) = fsel(:,:,3 ) + q t(:,:)171 fsel(:,:,3 ) = fsel(:,:,3 ) + qsr (:,:) + qns (:,:) 172 172 fsel(:,:,4 ) = fsel(:,:,4 ) + emp (:,:) 173 173 fsel(:,:,5 ) = fsel(:,:,5 ) + tb (:,:,1) - sst(:,:) … … 258 258 fsel(:,:,1 ) = utau(:,:) * umask(:,:,1) 259 259 fsel(:,:,2 ) = vtau(:,:) * vmask(:,:,1) 260 fsel(:,:,3 ) = qt (:,:) * tmask(:,:,1)260 fsel(:,:,3 ) = (qsr (:,:) + qnr (:,:)) * tmask(:,:,1) 261 261 fsel(:,:,4 ) = emp (:,:) * tmask(:,:,1) 262 262 fsel(:,:,5 ) = (tb (:,:,1) -sst(:,:)) *tmask(:,:,1) -
trunk/NEMO/OPA_SRC/TRA/trasbc.F90
r708 r714 132 132 #endif 133 133 IF( lk_vvl) THEN 134 zta = ro0cpr * ( qt(ji,jj) - qsr(ji,jj) ) * zse3t &! temperature : heat flux134 zta = ro0cpr * qns(ji,jj) * zse3t & ! temperature : heat flux 135 135 & - emp(ji,jj) * zsrau * tn(ji,jj,1) * zse3t ! & cooling/heating effet of EMP flux 136 136 zsa = 0.e0 ! No salinity concent./dilut. effect 137 137 ELSE 138 zta = ro0cpr * ( qt(ji,jj) - qsr(ji,jj)) * zse3t ! temperature : heat flux138 zta = ro0cpr * qns(ji,jj) * zse3t ! temperature : heat flux 139 139 zsa = emps(ji,jj) * zsrau * sn(ji,jj,1) * zse3t ! salinity : concent./dilut. effect 140 140 ENDIF -
trunk/NEMO/OPA_SRC/ZDF/zdfkpp.F90
r710 r714 459 459 zBosol(ji,jj) = grav * zthermal * qsr(ji,jj) 460 460 ! Non radiative surface buoyancy force 461 zBo (ji,jj) = grav * zthermal * ( qt(ji,jj) - qsr(ji,jj)) - grav * zhalin * emp(ji,jj)461 zBo (ji,jj) = grav * zthermal * qns(ji,jj) - grav * zhalin * emp(ji,jj) 462 462 ! Surface Temperature flux for non-local term 463 wt0(ji,jj) = - qt(ji,jj)* ro0cpr * tmask(ji,jj,1)463 wt0(ji,jj) = - ( qsr(ji,jj) + qns(ji,jj) )* ro0cpr * tmask(ji,jj,1) 464 464 ! Surface salinity flux for non-local term 465 465 ws0(ji,jj) = - ( emp(ji,jj) * sn(ji,jj,1) * rcs ) * tmask(ji,jj,1)
Note: See TracChangeset
for help on using the changeset viewer.