Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2026 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief Read pretabulated atom-relative RI-RS grids from data files.
10 : !> \par History
11 : !> 09.2026 created
12 : ! **************************************************************************************************
13 : MODULE gw_ri_rs_grid_from_file
14 : USE cp_files, ONLY: close_file,&
15 : open_file
16 : USE kinds, ONLY: default_path_length,&
17 : default_string_length,&
18 : dp
19 : USE particle_types, ONLY: particle_type
20 : USE post_scf_bandstructure_types, ONLY: post_scf_bandstructure_type,&
21 : rirs_grid_type
22 : #include "./base/base_uses.f90"
23 :
24 : IMPLICIT NONE
25 : PRIVATE
26 :
27 : PUBLIC :: read_ri_rs_grid_from_file
28 :
29 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'gw_ri_rs_grid_from_file'
30 :
31 : CONTAINS
32 :
33 : ! **************************************************************************************************
34 : !> \brief Read pretabulated atom-relative RI-RS grids from the CP2K data files.
35 : !> \param bs_env ...
36 : ! **************************************************************************************************
37 38 : SUBROUTINE read_ri_rs_grid_from_file(bs_env)
38 : TYPE(post_scf_bandstructure_type), POINTER :: bs_env
39 :
40 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_ri_rs_grid_from_file'
41 :
42 : CHARACTER(LEN=default_path_length) :: filepath
43 : CHARACTER(LEN=default_string_length) :: element_symbol
44 : INTEGER :: handle, iatom, ikind, natom, nkind
45 38 : TYPE(particle_type), DIMENSION(:), POINTER :: particle_set
46 38 : TYPE(rirs_grid_type), ALLOCATABLE :: kind_grids(:)
47 :
48 38 : CALL timeset(routineN, handle)
49 :
50 38 : CPASSERT(ASSOCIATED(bs_env%ri_rs%atomic_kind_set))
51 38 : nkind = SIZE(bs_env%ri_rs%atomic_kind_set)
52 38 : particle_set => bs_env%ri_rs%particle_set
53 38 : CPASSERT(ASSOCIATED(particle_set))
54 38 : natom = bs_env%n_atom
55 164 : ALLOCATE (kind_grids(nkind))
56 :
57 88 : DO ikind = 1, nkind
58 50 : element_symbol = TRIM(bs_env%ri_rs%atomic_kind_set(ikind)%element_symbol)
59 : CALL get_ri_rs_grid_filepath(element_symbol, bs_env%ri_rs%grid_select, &
60 50 : bs_env%ri_rs%grid_file_suffix, filepath)
61 50 : CALL read_ri_rs_grid_file(TRIM(filepath), kind_grids(ikind)%raw_points)
62 88 : kind_grids(ikind)%npts = SIZE(kind_grids(ikind)%raw_points, 2)
63 : END DO
64 :
65 38 : IF (ALLOCATED(bs_env%ri_rs%atomic_grids)) DEALLOCATE (bs_env%ri_rs%atomic_grids)
66 202 : ALLOCATE (bs_env%ri_rs%atomic_grids(natom))
67 126 : DO iatom = 1, natom
68 88 : ikind = particle_set(iatom)%atomic_kind%kind_number
69 126 : bs_env%ri_rs%atomic_grids(iatom) = kind_grids(ikind)
70 : END DO
71 88 : DEALLOCATE (kind_grids)
72 :
73 38 : CALL timestop(handle)
74 :
75 38 : END SUBROUTINE read_ri_rs_grid_from_file
76 :
77 : ! **************************************************************************************************
78 : !> \brief Construct the path of a pretabulated RI-RS grid file.
79 : !> \param element_symbol ...
80 : !> \param grid_select ...
81 : !> \param grid_file_suffix ...
82 : !> \param filepath ...
83 : ! **************************************************************************************************
84 50 : SUBROUTINE get_ri_rs_grid_filepath(element_symbol, grid_select, grid_file_suffix, filepath)
85 : CHARACTER(LEN=*), INTENT(IN) :: element_symbol
86 : INTEGER, INTENT(IN) :: grid_select
87 : CHARACTER(LEN=*), INTENT(IN) :: grid_file_suffix
88 : CHARACTER(LEN=*), INTENT(OUT) :: filepath
89 :
90 : CHARACTER(LEN=default_string_length) :: suffix
91 :
92 94 : SELECT CASE (grid_select)
93 : CASE (1)
94 44 : suffix = "_def2-tzvp-rs.ion"
95 : CASE (2)
96 4 : suffix = "_cc-pvtz-rs.ion"
97 : CASE (3)
98 2 : IF (LEN_TRIM(grid_file_suffix) > 0) THEN
99 2 : suffix = TRIM(grid_file_suffix)
100 : ELSE
101 0 : suffix = "_rirs.ion"
102 : END IF
103 : CASE DEFAULT
104 0 : CPABORT("Unknown grid_select (1=def2-TZVPP, 2=cc-pVTZ, 3=user-provided).")
105 : END SELECT
106 50 : filepath = "ri_rs_grid/"//TRIM(element_symbol)//TRIM(suffix)
107 50 : END SUBROUTINE get_ri_rs_grid_filepath
108 :
109 : ! **************************************************************************************************
110 : !> \brief Read Cartesian grid points from the existing CP2K RI-RS .ion format.
111 : !> \param filename Complete input filename; no suffix or element-name construction is performed.
112 : !> \param points Grid points in Bohr, indexed (alpha,l).
113 : ! **************************************************************************************************
114 50 : SUBROUTINE read_ri_rs_grid_file(filename, points)
115 : CHARACTER(LEN=*), INTENT(IN) :: filename
116 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :), &
117 : INTENT(OUT) :: points
118 :
119 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_ri_rs_grid_file'
120 :
121 : CHARACTER(LEN=default_path_length) :: line
122 : INTEGER :: handle, ierr, ipoint, iunit, npoints
123 : LOGICAL :: found_points
124 :
125 50 : CALL timeset(routineN, handle)
126 : CALL open_file(file_name=TRIM(filename), unit_number=iunit, &
127 50 : file_action='READ', file_status='OLD')
128 :
129 50 : CALL read_ri_rs_grid_header(iunit, filename, npoints)
130 :
131 50 : REWIND (iunit)
132 50 : found_points = .FALSE.
133 : DO
134 1792 : READ (iunit, '(A)', IOSTAT=ierr) line
135 1792 : IF (ierr /= 0) EXIT
136 1792 : IF (INDEX(line, '<grid_points>') > 0) THEN
137 : found_points = .TRUE.
138 : EXIT
139 : END IF
140 : END DO
141 50 : IF (.NOT. found_points) THEN
142 0 : CPABORT('RI-RS .ion file has no <grid_points> block: '//TRIM(filename))
143 : END IF
144 :
145 150 : ALLOCATE (points(3, npoints))
146 10174 : DO ipoint = 1, npoints
147 10124 : READ (iunit, *, IOSTAT=ierr) points(:, ipoint)
148 10174 : IF (ierr /= 0) CPABORT('Invalid grid point in RI-RS .ion file: '//TRIM(filename))
149 : END DO
150 50 : CALL close_file(unit_number=iunit)
151 50 : CALL timestop(handle)
152 100 : END SUBROUTINE read_ri_rs_grid_file
153 :
154 : ! **************************************************************************************************
155 : !> \brief Read and validate the point count in an open RI-RS grid file.
156 : !> \param iunit ...
157 : !> \param filename ...
158 : !> \param npoints ...
159 : ! **************************************************************************************************
160 50 : SUBROUTINE read_ri_rs_grid_header(iunit, filename, npoints)
161 : INTEGER, INTENT(IN) :: iunit
162 : CHARACTER(LEN=*), INTENT(IN) :: filename
163 : INTEGER, INTENT(OUT) :: npoints
164 :
165 : CHARACTER(LEN=default_path_length) :: line
166 : INTEGER :: colon, ierr
167 : LOGICAL :: found_size
168 :
169 50 : found_size = .FALSE.
170 50 : npoints = 0
171 : DO
172 1692 : READ (iunit, '(A)', IOSTAT=ierr) line
173 1692 : IF (ierr /= 0) EXIT
174 1692 : IF (INDEX(line, 'n points') > 0) THEN
175 50 : colon = INDEX(line, ':')
176 50 : IF (colon > 0) THEN
177 50 : READ (line(colon + 1:), *, IOSTAT=ierr) npoints
178 50 : found_size = ierr == 0 .AND. npoints > 0
179 : END IF
180 : EXIT
181 : END IF
182 : END DO
183 : IF (.NOT. found_size) THEN
184 0 : CPABORT('RI-RS .ion file has no valid n points field: '//TRIM(filename))
185 : END IF
186 50 : END SUBROUTINE read_ri_rs_grid_header
187 :
188 : END MODULE gw_ri_rs_grid_from_file
|
