Index: prepareQMX/qmx.in
===================================================================
--- prepareQMX/qmx.in	(revision 15)
+++ prepareQMX/qmx.in	(revision 16)
@@ -1,10 +0,0 @@
-high-level.program : TURBOMOLE
-high-level.options : opt
-low-level.program : VASP
-
-cluster.program : TURBOMOLE
-cluster.options : clusterX
-system.program : VASP
-
-
-job.program : QuasiNewton 
Index: prepareQMX/qmx.py
===================================================================
--- prepareQMX/qmx.py	(revision 15)
+++ prepareQMX/qmx.py	(revision 16)
@@ -1,19 +0,0 @@
-#!/usr/bin/env python
-from ase.embed import Embed
-from ase.io.vasp import read_vasp
-from ase.calculators.qmx import Qmx
-from ase.calculators.vasp import Vasp
-from ase.optimize import QuasiNewton
-
-high_level=Vasp(write_input=False)
-low_level=Vasp(write_input=False)
-qmx=Qmx(high_level, low_level)
-
-system=read_vasp('POSCAR')
-cluster=read_vasp('POSCAR')
-embed=Embed(system, cluster, cell_cluster="Auto")
-embed.embed()
-embed.set_calculator(qmx)
-
-job=QuasiNewton(embed, trajectory="embed.traj")
-job.run(fmax=0.01, steps=100)
Index: prepareQMX/qmxSTR.py
===================================================================
--- prepareQMX/qmxSTR.py	(revision 15)
+++ prepareQMX/qmxSTR.py	(revision 16)
@@ -28,6 +28,11 @@
         self.keywords['import']='ase.io.turbomole'
         self.keywords['class']='read_turbomole'
         self.keywords['class.options']='coord'
-
 #-------------------------------------------------------------------------------
-strDefinitions = [VASPStrDefinition(), TURBOMOLEStrDefintion()]
+class CLUSTERStrDefintion(Definition):
+    def __init__(self):
+        self.name='CLUSTER'
+        self.keywords = {}
+#-------------------------------------------------------------------------------
+strSystemDefinitions = [VASPStrDefinition(), TURBOMOLEStrDefintion()]
+strClusterDefinitions = [VASPStrDefinition(), TURBOMOLEStrDefintion(), CLUSTERStrDefintion()]
Index: prepareQMX/qmxWRITE.py
===================================================================
--- prepareQMX/qmxWRITE.py	(revision 15)
+++ prepareQMX/qmxWRITE.py	(revision 16)
@@ -25,7 +25,7 @@
                 arrImports.append(s)
 
     stream.write("\n")
-    for name in ['high-level', 'low-level', 'qmx']:
+    for name in ['high-level', 'low-level', 'hybrid']:
         for definition in definitions:
             if definition.system.lower() != name:
                 continue
@@ -37,7 +37,9 @@
                 stream.write(strSystem+'='+strClass+'('+strOptions+')\n')
 
     stream.write("\n")
-    for name in ['system', 'cluster', 'embed']:
+        
+    extraOptions = ''
+    for name in ['cluster', 'system', 'embed']:
         for definition in definitions:
             if definition.system.lower() != name:
                 continue
@@ -46,16 +48,20 @@
             strOptions=definition.getValue('class.options')
 
             if (strClass is not ''):
-                stream.write(strSystem+'='+strClass+"('"+strOptions+"')\n")
+                if name is 'embed':
+                    if strOptions is not '':
+                        strOptions += ', '
+                    strOptions += extraOptions
+                stream.write(strSystem+'='+strClass+"("+strOptions+")\n")
                 
+            if name is 'cluster' and definition.name is 'CLUSTER':
+                stream.write(strSystem+'=['+strOptions+']'+'\n')
+                extraOptions='cluster_pos=False'
+                
             if name is 'embed':
-                strMethod=definition.getValue('method')
-                strOptions=definition.getValue('method.options')
-                stream.write(strSystem+'.'+strMethod+'('+strOptions+')\n')
+                strOptions=definition.getValue('set_calculator')
+                stream.write(strSystem+'.set_calculator('+strOptions+')\n')
 
-                strOptions=definition.getValue('calculator')
-                stream.write(strSystem+'.set_calculator('+strOptions+')\n')
-        
     stream.write("\n")
     for definition in definitions:
         if (definition.system.lower() == 'job'):
