1C
2C     rt_tddft_input.F
3C
4C     Parses input deck for rt-tddft parameters.
5C
6C
7      subroutine rt_tddft_input (rtdb)
8      implicit none
9
10#include "rt_tddft.fh"
11#include "errquit.fh"
12#include "inp.fh"
13#include "rtdb.fh"
14#include "mafdecls.fh"
15#include "stdio.fh"
16
17
18C     == Inputs ==
19      integer, intent(in) :: rtdb
20
21
22C     == Parameters ==
23      character(*), parameter :: pname = "rt_tddft_input: "
24      integer, parameter      :: maxnum = 999999999 ! dummy # when using "*" as # prints
25
26C     == Varibles ==
27      logical done
28      character*255 test, curr_popt, loadopts, loadtarget
29      logical got_opt
30
31      double precision tmin, tmax, dt
32
33      logical prof, noprop, static, nodisk, matrix_checks
34C      logical dplot_do, subgs
35C      character*20 dplot_opts
36      integer nchecks, nprints, nrestarts
37C      integer nsnapshots
38      character*20 num_str
39      integer checklvl
40      character*16 field_name
41      character*255 tag_in, method
42      integer prop_method, exp_method
43      logical use_dmat
44
45
46C     (parameters that must be supplied--no defaults)
47      logical got_tmin, got_tmax, got_dt
48      integer nfields, nexcites
49
50
51
52C     == External ==
53      integer, external :: atoi
54
55
56C     (values which have no defaults and must be supplied)
57      got_tmin = .false.
58      got_tmax = .false.
59      got_dt = .false.
60
61      nfields = 0
62      nexcites = 0
63
64      if (.not.rtdb_put(rtdb,'rt_tddft:nfields',
65     $     mt_int,1,nfields))
66     $     call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
67
68      if (.not.rtdb_put(rtdb,'rt_tddft:nexcites',
69     $     mt_int,1,nexcites))
70     $     call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
71
72
73
74
75C
76C     Dump defaults into rtdb (will be replaced later if values found in
77C     input deck).
78C
79      call rt_tddft_input_put_defaults (rtdb)
80
81
82
83C
84C     Parse the input; we will put in rtdb later after checking.
85C
86      done = .false.
87      do while (.not. done)
88
89
90         if (.not. inp_read())
91     $        call errquit(pname//'Read failed input',0, INPUT_ERR)
92         if (.not. inp_a(test))
93     $        call errquit(pname//'Read failed keyword',0, INPUT_ERR)
94
95
96
97C
98C     Tag/title for the run.
99C
100         if (inp_compare(.false.,test,"tag")) then
101            if (.not.inp_a(tag_in)) then
102               call errquit (pname//"failed to parse tag",0,0)
103            endif
104
105            if (len_trim(tag_in) .gt. 24)
106     $           call errquit (pname//"tag too long, max length is 24",
107     $           0,0)
108
109            if (.not. rtdb_cput (rtdb, "rt_tddft:tag",
110     $           1, trim(tag_in)))
111     $           call errquit (pname//
112     $           "failed to put 'tag' target into rtdb",0,0)
113
114
115
116C
117C     Look for starting state to load.
118C
119         elseif (inp_compare(.false.,test,"load")) then
120            if (.not.inp_a(loadopts)) then
121               call errquit (pname//"failed to parse load option",0,0)
122            endif
123
124            if (loadopts .eq. "vectors") then
125               if (.not.inp_a(loadtarget)) then
126                  call errquit (pname//
127     $                 "failed to parse 'load vectors' target",0,0)
128               endif
129               if (.not. rtdb_cput (rtdb, "rt_tddft:init_movecs",
130     $              1, loadtarget)) call errquit (pname//
131     $              "failed to put 'init_movecs' target into rtdb",0,0)
132
133            elseif (loadopts .eq. "density") then
134               use_dmat=.true.
135               if (.not.rtdb_put(rtdb,'rt_tddft:use_dmat',mt_log,
136     $              1,use_dmat))
137     $            call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
138c              Read total (if CS) or alpha (if OS) dmat
139               if (.not.inp_a(loadtarget)) then
140                  call errquit (pname//
141     $                 "failed to parse 'load density' target",0,0)
142               endif
143               if (.not. rtdb_cput (rtdb, "rt_tddft:init_dmat1",
144     $              1, loadtarget)) call errquit (pname//
145     $              "failed to put 'init_density' target into rtdb",0,0)
146c              Read beta dmat if given (i.e. doing OS)
147               if (inp_a(loadtarget)) then
148                 if (.not. rtdb_cput (rtdb, "rt_tddft:init_dmat2",
149     $               1, loadtarget)) call errquit (pname//
150     $              "failed to put 'init_density' target into rtdb",0,0)
151               endif
152
153            elseif (loadopts .eq. "scf") then
154C     (no need to do anything--will use SCF vectors as starting point if nothing in rtdb)
155
156            elseif (loadopts .eq. "restart") then ! will look for *.rt_restart file
157               if (.not.rtdb_put (rtdb, "rt_tddft:restart",
158     $              mt_log, 1, .true.))
159     $              call errquit(pname//
160     $              'Write failed to rtdb',0,RTDB_ERR)
161
162            else
163               call errquit (pname//"invalid 'load' option: "
164     $              //loadopts,0,0)
165            endif
166
167
168C
169C     applied field sub-block
170C
171         elseif (inp_compare(.false.,test,'field')) then
172
173            if (.not. inp_a (field_name))
174     $           call errquit (pname//"failed to read field name",0,0)
175
176            nfields = nfields + 1
177
178            if (.not.rtdb_put(rtdb,'rt_tddft:nfields',mt_int,1,nfields))
179     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
180
181            call rt_tddft_input_field (rtdb, field_name, nfields)
182
183
184
185
186C
187C     visualization sub-block
188C
189         elseif (inp_compare(.false.,test,'visualization')) then
190            call rt_tddft_input_visualization (rtdb)
191
192
193C
194C     MO CAP sub-block
195C
196         elseif (inp_compare(.false.,test,'mocap')) then
197            call rt_tddft_input_mocap (rtdb)
198
199
200C
201C     Spatial CAP sub-block
202C
203         elseif (inp_compare(.false.,test,'spatialcap')) then
204            call rt_tddft_input_spatialcap (rtdb)
205
206
207C
208C     tmin
209C
210C     XXX HARDCODED TO -5*dt
211C
212c$$$         elseif (inp_compare(.false.,test,'tmin')) then
213c$$$            if (.not.inp_f(tmin)) call errquit (pname//
214c$$$     $           "tmin takes a float", 0, 0)
215c$$$            if (.not.rtdb_put(rtdb,'rt_tddft:tmin',mt_dbl,1,tmin))
216c$$$     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
217c$$$            got_tmin = .true.
218
219
220
221C
222C     tmax
223C
224         elseif (inp_compare(.false.,test,'tmax')) then
225            if (.not.inp_f(tmax)) call errquit (pname//
226     $           "tmax takes a float", 0, 0)
227            if (.not.rtdb_put(rtdb,'rt_tddft:tmax',mt_dbl,1,tmax))
228     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
229            got_tmax = .true.
230
231
232C
233C     dt
234C
235         elseif (inp_compare(.false.,test,'dt')) then
236            if (.not.inp_f(dt)) call errquit (pname//
237     $           "dt takes a float >= 0", 0, 0)
238            if (.not.rtdb_put(rtdb,'rt_tddft:dt',mt_dbl,1,dt))
239     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
240            got_dt = .true.
241
242
243C
244C     checklvl
245C
246         elseif (inp_compare(.false.,test,'checklvl')) then
247            if (.not.inp_i(checklvl)) call errquit (pname//
248     $           "checklvl takes a value of 1, 2, or 3", 0, 0)
249
250            if (.not.rtdb_put(rtdb,'rt_tddft:checklvl',
251     $           mt_int,1,checklvl))
252     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
253
254
255
256
257
258         elseif (inp_compare(.false.,test,'tolerances')) then
259            call rt_tddft_input_tolerances (rtdb)
260
261C
262C     nchecks
263C
264         elseif (inp_compare(.false.,test,'nchecks')) then
265            if (.not.inp_a(num_str))
266     $           call errquit (pname//
267     $           "nchecks takes an int >= 0 (or *)", 0, 0)
268
269            if ( trim(num_str) .eq. "*") then
270               nchecks = maxnum
271            else
272               nchecks = atoi (num_str)
273               if (nchecks < 0)
274     $              call errquit (pname//
275     $              "nchecks takes an int >= 0 (or *)", 0, 0)
276            endif
277
278            if (.not.rtdb_put(rtdb,'rt_tddft:nchecks',mt_int,1,nchecks))
279     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
280
281
282
283
284C
285C     nprints (note * means every step, uses large dummy number)
286C
287         elseif (inp_compare(.false.,test,'nprints')) then
288
289            if (.not.inp_a(num_str))
290     $           call errquit (pname//
291     $           "nprints takes an int >= 0 (or *)", 0, 0)
292
293            if ( trim(num_str) .eq. "*") then
294               nprints = maxnum
295            else
296               nprints = atoi (num_str)
297               if (nprints < 0)
298     $              call errquit (pname//
299     $              "nprints takes an int >= 0 (or *)", 0, 0)
300            endif
301
302            if (.not.rtdb_put(rtdb,'rt_tddft:nprints',mt_int,1,nprints))
303     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
304
305
306
307C
308C     nsnapshots
309C
310c$$$         elseif (inp_compare(.false.,test,'nsnapshots')) then
311c$$$            if (.not.inp_a(num_str))
312c$$$     $           call errquit (pname//
313c$$$     $           "nsnapshots takes an int >= 0 (or *)", 0, 0)
314c$$$
315c$$$            if ( trim(num_str) .eq. "*") then
316c$$$               nsnapshots = maxnum
317c$$$            else
318c$$$               nsnapshots = atoi (num_str)
319c$$$               if (nsnapshots < 0)
320c$$$     $              call errquit (pname//
321c$$$     $              "nsnapshots takes an int >= 0 (or *)", 0, 0)
322c$$$            endif
323c$$$
324c$$$            if (.not.rtdb_put(rtdb,'rt_tddft:nsnapshots',
325c$$$     $           mt_int,1,nsnapshots))
326c$$$     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
327
328
329C
330C     nrestarts
331C
332         elseif (inp_compare(.false.,test,'nrestarts')) then
333            if (.not.inp_a(num_str))
334     $           call errquit (pname//
335     $           "nrestarts takes an int >= 0 (or *)", 0, 0)
336
337            if ( trim(num_str) .eq. "*") then
338               nrestarts = maxnum
339            else
340               nrestarts = atoi (num_str)
341               if (nrestarts < 0)
342     $              call errquit (pname//
343     $              "nrestarts takes an int >= 0 (or *)", 0, 0)
344            endif
345
346            if (.not.rtdb_put(rtdb,'rt_tddft:nrestarts',
347     $           mt_int,1,nrestarts))
348     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
349
350
351
352
353C
354C     Propagator
355C
356         elseif (inp_compare(.false.,test,'propagator')) then
357
358            if (.not.inp_a(method)) then
359               call errquit (pname//
360     $              "failed to parse propagator method",0,0)
361            endif
362
363            if (method.eq."euler") then
364               prop_method = 1  ! euler
365
366            elseif (method.eq."rk4") then
367               prop_method = 2  ! 4th order runge-kutta
368
369            elseif (method.eq."magnus") then
370               prop_method = 3  ! 2nd order magnus
371
372            else
373               call errquit (pname//"invalid propagator: "//method)
374            endif
375
376            if (.not.rtdb_put(rtdb,'rt_tddft:prop_method',
377     $           mt_int,1,prop_method))
378     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
379
380
381
382C
383C     Exponentiation method
384C
385         elseif (inp_compare(.false.,test,'exp')) then
386
387            if (.not.inp_a(method)) then
388               call errquit (pname//
389     $              "failed to parse exponentiation method",0,0)
390            endif
391
392            if (method.eq."diag") then
393               exp_method = 2  ! diagonalization
394
395            elseif (method.eq."pseries") then
396               exp_method = 1  ! power series
397
398            elseif (method.eq."magnus") then
399               exp_method = 3   ! baker-campbell-hausdorff (disabled) !!since assumes e^X A e^-X
400               call errquit (pname//"BCH exp disabled"//method)
401
402            else
403               call errquit (pname//
404     $              "invalid exponentiation method: "//method)
405            endif
406
407            if (.not.rtdb_put(rtdb,'rt_tddft:exp_method',
408     $           mt_int,1,exp_method))
409     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
410
411
412
413
414C
415C     profiling
416C
417         elseif (inp_compare(.false.,test,'prof')) then
418            prof=.true.
419
420            if (.not.rtdb_put(rtdb,'rt_tddft:prof',mt_log,1,prof))
421     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
422
423
424C
425C     profiling
426C
427         elseif (inp_compare(.false.,test,'matrix_checks')) then
428            matrix_checks=.true.
429
430            if (.not.rtdb_put(rtdb,'rt_tddft:matrix_checks',
431     $           mt_log,1,matrix_checks))
432     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
433
434
435C
436C     Override DFT settings so that there is *absolutely* no I/O.  If
437C     "usedisk" is present, we will use DFT rtdb options, otherwise we
438C     disable everything (default).
439C
440         elseif (inp_compare(.false.,test,'usedisk')) then
441            call errquit (pname//"usedisk disabled",0,0)
442
443c$$$            nodisk=.false.
444c$$$            if (.not.rtdb_put(rtdb,'rt_tddft:nodisk',mt_log,1,nodisk))
445c$$$     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
446
447
448C
449C     post process using dplot? (i.e. make density plots)
450C
451c$$$         elseif (inp_compare(.false.,test,'dplot')) then
452c$$$            dplot_do=.true.
453c$$$
454c$$$            if (inp_a(dplot_opts)) then
455c$$$               if (trim(dplot_opts).eq."subgs") then
456c$$$                  subgs = .true.
457c$$$                  if (.not.rtdb_put(rtdb,'rt_tddft:subgs',
458c$$$     $                 mt_log,1,subgs))
459c$$$     $                 call errquit(pname//'Write failed to rtdb',
460c$$$     $                 0,RTDB_ERR)
461c$$$               else
462c$$$                  call errquit (pname//"invalid dplot option: "
463c$$$     $                 //trim(dplot_opts), 0, 0)
464c$$$               endif
465c$$$            endif
466c$$$
467c$$$            if (.not.rtdb_put(rtdb,'rt_tddft:dplot', mt_log,1,dplot_do))
468c$$$     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
469c$$$
470
471
472C
473C     noprop directive (dont propagate)
474C
475         elseif (inp_compare(.false.,test,'noprop')) then
476            noprop=.true.
477            if (.not.rtdb_put(rtdb,'rt_tddft:noprop',mt_log,1,noprop))
478     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
479
480
481C
482C     static directive (never rebuild Fock matrix)
483C
484         elseif (inp_compare(.false.,test,'static')) then
485            static = .true.
486            if (.not.rtdb_put(rtdb,'rt_tddft:static',mt_log,1,static))
487     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
488
489
490
491C
492C     Rules for applying fields to geometries
493C
494         elseif (inp_compare(.false.,test,"excite")) then
495            nexcites = nexcites + 1
496
497            if (.not.rtdb_put(rtdb,'rt_tddft:nexcites',
498     $           mt_int,1,nexcites))
499     $           call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
500
501            call rt_tddft_input_excite (rtdb, nexcites, nfields)
502
503
504C
505C     List of time-dependent system properties to print.
506C
507         elseif (inp_compare(.false.,test,'print')) then
508            call rt_tddft_input_print (rtdb)
509
510C
511C     end
512C
513         else if (inp_compare(.false.,test,'end')) then
514            done = .true.
515         else
516            call errquit(pname//'Unknown directive: '//trim(test),
517     $           0, INPUT_ERR)
518         endif
519      enddo  !end main parsing loop
520
521
522
523C      if (.not.got_tmin) call errquit(pname//"must supply tmin",0,0)
524      if (.not.got_tmax) call errquit(pname//"must supply tmax",0,0)
525      if (.not.got_dt) call errquit(pname//"must supply dt",0,0)
526
527
528C     HARDCODED TMIN
529C      tmin = -5.0*dt
530      tmin = 0d0
531
532
533      if (.not.rtdb_put(rtdb,'rt_tddft:tmin',mt_dbl,1,tmin))
534     $     call errquit(pname//'Write failed to rtdb',0,RTDB_ERR)
535
536
537C
538C     Check that parameters are valid and compatible with each other.
539C
540C      if (tmin.lt.0d0) call errquit(pname//"tmin must be > 0", 0, 0)
541      if (tmax.lt.0d0) call errquit(pname//"tmax must be > 0", 0, 0)
542      if (tmax.lt.tmin) call errquit(pname//"tmax must be > tmin", 0, 0)
543      if (dt.lt.0d0) call errquit(pname//"dt must be > 0", 0, 0)
544
545      end subroutine
546
547c $Id$
548