1C 2C rt_tddft_input_excite.F 3C 4C Input deck parser to rules about applying fields to geometries. 5C 6 subroutine rt_tddft_input_excite (rtdb, nexcites, nfields) 7 implicit none 8 9#include "errquit.fh" 10#include "mafdecls.fh" 11#include "stdio.fh" 12#include "global.fh" 13#include "msgids.fh" 14#include "util.fh" 15#include "rtdb.fh" 16#include "cdft.fh" 17#include "inp.fh" 18#include "geomP.fh" 19#include "geom.fh" 20#include "bas.fh" 21#include "sym.fh" 22#include "rt_tddft.fh" 23 24 25C == Inputs == 26 integer, intent(in) :: rtdb 27 integer, intent(in) :: nexcites 28 integer, intent(in) :: nfields 29 30C == Parameters == 31 character(len=*), parameter :: pname = "rt_tddft_input_excite: " 32 33 34C == Variables == 35 character*255 opt 36 type(rt_field_t) field 37 type(rt_excite_t) excite 38 integer igeom, ifield 39 logical found_geom, found_field 40 integer ao_indx 41 42 43 if (nexcites .gt. rt_max_excites) 44 $ call errquit (pname//"cannot exceed max excites", 0, 0) 45 46 47C 48C Get geom info from rtdb. Note: geomP.fh supplies ngeom_rtdb, 49C names_rtdb, max_geom_rtdb, ndipole 50C 51 if (.not.rtdb_get (rtdb, "geometry:ngeom", 52 $ mt_int, 1, ngeom_rtdb)) 53 $ call errquit (pname//"failed to get ngeom from rtdb",0,0) 54 55 if (ngeom_rtdb .lt. 1) 56 $ call errquit (pname//"invalid ngeom_rtdb < 1",0,0) 57 58c$$$ if (ngeom_rtdb .gt. nw_max_geoms) 59c$$$ $ call errquit (pname//"invalid ngeom_rtdb > nw_max_geoms",0,0) 60 61 if (ngeom_rtdb .gt. max_geom_rtdb) 62 $ call errquit(pname//"invalid ngeom_rtdb > max_geom_rtdb",0,0) 63 64 if (.not. rtdb_cget (rtdb,'geometry:names', 65 $ ngeom_rtdb,names_rtdb)) 66 $ call errquit (pname//"failed to read names from rtdb",0,0) 67 68 69 70C 71C Read in geom/atomic orbital to excite 72C 73c$$$ if (.not. inp_a (opt)) call errquit (pname// 74c$$$ $ "incorrect syntax--expecting: "// 75c$$$ $ "excite [geometry,ao] <target> with <field>",0,0) 76c$$$ 77c$$$ if (inp_compare(.false., opt, "geometry")) then !user want to excite a geometry 78c$$$ excite%ao_indx = -1 !excite all AO's 79c$$$ 80c$$$ if (.not. inp_a (opt)) call errquit (pname// 81c$$$ $ "expecting geometry tag",0,0) 82c$$$ 83c$$$ found_geom = .false. 84c$$$ do igeom = 1, ngeom_rtdb 85c$$$ if (names_rtdb(igeom) .eq. opt) then 86c$$$ if (found_geom) then 87c$$$ call errquit (pname// 88c$$$ $ "multiple matches to geom name: "//opt, 0, 0) 89c$$$ else 90c$$$ found_geom = .true. 91c$$$ excite%geom_indx = igeom 92c$$$ endif 93c$$$ endif 94c$$$ enddo 95c$$$ 96c$$$ if (.not. found_geom) 97c$$$ $ call errquit (pname// 98c$$$ $ "Specified geom to excite does not exist: "//opt, 0, 0) 99c$$$ 100c$$$ elseif (inp_compare(.false., opt, "ao")) then !user want to excite a particular atomic orbital 101c$$$ call errquit (pname//"ao excite not implemented yet",0,0) 102c$$$ excite%geom_indx = -1 !excite all geoms 103c$$$ 104c$$$ if (.not. inp_i (ao_indx)) call errquit (pname// 105c$$$ $ "ao takes a positive integer",0,0) 106c$$$ if (ao_indx < 1) call errquit (pname// 107c$$$ $ "ao takes a positive integer",0,0) 108c$$$ excite%ao_indx = ao_indx 109c$$$ 110c$$$ else 111c$$$ call errquit (pname//"invalid excitation target type: "// 112c$$$ $ trim(opt),0,0) 113c$$$ endif 114 115 116 if (.not. inp_a (opt)) call errquit (pname// 117 $ "expecting geometry tag",0,0) 118 119 found_geom = .false. 120 do igeom = 1, ngeom_rtdb 121 if (names_rtdb(igeom) .eq. opt) then 122 if (found_geom) then 123 call errquit (pname// 124 $ "multiple matches to geom name: "//opt, 0, 0) 125 else 126 found_geom = .true. 127 excite%geom_indx = igeom 128 endif 129 endif 130 enddo 131 132 if (.not. found_geom) 133 $ call errquit (pname// 134 $ "Specified geom to excite does not exist: "//opt, 0, 0) 135 136 excite%ao_indx = -1 !excite all AO's 137 138C 139C Then read in field to use. 140C 141 if (.not. inp_a (opt)) call errquit (pname// !get "with" 142 $ "incorrect syntax--expecting: "// 143 $ "excite <geom> with <field>",0,0) 144 145 if (opt .ne. "with") call errquit (pname// 146 $ "incorrect syntax--expecting: "// 147 $ "excite <geom> with <field>",0,0) 148 149 150 if (.not. inp_a (opt)) call errquit (pname// !get field 151 $ "incorrect syntax--expecting: "// 152 $ "excite <geom> with <field>",0,0) 153 154 155 found_field = .false. 156 do ifield = 1, nfields 157 call rt_tddft_field_rtdb_get (rtdb, ifield, field) 158 if (field%name .eq. opt) then 159 if (found_field) then 160 call errquit (pname// 161 $ "multiple matches to field name: "//opt, 0, 0) 162 else 163 found_field = .true. 164 excite%field_indx = ifield 165 endif 166 endif 167 enddo 168 169 if (.not. found_field) 170 $ call errquit (pname// 171 $ "Specified field does not exist: "//opt, 0, 0) 172 173 174C XXX FIX GEOM TO EXCITE IF -1 175 176 177c$$$C 178c$$$C Read in geometry to excite, check that it exists, and if so XXX 179c$$$C 180c$$$ if (.not. inp_a (opt)) call errquit (pname// !get geom 181c$$$ $ "incorrect syntax--expecting: excite <geom> with <field>") 182c$$$ 183c$$$ found_geom = .false. 184c$$$ do igeom = 1, ngeom_rtdb 185c$$$ if (names_rtdb(igeom) .eq. opt) then 186c$$$ if (found_geom) then 187c$$$ call errquit (pname// 188c$$$ $ "multiple matches to geom name: "//opt, 0, 0) 189c$$$ else 190c$$$ found_geom = .true. 191c$$$ excite%geom_indx = igeom 192c$$$ endif 193c$$$ endif 194c$$$ enddo 195c$$$ 196c$$$ if (.not. found_geom) 197c$$$ $ call errquit (pname// 198c$$$ $ "Specified geom to excite does not exist: "//opt, 0, 0) 199c$$$ 200c$$$ 201c$$$ 202c$$$C 203 204 205C 206C Check for additional options 207C 208c$$$ ao_indx = -99 209c$$$ if (inp_a (opt)) then 210c$$$ if (inp_compare (.false., opt, "ao")) then 211c$$$ if (.not. inp_i (ao_indx)) call errquit (pname// 212c$$$ $ "ao takes a positive integer",0,0) 213c$$$ if (ao_indx < 1) call errquit (pname// 214c$$$ $ "ao takes a positive integer",0,0) 215c$$$ excite%ao_indx = ao_indx 216c$$$ endif 217c$$$ endif 218 219 220C 221C Save excite rule to rtdb 222C 223 call rt_tddft_excite_rtdb_put (rtdb, nexcites, excite) 224 225 end subroutine 226 227 228C==================================================================== 229C 230C Generate entry name for field rtdb stuff (hack) 231C 232 subroutine rt_tddft_excite_rtdb_entry_name (i, name) 233 implicit none 234 235#include "errquit.fh" 236#include "rtdb.fh" 237#include "mafdecls.fh" 238#include "stdio.fh" 239#include "rt_tddft.fh" 240 241 242C == Inputs == 243 integer, intent(in) :: i 244 245 246C == Outputs == 247 character(len=*), intent(out) :: name !was 17 248 249 250C == Parameters == 251 character(len=*), parameter :: pname = 252 $ "rt_tddft_excite_rtdb_entry_name" 253 254 255C == Variables == 256 character*5 istring !note length 5 limit size of int 257 258 259 if ( (i .gt. rt_max_excites).or.(i .lt. 1) ) 260 $ call errquit(pname//"i must be between 1,rt_max_excites",0,0) 261 262 if (rt_max_fields .gt. 999) call errquit(pname// 263 $ "rt_max_excites too large; fix formatting", 0, 0) 264 265 write (istring, "(i0.5)") i 266 267 name = "rt_tddft:excite_"//trim(istring)//"_" 268 269 end subroutine 270 271 272C==================================================================== 273C 274C Load excite into rtbd. This is an ugly hack, but it's easier than 275C adding a custom struct to the rtdb routines. 276C 277 subroutine rt_tddft_excite_rtdb_put (rtdb, i, excite) 278 implicit none 279 280#include "rt_tddft.fh" 281#include "errquit.fh" 282#include "rtdb.fh" 283#include "mafdecls.fh" 284#include "stdio.fh" 285 286 287C == Inputs == 288 integer, intent(in) :: rtdb 289 integer, intent(in) :: i !index for the excite 290 type(rt_excite_t), intent(in) :: excite 291 292 293C == Parameters == 294 character(len=*), parameter :: pname="rt_tddft_excite_rtdb_put: " 295 296 297C == Variables == 298 character*32 basename 299 character*32 entry_name 300 301 if ( (i .gt. rt_max_excites).or.(i .lt. 1) ) 302 $ call errquit(pname//"i must be between 1,rt_max_excites",0,0) 303 304 call rt_tddft_excite_rtdb_entry_name (i, basename) 305 306 entry_name = trim(basename) // "geom_indx" 307 if (.not.rtdb_put(rtdb,entry_name,mt_int,1,excite%geom_indx)) 308 $ call errquit(pname//'Write failed to geom_indx rtdb', 309 $ 0,RTDB_ERR) 310 311 entry_name = trim(basename) // "field_indx" 312 if (.not.rtdb_put(rtdb,entry_name,mt_int,1,excite%field_indx)) 313 $ call errquit(pname//'Write failed to field_indx rtdb', 314 $ 0,RTDB_ERR) 315 316 entry_name = trim(basename) // "ao_indx" 317 if (.not.rtdb_put(rtdb,entry_name,mt_int,1,excite%ao_indx)) 318 $ call errquit(pname//'Write failed to ao_indx rtdb', 319 $ 0,RTDB_ERR) 320 321 end subroutine 322 323 324C==================================================================== 325C 326C Get excite from rtdb and put into struct 327C 328 subroutine rt_tddft_excite_rtdb_get (rtdb, i, excite) 329 implicit none 330 331#include "rt_tddft.fh" 332#include "errquit.fh" 333#include "rtdb.fh" 334#include "mafdecls.fh" 335#include "stdio.fh" 336 337 338C == Inputs == 339 integer, intent(in) :: rtdb 340 integer, intent(in) :: i !index for the excite 341 342 343C == Outputs == 344 type(rt_excite_t), intent(out) :: excite 345 346 347 348C == Parameters == 349 character(len=*), parameter :: pname="rt_tddft_excite_rtdb_get: " 350 351 352C == Variables == 353 character*32 basename 354 character*32 entry_name 355 356 357 if ( (i .gt. rt_max_excites).or.(i .lt. 1) ) 358 $ call errquit(pname//"i must be between 1,rt_max_excites",0,0) 359 360 call rt_tddft_excite_rtdb_entry_name (i, basename) 361 362 entry_name = trim(basename) // "geom_indx" 363 if (.not.rtdb_get(rtdb,entry_name,mt_int,1,excite%geom_indx)) 364 $ call errquit(pname//'Read failed for geom_indx rtdb', 365 $ 0,RTDB_ERR) 366 367 entry_name = trim(basename) // "field_indx" 368 if (.not.rtdb_get(rtdb,entry_name,mt_int,1,excite%field_indx)) 369 $ call errquit(pname//'Read failed for field_indx rtdb', 370 $ 0,RTDB_ERR) 371 372 entry_name = trim(basename) // "ao_indx" 373 if (.not.rtdb_get(rtdb,entry_name,mt_int,1,excite%ao_indx)) 374 $ call errquit(pname//'Read failed for ao_indx rtdb', 375 $ 0,RTDB_ERR) 376 377 end subroutine 378c $Id$ 379