Index: prepareQMX/qmxEMBED.py
===================================================================
--- prepareQMX/qmxEMBED.py	(revision 15)
+++ prepareQMX/qmxEMBED.py	(revision 16)
@@ -20,8 +20,7 @@
         self.keywords['import']='ase.embed'
         self.keywords['class']='Embed'
         self.keywords['class.options']='system, cluster, cell_cluster="Auto"'
-        self.keywords['method']='embed'
-        self.keywords['calculator']='qmx'
+        self.keywords['set_calculator']='qmx'
 
 #-------------------------------------------------------------------------------
 embedDefinitions = [EmbedDefinition()]
Index: prepareQMX/prepareQMX.py
===================================================================
--- prepareQMX/prepareQMX.py	(revision 15)
+++ prepareQMX/prepareQMX.py	(revision 16)
@@ -15,7 +15,7 @@
 #---------------------------------------------------------------------------------------------------
 from qmxEMBED import embedDefinitions
 from qmxJOB  import jobDefinitions
-from qmxSTR  import strDefinitions
+from qmxSTR  import strSystemDefinitions, strClusterDefinitions
 from qmxCALC import calcDefinitions, QmxCalcDefinition
 
 from qmxWRITE import writeData
@@ -121,9 +121,13 @@
         definitions.append(analyzeSection(calcDefinitions, section, lines_section))
         
     #structures (atoms)
-    if section == "cluster" or section == "system":
-        definitions.append(analyzeSection(strDefinitions, section, lines_section))
+    if section == "cluster":
+        definitions.append(analyzeSection(strClusterDefinitions, section, lines_section))
         
+    #structures (atoms)
+    if section == "system":
+        definitions.append(analyzeSection(strSystemDefinitions, section, lines_section))
+        
     #job
     if section == "job":
         definitions.append(analyzeSection(jobDefinitions, section, lines_section))
Index: prepareQMX/prepareQMX.GUI.py
===================================================================
--- prepareQMX/prepareQMX.GUI.py	(revision 15)
+++ prepareQMX/prepareQMX.GUI.py	(revision 16)
@@ -10,12 +10,13 @@
 #---------------------------------------------------------------------------------------------------
 from Tkinter import *
 from copy import deepcopy
+import tkFileDialog
 import StringIO, os
 
 #---------------------------------------------------------------------------------------------------
 from qmxEMBED import embedDefinitions
 from qmxJOB import jobDefinitions
-from qmxSTR import strDefinitions
+from qmxSTR import strClusterDefinitions, strSystemDefinitions
 from qmxCALC import calcDefinitions, QmxCalcDefinition
 
 from qmxWRITE import writeData
@@ -24,27 +25,24 @@
 #---------------------------------------------------------------------------------------------------
 #---------------------------------------------------------------------------------------------------
 class MyLabelFrame(LabelFrame):
-    def __init__(self, mainWnd, master, title, def_list, extras):
+    def __init__(self, mainWnd, master, title, name, def_list, fields):
         LabelFrame.__init__(self, master, text=title)
         self.mainWnd = mainWnd
         self.definitions=[]
-        self.title = title
 
         for definition in def_list:
             definition_new = deepcopy(definition)
-            definition_new.system = title
+            definition_new.system = name
             self.definitions.append(definition_new)
 
-        fields=['import', 'class', 'class.options']
-        if extras is not None:
-            for key in extras:
-                fields.append(key)
+        if fields is None:
+            fields = ['']
 
         self.entries=[]
         for key in fields:
             text=StringVar()
             text.set('')
-            self.entries.append((key, text, Entry(self, textvariable=text, width=15)))
+            self.entries.append((key, text, Entry(self, textvariable=text, width=25)))
             
     def grid(self, **arguments):
         LabelFrame.grid(self, arguments, sticky=(N+W+E))
@@ -54,18 +52,18 @@
     def doLayout(self):
         irow=0
         for (key, text, textfield) in self.entries:
-            Label(self, text=key, width=15, anchor=W).grid(row=irow)
-            textfield.grid(row=irow, column=1)
+            Label(self, text=key, width=15, anchor=W).grid(row=irow,  column=2)
+            textfield.grid(row=irow, column=3)
             textfield.bind("<KeyRelease>", self.update)
             irow+=1
         
         scrollbar=Scrollbar(self, orient=VERTICAL)
