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 system_oct_m 22 use accel_oct_m 23 use calc_mode_par_oct_m 24 use density_oct_m 25 use elf_oct_m 26 use energy_calc_oct_m 27 use geometry_oct_m 28 use global_oct_m 29 use grid_oct_m 30 use hamiltonian_elec_oct_m 31 use mesh_oct_m 32 use messages_oct_m 33 use modelmb_particles_oct_m 34 use mpi_oct_m 35 use multicomm_oct_m 36 use namespace_oct_m 37 use output_oct_m 38 use parser_oct_m 39 use poisson_oct_m 40 use profiling_oct_m 41 use space_oct_m 42 use simul_box_oct_m 43 use sort_oct_m 44 use states_abst_oct_m 45 use states_elec_oct_m 46 use states_elec_dim_oct_m 47 use v_ks_oct_m 48 use xc_oct_m 49 use xc_oep_oct_m 50 51 implicit none 52 53 private 54 public :: & 55 system_t, & 56 system_init, & 57 system_h_setup 58 59 type :: system_t 60 ! Components are public by default 61 type(space_t) :: space 62 type(geometry_t) :: geo 63 type(grid_t), pointer :: gr !< the mesh 64 type(states_elec_t), pointer :: st !< the states 65 type(v_ks_t) :: ks !< the Kohn-Sham potentials 66 type(output_t) :: outp !< the output 67 type(multicomm_t) :: mc !< index and domain communicators 68 type(namespace_t) :: namespace 69 type(hamiltonian_elec_t) :: hm 70 contains 71 final :: system_finalize 72 end type system_t 73 74contains 75 76 !---------------------------------------------------------- 77 function system_init(namespace) result(sys) 78 class(system_t), pointer :: sys 79 type(namespace_t), intent(in) :: namespace 80 81 type(profile_t), save :: prof 82 83 PUSH_SUB(system_init) 84 call profiling_in(prof,"SYSTEM_INIT") 85 86 SAFE_ALLOCATE(sys) 87 88 SAFE_ALLOCATE(sys%gr) 89 SAFE_ALLOCATE(sys%st) 90 91 sys%namespace = namespace 92 93 call messages_obsolete_variable(sys%namespace, 'SystemName') 94 95 call space_init(sys%space, sys%namespace) 96 97 call geometry_init(sys%geo, sys%namespace, sys%space) 98 call grid_init_stage_0(sys%gr, sys%namespace, sys%geo, sys%space) 99 call states_elec_init(sys%st, sys%namespace, sys%gr, sys%geo) 100 call sys%st%write_info(sys%namespace) 101 call grid_init_stage_1(sys%gr, sys%namespace, sys%geo) 102 ! if independent particles in N dimensions are being used, need to initialize them 103 ! after masses are set to 1 in grid_init_stage_1 -> derivatives_init 104 call modelmb_copy_masses (sys%st%modelmbparticles, sys%gr%der%masses) 105 106 call parallel_init() 107 108 call geometry_partition(sys%geo, sys%mc) 109 call kpoints_distribute(sys%st%d, sys%mc) 110 call states_elec_distribute_nodes(sys%st, sys%namespace, sys%mc) 111 call grid_init_stage_2(sys%gr, sys%namespace, sys%mc, sys%geo) 112 if(sys%st%symmetrize_density) call mesh_check_symmetries(sys%gr%mesh, sys%gr%sb) 113 114 call v_ks_nullify(sys%ks) 115 call output_init(sys%outp, sys%namespace, sys%gr%sb, sys%st, sys%st%nst, sys%ks) 116 call states_elec_densities_init(sys%st, sys%gr, sys%geo) 117 call states_elec_exec_init(sys%st, sys%namespace, sys%mc) 118 call elf_init(sys%namespace) 119 120 call v_ks_init(sys%ks, sys%namespace, sys%gr, sys%st, sys%geo, sys%mc) 121 122 call hamiltonian_elec_init(sys%hm, sys%namespace, sys%gr, sys%geo, sys%st, sys%ks%theory_level, & 123 sys%ks%xc, sys%mc, need_exchange = output_need_exchange(sys%outp) .or. sys%ks%oep%level /= XC_OEP_NONE) 124 125 if(poisson_is_multigrid(sys%hm%psolver)) call grid_create_multigrid(sys%gr, sys%namespace, sys%geo, sys%mc) 126 127 if (sys%hm%pcm%run_pcm .and. sys%mc%par_strategy /= P_STRATEGY_SERIAL .and. sys%mc%par_strategy /= P_STRATEGY_STATES) then 128 call messages_experimental('Parallel in domain calculations with PCM') 129 end if 130 131 call profiling_out(prof) 132 POP_SUB(system_init) 133 134 contains 135 136 ! --------------------------------------------------------- 137 subroutine parallel_init() 138 integer :: index_range(4) 139 140 PUSH_SUB(system_init.parallel_init) 141 142 ! store the ranges for these two indices (serves as initial guess 143 ! for parallelization strategy) 144 index_range(1) = sys%gr%mesh%np_global ! Number of points in mesh 145 index_range(2) = sys%st%nst ! Number of states 146 index_range(3) = sys%st%d%nik ! Number of k-points 147 index_range(4) = 100000 ! Some large number 148 149 ! create index and domain communicators 150 call multicomm_init(sys%mc, sys%namespace, mpi_world, calc_mode_par_parallel_mask(), calc_mode_par_default_parallel_mask(), & 151 mpi_world%size, index_range, (/ 5000, 1, 1, 1 /)) 152 153 POP_SUB(system_init.parallel_init) 154 end subroutine parallel_init 155 156 end function system_init 157 158 !---------------------------------------------------------- 159 subroutine system_finalize(sys) 160 type(system_t), intent(inout) :: sys 161 162 PUSH_SUB(system_finalize) 163 164 call hamiltonian_elec_end(sys%hm) 165 166 call multicomm_end(sys%mc) 167 168 call v_ks_end(sys%ks, sys%gr) 169 170 call output_end(sys%outp) 171 172 if(associated(sys%st)) then 173 call states_elec_end(sys%st) 174 SAFE_DEALLOCATE_P(sys%st) 175 end if 176 177 call geometry_end(sys%geo) 178 call simul_box_end(sys%gr%sb) 179 call grid_end(sys%gr) 180 181 call space_end(sys%space) 182 183 SAFE_DEALLOCATE_P(sys%gr) 184 185 POP_SUB(system_finalize) 186 end subroutine system_finalize 187 188 189 !---------------------------------------------------------- 190 subroutine system_h_setup(sys, calc_eigenval, calc_current) 191 type(system_t), intent(inout) :: sys 192 logical, optional, intent(in) :: calc_eigenval !< default is true 193 logical, optional, intent(in) :: calc_current !< default is true 194 195 integer, allocatable :: ind(:) 196 integer :: ist, ik 197 FLOAT, allocatable :: copy_occ(:) 198 logical :: calc_eigenval_ 199 logical :: calc_current_ 200 201 PUSH_SUB(system_h_setup) 202 203 calc_eigenval_ = optional_default(calc_eigenval, .true.) 204 calc_current_ = optional_default(calc_current, .true.) 205 call states_elec_fermi(sys%st, sys%namespace, sys%gr%mesh) 206 call density_calc(sys%st, sys%gr, sys%st%rho) 207 call v_ks_calc(sys%ks, sys%namespace, sys%hm, sys%st, sys%geo, calc_eigenval = calc_eigenval_, calc_current = calc_current_) ! get potentials 208 209 if(sys%st%restart_reorder_occs .and. .not. sys%st%fromScratch) then 210 message(1) = "Reordering occupations for restart." 211 call messages_info(1) 212 213 SAFE_ALLOCATE(ind(1:sys%st%nst)) 214 SAFE_ALLOCATE(copy_occ(1:sys%st%nst)) 215 216 do ik = 1, sys%st%d%nik 217 call sort(sys%st%eigenval(:, ik), ind) 218 copy_occ(1:sys%st%nst) = sys%st%occ(1:sys%st%nst, ik) 219 do ist = 1, sys%st%nst 220 sys%st%occ(ist, ik) = copy_occ(ind(ist)) 221 end do 222 end do 223 224 SAFE_DEALLOCATE_A(ind) 225 SAFE_DEALLOCATE_A(copy_occ) 226 end if 227 228 if(calc_eigenval_) call states_elec_fermi(sys%st, sys%namespace, sys%gr%mesh) ! occupations 229 call energy_calc_total(sys%namespace, sys%hm, sys%gr, sys%st) 230 231 POP_SUB(system_h_setup) 232 end subroutine system_h_setup 233 234end module system_oct_m 235 236!! Local Variables: 237!! mode: f90 238!! coding: utf-8 239!! End: 240