Changeset 5630 for branches/UKMO/dev_r5107_restart_func_and_date/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90
- Timestamp:
- 2015-07-23T18:05:51+02:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/UKMO/dev_r5107_restart_func_and_date/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90
r5500 r5630 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 … … 125 127 !! 126 128 INTEGER :: ji, jj, jk ! dummy loop indices 129 INTEGER :: jkbot ! 127 130 REAL(wp) :: zztmp, zztmpx, zztmpy ! 128 131 !! … … 148 151 CALL iom_put( "e3w" , fse3w_n(:,:,:) ) 149 152 ENDIF 153 154 CALL iom_put( "ssh" , sshn ) ! sea surface height 155 if( iom_use('ssh2') ) CALL iom_put( "ssh2", sshn(:,:) * sshn(:,:) ) ! square of sea surface height 150 156 151 157 CALL iom_put( "toce", tsn(:,:,:,jp_tem) ) ! 3D temperature … … 154 160 DO jj = 1, jpj 155 161 DO ji = 1, jpi 156 z2d(ji,jj) = tsn(ji,jj,MAX(mbathy(ji,jj),1),jp_tem) 162 jkbot = mbkt(ji,jj) 163 z2d(ji,jj) = tsn(ji,jj,jkbot,jp_tem) 157 164 END DO 158 165 END DO … … 165 172 DO jj = 1, jpj 166 173 DO ji = 1, jpi 167 z2d(ji,jj) = tsn(ji,jj,MAX(mbathy(ji,jj),1),jp_sal) 174 jkbot = mbkt(ji,jj) 175 z2d(ji,jj) = tsn(ji,jj,jkbot,jp_sal) 168 176 END DO 169 177 END DO 170 178 CALL iom_put( "sbs", z2d ) ! bottom salinity 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 ) 171 195 ENDIF 172 196 … … 176 200 DO jj = 1, jpj 177 201 DO ji = 1, jpi 178 z2d(ji,jj) = un(ji,jj,MAX(mbathy(ji,jj),1)) 202 jkbot = mbku(ji,jj) 203 z2d(ji,jj) = un(ji,jj,jkbot) 179 204 END DO 180 205 END DO 181 206 CALL iom_put( "sbu", z2d ) ! bottom i-current 182 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 183 213 184 214 CALL iom_put( "voce", vn(:,:,:) ) ! 3D j-current … … 187 217 DO jj = 1, jpj 188 218 DO ji = 1, jpi 189 z2d(ji,jj) = vn(ji,jj,MAX(mbathy(ji,jj),1)) 219 jkbot = mbkv(ji,jj) 220 z2d(ji,jj) = vn(ji,jj,jkbot) 190 221 END DO 191 222 END DO 192 223 CALL iom_put( "sbv", z2d ) ! bottom j-current 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 230 231 CALL iom_put( "woce", wn ) ! vertical velocity 232 IF( iom_use('w_masstr') .OR. iom_use('w_masstr2') ) THEN ! vertical mass transport & its square value 233 ! Caution: in the VVL case, it only correponds to the baroclinic mass transport. 234 z2d(:,:) = rau0 * e12t(:,:) 235 DO jk = 1, jpk 236 z3d(:,:,jk) = wn(:,:,jk) * z2d(:,:) 237 END DO 238 CALL iom_put( "w_masstr" , z3d ) 239 IF( iom_use('w_masstr2') ) CALL iom_put( "w_masstr2", z3d(:,:,:) * z3d(:,:,:) ) 193 240 ENDIF 194 241 … … 593 640 ENDIF 594 641 595 IF( .NOT. l k_cpl ) THEN642 IF( .NOT. ln_cpl ) THEN 596 643 CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping" , "W/m2" , & ! qrp 597 644 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) … … 602 649 ENDIF 603 650 604 IF( l k_cpl .AND. nn_ice <= 1 ) THEN651 IF( ln_cpl .AND. nn_ice <= 1 ) THEN 605 652 CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping" , "W/m2" , & ! qrp 606 653 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) … … 625 672 #endif 626 673 627 IF( l k_cpl .AND. nn_ice == 2 ) THEN674 IF( ln_cpl .AND. nn_ice == 2 ) THEN 628 675 CALL histdef( nid_T,"soicetem" , "Ice Surface Temperature" , "K" , & ! tn_ice 629 676 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) … … 780 827 ENDIF 781 828 782 IF( .NOT. l k_cpl ) THEN829 IF( .NOT. ln_cpl ) THEN 783 830 CALL histwrite( nid_T, "sohefldp", it, qrp , ndim_hT, ndex_hT ) ! heat flux damping 784 831 CALL histwrite( nid_T, "sowafldp", it, erp , ndim_hT, ndex_hT ) ! freshwater flux damping … … 786 833 CALL histwrite( nid_T, "sosafldp", it, zw2d , ndim_hT, ndex_hT ) ! salt flux damping 787 834 ENDIF 788 IF( l k_cpl .AND. nn_ice <= 1 ) THEN835 IF( ln_cpl .AND. nn_ice <= 1 ) THEN 789 836 CALL histwrite( nid_T, "sohefldp", it, qrp , ndim_hT, ndex_hT ) ! heat flux damping 790 837 CALL histwrite( nid_T, "sowafldp", it, erp , ndim_hT, ndex_hT ) ! freshwater flux damping … … 802 849 #endif 803 850 804 IF( l k_cpl .AND. nn_ice == 2 ) THEN851 IF( ln_cpl .AND. nn_ice == 2 ) THEN 805 852 CALL histwrite( nid_T, "soicetem", it, tn_ice(:,:,1) , ndim_hT, ndex_hT ) ! surf. ice temperature 806 853 CALL histwrite( nid_T, "soicealb", it, alb_ice(:,:,1), ndim_hT, ndex_hT ) ! ice albedo
Note: See TracChangeset
for help on using the changeset viewer.