-        self.listbox=Listbox(self, yscrollcommand=scrollbar.set, height=0, width=15, selectmode=SINGLE)
-        self.listbox.grid(row=0, column=2, rowspan=irow, columnspan=2, sticky=NS)
+        self.listbox=Listbox(self, yscrollcommand=scrollbar.set, height=5, width=15, selectmode=SINGLE, exportselection=0)
+        self.listbox.grid(row=0, column=0, rowspan=irow, sticky=NS)
         self.listbox.bind('<ButtonRelease-1>', self.select)
 
         scrollbar.config(command=self.listbox.yview)
-        scrollbar.grid(row=0, column=4, rowspan=irow, sticky=NS)
+        scrollbar.grid(row=0, column=1, rowspan=irow, sticky=NS)
        
         for item in self.definitions:
             self.listbox.insert(END, item.name)
@@ -93,80 +91,65 @@
     def __init__(self, master):
         Frame.__init__(self, master)
         
-        #--- Structures ---
-        bigFrame = LabelFrame(self, text='Structures')
-        irow=0; icol=0
-        def_list=embedDefinitions
-        extras=['method', 'method.options', 'calculator']
-        frame=MyLabelFrame(self, bigFrame, 'Embed', def_list, extras)
-        frame.doLayout()
-        frame.grid(row=irow, column=icol)
-        self.frames.append(frame)
-
-        icol+=1;
-        for label in ['System', 'Cluster']:
-            def_list=strDefinitions
-            frame=MyLabelFrame(self, bigFrame, label, def_list, None)
-            frame.doLayout()
-            frame.grid(row=irow, column=icol)
-            self.frames.append(frame)
-            icol+=1
-        bigFrame.grid(row=0, columnspan=3)
-
-        #--- Methods ---
-        bigFrame = LabelFrame(self, text='Methods')
-        irow=0; icol=0
-        def_list=[QmxCalcDefinition()]
-        frame=MyLabelFrame(self, bigFrame, 'Qmx', def_list, None)
-        frame.doLayout()
-        frame.grid(row=irow, column=icol)
-        self.frames.append(frame)
-
-        icol+=1
-        for label in ['Low-Level', 'High-Level']:
-            def_list=calcDefinitions
-            frame=MyLabelFrame(self, bigFrame, label, def_list, None)
-            frame.doLayout()
-            frame.grid(row=irow, column=icol)
-            self.frames.append(frame)
-            icol+=1
-        bigFrame.grid(row=1, columnspan=3)
+        fields=['import', 'class', 'class.options']
         
-        #--- General ---
-        bigFrame = LabelFrame(self, text='General')
-        extras=['method', 'method.options']
-        frame=MyLabelFrame(self, bigFrame, 'Job', jobDefinitions, extras)
-        frame.doLayout()
-        frame.grid()
-        self.frames.append(frame)
-        bigFrame.grid(row=2, sticky=NW)
-
-        frame=LabelFrame(self, text="PREVIEW")
-        scrollbar=Scrollbar(frame, orient=VERTICAL)
+        irow=0
+        for calcLabel, strLabel, strDefinitions in (
+                ("high-level",  "cluster", strClusterDefinitions), 
+                ("low-level",  "system", strSystemDefinitions)):
+            #column
+            icolLocal=0;
+            #build outFrame
+            outFrame = LabelFrame(self, text=calcLabel+" method - "+strLabel)
+            outFrame.grid(row=irow, columnspan=3,  sticky=W, pady=5)
+            irow+=1
+            
+            for title, def_list in ((calcLabel, calcDefinitions),  (strLabel, strDefinitions)):
+                # build frame
+                myFrame=MyLabelFrame(self, outFrame, title, title, def_list, fields)
+                myFrame.doLayout()
+                myFrame.grid(row=0, column=icolLocal, padx=5, pady=5)
+                # add frame for analyssis
+                self.frames.append(myFrame)
+                # move column
+                icolLocal += 1
+            
+        outFrame = LabelFrame(self, text="General")
+        outFrame.grid(row=irow, columnspan=3, pady=5); irow+=1
+        for title, name, def_list, extra_list,  irowLocal, icolLocal,  irspan in (
+            ("Optimization Method",  "job", jobDefinitions, fields + ['method', 'method.options'], 0, 0, 2), 
+            ("Embedding Method", "embed", embedDefinitions, fields + ['set_calculator'], 0, 1, 1), 
+            ("Hybrid Method", "hybrid", [QmxCalcDefinition()], fields, 1, 1, 1)):
+            # build frame
+            myFrame=MyLabelFrame(self, outFrame, title, name, def_list, extra_list)
+            myFrame.doLayout()
+            myFrame.grid(row=irowLocal, column=icolLocal,  rowspan=irspan, padx=5, pady=5)
+            # add frame for analyssis
+            self.frames.append(myFrame)
+            
+        outFrame=LabelFrame(self, text="PREVIEW")
+        scrollbar=Scrollbar(outFrame, orient=VERTICAL)
         scrollbar.set(0.0, 1.0)
