correct_cell_velocity Subroutine

private subroutine correct_cell_velocity(mesh, flow, bc, params, ustar, phi, local_u)

Updates cell-centered velocity using the pressure potential gradient.

Arguments

Type IntentOptional Attributes Name
type(mesh_t), intent(in) :: mesh
type(flow_mpi_t), intent(in) :: flow
type(bc_set_t), intent(in) :: bc
type(case_params_t), intent(in) :: params
real(kind=rk), intent(in) :: ustar(:,:)
real(kind=rk), intent(in) :: phi(:)
real(kind=rk), intent(out) :: local_u(:,:)

Calls

proc~~correct_cell_velocity~~CallsGraph proc~correct_cell_velocity mod_flow_projection::correct_cell_velocity proc~boundary_pressure_type mod_bc::boundary_pressure_type proc~correct_cell_velocity->proc~boundary_pressure_type proc~face_effective_neighbor mod_bc::face_effective_neighbor proc~correct_cell_velocity->proc~face_effective_neighbor proc~face_linear_scalar mod_flow_projection::face_linear_scalar proc~correct_cell_velocity->proc~face_linear_scalar proc~outward_normal mod_flow_projection::outward_normal proc~correct_cell_velocity->proc~outward_normal proc~is_periodic_face mod_bc::is_periodic_face proc~face_effective_neighbor->proc~is_periodic_face proc~face_neighbor_weight mod_flow_projection::face_neighbor_weight proc~face_linear_scalar->proc~face_neighbor_weight proc~face_neighbor_weight->proc~outward_normal proc~fatal_error mod_kinds::fatal_error proc~face_neighbor_weight->proc~fatal_error proc~patch_type_for_face mod_bc::patch_type_for_face proc~face_neighbor_weight->proc~patch_type_for_face proc~is_periodic_face->proc~patch_type_for_face

Called by

proc~~correct_cell_velocity~~CalledByGraph proc~correct_cell_velocity mod_flow_projection::correct_cell_velocity proc~advance_projection_step mod_flow_projection::advance_projection_step proc~advance_projection_step->proc~correct_cell_velocity program~lowmach_react_hex lowmach_react_hex program~lowmach_react_hex->proc~advance_projection_step

Source Code

   subroutine correct_cell_velocity(mesh, flow, bc, params, ustar, phi, local_u)
      type(mesh_t), intent(in) :: mesh
      type(flow_mpi_t), intent(in) :: flow
      type(bc_set_t), intent(in) :: bc
      type(case_params_t), intent(in) :: params
      real(rk), intent(in) :: ustar(:,:)
      real(rk), intent(in) :: phi(:)
      real(rk), intent(out) :: local_u(:,:)

      integer :: c, lf, f, nb
      real(rk) :: nvec(3), grad(3), phi_face

      local_u = zero

      do c = flow%first_cell, flow%last_cell
         grad = zero

         do lf = 1, mesh%ncell_faces(c)
            f = mesh%cell_faces(lf, c)
            nvec = outward_normal(mesh, f, c)

            nb = face_effective_neighbor(mesh, bc, f, c)

            if (nb > 0) then
               phi_face = face_linear_scalar(mesh, bc, f, c, nb, phi(c), phi(nb))
            else
               if (boundary_pressure_type(mesh, bc, f) == bc_dirichlet) then
                  phi_face = zero
               else
                  phi_face = phi(c)
               end if
            end if

            grad = grad + phi_face * nvec * mesh%faces(f)%area
         end do

         grad = grad / mesh%cells(c)%volume

         local_u(:, c) = ustar(:, c) - params%dt / params%rho * grad
      end do
   end subroutine correct_cell_velocity