ibmwallfun Subroutine

public subroutine ibmwallfun()

Uses

  • proc~~ibmwallfun~~UsesGraph proc~ibmwallfun ibmwallfun module~modfields modfields proc~ibmwallfun->module~modfields module~modglobal modglobal proc~ibmwallfun->module~modglobal module~modmpi modmpi proc~ibmwallfun->module~modmpi module~modstat_nc modstat_nc proc~ibmwallfun->module~modstat_nc module~modsubgriddata modsubgriddata proc~ibmwallfun->module~modsubgriddata decomp_2d decomp_2d module~modfields->decomp_2d mpi mpi module~modmpi->mpi module~modstat_nc->module~modmpi netcdf netcdf module~modstat_nc->netcdf

Arguments

None

Calls

proc~~ibmwallfun~~CallsGraph proc~ibmwallfun ibmwallfun interface~writestat_nc writestat_nc proc~ibmwallfun->interface~writestat_nc proc~diffc_corr diffc_corr proc~ibmwallfun->proc~diffc_corr proc~diffu_corr diffu_corr proc~ibmwallfun->proc~diffu_corr proc~diffv_corr diffv_corr proc~ibmwallfun->proc~diffv_corr proc~diffw_corr diffw_corr proc~ibmwallfun->proc~diffw_corr proc~wallfunheat wallfunheat proc~ibmwallfun->proc~wallfunheat proc~wallfunmom wallfunmom proc~ibmwallfun->proc~wallfunmom proc~writestat_1d_nc writestat_1D_nc proc~ibmwallfun->proc~writestat_1d_nc interface~writestat_nc->proc~writestat_1d_nc proc~writestat_2d_nc writestat_2D_nc interface~writestat_nc->proc~writestat_2d_nc proc~writestat_3d_nc writestat_3D_nc interface~writestat_nc->proc~writestat_3d_nc proc~writestat_3d_short_nc writestat_3D_short_nc interface~writestat_nc->proc~writestat_3d_short_nc proc~writestat_time_nc writestat_time_nc interface~writestat_nc->proc~writestat_time_nc zstart zstart proc~diffc_corr->zstart proc~diffu_corr->zstart proc~diffv_corr->zstart proc~diffw_corr->zstart mpi_allreduce mpi_allreduce proc~wallfunheat->mpi_allreduce proc~heat_transfer_coef_flux heat_transfer_coef_flux proc~wallfunheat->proc~heat_transfer_coef_flux proc~interp_velocity_c interp_velocity_c proc~wallfunheat->proc~interp_velocity_c proc~is_equal is_equal proc~wallfunheat->proc~is_equal proc~local_coords local_coords proc~wallfunheat->proc~local_coords proc~moist_flux moist_flux proc~wallfunheat->proc~moist_flux proc~trilinear_interp_var trilinear_interp_var proc~wallfunheat->proc~trilinear_interp_var proc~wallfunheat->zstart proc~wallfunmom->mpi_allreduce proc~alignment alignment proc~wallfunmom->proc~alignment proc~wallfunmom->proc~is_equal proc~wallfunmom->proc~local_coords proc~mom_transfer_coef_neutral mom_transfer_coef_neutral proc~wallfunmom->proc~mom_transfer_coef_neutral proc~mom_transfer_coef_stability mom_transfer_coef_stability proc~wallfunmom->proc~mom_transfer_coef_stability proc~wallfunmom->proc~trilinear_interp_var proc~wallfunmom->zstart nf90_inq_varid nf90_inq_varid proc~writestat_1d_nc->nf90_inq_varid nf90_put_var nf90_put_var proc~writestat_1d_nc->nf90_put_var nf90_sync nf90_sync proc~writestat_1d_nc->nf90_sync proc~alignment->proc~is_equal proc~local_coords->proc~is_equal proc~cross_product cross_product proc~local_coords->proc~cross_product proc~trilinear_interp_var->zstart proc~eval_corners eval_corners proc~trilinear_interp_var->proc~eval_corners proc~trilinear_interp trilinear_interp proc~trilinear_interp_var->proc~trilinear_interp proc~writestat_2d_nc->nf90_inq_varid proc~writestat_2d_nc->nf90_put_var proc~writestat_2d_nc->nf90_sync proc~writestat_3d_nc->nf90_inq_varid proc~writestat_3d_nc->nf90_put_var proc~writestat_3d_nc->nf90_sync proc~writestat_3d_short_nc->nf90_inq_varid proc~writestat_3d_short_nc->nf90_put_var proc~writestat_3d_short_nc->nf90_sync proc~writestat_time_nc->nf90_inq_varid proc~writestat_time_nc->nf90_put_var proc~writestat_time_nc->nf90_sync

Called by

proc~~ibmwallfun~~CalledByGraph proc~ibmwallfun ibmwallfun program~dalesurban DALESURBAN program~dalesurban->proc~ibmwallfun