-        self.text=Text(frame, yscrollcommand=scrollbar.set, width=60, heigh=10, state=DISABLED)
+        self.text=Text(outFrame, yscrollcommand=scrollbar.set, width=80, heigh=10, state=DISABLED)
         self.text.grid(sticky=(N+W+S+E))
         scrollbar.config(command=self.text.yview)
         scrollbar.grid(row=0, column=1, sticky=NS)        
-        frame.grid(row=2, column=1)
+        outFrame.grid(row=irow, column=0, sticky=W, pady=5);
 
-        frame = Frame(self)
-
-        #image = PhotoImage(file="exit.gif")
-        buttonExit=Button(frame, text="Quit!", command=self.quit)
-        #buttonExit.image = image
+        myFrame = Frame(self)
+        buttonExit=Button(myFrame, text="Quit!", command=self.quit)
         buttonExit.grid(row=1, column = 1, rowspan=3, sticky= W+E+N+S)
         
-        buttonSave=Button(frame, text='save settings', state=DISABLED)
+        buttonSave=Button(myFrame, text='save settings', state=DISABLED)
         buttonSave.grid(row=1, column = 2, sticky=(EW))
 
-        buttonScript=Button(frame, text='Write qmx.in input file', command=self.writeScript)
+        buttonScript=Button(myFrame, text='Write SCRIPT file', command=self.writeScript)
         buttonScript.grid(row=2, column = 2, sticky=EW)
 
-        buttonPython=Button(frame, text='Write qmx.py script file', command=self.writePython)
+        buttonPython=Button(myFrame, text='Write PYTHON file', command=self.writePython)
         buttonPython.grid(row=3, column = 2, sticky=EW)
 
-        frame.grid(row=2, column=2, sticky=SE)
+        myFrame.grid(row=irow, column=2, sticky=SE)
 
         #--- MainWindow ---
         self.grid(sticky=(W+N+S+E))
@@ -176,7 +159,10 @@
         quit()
     
     def writePython(self):
-        file = open("qmx.py","w")
+        fileName = tkFileDialog.asksaveasfile(mode='w')
+        if fileName is None:
+            return
+        file = open(fileName,"w")
         definitions=[]
         for frame in self.frames:
             definition = frame.getSelection()
@@ -189,7 +175,10 @@
         
         
     def writeScript(self):
-        file = open("qmx.in","w")
+        fileName = tkFileDialog.asksaveasfile(mode='w')
+        if fileName is None:
+            return
+        file = open(fileName,"w")
         for frame in self.frames:
             definition = frame.getSelection()
             if definition is None:
Index: prepareQMX/qmxCALC.py
===================================================================
--- prepareQMX/qmxCALC.py	(revision 15)
+++ prepareQMX/qmxCALC.py	(revision 16)
@@ -48,4 +48,87 @@
         self.keywords['class.options']='high_level, low_level'
 
 #-------------------------------------------------------------------------------
