Changeset 5463 for trunk/NEMOGCM
- Timestamp:
- 2015-06-22T16:57:34+02:00 (9 years ago)
- Location:
- trunk/NEMOGCM
- Files:
-
- 2 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/CONFIG/SHARED/field_def.xml
r5426 r5463 42 42 <field id="sssmin" long_name="min of sea surface salinity" field_ref="sss" operation="minimum" /> 43 43 <field id="sbs" long_name="sea bottom salinity" unit="1e-3" /> 44 45 <field id="taubot" long_name="bottom stress module" unit="N/m2" /> 44 46 45 47 <field id="ssh" long_name="sea surface height" standard_name="sea_surface_height_above_geoid" unit="m" /> … … 366 368 <field id="ssu" long_name="ocean surface current along i-axis" unit="m/s" /> 367 369 <field id="sbu" long_name="ocean bottom current along i-axis" unit="m/s" /> 370 <field id="ubar" long_name="ocean barotropic current along i-axis" unit="m/s" /> 368 371 <field id="uocetr_eff" long_name="Effective ocean transport along i-axis" standard_name="ocean_volume_x_transport" unit="m3/s" grid_ref="grid_U_3D" /> 369 372 <field id="uocet" long_name="ocean transport along i-axis times temperature (CRS)" unit="degC*m/s" grid_ref="grid_U_3D" /> … … 400 403 <field id="ssv" long_name="ocean surface current along j-axis" unit="m/s" /> 401 404 <field id="sbv" long_name="ocean bottom current along j-axis" unit="m/s" /> 405 <field id="vbar" long_name="ocean barotropic current along j-axis" unit="m/s" /> 402 406 <field id="vocetr_eff" long_name="Effective ocean transport along j-axis" standard_name="ocean_volume_y_transport" unit="m3/s" grid_ref="grid_V_3D" /> 403 407 <field id="vocet" long_name="ocean transport along j-axis times temperature (CRS)" unit="degC*m/s" grid_ref="grid_V_3D" /> -
trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90
r5461 r5463 46 46 USE iom 47 47 USE ioipsl 48 USE dynspg_oce, ONLY: un_adv, vn_adv ! barotropic velocities 49 48 50 #if defined key_lim2 49 51 USE limwri_2 … … 176 178 CALL iom_put( "sbs", z2d ) ! bottom salinity 177 179 ENDIF 180 181 IF ( iom_use("taubot") ) THEN ! bottom stress 182 z2d(:,:) = 0._wp 183 DO jj = 2, jpjm1 184 DO ji = fs_2, fs_jpim1 ! vector opt. 185 zztmpx = ( bfrua(ji ,jj) * un(ji ,jj,mbku(ji ,jj)) & 186 & + bfrua(ji-1,jj) * un(ji-1,jj,mbku(ji-1,jj)) ) 187 zztmpy = ( bfrva(ji, jj) * vn(ji,jj ,mbkv(ji,jj )) & 188 & + bfrva(ji,jj-1) * vn(ji,jj-1,mbkv(ji,jj-1)) ) 189 z2d(ji,jj) = rau0 * SQRT( zztmpx * zztmpx + zztmpy * zztmpy ) * tmask(ji,jj,1) 190 ! 191 ENDDO 192 ENDDO 193 CALL lbc_lnk( z2d, 'T', 1. ) 194 CALL iom_put( "taubot", z2d ) 195 ENDIF 178 196 179 197 CALL iom_put( "uoce", un(:,:,:) ) ! 3D i-current … … 188 206 CALL iom_put( "sbu", z2d ) ! bottom i-current 189 207 ENDIF 208 #if defined key_dynspg_ts 209 CALL iom_put( "ubar", un_adv(:,:) ) ! barotropic i-current 210 #else 211 CALL iom_put( "ubar", un_b(:,:) ) ! barotropic i-current 212 #endif 190 213 191 214 CALL iom_put( "voce", vn(:,:,:) ) ! 3D j-current … … 200 223 CALL iom_put( "sbv", z2d ) ! bottom j-current 201 224 ENDIF 225 #if defined key_dynspg_ts 226 CALL iom_put( "vbar", vn_adv(:,:) ) ! barotropic j-current 227 #else 228 CALL iom_put( "vbar", vn_b(:,:) ) ! barotropic j-current 229 #endif 202 230 203 231 CALL iom_put( "woce", wn ) ! vertical velocity
Note: See TracChangeset
for help on using the changeset viewer.