mg_fas_vcycle Subroutine

  subroutine mg_fas_vcycle(mg, highest_lvl, max_res, standalone)
    use mpi
    type(mg_t), intent(inout)       :: mg
    integer, intent(in), optional   :: highest_lvl
    !< Maximum level for V-cycle
    real(dp), intent(out), optional :: max_res
    !< Store max(abs(residual))
    !> Whether the V-cycle is called by itself (default: true)
    logical, intent(in), optional   :: standalone
    integer                         :: lvl, min_lvl, i, max_lvl, ierr
    real(dp)                        :: init_res, res
    logical                         :: is_standalone

    is_standalone = .true.
    if (present(standalone)) is_standalone = standalone

    call check_methods(mg)
    if (timer_smoother == -1) call mg_add_timers(mg)

    if (is_standalone) call mg_timer_start(mg%timers(timer_total_vcycle))

    if (mg%subtract_mean .and. .not. present(highest_lvl)) then
       ! Assume that this is a stand-alone call. For fully periodic solutions,
       ! ensure the mean source term is zero.
       call subtract_mean(mg, mg_irhs, .false.)
    end if

    min_lvl = mg%lowest_lvl
    max_lvl = mg%highest_lvl
    if (present(highest_lvl)) max_lvl = highest_lvl

    ! Ensure ghost cells are filled correctly
    if (is_standalone) then
       call mg_fill_ghost_cells_lvl(mg, max_lvl, mg_iphi)
    end if

    do lvl = max_lvl,  min_lvl+1, -1
       ! Downwards relaxation
       call smooth_boxes(mg, lvl, mg%n_cycle_down)

       ! Set rhs on coarse grid, restrict phi, and copy mg_iphi to mg_iold for the
       ! correction later
       call mg_timer_start(mg%timers(timer_update_coarse))
       call update_coarse(mg, lvl)
       call mg_timer_end(mg%timers(timer_update_coarse))
    end do

    call mg_timer_start(mg%timers(timer_coarse))
    if (.not. all(mg%boxes(mg%lvls(min_lvl)%ids)%rank == &
       mg%boxes(mg%lvls(min_lvl)%ids(1))%rank)) then
       error stop "Multiple CPUs for coarse grid (not implemented yet)"
    end if

    init_res = max_residual_lvl(mg, min_lvl)
    do i = 1, mg%max_coarse_cycles
       call smooth_boxes(mg, min_lvl, mg%n_cycle_up+mg%n_cycle_down)
       res = max_residual_lvl(mg, min_lvl)
       if (res < mg%residual_coarse_rel * init_res .or. res < &
          mg%residual_coarse_abs) exit
    end do
    call mg_timer_end(mg%timers(timer_coarse))

    ! Do the upwards part of the v-cycle in the tree
    do lvl = min_lvl+1, max_lvl
       ! Correct solution at this lvl using lvl-1 data
       ! phi = phi + prolong(phi_coarse - phi_old_coarse)
       call mg_timer_start(mg%timers(timer_correct))
       call correct_children(mg, lvl-1)

       ! Have to fill ghost cells after correction
       call mg_fill_ghost_cells_lvl(mg, lvl, mg_iphi)
       call mg_timer_end(mg%timers(timer_correct))

       ! Upwards relaxation
       call smooth_boxes(mg, lvl, mg%n_cycle_up)
    end do

    if (present(max_res)) then
       init_res = 0.0_dp
       do lvl = min_lvl, max_lvl
          res = max_residual_lvl(mg, lvl)
          init_res = max(res, init_res)
       end do
       call mpi_allreduce(init_res, max_res, 1, mpi_double, mpi_max, mg%comm,&
           ierr)
    end if

    ! Subtract mean(phi) from phi
    if (mg%subtract_mean) then
       call subtract_mean(mg, mg_iphi, .true.)
    end if

    if (is_standalone) call mg_timer_end(mg%timers(timer_total_vcycle))
  end subroutine mg_fas_vcycle