1!! Copyright (C) 2002-2006 M. Marques, A. Castro, A. Rubio, G. Bertsch 2!! 3!! This program is free software; you can redistribute it and/or modify 4!! it under the terms of the GNU General Public License as published by 5!! the Free Software Foundation; either version 2, or (at your option) 6!! any later version. 7!! 8!! This program is distributed in the hope that it will be useful, 9!! but WITHOUT ANY WARRANTY; without even the implied warranty of 10!! MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 11!! GNU General Public License for more details. 12!! 13!! You should have received a copy of the GNU General Public License 14!! along with this program; if not, write to the Free Software 15!! Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 16!! 02110-1301, USA. 17!! 18 19#include "global.h" 20 21module ground_state_oct_m 22 use calc_mode_par_oct_m 23 use global_oct_m 24 use grid_oct_m 25 use hamiltonian_elec_oct_m 26 use io_function_oct_m 27 use lcao_oct_m 28 use mesh_oct_m 29 use messages_oct_m 30 use mpi_oct_m 31 use multicomm_oct_m 32 use namespace_oct_m 33 use pcm_oct_m 34 use rdmft_oct_m 35 use restart_oct_m 36 use scf_oct_m 37 use states_abst_oct_m 38 use states_elec_oct_m 39 use states_elec_restart_oct_m 40 use system_oct_m 41 use v_ks_oct_m 42 43 implicit none 44 45 private 46 public :: & 47 ground_state_run_init, & 48 ground_state_run 49 50contains 51 52 subroutine ground_state_run_init() 53 54 PUSH_SUB(ground_state_run_init) 55 56 call calc_mode_par_set_parallelization(P_STRATEGY_STATES, default = .false.) 57#ifdef HAVE_SCALAPACK 58 call calc_mode_par_set_scalapack_compat() 59#endif 60 61 POP_SUB(ground_state_run_init) 62 end subroutine ground_state_run_init 63 64 ! --------------------------------------------------------- 65 subroutine ground_state_run(sys, fromScratch) 66 type(system_t), intent(inout) :: sys 67 logical, intent(inout) :: fromScratch 68 69 type(scf_t) :: scfv 70 type(restart_t) :: restart_load, restart_dump 71 integer :: ierr 72 type(rdm_t) :: rdm 73 74 PUSH_SUB(ground_state_run) 75 76 call messages_write('Info: Allocating ground state wave-functions') 77 call messages_info() 78 79 if(sys%st%parallel_in_states) then 80 call messages_experimental('State parallelization for ground state calculations') 81 end if 82 83 if (sys%hm%pcm%run_pcm) then 84 if (sys%hm%pcm%epsilon_infty /= sys%hm%pcm%epsilon_0 .and. sys%hm%pcm%tdlevel /= PCM_TD_EQ) then 85 message(1) = 'Non-equilbrium PCM is not active in a time-independent run.' 86 message(2) = 'You set epsilon_infty /= epsilon_0, but epsilon_infty is not relevant for CalculationMode = gs.' 87 message(3) = 'By definition, the ground state is in equilibrium with the solvent.' 88 message(4) = 'Therefore, the only relevant dielectric constant is the static one.' 89 message(5) = 'Nevertheless, the dynamical PCM response matrix is evaluated for benchamarking purposes.' 90 call messages_warning(5) 91 end if 92 end if 93 94 call states_elec_allocate_wfns(sys%st, sys%gr%mesh, packed=.true.) 95 96#ifdef HAVE_MPI 97 ! sometimes a deadlock can occur here (if some nodes can allocate and other cannot) 98 if(sys%st%dom_st_kpt_mpi_grp%comm > 0) call MPI_Barrier(sys%st%dom_st_kpt_mpi_grp%comm, mpi_err) 99#endif 100 call messages_write('Info: Ground-state allocation done.') 101 call messages_info() 102 103 if(.not. fromScratch) then 104 ! load wavefunctions 105 ! in RDMFT we need the full ground state 106 call restart_init(restart_load, sys%namespace, RESTART_GS, RESTART_TYPE_LOAD, sys%mc, ierr, & 107 mesh=sys%gr%mesh, exact = (sys%ks%theory_level == RDMFT)) 108 if(ierr == 0) & 109 call states_elec_load(restart_load, sys%namespace, sys%st, sys%gr, ierr) 110 111 if(ierr /= 0) then 112 call messages_write("Unable to read wavefunctions.") 113 call messages_new_line() 114 call messages_write("Starting from scratch!") 115 call messages_warning() 116 fromScratch = .true. 117 end if 118 end if 119 120 call write_canonicalized_xyz_file("exec", "initial_coordinates", sys%geo, sys%gr%mesh, sys%namespace) 121 122 if(sys%ks%theory_level /= RDMFT) then 123 call scf_init(scfv, sys%namespace, sys%gr, sys%geo, sys%st, sys%mc, sys%hm, sys%ks) 124 end if 125 126 if (fromScratch .and. sys%ks%theory_level /= RDMFT) then 127 call lcao_run(sys, lmm_r = scfv%lmm_r) 128 else 129 ! setup Hamiltonian 130 call messages_write('Info: Setting up Hamiltonian.') 131 call messages_info() 132 call system_h_setup(sys, calc_eigenval = .false., calc_current = .false.) 133 end if 134 135 call restart_init(restart_dump, sys%namespace, RESTART_GS, RESTART_TYPE_DUMP, sys%mc, ierr, mesh=sys%gr%mesh) 136 137 ! run self-consistency 138 call scf_state_info(sys%st) 139 140 if(sys%st%d%pack_states .and. hamiltonian_elec_apply_packed(sys%hm)) call sys%st%pack() 141 142 ! self-consistency for occupation numbers and natural orbitals in RDMFT 143 if(sys%ks%theory_level == RDMFT) then 144 call rdmft_init(rdm, sys%namespace, sys%gr, sys%st, sys%geo, sys%mc, fromScratch) 145 call scf_rdmft(rdm, sys%namespace, sys%gr, sys%geo, sys%st, sys%ks, sys%hm, sys%outp, restart_dump) 146 call rdmft_end(rdm, sys%gr) 147 else 148 if(.not. fromScratch) then 149 call scf_run(scfv, sys%namespace, sys%mc, sys%gr, sys%geo, sys%st, sys%ks, sys%hm, sys%outp, & 150 restart_load=restart_load, restart_dump=restart_dump) 151 call restart_end(restart_load) 152 else 153 call scf_run(scfv, sys%namespace, sys%mc, sys%gr, sys%geo, sys%st, sys%ks, sys%hm, sys%outp, & 154 restart_dump=restart_dump) 155 end if 156 157 call scf_end(scfv) 158 end if 159 160 call restart_end(restart_dump) 161 162 if(sys%st%d%pack_states .and. hamiltonian_elec_apply_packed(sys%hm)) call sys%st%unpack() 163 164 ! clean up 165 call states_elec_deallocate_wfns(sys%st) 166 167 POP_SUB(ground_state_run) 168 end subroutine ground_state_run 169 170end module ground_state_oct_m 171 172!! Local Variables: 173!! mode: f90 174!! coding: utf-8 175!! End: 176