-calcDefinitions = [VASPCalcDefinition(), TURBOMOLECalcDefinition(), MOPACCalcDefinition()]
+class LJCalcDefinition(Definition):
+    def __init__(self):
+        self.name='Lennard Jones'
+        self.system='LJ'
+        self.keywords = {}
+        self.keywords['import']='ase.calculators.lj'
+        self.keywords['class']='LennardJones'
+
+#-------------------------------------------------------------------------------
+class EMTCalcDefinition(Definition):
+    def __init__(self):
+        self.name='EMT'
+        self.system='EMT'
+        self.keywords = {}
+        self.keywords['import']='ase.calculators.emt'
+        self.keywords['class']='EMT'
+
+#-------------------------------------------------------------------------------
+class SiestaCalcDefinition(Definition):
+    def __init__(self):
+        self.name='Siesta'
+        self.system='Siesta'
+        self.keywords = {}
+        self.keywords['import']='ase.calculators.siesta'
+        self.keywords['class']='Siesta'
+
+#-------------------------------------------------------------------------------
+class DacapoCalcDefinition(Definition):
+    def __init__(self):
+        self.name='Dacapo'
+        self.system='Dacapo'
+        self.keywords = {}
+        self.keywords['import']='ase.calculators.dacapo'
+        self.keywords['class']='Dacapo'
+
+#-------------------------------------------------------------------------------
+class AimsCubeCalcDefinition(Definition):
+    def __init__(self):
+        self.name='AimsCube'
+        self.system='AimsCube'
+        self.keywords = {}
+        self.keywords['import']='ase.calculators.ase'
+        self.keywords['class']='AimsCube'
+
+#-------------------------------------------------------------------------------
+class AimsCalcDefinition(Definition):
+    def __init__(self):
+        self.name='Aims'
+        self.system='Aims'
+        self.keywords = {}
+        self.keywords['import']='ase.calculators.aims'
+        self.keywords['class']='Aims'
+
+#-------------------------------------------------------------------------------
+class ExcitingCalcDefinition(Definition):
+    def __init__(self):
+        self.name='Exciting'
+        self.system='Exciting'
+        self.keywords = {}
+        self.keywords['import']='ase.calculators.exciting'
+        self.keywords['class']='Exciting'
+
+#-------------------------------------------------------------------------------
+class DftbCalcDefinition(Definition):
+    def __init__(self):
+        self.name='Dftb'
+        self.system='Dftb'
+        self.keywords = {}
+        self.keywords['import']='ase.calculators.dftb'
+        self.keywords['class']='Dftb'
+
+#-------------------------------------------------------------------------------
+class GaussianCalcDefinition(Definition):
+    def __init__(self):
+        self.name='Gaussian'
+        self.system='Gaussian'
+        self.keywords = {}
+        self.keywords['import']='ase.calculators.gaussian'
+        self.keywords['class']='Gaussian'
+
+calcDefinitions = \
+  [ VASPCalcDefinition(), TURBOMOLECalcDefinition(), MOPACCalcDefinition(), LJCalcDefinition(), \
+    EMTCalcDefinition(), SiestaCalcDefinition(), DacapoCalcDefinition(), AimsCubeCalcDefinition(), \
+    AimsCalcDefinition(), ExcitingCalcDefinition(), DftbCalcDefinition(), GaussianCalcDefinition() ]
Index: ase/calculators/__init__.py
===================================================================
--- ase/calculators/__init__.py	(revision 15)
+++ ase/calculators/__init__.py	(revision 16)
@@ -16,12 +16,13 @@
 from ase.calculators.test import numeric_force, numeric_forces, TestPotential
 from ase.calculators.qmx import Qmx
 from ase.calculators.mopac import Mopac
+from ase.calculators.gaussian import Gaussian
 
 if _deprecate_things_from_ase_module:
     from ase.utils.deprecate import Deprecate
     _locals = locals()
     for name in ['LennardJones', 'EMT', 'Siesta', 'Dacapo', 'Vasp',
-                 'Aims', 'AimsCube', 'Turbomole', 'Exciting', 'Dftb', 'Qmx', 
+                 'Aims', 'AimsCube', 'Turbomole', 'Exciting', 'Dftb', 'Qmx', 'Gaussian', 
                  'SinglePointCalculator', 'numeric_force', 'numeric_forces',
                  'TestPotential']:
         obj = _locals[name]
Index: ase/calculators/qmx.py
===================================================================
--- ase/calculators/qmx.py	(revision 15)
+++ ase/calculators/qmx.py	(revision 16)
@@ -30,20 +30,26 @@
         self.print_forces = print_forces
 
     def get_energy_subsystem(self, path, calculator, atoms, force_consistent):
+        # writing output
+        line = "running energy in: " + path + "\n"
+        sys.stderr.write(line)
         # go to directory and calculate energies
