Octopus
scissor.F90
Go to the documentation of this file.
1!! Copyright (C) 2016 X. Andrade, N. Tancogne-Dejean
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., 59 Temple Place - Suite 330, Boston, MA
16!! 02111-1307, USA.
17
18#include "global.h"
19
20module scissor_oct_m
21 use batch_oct_m
23 use debug_oct_m
24 use global_oct_m
27 use mesh_oct_m
32 use parser_oct_m
33 use phase_oct_m
36 use space_oct_m
42
43 implicit none
44
45 private
46
47 public :: &
48 scissor_t, &
54
55 type scissor_t
56 private
57 logical, public :: apply = .false.
58 real(real64) :: gap
59 type(states_elec_t) :: gs_st
60 end type scissor_t
61
62 interface scissor_commute_r
64 end interface
65
66contains
67
68 subroutine scissor_init(this, namespace, space, st, mesh, kpoints, phase, gap, mc)
69 type(scissor_t), intent(inout) :: this
70 type(namespace_t), intent(in) :: namespace
71 class(space_t), intent(in) :: space
72 type(states_elec_t), intent(in) :: st
73 class(mesh_t), intent(in) :: mesh
74 type(kpoints_t), intent(in) :: kpoints
75 type(phase_t), intent(in) :: phase
76 real(real64), intent(in) :: gap
77 type(multicomm_t), intent(in) :: mc
78
79 type(restart_t) :: restart_gs
80 integer :: ierr
81 integer :: ist, ik
82 complex(real64), allocatable :: temp_state(:,:)
83
84 push_sub(scissor_init)
85
86 assert(.not. this%apply)
87 assert(.not. states_are_real(st))
88
89 call messages_print_with_emphasis(msg="TDScissor", namespace=namespace)
90
91 if (st%parallel_in_states) then
92 call messages_not_implemented("Scissor operator parallel in states", namespace=namespace)
93 end if
94 if (mesh%parallel_in_domains) then
95 call messages_not_implemented("Scissor operator parallel in domains", namespace=namespace)
96 end if
97
98 this%apply = .true.
99 this%gap = gap
100
101 write(message(1),'(a)') 'Start loading GS states.'
102 call messages_info(1, namespace=namespace)
103 !We need to load GS states and to store them in this%gs_st
104 call states_elec_copy(this%gs_st, st)
105
106 call restart_gs%init(namespace, restart_proj, restart_type_load, mc, ierr, mesh=mesh)
107 if (ierr /= 0) then
108 message(1) = "Unable to read states information."
109 call messages_fatal(1, namespace=namespace)
110 end if
111
112 call states_elec_load(restart_gs, namespace, space, this%gs_st, mesh, kpoints, ierr, label = ': gs for TDScissor')
113 if (ierr /= 0 .and. ierr /= (this%gs_st%st_end-this%gs_st%st_start+1)*this%gs_st%nik*this%gs_st%d%dim) then
114 message(1) = "Unable to read wavefunctions for TDScissor."
115 call messages_fatal(1, namespace=namespace)
116 end if
117 call restart_gs%end()
118
119 if (space%is_periodic() .and. .not. (kpoints_number(kpoints) == 1 .and. kpoints_point_is_gamma(kpoints, 1))) then
120
121 write(message(1),'(a)') 'Adding the phase for GS states.'
122 call messages_info(1, namespace=namespace)
123
124 safe_allocate(temp_state(1:mesh%np_part, 1:this%gs_st%d%dim))
125 ! We apply the phase to these states, as we need it for the projectors later
126 do ik=this%gs_st%d%kpt%start, this%gs_st%d%kpt%end
127
128 do ist = this%gs_st%st_start, this%gs_st%st_end
129 call states_elec_get_state(this%gs_st, mesh, ist, ik, temp_state)
130 call phase%apply_to_single(temp_state, mesh%np, this%gs_st%d%dim, ik, .false.)
131 call states_elec_set_state(this%gs_st, mesh, ist, ik,temp_state)
132 end do
133
134 end do
135
136 safe_deallocate_a(temp_state)
137 end if
138
139 call messages_print_with_emphasis(namespace=namespace)
140
141 pop_sub(scissor_init)
142 end subroutine scissor_init
143
144 subroutine scissor_end(this)
145 type(scissor_t), intent(inout) :: this
146
147 push_sub(scissor_end)
148
149 this%apply = .false.
150 call states_elec_end(this%gs_st)
151
152 pop_sub(scissor_end)
153 end subroutine scissor_end
155#include "undef.F90"
156#include "real.F90"
157#include "scissor_inc.F90"
158
159#include "undef.F90"
160#include "complex.F90"
161#include "scissor_inc.F90"
162end module scissor_oct_m
164!! Local Variables:
165!! mode: f90
166!! coding: utf-8
167!! End:
This module implements batches of mesh functions.
Definition: batch.F90:135
This module implements common operations on batches of mesh functions.
Definition: batch_ops.F90:118
logical pure function, public kpoints_point_is_gamma(this, ik)
Definition: kpoints.F90:1552
integer pure function, public kpoints_number(this)
Definition: kpoints.F90:1101
This module defines various routines, operating on mesh functions.
This module defines the meshes, which are used in Octopus.
Definition: mesh.F90:120
subroutine, public messages_print_with_emphasis(msg, iunit, namespace)
Definition: messages.F90:898
subroutine, public messages_not_implemented(feature, namespace)
Definition: messages.F90:1091
character(len=512), private msg
Definition: messages.F90:167
character(len=256), dimension(max_lines), public message
to be output by fatal, warning
Definition: messages.F90:162
subroutine, public messages_fatal(no_lines, only_root_writes, namespace)
Definition: messages.F90:410
subroutine, public messages_info(no_lines, iunit, debug_only, stress, all_nodes, namespace)
Definition: messages.F90:594
This module handles the communicators for the various parallelization strategies.
Definition: multicomm.F90:147
integer, parameter, public restart_proj
Definition: restart.F90:156
integer, parameter, public restart_type_load
Definition: restart.F90:183
subroutine zscissor_commute_r(this, mesh, ik, psi, gpsi)
Definition: scissor.F90:529
subroutine, public zscissor_apply(this, mesh, psib, hpsib)
Definition: scissor.F90:486
subroutine dscissor_commute_r(this, mesh, ik, psi, gpsi)
Definition: scissor.F90:360
subroutine, public dscissor_apply(this, mesh, psib, hpsib)
Definition: scissor.F90:317
subroutine, public scissor_end(this)
Definition: scissor.F90:240
subroutine, public scissor_init(this, namespace, space, st, mesh, kpoints, phase, gap, mc)
Definition: scissor.F90:164
pure logical function, public states_are_real(st)
This module handles spin dimensions of the states and the k-point distribution.
subroutine, public states_elec_end(st)
finalize the states_elec_t object
subroutine, public states_elec_copy(stout, stin, exclude_wfns, exclude_eigenval, special)
make a (selective) copy of a states_elec_t object
This module handles reading and writing restart information for the states_elec_t.
subroutine, public states_elec_load(restart, namespace, space, st, mesh, kpoints, ierr, iter, lr, lowest_missing, label, verbose, skip)
returns in ierr: <0 => Fatal error, or nothing read =0 => read all wavefunctions >0 => could only rea...
int true(void)