Source Code

   subroutine ibmwallfun
     use modglobal, only : libm, iwallmom, iwalltemp, xhat, yhat, zhat, ltempeq, lmoist, &
                           ib, ie, ih, ihc, jb, je, jh, jhc, kb, ke, kh, khc, nsv, totheatflux, totqflux, nfcts, rk3step, timee, nfcts, lwritefac, dt, dtfac, tfac, tnextfac
     use modfields, only : u0, v0, w0, thl0, qt0, sv0, up, vp, wp, thlp, qtp, svp, &
                           tau_x, tau_y, tau_z, thl_flux
     use modsubgriddata, only : ekm, ekh
     use modmpi, only : myid, comm3d, MPI_SUM, mpierr, MY_REAL
     use modstat_nc, only : writestat_nc, writestat_1D_nc, writestat_2D_nc

     real, allocatable :: rhs(:,:,:)
     integer n
     real :: thl_flux_sum, thl_flux_tot, mom_flux_sum, mom_flux_tot
     logical thl_flux_file_exists, mom_flux_file_exists

      if (.not. libm) return

      allocate(rhs(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh))

      if (iwallmom > 1) then
        rhs = up
        call wallfunmom(xhat, up, bound_info_u)
        tau_x(:,:,kb:ke+kh) = tau_x(:,:,kb:ke+kh) + (up - rhs)

        rhs = vp
        call wallfunmom(yhat, vp, bound_info_v)
        tau_y(:,:,kb:ke+kh) = tau_y(:,:,kb:ke+kh) + (vp - rhs)

        rhs = wp
        call wallfunmom(zhat, wp, bound_info_w)
        tau_z(:,:,kb:ke+kh) = tau_z(:,:,kb:ke+kh) + (wp - rhs)

        ! mom_flux_sum = sum(tau_x(ib:ie,jb:je,kb+1:ke) + tau_y(ib:ie,jb:je,kb+1:ke) + tau_z(ib:ie,jb:je,kb+1:ke))
        ! call MPI_ALLREDUCE(mom_flux_sum, mom_flux_tot, 1, MY_REAL, MPI_SUM, comm3d, mpierr)
        ! if (myid == 0) then
        !    if (rk3step == 3) then
        !         inquire(file="mom_flux.txt", exist=mom_flux_file_exists)
        !         if (mom_flux_file_exists) then
        !           open(12, file="mom_flux.txt", status="old", position="append", action="write")
        !         else
        !           open(12, file="mom_flux.txt", status="new", action="write")
        !         end if
        !         write(12, *) timee, -mom_flux_tot
        !         close(12)
        !    end if
        ! end if
      end if

      call diffu_corr
      call diffv_corr
      call diffw_corr

      if (ltempeq .or. lmoist .or. lwritefac) then
        rhs = thlp
        totheatflux = 0 ! Reset total heat flux to zero so we only account for that in this step.
        totqflux = 0
        call wallfunheat
        thl_flux(:,:,kb:ke+kh) = thl_flux(:,:,kb:ke+kh) + (thlp - rhs)
        if (ltempeq) call diffc_corr(thl0, thlp, ih, jh, kh)
        if (lmoist)  call diffc_corr(qt0, qtp, ih, jh, kh)

        ! thl_flux_sum = sum(thl_flux(ib:ie,jb:je,kb+1:ke))
        ! call MPI_ALLREDUCE(thl_flux_sum, thl_flux_tot, 1, MY_REAL, MPI_SUM, comm3d, mpierr)
        ! if (myid == 0) then
        !    if (rk3step == 3) then
        !         inquire(file="thl_flux.txt", exist=thl_flux_file_exists)
        !         if (thl_flux_file_exists) then
        !           open(12, file="thl_flux.txt", status="old", position="append", action="write")
        !         else
        !           open(12, file="thl_flux.txt", status="new", action="write")
        !         end if
        !         write(12, *) timee, thl_flux_tot
        !         close(12)
        !    end if
        ! end if
      end if

      do n = 1,nsv
        call diffc_corr(sv0(:,:,:,n), svp(:,:,:,n), ihc, jhc, khc)
      end do

      deallocate(rhs)

      if (lwritefac .and. rk3step==3) then

        if (myid == 0) then
            fac_tau_x_av = fac_tau_x_av + dt*fac_tau_x
            fac_tau_y_av = fac_tau_y_av + dt*fac_tau_y
            fac_tau_z_av = fac_tau_z_av + dt*fac_tau_z
            fac_pres_av = fac_pres_av + dt*fac_pres
            fac_htc_av = fac_htc_av + dt*fac_htc
            fac_cth_av = fac_cth_av + dt*fac_cth

            if (timee >= tnextfac) then
               tfac = timee - tfac
               allocate(varsfac(nfcts,nstatfac))
               varsfac(:,1) = fac_tau_x_av(1:nfcts)/tfac
               varsfac(:,2) = fac_tau_y_av(1:nfcts)/tfac
               varsfac(:,3) = fac_tau_z_av(1:nfcts)/tfac
               varsfac(:,4) = fac_pres_av(1:nfcts)/tfac
               varsfac(:,5) = fac_htc_av(1:nfcts)/tfac
               varsfac(:,6) = fac_cth_av(1:nfcts)/tfac
               call writestat_nc(ncidfac,1,tncstatfac,(/timee/),nrecfac,.true.)
               call writestat_1D_nc(ncidfac,nstatfac,ncstatfac,varsfac,nrecfac,nfcts)
               deallocate(varsfac)

               tfac = timee
               tnextfac = NINT((timee + dtfac))*1.0

               fac_tau_x_av = 0.
               fac_tau_y_av = 0.
               fac_tau_z_av = 0.
               fac_pres_av = 0.
               fac_htc_av = 0.
               fac_cth_av = 0.
            end if

        end if !myid

      end if

    end subroutine ibmwallfun