read_mopac_out_file: Read mopac output file
In conmolfields: Continuous molecular fields

Read mopac output file

1	read_mopac_out_file(fname, mol)

`fname`
`x`	mol

##---- Should be DIRECTLY executable !! ----
##-- ==>  Define data, use random,
##--	or do  help(data=index)  for the standard data sets.

## The function is currently defined as
read_mopac_out_file <- function(fname, mol) {
 
  # proceed to next line
  next_line <- function() {
    iline <<- iline + 1
	aline <<- arc_lines[iline]
  }
  
  get_num_heavy_atoms <- function(mol) {
	iheavy <- 0	
	for (atom in mol$atoms) {
	  if (atom$el != "H") iheavy <- iheavy + 1
	}
	return(iheavy)
  }
  
  arc_lines <- readLines(fname)
  nlines <- length(arc_lines)
  iline <- 0
  
  # Skip to superdelocalizabilities
  while (iline < nlines) {
    next_line()
    r <- regexpr("           SUPERDELOCALIZABILITIES", aline)
	if (r > 0) break
  }
  
  # Read several scalar values
  next_line()
  next_line()
  mulliken_electronegativity <- as.numeric(substr(aline, 34, 44))
  next_line()
  parr_pople_absolute_hardness <- as.numeric(substr(aline, 34, 44))
  next_line()
  schuurmann_mo_shift_alpha <- as.numeric(substr(aline, 34, 44))
  next_line()
  next_line()
  ehomo <- as.numeric(substr(aline, 34, 44))
  next_line()
  elumo <- as.numeric(substr(aline, 34, 44))
  for (i in 1:4) next_line()
  
  # Read arrays with superdelocalizabilities
  num_heavy_atoms <- get_num_heavy_atoms(mol)
  heavy_atom_index <- integer(num_heavy_atoms)
  Dn <- double(num_heavy_atoms)
  De <- double(num_heavy_atoms)
  qZ <- double(num_heavy_atoms)
  piS <- double(num_heavy_atoms)
  chomo <- double(num_heavy_atoms)
  clumo <- double(num_heavy_atoms)
  # Read Dn(r), De(r) and q(r)-Z(r)
  for (i in 1:num_heavy_atoms) {
    next_line()
    heavy_atom_index[i] <- as.integer(substr(aline, 5, 7))
	Dn[i] <- as.numeric(substr(aline, 11, 20))
	De[i] <- as.numeric(substr(aline, 24, 33))
	qZ[i] <- as.numeric(substr(aline, 37, 46))
  }
  for (i in 1:5) next_line()
  # Read piS
  for (i in 1:num_heavy_atoms) {
    next_line()
	piS[i] <- as.numeric(substr(aline, 11, 20))
  }
  for (i in 1:5) next_line()
  # Read homo and lumo
  for (i in 1:num_heavy_atoms) {
    next_line()
	chomo[i] <- as.numeric(substr(aline, 24, 33))
	clumo[i] <- as.numeric(substr(aline, 37, 46))
  }
  
  list(
    num_heavy_atoms = num_heavy_atoms,
	heavy_atom_index = heavy_atom_index,
	Dn = Dn,
	De = De,
	qZ = qZ,
	piS = piS,
	chomo = chomo,
	clumo = clumo,
	mulliken_electronegativity = mulliken_electronegativity,
	parr_pople_absolute_hardness = parr_pople_absolute_hardness,
	schuurmann_mo_shift_alpha = schuurmann_mo_shift_alpha,
	ehomo = ehomo,
	elumo = elumo
  )
}