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