source: trunk/src/run_ape.F90 @ 627

!! Copyright (C) 2004-2011 M. Oliveira, F. Nogueira
!!
!! This program is free software; you can redistribute it and/or modify
!! it under the terms of the GNU General Public License as published by
!! the Free Software Foundation; either version 2, or (at your option)
!! any later version.
!!
!! This program is distributed in the hope that it will be useful,
!! but WITHOUT ANY WARRANTY; without even the implied warranty of
!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
!! GNU General Public License for more details.
!!
!! You should have received a copy of the GNU General Public License
!! along with this program; if not, write to the Free Software
!! Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA
!! 02111-1307, USA.
!!
!! $Id$

#include "global.h"

module run_ape_m
  use global_m
  use oct_parser_m
  use utilities_m
  use messages_m
  use units_m
  use output_m
  use wave_equations_integrator_m
  use eigensolver_m
  use atom_m
  use numerical_tests_m
  implicit none


                    !---Global Variables---!

  !The stack
  integer :: stack(100) = 0
  integer :: n_inst     = 0

  !Modes
  integer, parameter :: AE         = 1,  &
                        PP_GEN     = 2,  &
                        PP_TEST    = 4,  &
                        PP_CONVERT = 8,  &
                        NUM_TESTS  = 16, &
                        XC_EVAL    = 32

  !Instructions
  integer, parameter :: ATOM_AE_CREATE  = 1,  &
                        ATOM_AE_SAVE    = 2,  &
                        ATOM_AE_LOAD    = 3,  &
                        ATOM_AE_XC_EVAL = 4,  &
                        ATOM_AE_END     = 5,  &
                        ATOM_PS_CREATE  = 11, &
                        ATOM_PS_SAVE    = 12, &
                        ATOM_PS_LOAD    = 13, &
                        ATOM_PS_TEST    = 14, &
                        ATOM_PS_END     = 15, &
                        ATOM_KB_END     = 16, &
                        NUMERICAL_TESTS = 30, &
                        LEAVE           = 100


                    !---Public/Private Statements---!

  private
  public :: run, AE, PP_GEN, PP_TEST, NUM_TESTS, XC_EVAL