-        print "running energy in: ", path
         os.chdir(path)
         atoms.set_calculator(calculator)
         energy = atoms.get_potential_energy()
+        # return path and result
         os.chdir("..")
         return energy
 
     def get_forces_subsystem(self, path, calculator, atoms):
+        # writing output
+        line = "running forces in: " + path + "\n"
+        sys.stderr.write(line)
         # go to directory and calculate forces
-        print "running forces in: ", path
         os.chdir(path)
         atoms.set_calculator(calculator)
         forces = atoms.get_forces()
+        # return path and result
         os.chdir("..")
         return forces
 
@@ -55,7 +61,7 @@
         # calculate energies
         energy = e_sys_lo - e_cl_lo + e_cl_hi
         # print energies
-        print "%20s = %15s - %15s + %15s" %("E(C:S)", "E(S-MM)", "E(C-MM)", "E(C-QM)")
+        print "%20s = %15s - %15s + %15s" %("E(C:S)", "E(S,LL)", "E(C,LL)", "E(C,HL)")
         print "%20f = %15f - %15f + %15f" %(energy, e_sys_lo, e_cl_lo, e_cl_hi)
         # set energies and return
         if force_consistent:
@@ -145,3 +151,6 @@
             print
         
         return f_sys_lo
+
+    def get_stress(self, atoms):
+    	return None
Index: ase/calculators/gaussian.py
===================================================================
--- ase/calculators/gaussian.py	(revision 0)
+++ ase/calculators/gaussian.py	(revision 16)
@@ -0,0 +1,251 @@
+"""
+This module is the Gaussian module for ASE,  based on the PyMD Gaussian Module by R. Bulo,  
+implemented by T. Kerber
+
+ @author:       Rosa Bulo
+ @organization: Vrije Universiteit Amsterdam
+ @contact:      bulo@few.vu.nl
+
+Torsten Kerber,  ENS LYON: 2011,  07,  11
+
+This work is supported by Award No. UK-C0017,  made by King Abdullah
+University of Science and Technology(KAUST)
+"""
+import os,  sys,  string
+
+import numpy as np
+from ase.units import Hartree,  Bohr
+from ase.calculators.general import Calculator
+import os, sys, re, math, commands, subprocess, time
+
+class Gaussian(Calculator):
+    """
+    Class for representing the settings of a ReaxffForceJob
+        
+    It contains the following variables:
+                        
+    basis:  
+        A string representing the basisset to be used
+
+    functional:
+        A string representing the DFT functional or wavefunction
+        method to be used.
+
+    """
+
+    def __init__(self, functional='BP86', basis='TZVP', atoms=None, inputfile=None):
+        self.functional = functional
+        self.basis = basis
+        self.restartdata = None
+        self.forces = None
+        self.atoms_old = None
+        if inputfile != None:
+            self = self.get_settings(inputfile)                
+        self.atoms = atoms
+
+    def set_functional(self, functional):
+        self.functional = functional
+
+    def set_basis(self, basis):
+        self.basis = basis
+
+    def setup_job(self,  type='force'):
+        ForceJob.setup_job(self, type=type)
+ 
+        os.chdir(self.dirname)
+        self.myfiles = GaussianFileManager()
+        # Add the restart restuls to the filemanager. The checksum will be empty('')
+        if self.restartdata != None:
+            self.r_m = GaussianResults(self.myfiles)
+            self.myfiles.add_results(self.r_m, dir=self.restartdata+'/')
+
+        os.chdir(self.dir)
+
+    def write_parfile(self, filename):
+        parfile = open(filename, 'w')
+        parfile.write(self.prm)
+        parfile.close()
+
+    def write_moldata(self,  type='force'):
+        out = commands.getoutput('mkdir SETUP')
+        out = commands.getoutput('ls SETUP/coord.*').split()
+        pdb = False
+        psf = False
+        prm = False
+
+        # Change labels variable if necessary
+        norderats = 0
+#        for el in self.labels:
+#            norderats += 1
+#            if norderats != self.atoms.pdb.get_number_of_atoms():
+#                self.labels = self.get_labels()
+
+        self.write_input(type)
+
+    def write_input(self, type):
+        input = self.get_input_block(type)
+        pyfile = open('ASE-gaussian.com', 'w')
+        pyfile.write(input)
+        pyfile.close()
+
+    def get_inp(self,  filename):
+        input = None
+        lis = commands.getoutput('ls %s'%(filename)).split('\n')
+        if not re.search('ls:', lis[0]):
+            infile = open(filename)
+            input = infile.readlines()
+            infile.close()
+        return input
+
+    def run(self, type='force'):
+        startrun = time.time()
+        # I have to write the input every timestep,  because that is where the coordinates are
+        # self.write_input(type)
+        self.write_moldata(type)
+
+        # Make sure that restart data is copied to the current directory
+#        if self.r_m != None:
+#            self.r_m.files.copy_job_result_files(self.r_m.fileid)
+
+        out = commands.getoutput('g09 ASE-gaussian.com')
+        endrun = time.time()
+        #print 'timings: ', endrun-startrun
+        # Write the output
+#        outputfile = open('ASE-gaussian.log')
+#        outfile = open('ASE-gaussian.out',  'w')
+#        for line in outputfile:
+#            outfile.write(line)
+#        outputfile.close()
+#        outfile.close()
+        # End write output
+        self.energy = self.read_energy()
+        if type == 'force':
+            self.forces = self.read_forces()
+        self.energy_zero = self.energy
+
+        # Change the results object to the new results
+        input = self.get_input_block(type)
+
+    def read_energy(self,  charges=True):
+        hartree2kcal = 627.5095
+
+        outfile = open('ASE-gaussian.log')
+        lines = outfile.readlines()
+        outfile.close()
+
+        chargelines = 0
+        if charges:
+            nats = len(self.atoms)
+            block = ''
+        for line in lines:
+            if re.search('SCF Done', line):
+                words = line.split()
+                energy = float(words[4]) * hartree2kcal
+            if charges:
+                if re.search(' Mulliken atomic charges:', line):
+                    chargelines += 1
+            if chargelines > 0 and chargelines <= nats+2:
+                chargelines += 1
+                block += line
+   
+        if charges:
+            chargefile = open('charge.out', 'a')
+            chargefile.write(block)
+            chargefile.close()
+
+        return energy
+
+    def read_forces(self):
+        factor = Hartree / Bohr
+        factor = 1.0
+
+        outfile = open('ASE-gaussian.log')
+        lines = outfile.readlines()
+        outfile.close()
+
+        outputforces = None
+        forces = None
+        nats = len(self.atoms)
+
+        for iline, line in enumerate(lines):
+            if re.search('Forces \(Ha', line):
+                forces = np.zeros((nats,  3),  float)
+                for iat in range(nats):
+                    forceline = lines[iline+iat+3]
+                    words = forceline.split()
+                    for idir,  word in enumerate(words[2:5]):
+                        forces[iat,  idir] = float(word) * factor
+                break
+        return forces
+
+    def get_input_block(self, type):
+        block = ''
+        block += '%chk=ASE-gaussian.chk\n'
+        block += '%Mem=256MB\n'
+        block += '%nproc=32\n'
+        block += 'MaxDisk=16gb\n'
+
+        block += '#'
+        block += '%s'%(self.functional)
+        block += '/'
+        block += '%s'%(self.basis)
+
+        if type == 'force':
+            block += 'Force '
+#        if self.r_m != None:
+#            block += 'Guess=Read'
+
+        block += '\n\nGaussian job\n\n'
+
+        charge = sum(self.atoms.get_charges())
+        block += '%i '%(charge)
+        block += '%i '%(1)
+        block += '\n'
+
+        coords = self.atoms.get_positions()
+        for i, at in enumerate(self.atoms.get_chemical_symbols()):
+            block += ' %2s'%(at)
+            coord = coords[i]
+            block += '  %16.5f %16.5f %16.5f'%(coord[0], coord[1], coord[2])
+            block += '\n'
+
+        if self.atoms.get_pbc().any(): 
+            cell = self.atoms.get_cell()
+           
+            for v in cell: 
+                block += 'TV %8.3f %8.3f %8.3f\n'%(v[0], v[1], v[2])
+        block += '\n'
+
+        return block
+
+
+    def get_settings(self,  filename=None):
+        settings = GaussianSettings()
+
+        if filename != None or filename == '':
+            input = self.get_inp(filename)
+        else:
+            settings.set_functional("")
+            return settings
+    
+        if input == None:
+            settings.set_functional("")
+            return settings
+
+        for i, line in enumerate(input):
+            if len(line.strip()) == 0:
+                continue
+            if line[0] == '#':
+                keyline = i
+                words = line[1:].split()
+            else:
+                settings.functional = words[0].split('/')[0]
+                settings.basis = words[0].split('/')[1]
+            break
+            return settings
+                
+    def update(self,  atoms):
+        if self.atoms_old != atoms:
+            self.atoms = atoms.copy()
+            self.atoms_old = atoms.copy()
+            self.run()
Index: ase/embed.py
===================================================================
--- ase/embed.py	(revision 15)
+++ ase/embed.py	(revision 16)
@@ -16,7 +16,7 @@
 
 class Embed(Atoms):
     #--- constructor of the Embed class ---
