Navigation:

Basin-Hopping

The basin-hopping algorithm is a stochastic simulation technique to efficiently find local minima on a potential energy surface (PES). The PES for a complex containing small molecules, such as water, can be relatively simple and will normally support only a few minima. For such systems finding this minimum is not a problem and any minimization technique will work. However the PESs for larger molecules, or larger numbers of molecules can support large numbers of minima. For such systems more robust methods are needed to find these minima.

References

See the innumerable papers by David Wales for details of how this algorithm works. Here are a select few:

Sample Orient command file

Data needed

In the example used here you will need to use the data for Model(1) of the pyridine dimer. Please pick up the data from that link. Make sure that the example provided works, and that the output agrees with the sample output provided.

warning

The data for Model(1) can be read by Orient 4.8 and earlier versions only.

Pyridine dimer : Orient 4.8 commands

warning This command file is suitable for Orient 4.8. It will not work for Orient 4.9.

pyridine_n2_basinhopping.ornt

Title pyr2 : Basin-hopping

UNITS BOHR

Parameters
      Sites     26 polarizable     26
      S-functions 100000
      Alphas 100000
      Parameter-sets 100000
      Pairs 100000
End

Types
      H1         Z     1
      H2         Z     1
      H3         Z     1
      H4         Z     1
      H5         Z     1
      N          Z     7
      C1         Z     6
      C2         Z     6
      C3         Z     6
      C4         Z     6
      C5         Z     6
End

Molecule  pyridine1 at  0.0 0.0 0.0 rotated by 0.0 about 0.0 0.0 1.0
  #include ./Model1/pyr.mom
End
Edit pyridine1
  #include ./Model1/pyr.axes
End
Polarizabilities for pyridine1
   ! Assumed that the pols are in the local-axes
   ! So they are read after axes are defined
   Read rank 1
     #include ./Model1/pyr.pol
   End
   Limit rank 1 for +++
      H1  H2  H3  H4  H5  N   C1  C2  C3  C4   +++
      C5 
End

Molecule  pyridine2 at  0.0 0.0 10.0 rotated by 0.0 about 0.0 0.0 1.0
  #include ./Model1/pyr.mom
End
Edit pyridine2
  #include ./Model1/pyr.axes
End
Polarizabilities for pyridine2
   ! Assumed that the pols are in the local-axes
   ! So they are read after axes are defined
   Read rank 1
     #include ./Model1/pyr.pol
   End
   Limit rank 1 for +++
      H1  H2  H3  H4  H5  N   C1  C2  C3  C4   +++
      C5 
End

Units Bohr Hartree

Pairs
   #include ./Model1/pyr2.pot
End

Units Bohr kJ/mol

Switch Induce On Iterate On
Options
  Induction Iterations 60  Convergence 1e-12
End

Options
  bfgs
  Iterations 1000
  convergence 1d-5
  show none
  plot none
end

Time

Basin-hopping
  ! restart save.cfg
  ! Verbose
  Allocate
    minima 100
  end
  Temperature 500
  Seed 777  ( optional )
  Steps
    number 500
    !  A serious calculation would need many more steps than this
    maxdisp 0.25
    maxrot 30 ( degrees )
    blocks rotations 50 translations 50 both 100
  end
  Difference energy 1e-5 MI 0.5
  Container radius 7.0
  Quenching
    reject -100
  end
  Files
    xyz    pyr2_model1_n2_run1.xyz
    minima pyr2_model1_n2_run1.geom
  end
end

Time

Finish

Run this job using (I have used orient-4.8 as the executable. You need to make sure that you use the correct executable):

This can take a while as the above example uses 500 steps in the search process. The output is also quite long, but we can use the **Cluster** program that is part of the CamCASP suite to analyse it and provide a sensible summary of the results.

=====Analysing the basin-hopping results using Cluster=====

The Cluster program that is part of the CamCASP suite of codes can be used to perform an analysis of the results of the Basin-Hopping simulations. Here is a sample code: <code | reorient-analyse.clt> Title Re-orient and analyse clusters

Global

End

Orient

End

Finish

</code> In this example, the outputs of two Basin-Hopping simulations is analysed. Cluster uses the files containing the **minima**. These will be the files names with the command:

in the Basin-Hopping commands in the Orient command file. Here these are RDX_run1.geom and RDX_run2.geom. The **APPEND** command tells Cluster to append the contents of the second file after those of the first. The **OFFSET 100** command tells it to change the numbers of the minima in the second file by 100. This ensures that we have unique numbers to all minima. Of course, if your files contain more than 100 minima each then a larger offset will be needed.

Run this using

It will run in a few seconds or so. The output will be fairly large, but the important parts are at the bottom and should be something like <code | partial output> Configuration statistics ! Config List name : CfgList-1 ! Configs read from file : RDX_run1.geom UNITS BOHR DEGREE KJ/MOL AMU BOHR^2 !


! INDEX ENERGY I_XX I_YY I_ZZ NUM-MOLS NUM-SIMILAR SYMMETRY

!


End

Finish

}}}

This listing is a summary of the minima that Cluster considers as unique. The fields include:

The clustering is done using the commands in the ANALYSE block:

  Analyse
    Sort
    Similarity  Moments  & Energy 
    Energy-Sigma 0.5 kJ/mol
    Moment-Sigma 40.0 AMU
    Similarity-cutoff 0.5
    Details
  End

Here the structure similarity is based on the moments of inertia **and** the energies. Each quantity is given a tolerance set by a standard deviation. This defines a similarity probability as a normal distribution with the specified width. Structures are deemed to be similar to a probability. The cutoff level is set by **Similarity-cutoff**. Here it is $0.5$ so structures with a similarity probability $p \ge 0.5$ are deemed to be similar.

The clustering will depend on the standard deviations chosen. This needs to be tuned to the problem of interest. It is quite likely that the **Energy-Sigma** is too large in the example above, and quite likely that the **Moment-Sigma** is too small. You need to experiment to decide.