The script is intended to be used from PyMol molecular visualization system and allows one to generate &QM_KIND, &LINK and &CELL parts of CP2K input easily. Given a molecular object name and QM part selection name (which can be defined using either approach available in PyMol) it generates output file with corresponding input groups.
To use the script:
1. Save it to some file
2. Start PyMol and issue the following command from it’s shell:
> run /path/to/
3. Select qm part of your system using either mechanism available ni pymol, the issue a command:
> gen_cp2k_qmmm filename, qm_part_selection_name, the_whole_system_object
The output will be saved to file «filename».
N.B.: The file will be overwritten!

The script:

# PyMol script to generate part of the &QMMM section of CP2K input.
#   gen_cp2k_qmmm filename, qm_part_selection_name, the_whole_system_object
#   N.B.: The output is written to file "filename". "Filename" is overwritten!
#   CP2K:
#   PyMol:

from pymol import cmd, stored

def gen_cp2k_qmmm (filename, qmsele, mmsele):
  f = open(filename, 'w')
  print >> f, "&QMMM"

# QMMM CELL setup
# n.b. 6A padding
# n.b. cell should be cubic in case of PSOLVER=WAVELET
  ext = cmd.get_extent(qmsele)
  cell = (ext[1][0]-ext[0][0] + 6, ext[1][1]-ext[0][1] + 6, ext[1][2]-ext[0][2] + 6)
  print >> f, "  &CELL\n    ABC %.2f %.2f %.2f\n  &END CELL" % cell

  cmd.iterate_state(1, qmsele, "stored.elems.append(elem)")
  kinds = set(stored.elems)
  for kind in kinds:
    print >> f, "  &QM_KIND", kind
    cmd.iterate_state(1, qmsele + ' & elem '+kind, "stored.ids.append(ID)")
    print >> f, "    MM_INDEX", 
    for id in stored.ids: 
      print >> f, id,
    print >> f, "\n  &END QM_KIND"

  stored.qm_ids = []
  cmd.iterate_state(1, qmsele, "stored.qm_ids.append(ID)")
  model = cmd.get_model(mmsele)
  for bond in
    [ia, ib] = bond.index
    a = model.atom[ia].id
    b = model.atom[ib].id
    link = 0
    if ((a in stored.qm_ids) and not (b in stored.qm_ids)):
      qmid = a
      mmid = b
      link = 1
    if ((b in stored.qm_ids) and not (a in stored.qm_ids)):
      mmid = a
      qmid = b
      link = 1
    if link: print >> f, "  &LINK\n    QM_INDEX %d\n    MM_INDEX %d\n    QM_KIND H\n  &END LINK" % (qmid, mmid)
  print >> f, "&END QMMM"

cmd.extend("gen_cp2k_qmmm", gen_cp2k_qmmm)