-    def __init__(self, system, cluster, cell_cluster = "Auto", cluster_pos = True):
+    def __init__(self, system, cluster, cell_cluster = "Auto", cluster_pos = True, linkType = 'H'):
         super(Embed, self).__init__()
         # define the atom map
         self.atom_map_sys_cl = []
@@ -39,6 +39,8 @@
         # set the cell of the system
         self.set_cell(system.get_cell())
         self.cell_cluster = cell_cluster
+        
+        self.embed(linkType)
         return
 
     #--- set the cluster ---
@@ -98,12 +100,12 @@
         return self.atoms_system
 
     #--- Embedding ---
-    def embed(self):
+    def embed(self, linkType):
         # is the cluster and the host system definied ?
         if self.atoms_cluster is None or self.atoms_system is None:
             return
         self.find_cluster()
-        self.set_linkatoms()
+        self.set_linkatoms(linkType)
         print "link atoms found: ", len(self.linkatoms)
         if self.cell_cluster == "System":
             self.atoms_cluster.set_cell(self.atoms_system.get_cell())
@@ -165,7 +167,7 @@
             if bSysOnly:
                 self.atom_map_sys_cl[iat_sys] = -1
 
-    def set_linkatoms(self, tol=15., linkAtom=None, debug=False):
+    def set_linkatoms(self, linkType, tol=15.):
         # local copies of xyz coordinates to avoid massive copying of xyz objects
         xyzs_cl=[]
         for atom_cl in self.atoms_cluster:
