1C
2C     rt_tddft_os_restart.F
3C
4C     Routines for dumping and reading in restart point info.
5C
6C     Open shell version.
7C
8      subroutine rt_tddft_os_restart_save (params, it, tt, g_zdens_ao)
9      implicit none
10
11#include "errquit.fh"
12#include "mafdecls.fh"
13#include "stdio.fh"
14#include "global.fh"
15#include "rtdb.fh"
16#include "util.fh"
17#include "cdft.fh"
18#include "matutils.fh"
19#include "rt_tddft.fh"
20
21C     == Inputs ==
22      type(rt_params_t), intent(in)     :: params
23      integer, intent(in)               :: it
24      double precision, intent(in)      :: tt
25      integer, intent(in)               :: g_zdens_ao(2)
26
27
28C     == Parameters ==
29      character(*), parameter :: pname = "rt_tddft_os_restart_save: "
30
31
32C     == Variables ==
33      integer n
34      integer g_densre_ao(2), g_densim_ao(2)
35      character(255) fname_re_alpha, fname_re_beta
36      character(255) fname_im_alpha, fname_im_beta
37      logical ok
38      double precision checksum
39      double precision elapsed
40
41
42      if (params%prof) call prof_start (elapsed)
43
44      if (params%ipol .ne. 2)
45     $     call errquit (pname//"only valid for ipol = 2",0,0)
46
47C
48C     First, store current time in rtdb.  Note we put both it and tt
49C     which is redundant, but we check that tt = tmin + (it-1)*dt when
50C     we restart.
51C
52      if (.not.rtdb_put(params%rtdb,'rt_tddft:rest_it', mt_int, 1, it))
53     $     call errquit(pname//'Write failed rest_it to rtdb',
54     $     0,RTDB_ERR)
55
56      if (.not.rtdb_put(params%rtdb,'rt_tddft:rest_tt', mt_dbl, 1, tt))
57     $     call errquit(pname//'Write failed rest_tt to rtdb',
58     $     0,RTDB_ERR)
59
60
61C
62C     Also store some other params (so we can check when we restart).
63C
64CXXX  [KAL]: ns_ao, ns_mo
65
66
67C
68C     Compute and store checksum for complex dens mat.  Note checksum is
69C     calculated using the alpha complex AO dens mat (beta part not
70C     used for checksum).
71C
72      checksum = zmat_checksum (g_zdens_ao(1))
73
74      if (.not.rtdb_put(params%rtdb,'rt_tddft:rest_checksum',
75     $     mt_dbl, 1, checksum))
76     $     call errquit(pname//'Write failed rest_checksum to rtdb',
77     $     0,RTDB_ERR)
78
79
80C
81C     Dump dens mat in AO basis to file.  Note we do real and imag parts
82C     separately.
83C
84      n = params%ns_ao   !alias for clarity
85
86      ok = .true.
87      ok = ok.and.ga_create(mt_dbl, n,n, "Re[Pa]", 0, 0, g_densre_ao(1))
88      ok = ok.and.ga_create(mt_dbl, n,n, "Re[Pb]", 0, 0, g_densre_ao(2))
89
90      ok = ok.and.ga_create(mt_dbl, n,n, "Im[Pa]", 0, 0, g_densim_ao(1))
91      ok = ok.and.ga_create(mt_dbl, n,n, "Im[Pb]", 0, 0, g_densim_ao(2))
92
93      if (.not.ok) call errquit(pname//"GA allocation failed", 0,GA_ERR)
94
95      call convert_z2d (g_zdens_ao(1), g_densre_ao(1), g_densim_ao(1))
96      call convert_z2d (g_zdens_ao(2), g_densre_ao(2), g_densim_ao(2))
97
98      call util_file_name ("densmat_ao_alpha_re",
99     $     .false., .false., fname_re_alpha)
100      call util_file_name ("densmat_ao_beta_re",
101     $     .false., .false., fname_re_beta)
102
103      call util_file_name ("densmat_ao_alpha_im",
104     $     .false., .false., fname_im_alpha)
105      call util_file_name ("densmat_ao_beta_im",
106     $     .false., .false., fname_im_beta)
107
108      if (.not. dmat_io_dump (g_densre_ao(1), fname_re_alpha))
109     $     call errquit (pname//"failed to dump densao_alpha_re",0,0)
110
111      if (.not. dmat_io_dump (g_densre_ao(2), fname_re_beta))
112     $     call errquit (pname//"failed to dump densao_beta_re",0,0)
113
114      if (.not. dmat_io_dump (g_densim_ao(1), fname_im_alpha))
115     $     call errquit (pname//"failed to dump densao_alpha_im",0,0)
116
117      if (.not. dmat_io_dump (g_densim_ao(2), fname_im_beta))
118     $     call errquit (pname//"failed to dump densao_beta_im",0,0)
119
120
121C
122C     Clean up.
123C
124      ok = .true.
125      ok = ok .and. ga_destroy(g_densre_ao(1))
126      ok = ok .and. ga_destroy(g_densre_ao(2))
127      ok = ok .and. ga_destroy(g_densim_ao(1))
128      ok = ok .and. ga_destroy(g_densim_ao(2))
129      if (.not.ok) call errquit(pname//"failed to destroy GAs", 0,0)
130
131      if (params%prof) call prof_end (elapsed, "Saving restart point")
132      end subroutine
133
134
135C====================================================================
136C
137C     Load info from previous restart point.
138C
139      subroutine rt_tddft_os_restart_load (params, it, tt, g_zdens_ao)
140      implicit none
141
142#include "errquit.fh"
143#include "mafdecls.fh"
144#include "stdio.fh"
145#include "global.fh"
146#include "rtdb.fh"
147#include "util.fh"
148#include "cdft.fh"
149#include "matutils.fh"
150#include "rt_tddft.fh"
151
152
153C     == Inputs ==
154      type(rt_params_t), intent(in) :: params
155
156
157C     == Outputs ==
158      integer, intent(in)           :: g_zdens_ao(2)
159
160
161C     == Parameters ==
162      character(*), parameter  :: pname = "rt_tddft_os_restart_load: "
163
164
165C     == Outputs ==
166      integer, intent(out)          :: it
167      double precision, intent(out) :: tt
168
169
170C     == Variables ==
171      double precision checksum, rest_checksum
172      integer n
173      integer g_densre_ao(2), g_densim_ao(2)
174      character(255) fname_re_alpha, fname_re_beta
175      character(255) fname_im_alpha, fname_im_beta
176      logical ok
177      double precision elapsed
178      double precision tdiff
179
180
181      call rt_tddft_os_confirm (params,'rt_tddft_os_restart.F')
182
183      if (params%prof) call prof_start (elapsed)
184
185C
186C     Read in previous time values.
187C
188      if (.not.rtdb_get (params%rtdb,"rt_tddft:rest_it",mt_int, 1, it))
189     $     call errquit (pname//"failed to read rest_it from rtdb",0,0)
190
191      if (.not.rtdb_get (params%rtdb,"rt_tddft:rest_tt",mt_dbl, 1, tt))
192     $     call errquit (pname//"failed to read rest_tt from rtdb",0,0)
193
194      tdiff = tt - (params%tmin + (it-1)*params%dt)
195      if ( abs(tdiff) > 1d-4 ) then
196         call errquit (pname//"inconsistent rest_it, rest_tt",0,0)
197      endif
198
199      if (.not.rtdb_get (params%rtdb,"rt_tddft:rest_checksum",
200     $     mt_dbl, 1, rest_checksum)) call errquit (pname//
201     $     "failed to read rest_checksum from rtdb",0,0)
202
203
204C
205C     Read in complex dens mat in AO basis.
206C
207      call util_file_name ("densmat_ao_alpha_re",
208     $     .false., .false., fname_re_alpha)
209
210      call util_file_name ("densmat_ao_beta_re",
211     $     .false., .false., fname_re_beta)
212
213      call util_file_name ("densmat_ao_alpha_im",
214     $     .false., .false., fname_im_alpha)
215
216      call util_file_name ("densmat_ao_beta_im",
217     $     .false., .false., fname_im_beta)
218
219      n = params%ns_ao   !alias for clarity
220
221      ok = .true.
222      ok = ok.and.ga_create(mt_dbl, n, n, "Re[P]", 0, 0, g_densre_ao(1))
223      ok = ok.and.ga_create(mt_dbl, n, n, "Re[P]", 0, 0, g_densre_ao(2))
224
225      ok = ok.and.ga_create(mt_dbl, n, n, "Im[P]", 0, 0, g_densim_ao(1))
226      ok = ok.and.ga_create(mt_dbl, n, n, "Im[P]", 0, 0, g_densim_ao(2))
227
228      if (.not. dmat_io_read (g_densre_ao(1), fname_re_alpha))
229     $     call errquit (pname//"failed to read densao_re_alpha",0,0)
230      if (.not. dmat_io_read (g_densre_ao(2), fname_re_beta))
231     $     call errquit (pname//"failed to read densao_re_beta",0,0)
232
233      if (.not. dmat_io_read (g_densim_ao(1), fname_im_alpha))
234     $     call errquit (pname//"failed to read densao_im_alpha",0,0)
235      if (.not. dmat_io_read (g_densim_ao(2), fname_im_beta))
236     $     call errquit (pname//"failed to read densao_im_beta",0,0)
237
238
239C
240C     Load real and im parts into complex GA and check that checksum is
241C     consistent.
242C
243      call convert_d2z (1d0, g_densre_ao(1), 1d0, g_densim_ao(1),
244     $     g_zdens_ao(1))
245      call convert_d2z (1d0, g_densre_ao(2), 1d0, g_densim_ao(2),
246     $     g_zdens_ao(2))
247
248
249C
250C     Recall, we did checksum using only alpha part of complex dens mat
251C     in AO basis.
252C
253      checksum = zmat_checksum (g_zdens_ao(1))
254
255      if ( abs(checksum - rest_checksum) > params%tol_zero)
256     $     call errquit (pname//
257     $     "bad checksum while importing density matrices",0,0)
258
259      ok = .true.
260      ok = ok .and. ga_destroy(g_densre_ao(1))
261      ok = ok .and. ga_destroy(g_densre_ao(2))
262      ok = ok .and. ga_destroy(g_densim_ao(1))
263      ok = ok .and. ga_destroy(g_densim_ao(2))
264      if (.not.ok) call errquit(pname//"failed to destroy GAs", 0,0)
265
266      if (params%prof) call prof_end (elapsed, "Loading restart point")
267      end subroutine
268c $Id$
269