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 - 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 2608 : 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 /= "")
81 :
82 2608 : ind = INDEX(file_name, ".", BACK=.TRUE.)
83 2608 : IF (.NOT. ind == 0) THEN
84 2608 : WRITE (result_file_name, *) file_name(1:ind - 1), ".", &
85 5216 : TRIM(ADJUSTL(extra))
86 : ELSE
87 0 : WRITE (result_file_name, *) TRIM(file_name), ".", extra
88 : END IF
89 2608 : result_file_name = TRIM(ADJUSTL(result_file_name))
90 2608 : CPASSERT(result_file_name /= "")
91 2608 : 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 /= "")
108 :
109 1427 : ind = INDEX(file_name, ".", BACK=.TRUE.)
110 1427 : IF (.NOT. ind == 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 /= "")
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 2864 : 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 2864 : CPASSERT(file_name /= "")
137 :
138 2864 : rval_to_string = ""
139 :
140 2864 : WRITE (rval_to_string, "(F16.2)") rvalue
141 2864 : ind = INDEX(file_name, ".", BACK=.TRUE.)
142 2864 : IF (.NOT. ind == 0) THEN
143 2864 : WRITE (result_file_name, *) file_name(1:ind - 1), "_T", &
144 5728 : TRIM(ADJUSTL(rval_to_string)), file_name(ind:LEN_TRIM(file_name))
145 : ELSE
146 0 : IF (LEN(file_name) == 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 2864 : result_file_name = TRIM(ADJUSTL(result_file_name))
154 2864 : CPASSERT(result_file_name /= "")
155 2864 : 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 /= "")
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 == 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) == 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 /= "")
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 12 : WRITE (f_unit) tmc_env%params%Temp(:), &
230 3 : tmc_env%m_env%gt_act%nr, &
231 84 : 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 15 : tmc_env%m_env%result_count(0:), &
237 33 : tmc_env%params%move_types%mv_weight, &
238 111 : tmc_env%params%move_types%acc_count, &
239 111 : tmc_env%params%move_types%mv_count, &
240 102 : tmc_env%params%move_types%subbox_acc_count, &
241 102 : tmc_env%params%move_types%subbox_count, &
242 39 : 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 : REAL(KIND=dp), PARAMETER :: eps_temp = 5.0E-03_dp
283 :
284 : INTEGER :: file_ptr, i, temp_size
285 : LOGICAL :: flag
286 2 : REAL(KIND=dp), ALLOCATABLE, DIMENSION(:) :: tmp_temp
287 : REAL(KIND=dp), DIMENSION(nr_mv_types) :: mv_weight_tmp
288 :
289 2 : CPASSERT(ASSOCIATED(tmc_env))
290 2 : CPASSERT(ASSOCIATED(tmc_env%m_env))
291 2 : CPASSERT(ASSOCIATED(tmc_env%params))
292 2 : CPASSERT(ASSOCIATED(tmc_env%m_env%gt_act))
293 :
294 2 : IF (file_name == tmc_default_restart_in_file_name) THEN
295 2 : INQUIRE (FILE=tmc_default_restart_in_file_name, EXIST=flag)
296 2 : CPASSERT(flag)
297 : CALL open_file(file_name=tmc_default_restart_in_file_name, file_status="OLD", &
298 2 : file_action="READ", unit_number=file_ptr)
299 2 : READ (file_ptr, *) file_name
300 2 : CALL close_file(unit_number=file_ptr)
301 : END IF
302 :
303 : CALL open_file(file_name=file_name, file_status="OLD", file_form="UNFORMATTED", &
304 2 : file_action="READ", unit_number=file_ptr)
305 2 : READ (file_ptr) temp_size
306 2 : IF (temp_size /= SIZE(tmc_env%params%Temp)) THEN
307 : CALL cp_abort(__LOCATION__, &
308 : "the actual specified temperatures does not "// &
309 0 : "fit in amount with the one from restart file ")
310 : END IF
311 6 : ALLOCATE (tmp_temp(temp_size))
312 2 : READ (file_ptr) tmp_temp(:), &
313 2 : tmc_env%m_env%gt_act%nr, &
314 56 : tmc_env%m_env%gt_act%rng_seed, &
315 2 : tmc_env%m_env%gt_act%rnd_nr, &
316 2 : tmc_env%m_env%gt_act%prob_acc, &
317 2 : tmc_env%m_env%gt_act%mv_conf, & !
318 2 : tmc_env%m_env%gt_act%mv_next_conf, & !
319 10 : tmc_env%m_env%result_count(0:), &
320 2 : mv_weight_tmp, & !
321 74 : tmc_env%params%move_types%acc_count, &
322 74 : tmc_env%params%move_types%mv_count, &
323 68 : tmc_env%params%move_types%subbox_acc_count, &
324 68 : tmc_env%params%move_types%subbox_count, & !
325 26 : tmc_env%params%cell%hmat, &
326 2 : job_counts, &
327 4 : timings
328 :
329 8 : IF (ANY(ABS(tmc_env%params%Temp(:) - tmp_temp(:)) >= eps_temp)) THEN
330 : CALL cp_abort(__LOCATION__, "the temperatures differ from the previous calculation. "// &
331 0 : "There were the following temperatures used:")
332 : END IF
333 22 : IF (ANY(mv_weight_tmp(:) /= tmc_env%params%move_types%mv_weight(:))) THEN
334 0 : CPWARN("The amount of mv types differs between the original and the restart run.")
335 : END IF
336 :
337 8 : DO i = 1, SIZE(tmc_env%params%Temp)
338 6 : tmc_env%m_env%gt_act%conf(i)%elem => tmc_env%m_env%result_list(i)%elem
339 6 : READ (file_ptr) tmc_env%m_env%result_list(i)%elem%nr, &
340 168 : tmc_env%m_env%result_list(i)%elem%rng_seed, &
341 384 : tmc_env%m_env%result_list(i)%elem%pos, &
342 384 : tmc_env%m_env%result_list(i)%elem%vel, &
343 24 : tmc_env%m_env%result_list(i)%elem%box_scale, &
344 6 : tmc_env%m_env%result_list(i)%elem%potential, &
345 6 : tmc_env%m_env%result_list(i)%elem%e_pot_approx, &
346 6 : tmc_env%m_env%result_list(i)%elem%ekin, &
347 6 : tmc_env%m_env%result_list(i)%elem%ekin_before_md, &
348 14 : tmc_env%m_env%result_list(i)%elem%temp_created
349 : END DO
350 2 : CALL close_file(unit_number=file_ptr)
351 2 : END SUBROUTINE read_restart_file
352 :
353 : !----------------------------------------------------------------------------
354 : ! printing configuration in file
355 : !----------------------------------------------------------------------------
356 :
357 : ! **************************************************************************************************
358 : !> \brief select the correct configuration to print out the
359 : !> (coordinates, forces, cell ...)
360 : !> \param result_list list of configurations for each temperature
361 : !> \param result_count list with number of Markov Chain number
362 : !> for each teperature (index 0 for global tree)
363 : !> \param conf_updated index of the updated (modified element)
364 : !> \param accepted acceptance flag
365 : !> \param tmc_params TMC environment parameters
366 : !> \author Mandes 02.2013
367 : ! **************************************************************************************************
368 9086 : SUBROUTINE write_result_list_element(result_list, result_count, conf_updated, &
369 : accepted, tmc_params)
370 : TYPE(elem_array_type), DIMENSION(:), POINTER :: result_list
371 : INTEGER, DIMENSION(:), POINTER :: result_count
372 : INTEGER :: conf_updated
373 : LOGICAL, INTENT(IN) :: accepted
374 : TYPE(tmc_param_type), POINTER :: tmc_params
375 :
376 : CHARACTER(LEN=*), PARAMETER :: routineN = 'write_result_list_element'
377 :
378 : CHARACTER(LEN=default_path_length) :: file_name
379 : INTEGER :: handle, i
380 :
381 4543 : file_name = ""
382 :
383 4543 : CPASSERT(ASSOCIATED(result_list))
384 4543 : CPASSERT(ASSOCIATED(result_count))
385 4543 : CPASSERT(ASSOCIATED(tmc_params))
386 4543 : CPASSERT(ASSOCIATED(tmc_params%Temp))
387 4543 : CPASSERT(conf_updated >= 0)
388 4543 : CPASSERT(conf_updated <= SIZE(tmc_params%Temp))
389 :
390 : ! start the timing
391 4543 : CALL timeset(routineN, handle)
392 :
393 4543 : IF (conf_updated == 0) THEN
394 : ! for debugging print every configuration of every temperature
395 0 : DO i = 1, SIZE(tmc_params%Temp)
396 0 : WRITE (file_name, *) "every_step_", TRIM(tmc_default_trajectory_file_name)
397 : CALL write_element_in_file(elem=result_list(i)%elem, &
398 : tmc_params=tmc_params, conf_nr=result_count(0), &
399 0 : file_name=expand_file_name_temp(file_name=file_name, rvalue=tmc_params%Temp(i)))
400 : END DO
401 : ELSE
402 4543 : IF ((.NOT. tmc_params%print_only_diff_conf) .OR. &
403 : (tmc_params%print_only_diff_conf .AND. accepted)) THEN
404 : CALL write_element_in_file(elem=result_list(conf_updated)%elem, &
405 : tmc_params=tmc_params, conf_nr=result_count(conf_updated), &
406 : file_name=expand_file_name_temp(file_name=TRIM(tmc_default_trajectory_file_name), &
407 1018 : rvalue=tmc_params%Temp(conf_updated)))
408 : END IF
409 : END IF
410 : ! end the timing
411 4543 : CALL timestop(handle)
412 4543 : END SUBROUTINE write_result_list_element
413 :
414 : ! **************************************************************************************************
415 : !> \brief writes the trajectory element in a file from sub tree element
416 : !> \param elem actual tree element to be printed out
417 : !> \param tmc_params TMC environment parameters
418 : !> \param temp_index ...
419 : !> \param file_name file name will be extended by type of file (pos, cell,...)
420 : !> \param conf_nr Markov chain element number
421 : !> \param conf_info whole header line
422 : !> \author Mandes 11.2012
423 : ! **************************************************************************************************
424 1018 : SUBROUTINE write_element_in_file(elem, tmc_params, temp_index, file_name, conf_nr, &
425 : conf_info)
426 : TYPE(tree_type), POINTER :: elem
427 : TYPE(tmc_param_type), POINTER :: tmc_params
428 : INTEGER, OPTIONAL :: temp_index
429 : CHARACTER(LEN=*), OPTIONAL :: file_name
430 : INTEGER, OPTIONAL :: conf_nr
431 : CHARACTER(LEN=*), OPTIONAL :: conf_info
432 :
433 : CHARACTER(LEN=*), PARAMETER :: routineN = 'write_element_in_file'
434 :
435 : CHARACTER(LEN=default_path_length) :: file_name_act, tmp_name
436 : CHARACTER(LEN=default_string_length) :: header
437 : INTEGER :: file_ptr, handle, i, nr_atoms
438 : LOGICAL :: file_exists, print_it
439 : REAL(KIND=dp) :: vol
440 : REAL(KIND=dp), DIMENSION(3, 3) :: hmat_scaled
441 :
442 1018 : file_name_act = ""
443 1018 : tmp_name = ""
444 1018 : header = ""
445 1018 : print_it = .TRUE.
446 :
447 0 : CPASSERT(ASSOCIATED(elem))
448 1018 : CPASSERT(ASSOCIATED(tmc_params))
449 1018 : CPASSERT(ASSOCIATED(tmc_params%atoms))
450 1018 : CPASSERT(PRESENT(conf_nr) .OR. PRESENT(conf_info))
451 :
452 1018 : IF (print_it) THEN
453 : ! start the timing
454 1018 : CALL timeset(routineN, handle)
455 :
456 : ! set default file name
457 1018 : IF (PRESENT(file_name)) THEN
458 1018 : CPASSERT(file_name /= "")
459 1018 : file_name_act = file_name
460 : ELSE
461 0 : CPASSERT(ASSOCIATED(tmc_params%Temp))
462 0 : CPASSERT(PRESENT(temp_index))
463 : file_name_act = expand_file_name_temp(file_name=tmc_default_trajectory_file_name, &
464 0 : rvalue=tmc_params%Temp(temp_index))
465 : END IF
466 :
467 1018 : nr_atoms = SIZE(elem%pos)/tmc_params%dim_per_elem
468 :
469 : ! set header (for coordinate or force file)
470 1018 : IF (tmc_params%print_trajectory .OR. tmc_params%print_forces) THEN
471 1018 : IF (PRESENT(conf_info)) THEN
472 0 : WRITE (header, *) TRIM(ADJUSTL(conf_info))
473 : ELSE
474 : !WRITE(header,FMT="(A,I8,A,F20.10)") " i = ", conf_nr,", E = ", elem%potential
475 1018 : WRITE (header, FMT="(A,I8,A,F20.10,F20.10,A,I8,I8)") "i =", conf_nr, " ,E =", &
476 2036 : elem%potential, elem%ekin, " st elem", elem%sub_tree_nr, elem%nr
477 : END IF
478 : END IF
479 :
480 : ! write the coordinates
481 1018 : IF (tmc_params%print_trajectory) THEN
482 1018 : tmp_name = expand_file_name_ending(file_name_act, "xyz")
483 : CALL open_file(file_name=tmp_name, file_status="UNKNOWN", &
484 : file_action="WRITE", file_position="APPEND", &
485 1018 : unit_number=file_ptr)
486 1018 : WRITE (file_ptr, FMT="(I8)") nr_atoms
487 1018 : WRITE (file_ptr, *) TRIM(header)
488 44643 : DO i = 1, SIZE(elem%pos), tmc_params%dim_per_elem
489 : WRITE (file_ptr, FMT="(A4,1X,1000F20.10)") &
490 43625 : TRIM(tmc_params%atoms((i - 1)/tmc_params%dim_per_elem + 1)%name), &
491 219143 : elem%pos(i:i + tmc_params%dim_per_elem - 1)*au2a
492 : END DO
493 1018 : CALL close_file(unit_number=file_ptr)
494 : END IF
495 :
496 : ! write the forces
497 1018 : IF (tmc_params%print_forces) THEN
498 331 : tmp_name = expand_file_name_ending(file_name_act, "frc")
499 : CALL open_file(file_name=tmp_name, file_status="UNKNOWN", &
500 : file_action="WRITE", file_position="APPEND", &
501 331 : unit_number=file_ptr)
502 331 : WRITE (file_ptr, FMT="(I8)") nr_atoms
503 331 : WRITE (file_ptr, *) TRIM(header)
504 7282 : DO i = 1, SIZE(elem%pos), tmc_params%dim_per_elem
505 : WRITE (file_ptr, FMT="(A4,1X,1000F20.10)") &
506 6951 : TRIM(tmc_params%atoms((i - 1)/tmc_params%dim_per_elem + 1)%name), &
507 35086 : elem%frc(i:i + tmc_params%dim_per_elem - 1)
508 : END DO
509 331 : CALL close_file(unit_number=file_ptr)
510 : END IF
511 :
512 : ! write the cell dipoles
513 1018 : IF (tmc_params%print_dipole) THEN
514 : CALL write_dipoles_in_file(file_name=file_name_act, &
515 0 : conf_nr=conf_nr, dip=elem%dipole)
516 : END IF
517 :
518 : ! write the cell file
519 1018 : IF (tmc_params%print_cell) THEN
520 392 : tmp_name = expand_file_name_ending(file_name_act, "cell")
521 : ! header
522 392 : INQUIRE (FILE=tmp_name, EXIST=file_exists) ! file_exists will be TRUE if the file exist
523 392 : IF (.NOT. file_exists) THEN
524 : CALL open_file(file_name=tmp_name, file_status="NEW", &
525 6 : file_action="WRITE", unit_number=file_ptr)
526 : WRITE (file_ptr, FMT='(A,9(7X,A2," [Angstrom]"),6X,A)') &
527 6 : "# MC step ", "Ax", "Ay", "Az", "Bx", "By", "Bz", "Cx", "Cy", "Cz", &
528 12 : "Volume [Angstrom^3]"
529 : ELSE
530 : CALL open_file(file_name=tmp_name, file_status="OLD", &
531 : file_action="WRITE", file_position="APPEND", &
532 386 : unit_number=file_ptr)
533 : END IF
534 : CALL get_scaled_cell(cell=tmc_params%cell, &
535 : box_scale=elem%box_scale, scaled_hmat=hmat_scaled, &
536 392 : vol=vol)
537 392 : WRITE (file_ptr, FMT="(I8,9(1X,F19.10),1X,F24.10)") conf_nr, &
538 5488 : hmat_scaled(:, :)*au2a, vol*au2a**3
539 : !TODO better cell output e.g. using cell_types routine
540 392 : CALL close_file(unit_number=file_ptr)
541 : END IF
542 :
543 : ! write the different energies
544 1018 : IF (tmc_params%print_energies) THEN
545 331 : tmp_name = expand_file_name_ending(file_name_act, "ener")
546 : ! header
547 331 : INQUIRE (FILE=tmp_name, EXIST=file_exists) ! file_exists will be TRUE if the file exist
548 331 : IF (.NOT. file_exists) THEN
549 : CALL open_file(file_name=tmp_name, file_status="NEW", &
550 3 : file_action="WRITE", unit_number=file_ptr)
551 : WRITE (file_ptr, FMT='(A,4A20)') &
552 3 : "# MC step ", " exact ", " approx ", " last SCF ", " kinetic "
553 : ELSE
554 : CALL open_file(file_name=tmp_name, file_status="OLD", &
555 : file_action="WRITE", file_position="APPEND", &
556 328 : unit_number=file_ptr)
557 : END IF
558 331 : WRITE (file_ptr, FMT="(I8,14F20.10)") conf_nr, elem%potential, elem%e_pot_approx, &
559 662 : elem%scf_energies(MOD(elem%scf_energies_count, 4) + 1), elem%ekin
560 331 : CALL close_file(unit_number=file_ptr)
561 : END IF
562 :
563 : ! end the timing
564 1018 : CALL timestop(handle)
565 : END IF
566 1018 : END SUBROUTINE write_element_in_file
567 :
568 : ! **************************************************************************************************
569 : !> \brief writes the cell dipoles in dipole trajectory file
570 : !> \param file_name ...
571 : !> \param conf_nr ...
572 : !> \param dip ...
573 : !> \param file_ext ...
574 : !> \param
575 : !> \author Mandes 11.2012
576 : ! **************************************************************************************************
577 500 : SUBROUTINE write_dipoles_in_file(file_name, conf_nr, dip, file_ext)
578 : CHARACTER(LEN=default_path_length) :: file_name
579 : INTEGER :: conf_nr
580 : REAL(KIND=dp), DIMENSION(:), POINTER :: dip
581 : CHARACTER(LEN=*), INTENT(in), OPTIONAL :: file_ext
582 :
583 : CHARACTER(LEN=default_path_length) :: file_name_tmp
584 : INTEGER :: file_ptr
585 : LOGICAL :: file_exists
586 :
587 500 : CPASSERT(ASSOCIATED(dip))
588 :
589 500 : IF (PRESENT(file_ext)) THEN
590 500 : CPASSERT(file_ext /= "")
591 500 : file_name_tmp = expand_file_name_ending(file_name, TRIM(file_ext))
592 : ELSE
593 0 : file_name_tmp = expand_file_name_ending(file_name, "dip")
594 : END IF
595 500 : INQUIRE (FILE=file_name_tmp, EXIST=file_exists)
596 500 : IF (.NOT. file_exists) THEN
597 : CALL open_file(file_name=file_name_tmp, file_status="NEW", &
598 3 : file_action="WRITE", unit_number=file_ptr)
599 3 : WRITE (file_ptr, FMT='(A8,10A20)') "# conf_nr", "dip_x [C Angstrom]", &
600 6 : "dip_y [C Angstrom]", "dip_z [C Angstrom]"
601 : ELSE
602 : CALL open_file(file_name=file_name_tmp, file_status="OLD", &
603 : file_action="WRITE", file_position="APPEND", &
604 497 : unit_number=file_ptr)
605 : END IF
606 2000 : WRITE (file_ptr, FMT="(I8,10F20.10)") conf_nr, dip(:)
607 500 : CALL close_file(unit_number=file_ptr)
608 500 : END SUBROUTINE write_dipoles_in_file
609 :
610 : !----------------------------------------------------------------------------
611 : ! read configuration from file
612 : !----------------------------------------------------------------------------
613 :
614 : ! **************************************************************************************************
615 : !> \brief read the trajectory element from a file from sub tree element
616 : !> \param elem actual tree element to be printed out
617 : !> \param tmc_ana TMC analysis environment parameters
618 : !> \param conf_nr Markov chain element number
619 : !> (input the old number and read only if conf nr from file is greater
620 : !> \param stat ...
621 : !> \author Mandes 03.2013
622 : ! **************************************************************************************************
623 2098 : SUBROUTINE read_element_from_file(elem, tmc_ana, conf_nr, stat)
624 : TYPE(tree_type), POINTER :: elem
625 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
626 : INTEGER :: conf_nr, stat
627 :
628 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_element_from_file'
629 :
630 : INTEGER :: conf_nr_old, handle, i_tmp
631 : LOGICAL :: files_conf_missmatch
632 :
633 1049 : stat = TMC_STATUS_OK
634 1049 : conf_nr_old = conf_nr
635 1049 : files_conf_missmatch = .FALSE.
636 :
637 1049 : CPASSERT(ASSOCIATED(elem))
638 1049 : CPASSERT(ASSOCIATED(tmc_ana))
639 1049 : CPASSERT(ASSOCIATED(tmc_ana%atoms))
640 :
641 : ! start the timing
642 1049 : CALL timeset(routineN, handle)
643 :
644 : ! read the coordinates
645 1049 : IF (tmc_ana%id_traj > 0) THEN
646 1049 : i_tmp = conf_nr_old
647 : CALL read_pos_from_file(elem=elem, tmc_ana=tmc_ana, stat=stat, &
648 1049 : conf_nr=i_tmp)
649 1049 : IF (stat == TMC_STATUS_WAIT_FOR_NEW_TASK) THEN
650 : CALL cp_warn(__LOCATION__, &
651 : 'end of position file reached at line '// &
652 : cp_to_string(REAL(tmc_ana%lc_traj, KIND=dp))//", last element "// &
653 18 : cp_to_string(tmc_ana%last_elem%nr))
654 : ELSE
655 1031 : CPASSERT(i_tmp > conf_nr_old)
656 1031 : conf_nr = i_tmp
657 1031 : elem%nr = i_tmp
658 : END IF
659 : END IF
660 :
661 : ! read the forces
662 : ! TODO if necessary
663 :
664 : ! read the dipoles file
665 1049 : IF (tmc_ana%id_dip > 0 .AND. stat == TMC_STATUS_OK) THEN
666 0 : i_tmp = conf_nr_old
667 : search_conf_dip: DO
668 : CALL read_dipole_from_file(elem=elem, tmc_ana=tmc_ana, stat=stat, &
669 0 : conf_nr=i_tmp)
670 0 : IF (stat == TMC_STATUS_WAIT_FOR_NEW_TASK) THEN
671 : CALL cp_warn(__LOCATION__, &
672 : 'end of dipole file reached at line'// &
673 0 : cp_to_string(REAL(tmc_ana%lc_dip, KIND=dp)))
674 0 : EXIT search_conf_dip
675 : END IF
676 : ! check consitence with pos file
677 0 : IF (tmc_ana%id_traj > 0) THEN
678 0 : IF (i_tmp == conf_nr) THEN
679 : files_conf_missmatch = .FALSE.
680 : EXIT search_conf_dip
681 : ELSE
682 : ! the configuration numbering differ from the position file,
683 : ! but we keep on searching for the correct configuration
684 : files_conf_missmatch = .TRUE.
685 : END IF
686 : ! if no pos file, just take the next conf
687 0 : ELSE IF (i_tmp > conf_nr_old) THEN
688 0 : conf_nr = i_tmp
689 0 : elem%nr = i_tmp
690 0 : EXIT search_conf_dip
691 : END IF
692 : END DO search_conf_dip
693 : END IF
694 :
695 : ! read the cell file
696 1049 : IF (tmc_ana%id_cell > 0 .AND. stat == TMC_STATUS_OK) THEN
697 : search_conf_cell: DO
698 : CALL read_cell_from_file(elem=elem, tmc_ana=tmc_ana, stat=stat, &
699 1206 : conf_nr=i_tmp)
700 1206 : IF (stat == TMC_STATUS_WAIT_FOR_NEW_TASK) THEN
701 : CALL cp_warn(__LOCATION__, &
702 : 'end of cell file reached at line at line'// &
703 0 : cp_to_string(REAL(tmc_ana%lc_cell, KIND=dp)))
704 0 : EXIT search_conf_cell
705 : END IF
706 : ! check consitence with pos file
707 1206 : IF (tmc_ana%id_traj > 0) THEN
708 1206 : IF (i_tmp == conf_nr) THEN
709 : files_conf_missmatch = .FALSE.
710 : EXIT search_conf_cell
711 : ELSE
712 : ! the configuration numbering differ from the position file,
713 : ! but we keep on searching for the correct configuration
714 : files_conf_missmatch = .TRUE.
715 : END IF
716 : ! if no pos file, just take the next conf
717 0 : ELSE IF (i_tmp > conf_nr_old) THEN
718 0 : conf_nr = i_tmp
719 0 : elem%nr = i_tmp
720 0 : EXIT search_conf_cell
721 : END IF
722 : END DO search_conf_cell
723 :
724 : END IF
725 :
726 : ! write the different energies
727 : ! TODO if necessary
728 :
729 1049 : IF (files_conf_missmatch) THEN
730 : CALL cp_warn(__LOCATION__, &
731 : 'there is a missmatch in the configuration numbering. '// &
732 : "Read number of lines (pos|cell|dip)"// &
733 : cp_to_string(tmc_ana%lc_traj)//"|"// &
734 : cp_to_string(tmc_ana%lc_cell)//"|"// &
735 0 : cp_to_string(tmc_ana%lc_dip))
736 : END IF
737 :
738 : ! end the timing
739 1049 : CALL timestop(handle)
740 1049 : END SUBROUTINE read_element_from_file
741 :
742 : ! **************************************************************************************************
743 : !> \brief search for the next configurational position in file
744 : !> \param elem actual tree element to be read
745 : !> \param tmc_ana ...
746 : !> \param stat ...
747 : !> \param conf_nr Markov chain element number
748 : !> (input the old number and read only if conf nr from file is greater
749 : !> \param header_info ...
750 : !> \author Mandes 03.2013
751 : ! **************************************************************************************************
752 2098 : SUBROUTINE read_pos_from_file(elem, tmc_ana, stat, conf_nr, header_info)
753 : TYPE(tree_type), POINTER :: elem
754 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
755 : INTEGER :: stat, conf_nr
756 : CHARACTER(LEN=*), OPTIONAL :: header_info
757 :
758 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_pos_from_file'
759 :
760 : CHARACTER(LEN=default_string_length) :: c_tmp
761 : INTEGER :: handle, i, i_tmp, status
762 :
763 1049 : stat = TMC_STATUS_FAILED
764 :
765 0 : CPASSERT(ASSOCIATED(elem))
766 1049 : CPASSERT(ASSOCIATED(elem%pos))
767 1049 : CPASSERT(ASSOCIATED(tmc_ana))
768 1049 : CPASSERT(tmc_ana%id_traj > 0)
769 :
770 : ! start the timing
771 1049 : CALL timeset(routineN, handle)
772 :
773 : search_next_conf: DO
774 6105 : c_tmp(:) = " "
775 6105 : tmc_ana%lc_traj = tmc_ana%lc_traj + 1
776 6105 : READ (tmc_ana%id_traj, '(A)', IOSTAT=status) c_tmp(:)
777 6105 : IF (status > 0) THEN
778 : CALL cp_abort(__LOCATION__, &
779 : "configuration header read error at line: "// &
780 0 : cp_to_string(tmc_ana%lc_traj)//": "//c_tmp)
781 : END IF
782 6105 : IF (status < 0) THEN ! end of file reached
783 18 : stat = TMC_STATUS_WAIT_FOR_NEW_TASK
784 18 : EXIT search_next_conf
785 : END IF
786 6087 : IF (INDEX(c_tmp, "=") > 0) THEN
787 1206 : READ (c_tmp(INDEX(c_tmp, "=") + 1:), *, IOSTAT=status) i_tmp ! read the configuration number
788 1206 : IF (status /= 0) THEN
789 : CALL cp_abort(__LOCATION__, &
790 : "configuration header read error (for conf nr) at line: "// &
791 0 : cp_to_string(tmc_ana%lc_traj))
792 : END IF
793 1206 : IF (i_tmp > conf_nr) THEN
794 : ! TODO we could also read the energy ...
795 1031 : conf_nr = i_tmp
796 1031 : IF (PRESENT(header_info)) header_info = c_tmp
797 1031 : stat = TMC_STATUS_OK
798 1031 : EXIT search_next_conf
799 : END IF
800 : END IF
801 : END DO search_next_conf
802 :
803 1049 : IF (stat == TMC_STATUS_OK) THEN
804 22682 : pos_loop: DO i = 1, SIZE(elem%pos), tmc_ana%dim_per_elem
805 21651 : tmc_ana%lc_traj = tmc_ana%lc_traj + 1
806 : READ (tmc_ana%id_traj, FMT="(A4,1X,1000F20.10)", IOSTAT=status) &
807 86604 : c_tmp, elem%pos(i:i + tmc_ana%dim_per_elem - 1)
808 22682 : IF (status /= 0) THEN
809 : CALL cp_abort(__LOCATION__, &
810 : "configuration pos read error at line: "// &
811 0 : cp_to_string(tmc_ana%lc_traj))
812 : END IF
813 : END DO pos_loop
814 65984 : elem%pos(:) = elem%pos(:)/au2a
815 : END IF
816 :
817 : ! end the timing
818 1049 : CALL timestop(handle)
819 1049 : END SUBROUTINE read_pos_from_file
820 :
821 : ! **************************************************************************************************
822 : !> \brief search for the dipole entry
823 : !> \param elem actual tree element to be read
824 : !> \param tmc_ana ...
825 : !> \param stat ...
826 : !> \param conf_nr Markov chain element number
827 : !> (input the old number and read only if conf nr from file is greater
828 : !> \author Mandes 03.2013
829 : ! **************************************************************************************************
830 0 : SUBROUTINE read_dipole_from_file(elem, tmc_ana, stat, conf_nr)
831 : TYPE(tree_type), POINTER :: elem
832 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
833 : INTEGER :: stat, conf_nr
834 :
835 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_dipole_from_file'
836 :
837 : CHARACTER(LEN=250) :: c_tmp
838 : INTEGER :: handle, status
839 :
840 0 : stat = TMC_STATUS_FAILED
841 :
842 0 : CPASSERT(ASSOCIATED(elem))
843 0 : CPASSERT(ASSOCIATED(elem%dipole))
844 0 : CPASSERT(ASSOCIATED(tmc_ana))
845 0 : CPASSERT(tmc_ana%id_dip > 0)
846 :
847 : ! start the timing
848 0 : CALL timeset(routineN, handle)
849 0 : tmc_ana%lc_dip = tmc_ana%lc_dip + 1
850 0 : READ (tmc_ana%id_dip, FMT="(A)", IOSTAT=status) c_tmp
851 0 : IF (status == 0) THEN
852 : ! skip the initial line (header)
853 0 : IF (INDEX(c_tmp, "#") > 0) THEN
854 0 : tmc_ana%lc_dip = tmc_ana%lc_dip + 1
855 0 : READ (tmc_ana%id_dip, FMT="(A)", IOSTAT=status) c_tmp
856 : END IF
857 : END IF
858 0 : IF (status == 0) THEN
859 : READ (c_tmp, FMT="(I8,10F20.10)", IOSTAT=status) &
860 0 : conf_nr, elem%dipole(:)
861 : END IF
862 0 : IF (status == 0) THEN ! success
863 0 : stat = TMC_STATUS_OK
864 0 : ELSE IF (status < 0) THEN ! end of file reached
865 0 : stat = TMC_STATUS_WAIT_FOR_NEW_TASK
866 : ELSE
867 : IF (status /= 0) THEN
868 0 : CPWARN("configuration dipole read error at line: "//cp_to_string(tmc_ana%lc_dip))
869 : END IF
870 0 : stat = TMC_STATUS_FAILED
871 : END IF
872 :
873 : ! end the timing
874 0 : CALL timestop(handle)
875 0 : END SUBROUTINE read_dipole_from_file
876 :
877 : ! **************************************************************************************************
878 : !> \brief search for the cell entry
879 : !> \param elem actual tree element to be read
880 : !> \param tmc_ana ...
881 : !> \param stat ...
882 : !> \param conf_nr Markov chain element number
883 : !> (input the old number and read only if conf nr from file is greater
884 : !> \author Mandes 03.2013
885 : ! **************************************************************************************************
886 2412 : SUBROUTINE read_cell_from_file(elem, tmc_ana, stat, conf_nr)
887 : TYPE(tree_type), POINTER :: elem
888 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
889 : INTEGER :: stat, conf_nr
890 :
891 : CHARACTER(LEN=*), PARAMETER :: routineN = 'read_cell_from_file'
892 :
893 : CHARACTER(LEN=250) :: c_tmp
894 : INTEGER :: handle, status
895 : REAL(KIND=dp) :: r_tmp
896 : REAL(KIND=dp), DIMENSION(3, 3) :: hmat
897 :
898 1206 : stat = TMC_STATUS_FAILED
899 :
900 1206 : CPASSERT(ASSOCIATED(elem))
901 1206 : CPASSERT(ASSOCIATED(tmc_ana))
902 1206 : CPASSERT(ASSOCIATED(tmc_ana%cell))
903 1206 : CPASSERT(tmc_ana%id_cell > 0)
904 :
905 : ! start the timing
906 1206 : CALL timeset(routineN, handle)
907 :
908 1206 : tmc_ana%lc_cell = tmc_ana%lc_cell + 1
909 1206 : READ (tmc_ana%id_cell, FMT="(A)", IOSTAT=status) c_tmp
910 1206 : IF (status == 0) THEN
911 : ! skip the initial line (header)
912 1206 : IF (INDEX(c_tmp, "#") > 0) THEN
913 18 : tmc_ana%lc_cell = tmc_ana%lc_cell + 1
914 18 : READ (tmc_ana%id_cell, FMT="(A)", IOSTAT=status) c_tmp
915 : END IF
916 : END IF
917 1206 : IF (status == 0) THEN
918 1206 : READ (c_tmp, FMT="(I8,9(1X,F19.10),1X,F24.10)", IOSTAT=status) conf_nr, &
919 2412 : hmat(:, :), r_tmp
920 : END IF
921 1206 : IF (status < 0) THEN ! end of file reached
922 0 : stat = TMC_STATUS_WAIT_FOR_NEW_TASK
923 1206 : ELSE IF (status > 0) THEN
924 : IF (status /= 0) THEN
925 0 : CPABORT("configuration cell read error at line: "//cp_to_string(tmc_ana%lc_cell))
926 : END IF
927 0 : stat = TMC_STATUS_FAILED
928 : ELSE
929 1206 : IF (elem%nr < 0) elem%nr = conf_nr
930 15678 : hmat(:, :) = hmat(:, :)/au2a
931 : ! get the box scaling
932 : CALL get_cell_scaling(cell=tmc_ana%cell, scaled_hmat=hmat, &
933 1206 : box_scale=elem%box_scale)
934 1206 : stat = TMC_STATUS_OK
935 : END IF
936 : ! end the timing
937 1206 : CALL timestop(handle)
938 1206 : END SUBROUTINE read_cell_from_file
939 :
940 : !----------------------------------------------------------------------------
941 : ! get the configurations from file and calc
942 : !----------------------------------------------------------------------------
943 :
944 : ! **************************************************************************************************
945 : !> \brief opens the files for reading configurations data to analyze
946 : !> \param tmc_ana ...
947 : !> \param stat ...
948 : !> \param dir_ind ...
949 : !> \param
950 : !> \author Mandes 02.2013
951 : ! **************************************************************************************************
952 36 : SUBROUTINE analyse_files_open(tmc_ana, stat, dir_ind)
953 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
954 : INTEGER :: stat
955 : INTEGER, OPTIONAL :: dir_ind
956 :
957 : CHARACTER(LEN=*), PARAMETER :: routineN = 'analyse_files_open'
958 :
959 : CHARACTER(LEN=default_path_length) :: dir_name, file_name_act, file_name_temp
960 : INTEGER :: handle
961 : LOGICAL :: file_exists
962 :
963 18 : CPASSERT(ASSOCIATED(tmc_ana))
964 :
965 18 : stat = TMC_STATUS_WAIT_FOR_NEW_TASK
966 :
967 : ! start the timing
968 18 : CALL timeset(routineN, handle)
969 :
970 18 : IF (PRESENT(dir_ind)) THEN
971 18 : CPASSERT(ASSOCIATED(tmc_ana%dirs))
972 18 : CPASSERT(dir_ind > 0)
973 18 : CPASSERT(dir_ind <= SIZE(tmc_ana%dirs))
974 :
975 18 : IF (INDEX(tmc_ana%dirs(dir_ind), "/", BACK=.TRUE.) == &
976 : LEN_TRIM(tmc_ana%dirs(dir_ind))) THEN
977 18 : dir_name = TRIM(tmc_ana%dirs(dir_ind))
978 : ELSE
979 0 : dir_name = TRIM(tmc_ana%dirs(dir_ind))//"/"
980 : END IF
981 : ELSE
982 0 : dir_name = "./"
983 : END IF
984 :
985 : ! open the files
986 : file_name_temp = expand_file_name_temp( &
987 : file_name=tmc_default_trajectory_file_name, &
988 18 : rvalue=tmc_ana%temperature)
989 : ! position file
990 18 : IF (tmc_ana%costum_pos_file_name /= "") THEN
991 0 : file_name_act = TRIM(dir_name)//tmc_ana%costum_pos_file_name
992 : ELSE
993 : file_name_act = TRIM(dir_name)// &
994 18 : expand_file_name_ending(file_name_temp, "xyz")
995 : END IF
996 18 : INQUIRE (FILE=file_name_act, EXIST=file_exists)
997 18 : IF (file_exists) THEN
998 : CALL open_file(file_name=file_name_act, file_status="OLD", &
999 18 : file_action="READ", unit_number=tmc_ana%id_traj)
1000 18 : WRITE (tmc_ana%io_unit, FMT='(T2,A,"| ",A,T41,A40)') "TMC_ANA", &
1001 36 : "read xyz file", TRIM(file_name_act)
1002 : END IF
1003 :
1004 : ! cell file
1005 18 : IF (tmc_ana%costum_cell_file_name /= "") THEN
1006 0 : file_name_act = TRIM(dir_name)//tmc_ana%costum_cell_file_name
1007 : ELSE
1008 : file_name_act = TRIM(dir_name)// &
1009 18 : expand_file_name_ending(file_name_temp, "cell")
1010 : END IF
1011 18 : INQUIRE (FILE=file_name_act, EXIST=file_exists)
1012 18 : IF (file_exists) THEN
1013 : CALL open_file(file_name=file_name_act, file_status="OLD", &
1014 18 : file_action="READ", unit_number=tmc_ana%id_cell)
1015 18 : WRITE (tmc_ana%io_unit, FMT='(T2,A,"| ",A,T41,A40)') "TMC_ANA", &
1016 36 : "read cell file", TRIM(file_name_act)
1017 : END IF
1018 :
1019 : ! dipole file
1020 18 : IF (tmc_ana%costum_dip_file_name /= "") THEN
1021 18 : file_name_act = TRIM(dir_name)//tmc_ana%costum_dip_file_name
1022 : ELSE
1023 : file_name_act = TRIM(dir_name)// &
1024 0 : expand_file_name_ending(file_name_temp, "dip")
1025 : END IF
1026 18 : INQUIRE (FILE=file_name_act, EXIST=file_exists)
1027 18 : IF (file_exists) THEN
1028 : CALL open_file(file_name=file_name_act, file_status="OLD", &
1029 0 : file_action="READ", unit_number=tmc_ana%id_dip)
1030 0 : WRITE (tmc_ana%io_unit, FMT='(T2,A,"| ",A,T41,A40)') "TMC_ANA", &
1031 0 : "read dip file", TRIM(file_name_act)
1032 : END IF
1033 :
1034 18 : IF (tmc_ana%id_traj > 0 .OR. tmc_ana%id_cell > 0 .OR. &
1035 : tmc_ana%id_dip > 0) THEN
1036 18 : stat = TMC_STATUS_OK
1037 : ELSE
1038 : CALL cp_warn(__LOCATION__, &
1039 : "There is no file to open for temperature "//cp_to_string(tmc_ana%temperature)// &
1040 0 : "K in directory "//TRIM(dir_name))
1041 : END IF
1042 : ! end the timing
1043 18 : CALL timestop(handle)
1044 18 : END SUBROUTINE analyse_files_open
1045 :
1046 : ! **************************************************************************************************
1047 : !> \brief close the files for reading configurations data to analyze
1048 : !> \param tmc_ana ...
1049 : !> \param
1050 : !> \author Mandes 02.2013
1051 : ! **************************************************************************************************
1052 36 : SUBROUTINE analyse_files_close(tmc_ana)
1053 : TYPE(tmc_analysis_env), POINTER :: tmc_ana
1054 :
1055 : CHARACTER(LEN=*), PARAMETER :: routineN = 'analyse_files_close'
1056 :
1057 : INTEGER :: handle
1058 :
1059 18 : CPASSERT(ASSOCIATED(tmc_ana))
1060 :
1061 : ! start the timing
1062 18 : CALL timeset(routineN, handle)
1063 :
1064 : ! position file
1065 18 : IF (tmc_ana%id_traj > 0) CALL close_file(unit_number=tmc_ana%id_traj)
1066 :
1067 : ! cell file
1068 18 : IF (tmc_ana%id_cell > 0) CALL close_file(unit_number=tmc_ana%id_cell)
1069 :
1070 : ! dipole file
1071 18 : IF (tmc_ana%id_dip > 0) CALL close_file(unit_number=tmc_ana%id_dip)
1072 :
1073 : ! end the timing
1074 18 : CALL timestop(handle)
1075 18 : END SUBROUTINE analyse_files_close
1076 :
1077 : END MODULE tmc_file_io
|