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 50 : CALL get_ri_rs_grid_filepath(element_symbol, bs_env%ri_rs%grid_select, filepath)
60 50 : CALL read_ri_rs_grid_file(TRIM(filepath), kind_grids(ikind)%raw_points)
61 88 : kind_grids(ikind)%npts = SIZE(kind_grids(ikind)%raw_points, 2)
62 : END DO
63 :
64 38 : IF (ALLOCATED(bs_env%ri_rs%atomic_grids)) DEALLOCATE (bs_env%ri_rs%atomic_grids)
65 202 : ALLOCATE (bs_env%ri_rs%atomic_grids(natom))
66 126 : DO iatom = 1, natom
67 88 : ikind = particle_set(iatom)%atomic_kind%kind_number
68 126 : bs_env%ri_rs%atomic_grids(iatom) = kind_grids(ikind)
69 : END DO
70 88 : DEALLOCATE (kind_grids)
71 :
72 38 : CALL timestop(handle)
73 :
74 38 : END SUBROUTINE read_ri_rs_grid_from_file
75 :
76 : ! **************************************************************************************************
77 : !> \brief Construct the path of a pretabulated RI-RS grid file.
78 : !> \param element_symbol ...
79 : !> \param grid_select ...
80 : !> \param filepath ...
81 : ! **************************************************************************************************
82 50 : SUBROUTINE get_ri_rs_grid_filepath(element_symbol, grid_select, filepath)
83 : CHARACTER(LEN=*), INTENT(IN) :: element_symbol
84 : INTEGER, INTENT(IN) :: grid_select
85 : CHARACTER(LEN=*), INTENT(OUT) :: filepath
86 :
87 : CHARACTER(LEN=default_string_length) :: suffix
88 :
89 94 : SELECT CASE (grid_select)
90 : CASE (1)
91 44 : suffix = "_def2-tzvp-rs.ion"
92 : CASE (2)
93 4 : suffix = "_cc-pvtz-rs.ion"
94 : CASE (3)
95 2 : suffix = "_rirs.ion"
96 : CASE DEFAULT
97 0 : CPABORT("Unknown grid_select (1=def2-TZVPP, 2=cc-pVTZ, 3=user-provided).")
98 : END SELECT
99 50 : filepath = "ri_rs_grid/"//TRIM(element_symbol)//TRIM(suffix)
100 50 : END SUBROUTINE get_ri_rs_grid_filepath
101 :
102 : ! **************************************************************************************************
103 : !> \brief Read Cartesian grid points from the existing CP2K RI-RS .ion format.
104 : !> \param filename Complete input filename; no suffix or element-name construction is performed.
105 : !> \param points Grid points in Bohr, indexed (alpha,l).
106 : ! **************************************************************************************************
107 50 : SUBROUTINE read_ri_rs_grid_file(filename, points)
108 : CHARACTER(LEN=*), INTENT(IN) :: filename
109 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:, :), &
110 : INTENT(OUT) :: points
111 :
112 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_ri_rs_grid_file'
113 :
114 : CHARACTER(LEN=default_path_length) :: line
115 : INTEGER :: handle, ierr, ipoint, iunit, npoints
116 : LOGICAL :: found_points
117 :
118 50 : CALL timeset(routineN, handle)
119 : CALL open_file(file_name=TRIM(filename), unit_number=iunit, &
120 50 : file_action='READ', file_status='OLD')
121 :
122 50 : CALL read_ri_rs_grid_header(iunit, filename, npoints)
123 :
124 50 : REWIND (iunit)
125 50 : found_points = .FALSE.
126 : DO
127 1792 : READ (iunit, '(A)', IOSTAT=ierr) line
128 1792 : IF (ierr /= 0) EXIT
129 1792 : IF (INDEX(line, '<grid_points>') > 0) THEN
130 : found_points = .TRUE.
131 : EXIT
132 : END IF
133 : END DO
134 50 : IF (.NOT. found_points) THEN
135 0 : CPABORT('RI-RS .ion file has no <grid_points> block: '//TRIM(filename))
136 : END IF
137 :
138 150 : ALLOCATE (points(3, npoints))
139 10174 : DO ipoint = 1, npoints
140 10124 : READ (iunit, *, IOSTAT=ierr) points(:, ipoint)
141 10174 : IF (ierr /= 0) CPABORT('Invalid grid point in RI-RS .ion file: '//TRIM(filename))
142 : END DO
143 50 : CALL close_file(unit_number=iunit)
144 50 : CALL timestop(handle)
145 100 : END SUBROUTINE read_ri_rs_grid_file
146 :
147 : ! **************************************************************************************************
148 : !> \brief Read and validate the point count in an open RI-RS grid file.
149 : !> \param iunit ...
150 : !> \param filename ...
151 : !> \param npoints ...
152 : ! **************************************************************************************************
153 50 : SUBROUTINE read_ri_rs_grid_header(iunit, filename, npoints)
154 : INTEGER, INTENT(IN) :: iunit
155 : CHARACTER(LEN=*), INTENT(IN) :: filename
156 : INTEGER, INTENT(OUT) :: npoints
157 :
158 : CHARACTER(LEN=default_path_length) :: line
159 : INTEGER :: colon, ierr
160 : LOGICAL :: found_size
161 :
162 50 : found_size = .FALSE.
163 50 : npoints = 0
164 : DO
165 1692 : READ (iunit, '(A)', IOSTAT=ierr) line
166 1692 : IF (ierr /= 0) EXIT
167 1692 : IF (INDEX(line, 'n points') > 0) THEN
168 50 : colon = INDEX(line, ':')
169 50 : IF (colon > 0) THEN
170 50 : READ (line(colon + 1:), *, IOSTAT=ierr) npoints
171 50 : found_size = ierr == 0 .AND. npoints > 0
172 : END IF
173 : EXIT
174 : END IF
175 : END DO
176 : IF (.NOT. found_size) THEN
177 0 : CPABORT('RI-RS .ion file has no valid n points field: '//TRIM(filename))
178 : END IF
179 50 : END SUBROUTINE read_ri_rs_grid_header
180 :
181 : END MODULE gw_ri_rs_grid_from_file
|