Fix for wave forcing output (#126)

* fix wave forcing output

* remove noutput

src/Core/schism_glbl.F90

@@ -475,6 +475,7 @@ module schism_glbl
                                   &elev_nudge(:),uv_nudge(:),fluxprc(:),fluxevp(:), &
                                   &dav(:,:),elevmax(:),dav_max(:,:),dav_maxmag(:), & 
                                   &etaic(:),diffmax(:),diffmin(:),dfq1(:,:),dfq2(:,:)
+  real(rkind),save,allocatable,target :: rsxx(:), rsxy(:), rsyy(:)
 
   !(2,npa). ocean-ice stress (junk if no ice) [m^2/s/s]
   real(rkind),save,allocatable :: tau_oi(:,:)

src/Hydro/misc_subs.F90

@@ -6178,27 +6178,18 @@ end subroutine smooth_2dvar
 !     Compute wave force using Longuet-Higgins Stewart formulation
       subroutine compute_wave_force_lon(RSXX0,RSXY0,RSYY0)
       use schism_glbl, only : rkind,nsa,np,npa,nvrt,rho0,idry,idry_s,dp,dps,hmin_radstress, &
-     &WWAVE_FORCE,errmsg,it_main,time_stamp,ipgl,id_out_ww3
+     &WWAVE_FORCE,errmsg,it_main,time_stamp,ipgl,id_out_ww3, rsxx, rsxy, rsyy
       use schism_msgp
       use schism_io, only: writeout_nc
       implicit none
 !TODO: change to intent(in)
       REAL(rkind), intent(in) :: RSXX0(np),RSXY0(np),RSYY0(np) !from WW3, [N/m]
 
-      REAL(rkind) :: RSXX(npa),RSXY(npa),RSYY(npa) !from WW3, [N/m]
       !REAL(rkind), allocatable :: DSXX3D(:,:,:),DSXY3D(:,:,:),DSYY3D(:,:,:)
       REAL(rkind) :: DSXX3D(2,NVRT,nsa),DSXY3D(2,NVRT,nsa),DSYY3D(2,NVRT,nsa)
       integer :: IS,i
       REAL(rkind) :: HTOT,sum1,sum2,sum3,tmp
 
-!      allocate(DSXX3D(2,NVRT,nsa), DSYY3D(2,NVRT,nsa),DSXY3D(2,NVRT,nsa),stat=i)
-!      if(i/=0) call parallel_abort('compute_wave_force_lon, alloc')
-
-      !Output for check
-      call writeout_nc(id_out_ww3(1),'RSXX',1,1,np,RSXX0)
-      call writeout_nc(id_out_ww3(2),'RSXY',1,1,np,RSXY0)
-      call writeout_nc(id_out_ww3(3),'RSYY',1,1,np,RSYY0)
-
       !Check
       sum1=sum(RSXX0)
       sum2=sum(RSXY0)

src/Hydro/schism_init.F90

@@ -1423,7 +1423,8 @@ subroutine schism_init(iorder,indir,iths,ntime)
          &  iwater_type(npa),rho_mean(nvrt,nea),erho(nvrt,nea),& 
          & surf_t1(npa),surf_t2(npa),surf_t(npa),etaic(npa),sav_alpha(npa), &
          & sav_h(npa),sav_nv(npa),sav_di(npa),sav_cd(npa), &
-         & wwave_force(2,nvrt,nsa),btaun(npa),stat=istat)
+         & wwave_force(2,nvrt,nsa),btaun(npa), &
+         & rsxx(npa), rsxy(npa), rsyy(npa), stat=istat)
       if(istat/=0) call parallel_abort('INIT: other allocation failure')
 
 !     Tracers
@@ -1770,6 +1771,7 @@ subroutine schism_init(iorder,indir,iths,ntime)
       dfq1=0.d0; dfq2=0.d0 !for hotstart
       fluxevp=0.d0; fluxprc=0.d0
       prec_rain=0.d0; prec_snow=0.d0
+      rsxx=0.d0; rsxy=0.d0; rsyy=0.d0
 
 !     Fort.12 flags
 !      ifort12=0

src/Hydro/schism_step.F90

@@ -9183,6 +9183,20 @@ subroutine schism_step(it)
 
 #endif
 
+#if defined USE_WW3
+        ! Eastward wave radiation stress
+        icount=icount+1
+        call writeout_nc(id_out_ww3(1),'rsxx',1,1,npa,rsxx)
+
+        ! Eastward northward wave radiation stress
+        icount=icount+1
+        call writeout_nc(id_out_ww3(2),'rsxy',1,1,npa,rsxy)
+
+        ! Northward wave radiation stress
+        icount=icount+1
+        call writeout_nc(id_out_ww3(3),'rsyy',1,1,npa,rsyy)
+#endif
+
 #ifdef USE_MARSH
         if(iof_marsh(1)==1) call writeout_nc(id_out_var(noutput+5), &
      &'marsh_flag',4,1,nea,dble(imarsh))