contains

  subroutine run()
    !-----------------------------------------------------------------------!
    ! This subroutine reads the run mode and puts the corresponding tasks   !
    ! in a stack. It then performs all the tasks that are in the stack by   !
    ! calling the appropriate subroutines.                                  !
    !                                                                       !
    !-----------------------------------------------------------------------!
    integer             :: calcmode, mode_test, inst, message_length
    type(integrator_t)  :: integrator_sp, integrator_dp
    type(eigensolver_t) :: eigensolver
    type(atom_t)        :: ae_atom, ps_atom, kb_atom

    call push_sub("run")

    !Read calculation mode
    call oct_parse_int('CalculationMode', AE + PP_GEN, calcmode)

    !Write what the code is going to do and check the input
    message(1) = "Calculation Type:"
    message_length = 1
    mode_test = calcmode
    if (iand(calcmode, AE) /= 0) then
      message_length = message_length + 1
      message(message_length) = "  Atomic Calculation"
      mode_test = mode_test - iand(AE, calcmode)
    end if
    if (iand(calcmode, PP_GEN) /= 0) then
      message_length = message_length + 1
      message(message_length) = "  PseudoPotential Generation"
      mode_test = mode_test - iand(PP_GEN, calcmode)
    end if
    if (iand(calcmode, PP_TEST) /= 0) then
      message_length = message_length + 1
      message(message_length) = "  PseudoPotential Test"
      mode_test = mode_test - iand(PP_TEST, calcmode)
    end if
    if (iand(calcmode, PP_CONVERT) /= 0) then
      message_length = message_length + 1
      message(message_length) = "  PseudoPotential File Conversion"
      mode_test = mode_test - iand(PP_CONVERT, calcmode)
    end if
    if (iand(calcmode, NUM_TESTS) /= 0) then
      message_length = message_length + 1
      message(message_length) = "  Numerical Tests"
      mode_test = mode_test - iand(NUM_TESTS, calcmode)
    end if
    if (iand(calcmode, XC_EVAL) /= 0) then
      message_length = message_length + 1
      message(message_length) = "  Exchange-Correlation Evaluation"
      mode_test = mode_test - iand(XC_EVAL, calcmode)
    end if
    if (mode_test /= 0) then
      write(message(1), '(A,I2,A)') "Input: '", calcmode, "' is not a valid CalculationMode"
      call write_fatal(1)
    else
      call write_info(message_length)
    end if

    !Initialize stack
    !Remember this is a stack: first in last out!
    call push(LEAVE) !Last thing to do is always to exit

    !Are we going to convert the pseudopotential data?
    if (iand(calcmode, PP_CONVERT) /= 0) then
      call push(ATOM_PS_END)
      call push(ATOM_PS_LOAD)
    end if

    !Are we going to test the pseudopotentials?
    if (iand(calcmode, PP_TEST) /= 0) then
      call push(ATOM_PS_END)
      call push(ATOM_PS_TEST)
      call push(ATOM_PS_LOAD)
    end if

    !Are we going to do a pseudopotential generation?
    if (iand(calcmode, PP_GEN) /= 0) then
      call push(ATOM_AE_END)
      call push(ATOM_PS_END)
      call push(ATOM_KB_END)
      call push(ATOM_PS_SAVE)
      call push(ATOM_PS_CREATE)
      call push(ATOM_AE_LOAD)
    end if

    !Are we going to do an all-electron calculation?
    if (iand(calcmode, AE) /= 0) then
      call push(ATOM_AE_END)
      call push(ATOM_AE_SAVE)
      call push(ATOM_AE_CREATE)
    end if

    !Are we going to perform numerical tests?
    if (iand(calcmode, NUM_TESTS) /= 0) then
      call push(NUMERICAL_TESTS)
    end if

    !Are we going to evaluate the xc potential and energy?
    if (iand(calcmode, XC_EVAL) /= 0) then
      call push(ATOM_AE_END)
      call push(ATOM_AE_XC_EVAL)
      call push(ATOM_AE_LOAD)
    end if

    !Initialize things needed for all modes
    call units_init()
    call eigensolver_init(eigensolver)

    !Main loop: do everything that is in the stack
    do
      call pop(inst)
      select case (inst)
      case (ATOM_AE_CREATE)
        message(1) = ""
        message(2) = ""
        message(3) = str_center("-- All Electron Calculation --", 70)
        call write_info(3)

        call output_init("ae")
        call atom_null(ae_atom)
        call atom_create_ae(ae_atom, eigensolver)
        call atom_output(ae_atom, "ae")
        call output_end("ae")

      case (ATOM_AE_SAVE)
        call atom_save(ae_atom, "ae")

      case (ATOM_AE_LOAD)
        call atom_null(ae_atom)
        call atom_load(ae_atom, "ae")

      case (ATOM_AE_XC_EVAL)
        message(1) = ""
        message(2) = ""
        message(3) = str_center("-- Exchage-Correlation Evaluation --", 70)
        call output_init("xc")
        call atom_xc_eval(ae_atom)
        call output_end("xc")

      case (ATOM_AE_END)
        call atom_end(ae_atom)

      case (ATOM_PS_CREATE)
        message(1) = ""
        message(2) = ""
        message(3) = str_center("-- Pseudopotential Generation --", 70)
        call write_info(3)

        call output_init("pp")
        call output_init("kb")
        call atom_null(ps_atom)
        call atom_null(kb_atom)
        call atom_create_ps(ae_atom, ps_atom, kb_atom, eigensolver)
        call atom_output(ps_atom, "pp")
        call atom_output(kb_atom, "kb")
        call output_end("pp")
        call output_end("kb")

      case (ATOM_PS_SAVE)
        call atom_save(ps_atom, "pp")

      case (ATOM_PS_LOAD)
        call atom_null(ps_atom)
        call atom_load(ps_atom, "pp")

      case (ATOM_PS_TEST)
        message(1) = ""
        message(2) = ""
        message(3) = str_center("-- Pseudopotential Testing --", 70)
        call write_info(3)

        call output_init("tests")
        call atom_test(ps_atom, eigensolver)
        call output_end("tests")
      case (ATOM_PS_END)
        call atom_end(ps_atom)
      case (NUMERICAL_TESTS)
        message(1) = ""
        message(2) = ""
        message(3) = str_center("-- Numerical Tests --", 70)
        call output_init("nt")
        call numerical_tests_run()
        call output_end("nt")
      case (LEAVE)
        exit
      end select
    end do

    !End things needed for all modes
    call eigensolver_end(eigensolver)

    call pop_sub()
  contains

    subroutine push(inst)
      integer, intent(in) :: inst

      if (n_inst == 100) then
        message(1) = "Too many instructions in stack"
        call write_fatal(1)
      end if
      n_inst = n_inst + 1
      stack(n_inst) = inst

    end subroutine push

    subroutine pop(inst)
      integer, intent(out) :: inst

      if (n_inst == 0) then
        message(1) = "Not enough instructions in stack"
        call write_fatal(1)
      end if
      inst = stack(n_inst)
      stack(n_inst) = 0
      n_inst = n_inst - 1

    end subroutine pop

  end subroutine run

end module run_ape_m