@@ -173,9 +175,6 @@
         xyzs_sys=[]
         for atom_sys in self.atoms_system:
             xyzs_sys.append(atom_sys.get_position())
-        # set the standard link atom
-        if linkAtom is None:
-            linkAtom ='H'
         # number of atoms in the cluster and the system
         nat_cl=len(self.atoms_cluster)
         nat_sys=len(self.atoms_system)
@@ -235,10 +234,6 @@
                     r = np.sqrt(np.dot(xyz_diff, xyz_diff))
                     # ratio of the distance to the sum of covalent radius
                     f = r / (r_cl + r_sys)
-                    if debug:
-                        print "Covalent radii = ",r_cl, r_sys
-                        print "Distance ", f
-                        print "tol = ",(1+tol/100.),(1-tol/100.),(1-2*tol/100.)
                     if f <= (1+tol/100.) and f >= (1-2*tol/100.):
                         s = cell_L, self.atom_map_cl_sys[iat_cl], iat_sys, r_cl
                         bonds.append(s)
@@ -246,7 +241,7 @@
                     if f <= (1-2*tol/100.):
                         raise RuntimeError("QMX: The cluster atom", iat_cl, " and the system atom", iat_sys, "came too close")
 
-        r_h = covalent_radii[atomic_numbers[linkAtom]]
+        r_h = covalent_radii[atomic_numbers[linkType]]
         for bond in bonds:
             cell_L, iat_cl_sys, iat_sys, r_cl = bond
             # assign the tags for the border atoms
@@ -259,7 +254,7 @@
             # determine position of the link atom
             xyz_diff += xyzs_sys[iat_cl_sys]
             # create link atom
-            atom = Atom(symbol=linkAtom, position=xyz_diff, tag=50)
+            atom = Atom(symbol=linkType, position=xyz_diff, tag=50)
             # add atom to cluster
             self.atoms_cluster.append(atom)
             # add atom to the linkatoms
