Line data Source code
1 : !--------------------------------------------------------------------------------------------------!
2 : ! CP2K: A general program to perform molecular dynamics simulations !
3 : ! Copyright 2000-2024 CP2K developers group <https://cp2k.org> !
4 : ! !
5 : ! SPDX-License-Identifier: GPL-2.0-or-later !
6 : !--------------------------------------------------------------------------------------------------!
7 :
8 : ! **************************************************************************************************
9 : !> \brief - writing and printing the files, trajectory (pos, cell, dipoles) as
10 : !> well as restart files
11 : !> - usually just the Markov Chain elements are regarded, the elements
12 : !> beside this trajectory are neglected
13 : !> - futrthermore (by option) just the accepted configurations
14 : !> are print out to reduce the file sizes
15 : !> \par History
16 : !> 12.2012 created [Mandes Schoenherr]
17 : !> \author Mandes
18 : ! **************************************************************************************************
19 :
20 : MODULE tmc_file_io
21 : USE cp_files, ONLY: close_file,&
22 : open_file
23 : USE cp_log_handling, ONLY: cp_to_string
24 : USE kinds, ONLY: default_path_length,&
25 : default_string_length,&
26 : dp
27 : USE physcon, ONLY: au2a => angstrom
28 : USE tmc_analysis_types, ONLY: tmc_analysis_env
29 : USE tmc_calculations, ONLY: get_cell_scaling,&
30 : get_scaled_cell
31 : USE tmc_move_types, ONLY: nr_mv_types
32 : USE tmc_stati, ONLY: TMC_STATUS_FAILED,&
33 : TMC_STATUS_OK,&
34 : TMC_STATUS_WAIT_FOR_NEW_TASK,&
35 : tmc_default_restart_in_file_name,&
36 : tmc_default_restart_out_file_name,&
37 : tmc_default_trajectory_file_name
38 : USE tmc_tree_types, ONLY: elem_array_type,&
39 : tree_type
40 : USE tmc_types, ONLY: tmc_env_type,&
41 : tmc_param_type
42 : #include "../base/base_uses.f90"
43 :
44 : IMPLICIT NONE
45 :
46 : PRIVATE
47 :
48 : CHARACTER(len=*), PARAMETER, PRIVATE :: moduleN = 'tmc_file_io'
49 :
50 : ! filename manipulation
51 : PUBLIC :: expand_file_name_char, expand_file_name_temp, expand_file_name_int
52 : ! read/write restart file
53 : PUBLIC :: print_restart_file, read_restart_file
54 : ! write the configuration
55 : PUBLIC :: write_result_list_element
56 : PUBLIC :: write_element_in_file
57 : PUBLIC :: write_dipoles_in_file
58 : ! analysis read
59 : PUBLIC :: analyse_files_open, read_element_from_file, analyse_files_close
60 :
61 : CONTAINS
62 :
63 : !------------------------------------------------------------------------------
64 : ! routines for manipulating the file name
65 : !------------------------------------------------------------------------------
66 : ! **************************************************************************************************
67 : !> \brief placing a character string at the end of a file name
68 : !> (instead of the ending)
69 : !> \param file_name original file name
70 : !> \param extra string to be added before the file extension
71 : !> \return the new filename
72 : !> \author Mandes 11.2012
73 : ! **************************************************************************************************
74 2615 : FUNCTION expand_file_name_ending(file_name, extra) RESULT(result_file_name)
75 : CHARACTER(LEN=*) :: file_name, extra
76 : CHARACTER(LEN=default_path_length) :: result_file_name
77 :
78 : INTEGER :: ind
79 :
80 0 : CPASSERT(file_name .NE. "")
81 :
82 2615 : ind = INDEX(file_name, ".", BACK=.TRUE.)
83 2615 : IF (.NOT. ind .EQ. 0) THEN
84 2615 : WRITE (result_file_name, *) file_name(1:ind - 1), ".", &
85 5230 : TRIM(ADJUSTL(extra))
86 : ELSE
87 0 : WRITE (result_file_name, *) TRIM(file_name), ".", extra
88 : END IF
89 2615 : result_file_name = TRIM(ADJUSTL(result_file_name))
90 2615 : CPASSERT(result_file_name .NE. "")
91 2615 : END FUNCTION expand_file_name_ending
92 :
93 : ! **************************************************************************************************
94 : !> \brief placing a character string at the end of a file name
95 : !> (before the file extension)
96 : !> \param file_name original file name
97 : !> \param extra string to be added before the file extension
98 : !> \return the new filename
99 : !> \author Mandes 11.2012
100 : ! **************************************************************************************************
101 1427 : FUNCTION expand_file_name_char(file_name, extra) RESULT(result_file_name)
102 : CHARACTER(LEN=*) :: file_name, extra
103 : CHARACTER(LEN=default_path_length) :: result_file_name
104 :
105 : INTEGER :: ind
106 :
107 0 : CPASSERT(file_name .NE. "")
108 :
109 1427 : ind = INDEX(file_name, ".", BACK=.TRUE.)
110 1427 : IF (.NOT. ind .EQ. 0) THEN
111 1427 : WRITE (result_file_name, *) file_name(1:ind - 1), "_", &
112 2854 : TRIM(ADJUSTL(extra)), file_name(ind:LEN_TRIM(file_name))
113 : ELSE
114 0 : WRITE (result_file_name, *) TRIM(file_name), "_", extra
115 : END IF
116 1427 : result_file_name = TRIM(ADJUSTL(result_file_name))
117 1427 : CPASSERT(result_file_name .NE. "")
118 1427 : END FUNCTION expand_file_name_char
119 :
120 : ! **************************************************************************************************
121 : !> \brief placing the temperature at the end of a file name
122 : !> (before the file extension)
123 : !> \param file_name original file name
124 : !> \param rvalue temperature to be added
125 : !> \return the new filename
126 : !> \author Mandes 11.2012
127 : ! **************************************************************************************************
128 2871 : FUNCTION expand_file_name_temp(file_name, rvalue) RESULT(result_file_name)
129 : CHARACTER(LEN=*) :: file_name
130 : REAL(KIND=dp) :: rvalue
131 : CHARACTER(LEN=default_path_length) :: result_file_name
132 :
133 : CHARACTER(LEN=18) :: rval_to_string
134 : INTEGER :: ind
135 :
136 2871 : CPASSERT(file_name .NE. "")
137 :
138 2871 : rval_to_string = ""
139 :
140 2871 : WRITE (rval_to_string, "(F16.2)") rvalue
141 2871 : ind = INDEX(file_name, ".", BACK=.TRUE.)
142 2871 : IF (.NOT. ind .EQ. 0) THEN
143 2871 : WRITE (result_file_name, *) file_name(1:ind - 1), "_T", &
144 5742 : TRIM(ADJUSTL(rval_to_string)), file_name(ind:LEN_TRIM(file_name))
145 : ELSE
146 0 : IF (LEN(file_name) .EQ. 0) THEN
147 0 : WRITE (result_file_name, *) TRIM(file_name), "T", TRIM(ADJUSTL(rval_to_string)), &
148 0 : file_name(ind:LEN_TRIM(file_name))
149 : ELSE
150 0 : WRITE (result_file_name, *) TRIM(file_name), "_T", TRIM(ADJUSTL(rval_to_string))
151 : END IF
152 : END IF
153 2871 : result_file_name = TRIM(ADJUSTL(result_file_name))
154 2871 : CPASSERT(result_file_name .NE. "")
155 2871 : END FUNCTION expand_file_name_temp
156 :
157 : ! **************************************************************************************************
158 : !> \brief placing an integer at the end of a file name
159 : !> (before the file extension)
160 : !> \param file_name original file name
161 : !> \param ivalue number to be added
162 : !> \return the new filename
163 : !> \author Mandes 11.2012
164 : ! **************************************************************************************************
165 19 : FUNCTION expand_file_name_int(file_name, ivalue) RESULT(result_file_name)
166 : CHARACTER(LEN=*) :: file_name
167 : INTEGER :: ivalue
168 : CHARACTER(LEN=default_path_length) :: result_file_name
169 :
170 : CHARACTER(LEN=18) :: rval_to_string
171 : INTEGER :: ind
172 :
173 19 : CPASSERT(file_name .NE. "")
174 :
175 19 : rval_to_string = ""
176 :
177 19 : WRITE (rval_to_string, *) ivalue
178 19 : ind = INDEX(file_name, ".", BACK=.TRUE.)
179 19 : IF (.NOT. ind .EQ. 0) THEN
180 19 : WRITE (result_file_name, *) file_name(1:ind - 1), "_", &
181 38 : TRIM(ADJUSTL(rval_to_string)), file_name(ind:LEN_TRIM(file_name))
182 : ELSE
183 0 : IF (LEN(file_name) .EQ. 0) THEN
184 0 : WRITE (result_file_name, *) TRIM(file_name), "", TRIM(ADJUSTL(rval_to_string)), &
185 0 : file_name(ind:LEN_TRIM(file_name))
186 : ELSE
187 0 : WRITE (result_file_name, *) TRIM(file_name), "_", TRIM(ADJUSTL(rval_to_string)), &
188 0 : file_name(ind:LEN_TRIM(file_name))
189 : END IF
190 : END IF
191 19 : result_file_name = TRIM(ADJUSTL(result_file_name))
192 19 : CPASSERT(result_file_name .NE. "")
193 19 : END FUNCTION expand_file_name_int
194 :
195 : !------------------------------------------------------------------------------
196 : ! routines for reading and writing RESTART file
197 : !------------------------------------------------------------------------------
198 : ! **************************************************************************************************
199 : !> \brief prints out the TMC restart files with all last configurations and
200 : !> counters etc.
201 : !> \param tmc_env the tmc environment, storing result lists and counters an in
202 : !> temperatures
203 : !> \param job_counts the counters for counting the submitted different job types
204 : !> \param timings ...
205 : !> \author Mandes 11.2012
206 : ! **************************************************************************************************
207 3 : SUBROUTINE print_restart_file(tmc_env, job_counts, timings)
208 : TYPE(tmc_env_type), POINTER :: tmc_env
209 : INTEGER, DIMENSION(:) :: job_counts
210 : REAL(KIND=dp), DIMENSION(4) :: timings
211 :
212 : CHARACTER(LEN=default_path_length) :: c_tmp, file_name
213 : INTEGER :: f_unit, i
214 :
215 3 : c_tmp = ""
216 3 : CPASSERT(ASSOCIATED(tmc_env))
217 3 : CPASSERT(ASSOCIATED(tmc_env%m_env))
218 3 : CPASSERT(ASSOCIATED(tmc_env%params))
219 3 : CPASSERT(ASSOCIATED(tmc_env%m_env%gt_act))
220 :
221 3 : WRITE (c_tmp, FMT='(I9.9)') tmc_env%m_env%result_count(0)
222 : file_name = TRIM(expand_file_name_char( &
223 : file_name=tmc_default_restart_out_file_name, &
224 3 : extra=c_tmp))
225 : CALL open_file(file_name=file_name, file_status="REPLACE", &
226 : file_action="WRITE", file_form="UNFORMATTED", &
227 3 : unit_number=f_unit)
228 3 : WRITE (f_unit) SIZE(tmc_env%params%Temp)
229 3 : WRITE (f_unit) tmc_env%params%Temp(:), &
230 3 : tmc_env%m_env%gt_act%nr, &
231 3 : tmc_env%m_env%gt_act%rng_seed, &
232 3 : tmc_env%m_env%gt_act%rnd_nr, &
233 3 : tmc_env%m_env%gt_act%prob_acc, &
234 3 : tmc_env%m_env%gt_act%mv_conf, &
235 3 : tmc_env%m_env%gt_act%mv_next_conf, &
236 3 : tmc_env%m_env%result_count(0:), &
237 3 : tmc_env%params%move_types%mv_weight, &
238 3 : tmc_env%params%move_types%acc_count, &
239 3 : tmc_env%params%move_types%mv_count, &
240 3 : tmc_env%params%move_types%subbox_acc_count, &
241 3 : tmc_env%params%move_types%subbox_count, &
242 3 : tmc_env%params%cell%hmat, &
243 3 : job_counts, &
244 6 : timings
245 12 : DO i = 1, SIZE(tmc_env%params%Temp)
246 9 : WRITE (f_unit) tmc_env%m_env%result_list(i)%elem%nr, &
247 252 : tmc_env%m_env%result_list(i)%elem%rng_seed, &
248 576 : tmc_env%m_env%result_list(i)%elem%pos, &
249 576 : tmc_env%m_env%result_list(i)%elem%vel, &
250 36 : tmc_env%m_env%result_list(i)%elem%box_scale, &
251 9 : tmc_env%m_env%result_list(i)%elem%potential, &
252 9 : tmc_env%m_env%result_list(i)%elem%e_pot_approx, &
253 9 : tmc_env%m_env%result_list(i)%elem%ekin, &
254 9 : tmc_env%m_env%result_list(i)%elem%ekin_before_md, &
255 21 : tmc_env%m_env%result_list(i)%elem%temp_created
256 : END DO
257 3 : CALL close_file(unit_number=f_unit)
258 : ! write the file, where the restart file name is written in
259 : CALL open_file(file_name=tmc_default_restart_in_file_name, &
260 : file_action="WRITE", file_status="REPLACE", &
261 3 : unit_number=f_unit)
262 3 : WRITE (f_unit, *) TRIM(file_name)
263 3 : CALL close_file(unit_number=f_unit)
264 3 : END SUBROUTINE print_restart_file
265 :
266 : ! **************************************************************************************************
267 : !> \brief reads the TMC restart file with all last configurations and
268 : !> counters etc.
269 : !> \param tmc_env the tmc environment, storing result lists and counters an in
270 : !> temperatures
271 : !> \param job_counts the counters for counting the submitted different job types
272 : !> \param timings ...
273 : !> \param file_name the restart file name
274 : !> \author Mandes 11.2012
275 : ! **************************************************************************************************
276 2 : SUBROUTINE read_restart_file(tmc_env, job_counts, timings, file_name)
277 : TYPE(tmc_env_type), POINTER :: tmc_env
278 : INTEGER, DIMENSION(:) :: job_counts
279 : REAL(KIND=dp), DIMENSION(4) :: timings
280 : CHARACTER(LEN=*) :: file_name
281 :
282 : INTEGER :: file_ptr, i, temp_size
283 : LOGICAL :: flag
284 2 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: tmp_temp
285 : REAL(KIND=dp), DIMENSION(nr_mv_types) :: mv_weight_tmp
286 :
287 2 : CPASSERT(ASSOCIATED(tmc_env))
288 2 : CPASSERT(ASSOCIATED(tmc_env%m_env))
289 2 : CPASSERT(ASSOCIATED(tmc_env%params))
290 2 : CPASSERT(ASSOCIATED(tmc_env%m_env%gt_act))
291 :
292 2 : IF (file_name .EQ. tmc_default_restart_in_file_name) THEN
293 2 : INQUIRE (FILE=tmc_default_restart_in_file_name, EXIST=flag)
294 2 : CPASSERT(flag)
295 : CALL open_file(file_name=tmc_default_restart_in_file_name, file_status="OLD", &
296 2 : file_action="READ", unit_number=file_ptr)
297 2 : READ (file_ptr, *) file_name
298 2 : CALL close_file(unit_number=file_ptr)
299 : END IF
300 :
301 : CALL open_file(file_name=file_name, file_status="OLD", file_form="UNFORMATTED", &
302 2 : file_action="READ", unit_number=file_ptr)
303 2 : READ (file_ptr) temp_size
304 2 : IF (temp_size .NE. SIZE(tmc_env%params%Temp)) &
305 : CALL cp_abort(__LOCATION__, &
306 : "the actual specified temperatures does not "// &
307 0 : "fit in amount with the one from restart file ")
308 6 : ALLOCATE (tmp_temp(temp_size))
309 2 : READ (file_ptr) tmp_temp(:), &
310 2 : tmc_env%m_env%gt_act%nr, &
311 2 : tmc_env%m_env%gt_act%rng_seed, &
312 2 : tmc_env%m_env%gt_act%rnd_nr, &
313 2 : tmc_env%m_env%gt_act%prob_acc, &
314 2 : tmc_env%m_env%gt_act%mv_conf, & !
315 2 : tmc_env%m_env%gt_act%mv_next_conf, & !
316 2 : tmc_env%m_env%result_count(0:), &
317 2 : mv_weight_tmp, & !
318 2 : tmc_env%params%move_types%acc_count, &
319 2 : tmc_env%params%move_types%mv_count, &
320 2 : tmc_env%params%move_types%subbox_acc_count, &
321 2 : tmc_env%params%move_types%subbox_count, & !
322 2 : tmc_env%params%cell%hmat, &
323 2 : job_counts, &
324 4 : timings
325 :
326 8 : IF (ANY(ABS(tmc_env%params%Temp(:) - tmp_temp(:)) .GE. 0.005)) &
327 : CALL cp_abort(__LOCATION__, "the temperatures differ from the previous calculation. "// &
328 0 : "There were the following temperatures used:")
329 22 : IF (ANY(mv_weight_tmp(:) .NE. tmc_env%params%move_types%mv_weight(:))) THEN
330 0 : CPWARN("The amount of mv types differs between the original and the restart run.")
331 : END IF
332 :
333 8 : DO i = 1, SIZE(tmc_env%params%Temp)
334 6 : tmc_env%m_env%gt_act%conf(i)%elem => tmc_env%m_env%result_list(i)%elem
335 6 : READ (file_ptr) tmc_env%m_env%result_list(i)%elem%nr, &
336 168 : tmc_env%m_env%result_list(i)%elem%rng_seed, &
337 384 : tmc_env%m_env%result_list(i)%elem%pos, &
338 384 : tmc_env%m_env%result_list(i)%elem%vel, &
339 24 : tmc_env%m_env%result_list(i)%elem%box_scale, &
340 6 : tmc_env%m_env%result_list(i)%elem%potential, &
341 6 : tmc_env%m_env%result_list(i)%elem%e_pot_approx, &
342 6 : tmc_env%m_env%result_list(i)%elem%ekin, &
343 6 : tmc_env%m_env%result_list(i)%elem%ekin_before_md, &
344 14 : tmc_env%m_env%result_list(i)%elem%temp_created
345 : END DO
346 2 : CALL close_file(unit_number=file_ptr)
347 2 : END SUBROUTINE read_restart_file
348 :
349 : !----------------------------------------------------------------------------
350 : ! printing configuration in file
351 : !----------------------------------------------------------------------------
352 :
353 : ! **************************************************************************************************
354 : !> \brief select the correct configuration to print out the
355 : !> (coordinates, forces, cell ...)
356 : !> \param result_list list of configurations for each temperature
357 : !> \param result_count list with number of Markov Chain number
358 : !> for each teperature (index 0 for global tree)
359 : !> \param conf_updated index of the updated (modified element)
360 : !> \param accepted acceptance flag
361 : !> \param tmc_params TMC environment parameters
362 : !> \author Mandes 02.2013
363 : ! **************************************************************************************************
364 9350 : SUBROUTINE write_result_list_element(result_list, result_count, conf_updated, &
365 : accepted, tmc_params)
366 : TYPE(elem_array_type), DIMENSION(:), POINTER :: result_list
367 : INTEGER, DIMENSION(:), POINTER :: result_count
368 : INTEGER :: conf_updated
369 : LOGICAL, INTENT(IN) :: accepted
370 : TYPE(tmc_param_type), POINTER :: tmc_params
371 :
372 : CHARACTER(LEN=*), PARAMETER :: routineN = 'write_result_list_element'
373 :
374 : CHARACTER(LEN=default_path_length) :: file_name
375 : INTEGER :: handle, i
376 :
377 4675 : file_name = ""
378 :
379 4675 : CPASSERT(ASSOCIATED(result_list))
380 4675 : CPASSERT(ASSOCIATED(result_count))
381 4675 : CPASSERT(ASSOCIATED(tmc_params))
382 4675 : CPASSERT(ASSOCIATED(tmc_params%Temp))
383 4675 : CPASSERT(conf_updated .GE. 0)
384 4675 : CPASSERT(conf_updated .LE. SIZE(tmc_params%Temp))
385 :
386 : ! start the timing
387 4675 : CALL timeset(routineN, handle)
388 :
389 4675 : IF (conf_updated .EQ. 0) THEN
390 : ! for debugging print every configuration of every temperature
391 0 : DO i = 1, SIZE(tmc_params%Temp)
392 0 : WRITE (file_name, *) "every_step_", TRIM(tmc_default_trajectory_file_name)
393 : CALL write_element_in_file(elem=result_list(i)%elem, &
394 : tmc_params=tmc_params, conf_nr=result_count(0), &
395 0 : file_name=expand_file_name_temp(file_name=file_name, rvalue=tmc_params%Temp(i)))
396 : END DO
397 : ELSE
398 4675 : IF ((.NOT. tmc_params%print_only_diff_conf) .OR. &
399 : (tmc_params%print_only_diff_conf .AND. accepted)) THEN
400 : CALL write_element_in_file(elem=result_list(conf_updated)%elem, &
401 : tmc_params=tmc_params, conf_nr=result_count(conf_updated), &
402 : file_name=expand_file_name_temp(file_name=TRIM(tmc_default_trajectory_file_name), &
403 1025 : rvalue=tmc_params%Temp(conf_updated)))
404 : END IF
405 : END IF
406 : ! end the timing
407 4675 : CALL timestop(handle)
408 4675 : END SUBROUTINE write_result_list_element
409 :
410 : ! **************************************************************************************************
411 : !> \brief writes the trajectory element in a file from sub tree element
412 : !> \param elem actual tree element to be printed out
413 : !> \param tmc_params TMC environment parameters
414 : !> \param temp_index ...
415 : !> \param file_name file name will be extended by type of file (pos, cell,...)
416 : !> \param conf_nr Markov chain element number
417 : !> \param conf_info whole header line
418 : !> \author Mandes 11.2012
419 : ! **************************************************************************************************
420 1025 : SUBROUTINE write_element_in_file(elem, tmc_params, temp_index, file_name, conf_nr, &
421 : conf_info)
422 : TYPE(tree_type), POINTER :: elem
423 : TYPE(tmc_param_type), POINTER :: tmc_params
424 : INTEGER, OPTIONAL :: temp_index
425 : CHARACTER(LEN=*), OPTIONAL :: file_name
426 : INTEGER, OPTIONAL :: conf_nr
427 : CHARACTER(LEN=*), OPTIONAL :: conf_info
428 :
429 : CHARACTER(LEN=*), PARAMETER :: routineN = 'write_element_in_file'
430 :
431 : CHARACTER(LEN=default_path_length) :: file_name_act, tmp_name
432 : CHARACTER(LEN=default_string_length) :: header
433 : INTEGER :: file_ptr, handle, i, nr_atoms
434 : LOGICAL :: file_exists, print_it
435 : REAL(KIND=dp) :: vol
436 : REAL(KIND=dp), DIMENSION(3, 3) :: hmat_scaled
437 :
438 1025 : file_name_act = ""
439 1025 : tmp_name = ""
440 1025 : header = ""
441 1025 : print_it = .TRUE.
442 :
443 0 : CPASSERT(ASSOCIATED(elem))
444 1025 : CPASSERT(ASSOCIATED(tmc_params))
445 1025 : CPASSERT(ASSOCIATED(tmc_params%atoms))
446 1025 : CPASSERT(PRESENT(conf_nr) .OR. PRESENT(conf_info))
447 :
448 1025 : IF (print_it) THEN
449 : ! start the timing
450 1025 : CALL timeset(routineN, handle)
451 :
452 : ! set default file name
453 1025 : IF (PRESENT(file_name)) THEN
454 1025 : CPASSERT(file_name .NE. "")
455 1025 : file_name_act = file_name
456 : ELSE
457 0 : CPASSERT(ASSOCIATED(tmc_params%Temp))
458 0 : CPASSERT(PRESENT(temp_index))
459 : file_name_act = expand_file_name_temp(file_name=tmc_default_trajectory_file_name, &
460 0 : rvalue=tmc_params%Temp(temp_index))
461 : END IF
462 :
463 1025 : nr_atoms = SIZE(elem%pos)/tmc_params%dim_per_elem
464 :
465 : ! set header (for coordinate or force file)
466 1025 : IF (tmc_params%print_trajectory .OR. tmc_params%print_forces) THEN
467 1025 : IF (PRESENT(conf_info)) THEN
468 0 : WRITE (header, *) TRIM(ADJUSTL(conf_info))
469 : ELSE
470 : !WRITE(header,FMT="(A,I8,A,F20.10)") " i = ", conf_nr,", E = ", elem%potential
471 1025 : WRITE (header, FMT="(A,I8,A,F20.10,F20.10,A,I8,I8)") "i =", conf_nr, " ,E =", &
472 2050 : elem%potential, elem%ekin, " st elem", elem%sub_tree_nr, elem%nr
473 : END IF
474 : END IF
475 :
476 : ! write the coordinates
477 1025 : IF (tmc_params%print_trajectory) THEN
478 1025 : tmp_name = expand_file_name_ending(file_name_act, "xyz")
479 : CALL open_file(file_name=tmp_name, file_status="UNKNOWN", &
480 : file_action="WRITE", file_position="APPEND", &
481 1025 : unit_number=file_ptr)
482 1025 : WRITE (file_ptr, FMT="(I8)") nr_atoms
483 1025 : WRITE (file_ptr, *) TRIM(header)
484 44797 : DO i = 1, SIZE(elem%pos), tmc_params%dim_per_elem
485 : WRITE (file_ptr, FMT="(A4,1X,1000F20.10)") &
486 43772 : TRIM(tmc_params%atoms((i - 1)/tmc_params%dim_per_elem + 1)%name), &
487 219885 : elem%pos(i:i + tmc_params%dim_per_elem - 1)*au2a
488 : END DO
489 1025 : CALL close_file(unit_number=file_ptr)
490 : END IF
491 :
492 : ! write the forces
493 1025 : IF (tmc_params%print_forces) THEN
494 331 : tmp_name = expand_file_name_ending(file_name_act, "frc")
495 : CALL open_file(file_name=tmp_name, file_status="UNKNOWN", &
496 : file_action="WRITE", file_position="APPEND", &
497 331 : unit_number=file_ptr)
498 331 : WRITE (file_ptr, FMT="(I8)") nr_atoms
499 331 : WRITE (file_ptr, *) TRIM(header)
500 7282 : DO i = 1, SIZE(elem%pos), tmc_params%dim_per_elem
501 : WRITE (file_ptr, FMT="(A4,1X,1000F20.10)") &
502 6951 : TRIM(tmc_params%atoms((i - 1)/tmc_params%dim_per_elem + 1)%name), &
503 14233 : elem%frc(i:i + tmc_params%dim_per_elem - 1)
504 : END DO
505 331 : CALL close_file(unit_number=file_ptr)
506 : END IF
507 :
508 : ! write the cell dipoles
509 1025 : IF (tmc_params%print_dipole) THEN
510 : CALL write_dipoles_in_file(file_name=file_name_act, &
511 0 : conf_nr=conf_nr, dip=elem%dipole)
512 : END IF
513 :
514 : ! write the cell file
515 1025 : IF (tmc_params%print_cell) THEN
516 392 : tmp_name = expand_file_name_ending(file_name_act, "cell")
517 : ! header
518 392 : INQUIRE (FILE=tmp_name, EXIST=file_exists) ! file_exists will be TRUE if the file exist
519 392 : IF (.NOT. file_exists) THEN
520 : CALL open_file(file_name=tmp_name, file_status="NEW", &
521 6 : file_action="WRITE", unit_number=file_ptr)
522 : WRITE (file_ptr, FMT='(A,9(7X,A2," [Angstrom]"),6X,A)') &
523 6 : "# MC step ", "Ax", "Ay", "Az", "Bx", "By", "Bz", "Cx", "Cy", "Cz", &
524 12 : "Volume [Angstrom^3]"
525 : ELSE
526 : CALL open_file(file_name=tmp_name, file_status="OLD", &
527 : file_action="WRITE", file_position="APPEND", &
528 386 : unit_number=file_ptr)
529 : END IF
530 : CALL get_scaled_cell(cell=tmc_params%cell, &
531 : box_scale=elem%box_scale, scaled_hmat=hmat_scaled, &
532 392 : vol=vol)
533 392 : WRITE (file_ptr, FMT="(I8,9(1X,F19.10),1X,F24.10)") conf_nr, &
534 5488 : hmat_scaled(:, :)*au2a, vol*au2a**3
535 : !TODO better cell output e.g. using cell_types routine
536 392 : CALL close_file(unit_number=file_ptr)
537 : END IF
538 :
539 : ! write the different energies
540 1025 : IF (tmc_params%print_energies) THEN
541 331 : tmp_name = expand_file_name_ending(file_name_act, "ener")
542 : ! header
543 331 : INQUIRE (FILE=tmp_name, EXIST=file_exists) ! file_exists will be TRUE if the file exist
544 331 : IF (.NOT. file_exists) THEN
545 : CALL open_file(file_name=tmp_name, file_status="NEW", &
546 3 : file_action="WRITE", unit_number=file_ptr)
547 : WRITE (file_ptr, FMT='(A,4A20)') &
548 3 : "# MC step ", " exact ", " approx ", " last SCF ", " kinetic "
549 : ELSE
550 : CALL open_file(file_name=tmp_name, file_status="OLD", &
551 : file_action="WRITE", file_position="APPEND", &
552 328 : unit_number=file_ptr)
553 : END IF
554 331 : WRITE (file_ptr, FMT="(I8,14F20.10)") conf_nr, elem%potential, elem%e_pot_approx, &
555 662 : elem%scf_energies(MOD(elem%scf_energies_count, 4) + 1), elem%ekin
556 331 : CALL close_file(unit_number=file_ptr)
557 : END IF
558 :
559 : ! end the timing
560 1025 : CALL timestop(handle)
561 : END IF
562 1025 : END SUBROUTINE write_element_in_file
563 :
564 : ! **************************************************************************************************
565 : !> \brief writes the cell dipoles in dipole trajectory file
566 : !> \param file_name ...
567 : !> \param conf_nr ...
568 : !> \param dip ...
569 : !> \param file_ext ...
570 : !> \param
571 : !> \author Mandes 11.2012
572 : ! **************************************************************************************************
573 500 : SUBROUTINE write_dipoles_in_file(file_name, conf_nr, dip, file_ext)
574 : CHARACTER(LEN=default_path_length) :: file_name
575 : INTEGER :: conf_nr
576 : REAL(KIND=dp), DIMENSION(:), POINTER :: dip
577 : CHARACTER(LEN=*), INTENT(in), OPTIONAL :: file_ext
578 :
579 : CHARACTER(LEN=default_path_length) :: file_name_tmp
580 : INTEGER :: file_ptr
581 : LOGICAL :: file_exists
582 :
583 500 : CPASSERT(ASSOCIATED(dip))
584 :
585 500 : IF (PRESENT(file_ext)) THEN
586 500 : CPASSERT(file_ext .NE. "")
587 500 : file_name_tmp = expand_file_name_ending(file_name, TRIM(file_ext))
588 : ELSE
589 0 : file_name_tmp = expand_file_name_ending(file_name, "dip")
590 : END IF
591 500 : INQUIRE (FILE=file_name_tmp, EXIST=file_exists)
592 500 : IF (.NOT. file_exists) THEN
593 : CALL open_file(file_name=file_name_tmp, file_status="NEW", &
594 3 : file_action="WRITE", unit_number=file_ptr)
595 3 : WRITE (file_ptr, FMT='(A8,10A20)') "# conf_nr", "dip_x [C Angstrom]", &
596 6 : "dip_y [C Angstrom]", "dip_z [C Angstrom]"
597 : ELSE
598 : CALL open_file(file_name=file_name_tmp, file_status="OLD", &
599 : file_action="WRITE", file_position="APPEND", &
600 497 : unit_number=file_ptr)
601 : END IF
602 500 : WRITE (file_ptr, FMT="(I8,10F20.10)") conf_nr, dip(:)
603 500 : CALL close_file(unit_number=file_ptr)
604 500 : END SUBROUTINE write_dipoles_in_file
605 :
606 : !----------------------------------------------------------------------------
607 : ! read configuration from file
608 : !----------------------------------------------------------------------------
609 :
610 : ! **************************************************************************************************
611 : !> \brief read the trajectory element from a file from sub tree element
612 : !> \param elem actual tree element to be printed out
613 : !> \param tmc_ana TMC analysis environment parameters
614 : !> \param conf_nr Markov chain element number
615 : !> (input the old number and read only if conf nr from file is greater
616 : !> \param stat ...
617 : !> \author Mandes 03.2013
618 : ! **************************************************************************************************
619 2098 : SUBROUTINE read_element_from_file(elem, tmc_ana, conf_nr, stat)
620 : TYPE(tree_type), POINTER :: elem
621 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
622 : INTEGER :: conf_nr, stat
623 :
624 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_element_from_file'
625 :
626 : INTEGER :: conf_nr_old, handle, i_tmp
627 : LOGICAL :: files_conf_missmatch
628 :
629 1049 : stat = TMC_STATUS_OK
630 1049 : conf_nr_old = conf_nr
631 1049 : files_conf_missmatch = .FALSE.
632 :
633 1049 : CPASSERT(ASSOCIATED(elem))
634 1049 : CPASSERT(ASSOCIATED(tmc_ana))
635 1049 : CPASSERT(ASSOCIATED(tmc_ana%atoms))
636 :
637 : ! start the timing
638 1049 : CALL timeset(routineN, handle)
639 :
640 : ! read the coordinates
641 1049 : IF (tmc_ana%id_traj .GT. 0) THEN
642 1049 : i_tmp = conf_nr_old
643 : CALL read_pos_from_file(elem=elem, tmc_ana=tmc_ana, stat=stat, &
644 1049 : conf_nr=i_tmp)
645 1049 : IF (stat .EQ. TMC_STATUS_WAIT_FOR_NEW_TASK) THEN
646 : CALL cp_warn(__LOCATION__, &
647 : 'end of position file reached at line '// &
648 : cp_to_string(REAL(tmc_ana%lc_traj, KIND=dp))//", last element "// &
649 18 : cp_to_string(tmc_ana%last_elem%nr))
650 : ELSE
651 1031 : CPASSERT(i_tmp .GT. conf_nr_old)
652 1031 : conf_nr = i_tmp
653 1031 : elem%nr = i_tmp
654 : END IF
655 : END IF
656 :
657 : ! read the forces
658 : ! TODO if necessary
659 :
660 : ! read the dipoles file
661 1049 : IF (tmc_ana%id_dip .GT. 0 .AND. stat .EQ. TMC_STATUS_OK) THEN
662 0 : i_tmp = conf_nr_old
663 : search_conf_dip: DO
664 : CALL read_dipole_from_file(elem=elem, tmc_ana=tmc_ana, stat=stat, &
665 0 : conf_nr=i_tmp)
666 0 : IF (stat .EQ. TMC_STATUS_WAIT_FOR_NEW_TASK) THEN
667 : CALL cp_warn(__LOCATION__, &
668 : 'end of dipole file reached at line'// &
669 0 : cp_to_string(REAL(tmc_ana%lc_dip, KIND=dp)))
670 0 : EXIT search_conf_dip
671 : END IF
672 : ! check consitence with pos file
673 0 : IF (tmc_ana%id_traj .GT. 0) THEN
674 0 : IF (i_tmp .EQ. conf_nr) THEN
675 : files_conf_missmatch = .FALSE.
676 : EXIT search_conf_dip
677 : ELSE
678 : ! the configuration numbering differ from the position file,
679 : ! but we keep on searching for the correct configuration
680 : files_conf_missmatch = .TRUE.
681 : END IF
682 : ! if no pos file, just take the next conf
683 0 : ELSE IF (i_tmp .GT. conf_nr_old) THEN
684 0 : conf_nr = i_tmp
685 0 : elem%nr = i_tmp
686 0 : EXIT search_conf_dip
687 : END IF
688 : END DO search_conf_dip
689 : END IF
690 :
691 : ! read the cell file
692 1049 : IF (tmc_ana%id_cell .GT. 0 .AND. stat .EQ. TMC_STATUS_OK) THEN
693 : search_conf_cell: DO
694 : CALL read_cell_from_file(elem=elem, tmc_ana=tmc_ana, stat=stat, &
695 1206 : conf_nr=i_tmp)
696 1206 : IF (stat .EQ. TMC_STATUS_WAIT_FOR_NEW_TASK) THEN
697 : CALL cp_warn(__LOCATION__, &
698 : 'end of cell file reached at line at line'// &
699 0 : cp_to_string(REAL(tmc_ana%lc_cell, KIND=dp)))
700 0 : EXIT search_conf_cell
701 : END IF
702 : ! check consitence with pos file
703 1206 : IF (tmc_ana%id_traj .GT. 0) THEN
704 1206 : IF (i_tmp .EQ. conf_nr) THEN
705 : files_conf_missmatch = .FALSE.
706 : EXIT search_conf_cell
707 : ELSE
708 : ! the configuration numbering differ from the position file,
709 : ! but we keep on searching for the correct configuration
710 : files_conf_missmatch = .TRUE.
711 : END IF
712 : ! if no pos file, just take the next conf
713 0 : ELSE IF (i_tmp .GT. conf_nr_old) THEN
714 0 : conf_nr = i_tmp
715 0 : elem%nr = i_tmp
716 0 : EXIT search_conf_cell
717 : END IF
718 : END DO search_conf_cell
719 :
720 : END IF
721 :
722 : ! write the different energies
723 : ! TODO if necessary
724 :
725 1049 : IF (files_conf_missmatch) &
726 : CALL cp_warn(__LOCATION__, &
727 : 'there is a missmatch in the configuration numbering. '// &
728 : "Read number of lines (pos|cell|dip)"// &
729 : cp_to_string(tmc_ana%lc_traj)//"|"// &
730 : cp_to_string(tmc_ana%lc_cell)//"|"// &
731 0 : cp_to_string(tmc_ana%lc_dip))
732 :
733 : ! end the timing
734 1049 : CALL timestop(handle)
735 1049 : END SUBROUTINE read_element_from_file
736 :
737 : ! **************************************************************************************************
738 : !> \brief search for the next configurational position in file
739 : !> \param elem actual tree element to be read
740 : !> \param tmc_ana ...
741 : !> \param stat ...
742 : !> \param conf_nr Markov chain element number
743 : !> (input the old number and read only if conf nr from file is greater
744 : !> \param header_info ...
745 : !> \author Mandes 03.2013
746 : ! **************************************************************************************************
747 2098 : SUBROUTINE read_pos_from_file(elem, tmc_ana, stat, conf_nr, header_info)
748 : TYPE(tree_type), POINTER :: elem
749 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
750 : INTEGER :: stat, conf_nr
751 : CHARACTER(LEN=*), OPTIONAL :: header_info
752 :
753 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_pos_from_file'
754 :
755 : CHARACTER(LEN=default_string_length) :: c_tmp
756 : INTEGER :: handle, i, i_tmp, status
757 :
758 1049 : stat = TMC_STATUS_FAILED
759 :
760 0 : CPASSERT(ASSOCIATED(elem))
761 1049 : CPASSERT(ASSOCIATED(elem%pos))
762 1049 : CPASSERT(ASSOCIATED(tmc_ana))
763 1049 : CPASSERT(tmc_ana%id_traj .GT. 0)
764 :
765 : ! start the timing
766 1049 : CALL timeset(routineN, handle)
767 :
768 : search_next_conf: DO
769 6105 : c_tmp(:) = " "
770 6105 : tmc_ana%lc_traj = tmc_ana%lc_traj + 1
771 6105 : READ (tmc_ana%id_traj, '(A)', IOSTAT=status) c_tmp(:)
772 6105 : IF (status .GT. 0) &
773 : CALL cp_abort(__LOCATION__, &
774 : "configuration header read error at line: "// &
775 0 : cp_to_string(tmc_ana%lc_traj)//": "//c_tmp)
776 6105 : IF (status .LT. 0) THEN ! end of file reached
777 18 : stat = TMC_STATUS_WAIT_FOR_NEW_TASK
778 18 : EXIT search_next_conf
779 : END IF
780 6087 : IF (INDEX(c_tmp, "=") .GT. 0) THEN
781 1206 : READ (c_tmp(INDEX(c_tmp, "=") + 1:), *, IOSTAT=status) i_tmp ! read the configuration number
782 1206 : IF (status .NE. 0) &
783 : CALL cp_abort(__LOCATION__, &
784 : "configuration header read error (for conf nr) at line: "// &
785 0 : cp_to_string(tmc_ana%lc_traj))
786 1206 : IF (i_tmp .GT. conf_nr) THEN
787 : ! TODO we could also read the energy ...
788 1031 : conf_nr = i_tmp
789 1031 : IF (PRESENT(header_info)) header_info = c_tmp
790 1031 : stat = TMC_STATUS_OK
791 1031 : EXIT search_next_conf
792 : END IF
793 : END IF
794 : END DO search_next_conf
795 :
796 1049 : IF (stat .EQ. TMC_STATUS_OK) THEN
797 22682 : pos_loop: DO i = 1, SIZE(elem%pos), tmc_ana%dim_per_elem
798 21651 : tmc_ana%lc_traj = tmc_ana%lc_traj + 1
799 : READ (tmc_ana%id_traj, FMT="(A4,1X,1000F20.10)", IOSTAT=status) &
800 21651 : c_tmp, elem%pos(i:i + tmc_ana%dim_per_elem - 1)
801 22682 : IF (status .NE. 0) THEN
802 : CALL cp_abort(__LOCATION__, &
803 : "configuration pos read error at line: "// &
804 0 : cp_to_string(tmc_ana%lc_traj))
805 : END IF
806 : END DO pos_loop
807 65984 : elem%pos(:) = elem%pos(:)/au2a
808 : END IF
809 :
810 : ! end the timing
811 1049 : CALL timestop(handle)
812 1049 : END SUBROUTINE read_pos_from_file
813 :
814 : ! **************************************************************************************************
815 : !> \brief search for the dipole entry
816 : !> \param elem actual tree element to be read
817 : !> \param tmc_ana ...
818 : !> \param stat ...
819 : !> \param conf_nr Markov chain element number
820 : !> (input the old number and read only if conf nr from file is greater
821 : !> \author Mandes 03.2013
822 : ! **************************************************************************************************
823 0 : SUBROUTINE read_dipole_from_file(elem, tmc_ana, stat, conf_nr)
824 : TYPE(tree_type), POINTER :: elem
825 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
826 : INTEGER :: stat, conf_nr
827 :
828 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_dipole_from_file'
829 :
830 : CHARACTER(LEN=250) :: c_tmp
831 : INTEGER :: handle, status
832 :
833 0 : stat = TMC_STATUS_FAILED
834 :
835 0 : CPASSERT(ASSOCIATED(elem))
836 0 : CPASSERT(ASSOCIATED(elem%dipole))
837 0 : CPASSERT(ASSOCIATED(tmc_ana))
838 0 : CPASSERT(tmc_ana%id_dip .GT. 0)
839 :
840 : ! start the timing
841 0 : CALL timeset(routineN, handle)
842 0 : tmc_ana%lc_dip = tmc_ana%lc_dip + 1
843 0 : READ (tmc_ana%id_dip, FMT="(A)", IOSTAT=status) c_tmp
844 0 : IF (status .EQ. 0) THEN
845 : ! skip the initial line (header)
846 0 : IF (INDEX(c_tmp, "#") .GT. 0) THEN
847 0 : tmc_ana%lc_dip = tmc_ana%lc_dip + 1
848 0 : READ (tmc_ana%id_dip, FMT="(A)", IOSTAT=status) c_tmp
849 : END IF
850 : END IF
851 0 : IF (status .EQ. 0) THEN
852 : READ (c_tmp, FMT="(I8,10F20.10)", IOSTAT=status) &
853 0 : conf_nr, elem%dipole(:)
854 : END IF
855 0 : IF (status .EQ. 0) THEN ! success
856 0 : stat = TMC_STATUS_OK
857 0 : ELSE IF (status .LT. 0) THEN ! end of file reached
858 0 : stat = TMC_STATUS_WAIT_FOR_NEW_TASK
859 : ELSE
860 : IF (status .NE. 0) THEN
861 0 : CPWARN("configuration dipole read error at line: "//cp_to_string(tmc_ana%lc_dip))
862 : END IF
863 0 : stat = TMC_STATUS_FAILED
864 : END IF
865 :
866 : ! end the timing
867 0 : CALL timestop(handle)
868 0 : END SUBROUTINE read_dipole_from_file
869 :
870 : ! **************************************************************************************************
871 : !> \brief search for the cell entry
872 : !> \param elem actual tree element to be read
873 : !> \param tmc_ana ...
874 : !> \param stat ...
875 : !> \param conf_nr Markov chain element number
876 : !> (input the old number and read only if conf nr from file is greater
877 : !> \author Mandes 03.2013
878 : ! **************************************************************************************************
879 2412 : SUBROUTINE read_cell_from_file(elem, tmc_ana, stat, conf_nr)
880 : TYPE(tree_type), POINTER :: elem
881 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
882 : INTEGER :: stat, conf_nr
883 :
884 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_cell_from_file'
885 :
886 : CHARACTER(LEN=250) :: c_tmp
887 : INTEGER :: handle, status
888 : REAL(KIND=dp) :: r_tmp
889 : REAL(KIND=dp), DIMENSION(3, 3) :: hmat
890 :
891 1206 : stat = TMC_STATUS_FAILED
892 :
893 1206 : CPASSERT(ASSOCIATED(elem))
894 1206 : CPASSERT(ASSOCIATED(tmc_ana))
895 1206 : CPASSERT(ASSOCIATED(tmc_ana%cell))
896 1206 : CPASSERT(tmc_ana%id_cell .GT. 0)
897 :
898 : ! start the timing
899 1206 : CALL timeset(routineN, handle)
900 :
901 1206 : tmc_ana%lc_cell = tmc_ana%lc_cell + 1
902 1206 : READ (tmc_ana%id_cell, FMT="(A)", IOSTAT=status) c_tmp
903 1206 : IF (status .EQ. 0) THEN
904 : ! skip the initial line (header)
905 1206 : IF (INDEX(c_tmp, "#") .GT. 0) THEN
906 18 : tmc_ana%lc_cell = tmc_ana%lc_cell + 1
907 18 : READ (tmc_ana%id_cell, FMT="(A)", IOSTAT=status) c_tmp
908 : END IF
909 : END IF
910 1206 : IF (status .EQ. 0) THEN
911 1206 : READ (c_tmp, FMT="(I8,9(1X,F19.10),1X,F24.10)", IOSTAT=status) conf_nr, &
912 2412 : hmat(:, :), r_tmp
913 : END IF
914 1206 : IF (status .LT. 0) THEN ! end of file reached
915 0 : stat = TMC_STATUS_WAIT_FOR_NEW_TASK
916 1206 : ELSE IF (status .GT. 0) THEN
917 : IF (status .NE. 0) &
918 0 : CPABORT("configuration cell read error at line: "//cp_to_string(tmc_ana%lc_cell))
919 0 : stat = TMC_STATUS_FAILED
920 : ELSE
921 1206 : IF (elem%nr .LT. 0) elem%nr = conf_nr
922 15678 : hmat(:, :) = hmat(:, :)/au2a
923 : ! get the box scaling
924 : CALL get_cell_scaling(cell=tmc_ana%cell, scaled_hmat=hmat, &
925 1206 : box_scale=elem%box_scale)
926 1206 : stat = TMC_STATUS_OK
927 : END IF
928 : ! end the timing
929 1206 : CALL timestop(handle)
930 1206 : END SUBROUTINE read_cell_from_file
931 :
932 : !----------------------------------------------------------------------------
933 : ! get the configurations from file and calc
934 : !----------------------------------------------------------------------------
935 :
936 : ! **************************************************************************************************
937 : !> \brief opens the files for reading configurations data to analyze
938 : !> \param tmc_ana ...
939 : !> \param stat ...
940 : !> \param dir_ind ...
941 : !> \param
942 : !> \author Mandes 02.2013
943 : ! **************************************************************************************************
944 36 : SUBROUTINE analyse_files_open(tmc_ana, stat, dir_ind)
945 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
946 : INTEGER :: stat
947 : INTEGER, OPTIONAL :: dir_ind
948 :
949 : CHARACTER(LEN=*), PARAMETER :: routineN = 'analyse_files_open'
950 :
951 : CHARACTER(LEN=default_path_length) :: dir_name, file_name_act, file_name_temp
952 : INTEGER :: handle
953 : LOGICAL :: file_exists
954 :
955 18 : CPASSERT(ASSOCIATED(tmc_ana))
956 :
957 18 : stat = TMC_STATUS_WAIT_FOR_NEW_TASK
958 :
959 : ! start the timing
960 18 : CALL timeset(routineN, handle)
961 :
962 18 : IF (PRESENT(dir_ind)) THEN
963 18 : CPASSERT(ASSOCIATED(tmc_ana%dirs))
964 18 : CPASSERT(dir_ind .GT. 0)
965 18 : CPASSERT(dir_ind .LE. SIZE(tmc_ana%dirs))
966 :
967 18 : IF (INDEX(tmc_ana%dirs(dir_ind), "/", BACK=.TRUE.) .EQ. &
968 : LEN_TRIM(tmc_ana%dirs(dir_ind))) THEN
969 18 : dir_name = TRIM(tmc_ana%dirs(dir_ind))
970 : ELSE
971 0 : dir_name = TRIM(tmc_ana%dirs(dir_ind))//"/"
972 : END IF
973 : ELSE
974 0 : dir_name = "./"
975 : END IF
976 :
977 : ! open the files
978 : file_name_temp = expand_file_name_temp( &
979 : file_name=tmc_default_trajectory_file_name, &
980 18 : rvalue=tmc_ana%temperature)
981 : ! position file
982 18 : IF (tmc_ana%costum_pos_file_name .NE. "") THEN
983 0 : file_name_act = TRIM(dir_name)//tmc_ana%costum_pos_file_name
984 : ELSE
985 : file_name_act = TRIM(dir_name)// &
986 18 : expand_file_name_ending(file_name_temp, "xyz")
987 : END IF
988 18 : INQUIRE (FILE=file_name_act, EXIST=file_exists)
989 18 : IF (file_exists) THEN
990 : CALL open_file(file_name=file_name_act, file_status="OLD", &
991 18 : file_action="READ", unit_number=tmc_ana%id_traj)
992 18 : WRITE (tmc_ana%io_unit, FMT='(T2,A,"| ",A,T41,A40)') "TMC_ANA", &
993 36 : "read xyz file", TRIM(file_name_act)
994 : END IF
995 :
996 : ! cell file
997 18 : IF (tmc_ana%costum_cell_file_name .NE. "") THEN
998 0 : file_name_act = TRIM(dir_name)//tmc_ana%costum_cell_file_name
999 : ELSE
1000 : file_name_act = TRIM(dir_name)// &
1001 18 : expand_file_name_ending(file_name_temp, "cell")
1002 : END IF
1003 18 : INQUIRE (FILE=file_name_act, EXIST=file_exists)
1004 18 : IF (file_exists) THEN
1005 : CALL open_file(file_name=file_name_act, file_status="OLD", &
1006 18 : file_action="READ", unit_number=tmc_ana%id_cell)
1007 18 : WRITE (tmc_ana%io_unit, FMT='(T2,A,"| ",A,T41,A40)') "TMC_ANA", &
1008 36 : "read cell file", TRIM(file_name_act)
1009 : END IF
1010 :
1011 : ! dipole file
1012 18 : IF (tmc_ana%costum_dip_file_name .NE. "") THEN
1013 18 : file_name_act = TRIM(dir_name)//tmc_ana%costum_dip_file_name
1014 : ELSE
1015 : file_name_act = TRIM(dir_name)// &
1016 0 : expand_file_name_ending(file_name_temp, "dip")
1017 : END IF
1018 18 : INQUIRE (FILE=file_name_act, EXIST=file_exists)
1019 18 : IF (file_exists) THEN
1020 : CALL open_file(file_name=file_name_act, file_status="OLD", &
1021 0 : file_action="READ", unit_number=tmc_ana%id_dip)
1022 0 : WRITE (tmc_ana%io_unit, FMT='(T2,A,"| ",A,T41,A40)') "TMC_ANA", &
1023 0 : "read dip file", TRIM(file_name_act)
1024 : END IF
1025 :
1026 18 : IF (tmc_ana%id_traj .GT. 0 .OR. tmc_ana%id_cell .GT. 0 .OR. &
1027 : tmc_ana%id_dip .GT. 0) THEN
1028 18 : stat = TMC_STATUS_OK
1029 : ELSE
1030 : CALL cp_warn(__LOCATION__, &
1031 : "There is no file to open for temperature "//cp_to_string(tmc_ana%temperature)// &
1032 0 : "K in directory "//TRIM(dir_name))
1033 : END IF
1034 : ! end the timing
1035 18 : CALL timestop(handle)
1036 18 : END SUBROUTINE analyse_files_open
1037 :
1038 : ! **************************************************************************************************
1039 : !> \brief close the files for reading configurations data to analyze
1040 : !> \param tmc_ana ...
1041 : !> \param
1042 : !> \author Mandes 02.2013
1043 : ! **************************************************************************************************
1044 36 : SUBROUTINE analyse_files_close(tmc_ana)
1045 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
1046 :
1047 : CHARACTER(LEN=*), PARAMETER :: routineN = 'analyse_files_close'
1048 :
1049 : INTEGER :: handle
1050 :
1051 18 : CPASSERT(ASSOCIATED(tmc_ana))
1052 :
1053 : ! start the timing
1054 18 : CALL timeset(routineN, handle)
1055 :
1056 : ! position file
1057 18 : IF (tmc_ana%id_traj .GT. 0) CALL close_file(unit_number=tmc_ana%id_traj)
1058 :
1059 : ! cell file
1060 18 : IF (tmc_ana%id_cell .GT. 0) CALL close_file(unit_number=tmc_ana%id_cell)
1061 :
1062 : ! dipole file
1063 18 : IF (tmc_ana%id_dip .GT. 0) CALL close_file(unit_number=tmc_ana%id_dip)
1064 :
1065 : ! end the timing
1066 18 : CALL timestop(handle)
1067 18 : END SUBROUTINE analyse_files_close
1068 :
1069 : END MODULE tmc_file_io
|
