Changeset 7646 for trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90
- Timestamp:
- 2017-02-06T10:25:03+01:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMOGCM/NEMO/OPA_SRC/DIA/diawri.F90
r6387 r7646 41 41 USE zdfddm ! vertical physics: double diffusion 42 42 USE diahth ! thermocline diagnostics 43 USE wet_dry ! wetting and drying 44 USE sbcwave ! wave parameters 43 45 ! 44 46 USE lbclnk ! ocean lateral boundary conditions (or mpp link) … … 153 155 154 156 CALL iom_put( "ssh" , sshn ) ! sea surface height 157 IF( iom_use("wetdep") ) & ! wet depth 158 CALL iom_put( "wetdep" , ht_wd(:,:) + sshn(:,:) ) 155 159 156 160 CALL iom_put( "toce", tsn(:,:,:,jp_tem) ) ! 3D temperature … … 302 306 CALL iom_put( "hdiv", hdivn ) ! Horizontal divergence 303 307 ! 304 IF( iom_use("u_masstr") .OR. iom_use("u_ heattr") .OR. iom_use("u_salttr") ) THEN308 IF( iom_use("u_masstr") .OR. iom_use("u_masstr_vint") .OR. iom_use("u_heattr") .OR. iom_use("u_salttr") ) THEN 305 309 z3d(:,:,jpk) = 0.e0 310 z2d(:,:) = 0.e0 306 311 DO jk = 1, jpkm1 307 312 z3d(:,:,jk) = rau0 * un(:,:,jk) * e2u(:,:) * e3u_n(:,:,jk) * umask(:,:,jk) 313 z2d(:,:) = z2d(:,:) + z3d(:,:,jk) 308 314 END DO 309 315 CALL iom_put( "u_masstr", z3d ) ! mass transport in i-direction 316 CALL iom_put( "u_masstr_vint", z2d ) ! mass transport in i-direction vertical sum 310 317 ENDIF 311 318 … … 370 377 CALL iom_put( "v_salttr", 0.5 * z2d ) ! heat transport in j-direction 371 378 ENDIF 379 380 ! Vertical integral of temperature 381 IF( iom_use("tosmint") ) THEN 382 z2d(:,:)=0._wp 383 DO jk = 1, jpkm1 384 DO jj = 2, jpjm1 385 DO ji = fs_2, fs_jpim1 ! vector opt. 386 z2d(ji,jj) = z2d(ji,jj) + rau0 * e3t_n(ji,jj,jk) * tsn(ji,jj,jk,jp_tem) 387 END DO 388 END DO 389 END DO 390 CALL lbc_lnk( z2d, 'T', -1. ) 391 CALL iom_put( "tosmint", z2d ) 392 ENDIF 393 394 ! Vertical integral of salinity 395 IF( iom_use("somint") ) THEN 396 z2d(:,:)=0._wp 397 DO jk = 1, jpkm1 398 DO jj = 2, jpjm1 399 DO ji = fs_2, fs_jpim1 ! vector opt. 400 z2d(ji,jj) = z2d(ji,jj) + rau0 * e3t_n(ji,jj,jk) * tsn(ji,jj,jk,jp_sal) 401 END DO 402 END DO 403 END DO 404 CALL lbc_lnk( z2d, 'T', -1. ) 405 CALL iom_put( "somint", z2d ) 406 ENDIF 407 408 CALL iom_put( "bn2", rn2 ) !Brunt-Vaisala buoyancy frequency (N^2) 372 409 ! 373 410 CALL wrk_dealloc( jpi , jpj , z2d ) … … 666 703 CALL histdef( nid_T, "so28chgt", "Depth of 28C isotherm" , "m" , & ! hd28 667 704 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 668 CALL histdef( nid_T, "sohtc300", "Heat content 300 m" , " W", & ! htc3705 CALL histdef( nid_T, "sohtc300", "Heat content 300 m" , "J/m2" , & ! htc3 669 706 & jpi, jpj, nh_T, 1 , 1, 1 , -99 , 32, clop, zsto, zout ) 670 707 #endif … … 682 719 CALL histdef( nid_U, "vozocrtx", "Zonal Current" , "m/s" , & ! un 683 720 & jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) 721 IF( ln_wave .AND. ln_sdw) THEN 722 CALL histdef( nid_U, "sdzocrtx", "Stokes Drift Zonal Current" , "m/s" , & ! usd 723 & jpi, jpj, nh_U, ipk, 1, ipk, nz_U, 32, clop, zsto, zout ) 724 ENDIF 684 725 ! !!! nid_U : 2D 685 726 CALL histdef( nid_U, "sozotaux", "Wind Stress along i-axis" , "N/m2" , & ! utau … … 691 732 CALL histdef( nid_V, "vomecrty", "Meridional Current" , "m/s" , & ! vn 692 733 & jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) 734 IF( ln_wave .AND. ln_sdw) THEN 735 CALL histdef( nid_V, "sdmecrty", "Stokes Drift Meridional Current" , "m/s" , & ! vsd 736 & jpi, jpj, nh_V, ipk, 1, ipk, nz_V, 32, clop, zsto, zout ) 737 ENDIF 693 738 ! !!! nid_V : 2D 694 739 CALL histdef( nid_V, "sometauy", "Wind Stress along j-axis" , "N/m2" , & ! vtau … … 707 752 IF( lk_zdfddm ) THEN 708 753 CALL histdef( nid_W,"voddmavs","Salt Vertical Eddy Diffusivity" , "m2/s" , & ! avs 754 & jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout ) 755 ENDIF 756 757 IF( ln_wave .AND. ln_sdw) THEN 758 CALL histdef( nid_W, "sdvecrtz", "Stokes Drift Vertical Current" , "m/s" , & ! wsd 709 759 & jpi, jpj, nh_W, ipk, 1, ipk, nz_W, 32, clop, zsto, zout ) 710 760 ENDIF … … 829 879 ENDIF 830 880 881 IF( ln_wave .AND. ln_sdw ) THEN 882 CALL histwrite( nid_U, "sdzocrtx", it, usd , ndim_U , ndex_U ) ! i-StokesDrift-current 883 CALL histwrite( nid_V, "sdmecrty", it, vsd , ndim_V , ndex_V ) ! j-StokesDrift-current 884 CALL histwrite( nid_W, "sdvecrtz", it, wsd , ndim_T , ndex_T ) ! StokesDrift vert. current 885 ENDIF 886 831 887 ! 3. Close all files 832 888 ! --------------------------------------- … … 912 968 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 913 969 ! 914 CALL histdef( id_i, "ahtu" , "u-eddy diffusivity" , "m2/s" , & ! zonal current 915 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 916 CALL histdef( id_i, "ahtv" , "v-eddy diffusivity" , "m2/s" , & ! meridonal current 917 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 918 CALL histdef( id_i, "ahmt" , "t-eddy viscosity" , "m2/s" , & ! zonal current 919 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 920 CALL histdef( id_i, "ahmf" , "f-eddy viscosity" , "m2/s" , & ! meridonal current 921 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 970 IF( ALLOCATED(ahtu) ) THEN 971 CALL histdef( id_i, "ahtu" , "u-eddy diffusivity" , "m2/s" , & ! zonal current 972 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 973 CALL histdef( id_i, "ahtv" , "v-eddy diffusivity" , "m2/s" , & ! meridonal current 974 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 975 ENDIF 976 IF( ALLOCATED(ahmt) ) THEN 977 CALL histdef( id_i, "ahmt" , "t-eddy viscosity" , "m2/s" , & ! zonal current 978 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 979 CALL histdef( id_i, "ahmf" , "f-eddy viscosity" , "m2/s" , & ! meridonal current 980 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 981 ENDIF 922 982 ! 923 983 CALL histdef( id_i, "sowaflup", "Net Upward Water Flux" , "Kg/m2/S", & ! net freshwater … … 939 999 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 940 1000 ENDIF 1001 ! 1002 IF( ln_wave .AND. ln_sdw ) THEN 1003 CALL histdef( id_i, "sdzocrtx", "Stokes Drift Zonal" , "m/s" , & ! StokesDrift zonal current 1004 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 1005 CALL histdef( id_i, "sdmecrty", "Stokes Drift Merid" , "m/s" , & ! StokesDrift meridonal current 1006 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 1007 CALL histdef( id_i, "sdvecrtz", "Stokes Drift Vert" , "m/s" , & ! StokesDrift vertical current 1008 & jpi, jpj, nh_i, jpk, 1, jpk, nz_i, 32, clop, zsto, zout ) 1009 ENDIF 941 1010 942 1011 #if defined key_lim2 … … 963 1032 CALL histwrite( id_i, "vovecrtz", kt, wn , jpi*jpj*jpk, idex ) ! now k-velocity 964 1033 ! 965 CALL histwrite( id_i, "ahtu" , kt, ahtu , jpi*jpj*jpk, idex ) ! aht at u-point 966 CALL histwrite( id_i, "ahtv" , kt, ahtv , jpi*jpj*jpk, idex ) ! - at v-point 967 CALL histwrite( id_i, "ahmt" , kt, ahmt , jpi*jpj*jpk, idex ) ! ahm at t-point 968 CALL histwrite( id_i, "ahmf" , kt, ahmf , jpi*jpj*jpk, idex ) ! - at f-point 969 ! 970 CALL histwrite( id_i, "sowaflup", kt, emp-rnf , jpi*jpj , idex ) ! freshwater budget 1034 IF( ALLOCATED(ahtu) ) THEN 1035 CALL histwrite( id_i, "ahtu" , kt, ahtu , jpi*jpj*jpk, idex ) ! aht at u-point 1036 CALL histwrite( id_i, "ahtv" , kt, ahtv , jpi*jpj*jpk, idex ) ! - at v-point 1037 ENDIF 1038 IF( ALLOCATED(ahmt) ) THEN 1039 CALL histwrite( id_i, "ahmt" , kt, ahmt , jpi*jpj*jpk, idex ) ! ahm at t-point 1040 CALL histwrite( id_i, "ahmf" , kt, ahmf , jpi*jpj*jpk, idex ) ! - at f-point 1041 ENDIF 1042 ! 1043 CALL histwrite( id_i, "sowaflup", kt, emp - rnf , jpi*jpj , idex ) ! freshwater budget 971 1044 CALL histwrite( id_i, "sohefldo", kt, qsr + qns , jpi*jpj , idex ) ! total heat flux 972 1045 CALL histwrite( id_i, "soshfldo", kt, qsr , jpi*jpj , idex ) ! solar heat flux … … 979 1052 CALL histwrite( id_i, "vovvle3t", kt, e3t_n (:,:,:) , jpi*jpj*jpk, idex )! T-cell thickness 980 1053 END IF 1054 1055 IF( ln_wave .AND. ln_sdw ) THEN 1056 CALL histwrite( id_i, "sdzocrtx", kt, usd , jpi*jpj*jpk, idex) ! now StokesDrift i-velocity 1057 CALL histwrite( id_i, "sdmecrty", kt, vsd , jpi*jpj*jpk, idex) ! now StokesDrift j-velocity 1058 CALL histwrite( id_i, "sdvecrtz", kt, wsd , jpi*jpj*jpk, idex) ! now StokesDrift k-velocity 1059 ENDIF 1060 981 1061 ! 3. Close the file 982 1062 ! -----------------
Note: See TracChangeset
for help on using the changeset viewer.