Gamess (US) frequently asked questions. Part 7: How to distinguish alpha from beta orbitals in the $VEC deck

Each line in a $VEC group contains the coefficients of five basis functions for a given orbital. These are formatted in a special way, with seven numbers in each line. These numbers are:

1st) the number of the orbital to which the coefficients belong (written with at most two characters, so that 1 means orbital 1, .. , 99 means orbital 99, 00 means orbital 100) . This number is repeated in the beginning of every line, until all coefficients for that orbital have been written

2nd) this number tells the program how to assign the coefficients to the basis functions. "1" means that the coefficients are for basis functions 1-5, "2" means that the coefficients are for basis functions 5-10, etc. In general , that number "n" directs the program to assign the five coefficients present in the line to basis functions 5*(n-1)+1 to 5*n.

3rd to 7th) coefficients of five basis functions

BETA orbitals are punched as a group immediately after all ALPHA orbitals.

This format entails that in molecules with more than 100 orbitals the $VEC group contains several blocks with the same 1st number. For example, in a molecule with 200 orbitals, alpha orbital 27 is described by the first block of lines beginning with "27", and alpha orbital 127 is described by the SECOND block of lines beginning with "27".

I usually find the beginning of the BETA orbitals by repeating a search for the string " 1 1" : if that string is preceded by a block beginning with "00 1", it usually refers to orbitals 101, or 201, etc. (the exception being those systems with exactly 100, 200, etc. orbitals). If string " 1 1" is NOT preceded by a block beginning with "00 1", you are sure to have found the beginnning of the BETA orbitals

Gamess (US) frequently asked questions part 6: Obtaining proper SCF convergence (Anti-)ferromagnetic coupled Fe-S clusters

Obtaining SCF convergence of FeS clusters is a very demanding task.

The problem in FeS clusters is the arrangement of spins on the Fe atoms: if you have a cluster with 4 Fe atoms, each of them with 5 up-spins, and a total spin of zero, the arrangement of spins on the atoms could be

• Fe1 and Fe2  up-spin, Fe3 and Fe4 down-spin; or
• Fe1 and Fe4  up-spin, Fe2 and Fe3 down-spin; or
• Fe1 and Fe3  up-spin, Fe2 and Fe4 down-spin;

The problem is compounded if you have a mixture of Fe²⁺ and Fe³⁺, which may lead to 12 (or more) different spin arrangements, depending on the number of Fe²⁺ atoms. However, if you have a good guess SCF for one instance instance, you may simply substitute the coordinates of Fe2 with those of Fe4 to get a comparably good guess for the second instance, and so forth... This is the approach suggested by Greco, Fantucci, Ryde, de Gioia (2011) Int. J. Quantum Chem. 111, 3949-3960. Obtaining the guess for one of the instances is in itself quite difficult, and I usually follow the approach outlined by Szilagyi, R. K. and Winslow, M. A. (2006) J. Comput. Chem., 27: 1385–1397  .

It goes like this:

- obtain orbitals for bare Fe²⁺, Fe³⁺, S^2-, and isolated ligands, with proper spins on the Fe atoms (5/2 for Fe³⁺, 2 for Fe²⁺)

- Manually split the "alpha/up" and "beta/down" portions of the resulting  $VEC groups. For example, assuming you have a system with three Fe atoms (two Fe²⁺ and one Fe³⁺) with total spin S=5/2 and the $VEC groups for bare Fe²⁺ and bare Fe³⁺, you should cut the $VEC groups of Fe²⁺ and Fe³⁺ as:

$VEC  for the alpha (up) electrons of Fe²⁺   (let's call it "Fe2+_5_d_electrons")
$VEC  for the alpha (up) electrons of Fe³⁺   (let's call it "Fe3+_5_d_electrons")
$VEC  for the beta (down) electrons of Fe²⁺   (let's call it "Fe2+_1_d_electron")
$VEC  for the beta (down) electrons of Fe³⁺   (let's call it "Fe3+_0_d_electrons")

The total spin S=5/2 in this sample problem implies that  both Fe²⁺ atoms spins should annull each other, i.e., one Fe²⁺ is mostly "up" and the other is mostly "down". Building the new guess for the "up" electrons should therefore include:

"Fe2+_5_d_electrons" for one of the  Fe²⁺ ions,
"Fe2+_1_d_electrons" for the other  Fe²⁺,
"Fe3+_5_d_electrons" for the Fe³⁺

Building the new guess for the "down" electrons should  include:
"Fe2+_1_d_electrons" for the FIRST Fe²⁺ ions,
"Fe2+_5_d_electrons" for the other Fe²⁺,
"Fe3+_0_d_electrons" for the Fe³⁺

- combine the orbitals using the small utility called combo, which you may obtain from Alex Granovsky's Firefly website.

- Manually paste the "alpha" and "beta" guesses  into a single $vec group, which would be the proper guess.

- cross all your fingers and toes, and expect it to converge into the proper state. If it does not converge, change convergers (SOSCF=.T. DIIS=.F.), onset of SOSCF (SOGTOL=1e-3) , etc.

- after SCF optimization using this guess, manually scramble the ordering of Fe atoms in your input, to ascertain whether a lower energy solution can be obtained with a different spin distribution.



Good Luck!

Gamess (US) frequently asked questions Part 5: "THE VIBRATIONAL ANALYSIS IS NOT VALID"

Gamess (US) and Firefly by default assume geometric convergence has been achieved when the maximum gradient is below 1e-4 and the RMS gradient is smaller  than 1/3 of the maximum gradient.  This convergence criterion may be changed by the user with


 $STATPT OPTTOL=<your desired convergence criterion> $END


It is well known that the vibrational analysis is strictly valid mathematically when the Hessian is computed in true stationary points (i.e when the gradient is exactly equal to zero). If the maximum gradient is sufficiently close to zero, the vibrational analysis (although not absolutely correct) is still close enough to the "true" solution for all practical purposes.



This introduction brings us to today's FAQ. A recurring question in both the Gamess-US list and the Firefly forums concerns the message often printed by the program after a vibrational analysis:

*THIS IS NOT A STATIONARY POINT ON THE MOLECULAR PES THE VIBRATIONAL ANALYSIS IS NOT VALID*

This message arises from the way gradients are analyzed by Gamess: gradients are originally computed in one set of coordinates (cartesian coordinates, I believe) , and then transformed into the  coordinate system specified by the user. Optimizations stop when the "transformed gradient" lies below OPTTOL, but Gamess uses the original, non-transformed, gradient to decide whether to consider the geometry as a stationary point on the molecular PES.  Therefore, if  the geometry is converged, the scary message in capital letters above may be safely disregarded. When in doubt, simply decrease your OPTTOL value, continue the optimization and re-compute the hessian.

Gamess (US) frequently asked questions Part 4: The rungms script

Third guest post by Kirill Berezovsky (Petrozadovsk State University).

Gamess (US) is run though a script provided in the installation package. This script MUST be adapted by the user. The rungms-script below has been adapted to be used with the Gamess (US) installation described in earlier posts. Happy computing :-)

RUNGMS script

#!/bin/csh

set SCR=/scr/$USER

set USERSCR=~$USER/scr

set GMSPATH=/usr/local/gamess

set JOB=$1

set VERNO=$2

if ($VERNO == cpu) set TARGET=mpi

if ($VERNO == gpu) set TARGET=ga

if (null$VERNO == null) set VERNO=cpu

set master=`hostname`

printf "\n * Started at:           `date`"

set DSK=`df -m $SCR | awk 'NR==2{print$4}'`

printf "\n * Available disk space: $DSK MB"

printf "\n * Temporary files in:   $SCR and in $USERSCR\n\n"

limit stacksize 8192

if ($JOB:r.inp == $JOB) set JOB=$JOB:r

if (-e $JOB.inp) then

   cp  $JOB.inp  $SCR/$JOB.F05

else

   echo "Input file $JOB.inp not found"

   exit 4

endif

source $GMSPATH/gms-files.csh

if (-e $HOME/.gmsrc) source $HOME/.gmsrc

set ngddi=`grep -i '^ \$GDDI' $SCR/$JOB.F05 | grep -iv 'NGROUP=0 ' | wc -l`

if ($ngddi > 0) then

   set GDDIjob=true

   echo "This is a GDDI run, keeping various output files on local disks"

   set echo

   setenv  OUTPUT $SCR/$JOB.F06

   setenv   PUNCH $SCR/$JOB.F07

   unset echo

else

   set GDDIjob=false

endif

if ((-e $PUNCH) || (-e $MAKEFP) || (-e $TRAJECT) || (-e $RESTART) ) then

   echo "Please save, rename, or erase these files from a previous run:"

   echo "     $PUNCH,"

   echo "     $TRAJECT,"

   echo "     $RESTART, and/or"

   echo "     $MAKEFP,"

   echo "and then resubmit this computation."

   exit 4

endif

#----------------------------------------------------------------------

if ($TARGET == mpi) then

 set NCPUS=`grep cores /proc/cpuinfo | wc -l`

 echo " * CPU cores:   $NCPUS "

 echo " * GPU devices: not used"

 echo " "

 setenv TRAJECT $USERSCR/$JOB.trj

 setenv RESTART $USERSCR/$JOB.rst

 setenv INPUT $SCR/$JOB.F05

 setenv PUNCH $USERSCR/$JOB.dat

 if ( -e $TRAJECT ) rm $TRAJECT

 if ( -e  $PUNCH ) rm $PUNCH

 if ( -e  $RESTART ) rm $RESTART

 setenv LD_LIBRARY_PATH /opt/intel/impi/4.0.2.003/intel64/lib:$LD_LIBRARY_PATH

 set path= ( /opt/intel/impi/4.0.2.003/intel64/bin $path )

mpdboot

 mpiexec -n $NCPUS $GMSPATH/gamess.$VERNO.x

mpdallexit

 cp $PUNCH .

endif

#----------------------------------------------------------------------

if ($TARGET == ga) then

 set PPN=1

 set NCPUS=1

 @ NPROCS = $NCPUS

 setenv HOSTFILE $SCR/$JOB.nodes.mpd

 if (-e $HOSTFILE) rm $HOSTFILE

 touch $HOSTFILE

  

 echo `hostname` >> $HOSTFILE

 set NNODES=1

 setenv PROCFILE $SCR/$JOB.processes.mpd

 if (-e $PROCFILE) rm $PROCFILE

 touch $PROCFILE

 echo "-n $NPROCS -host `hostname` $GMSPATH/gamess.$VERNO.x" >> $PROCFILE

 set path=(/opt/intel/impi/4.0.2.003/intel64/bin $path)

   setenv I_MPI_WAIT_MODE enable

   setenv I_MPI_PIN disable

   setenv I_MPI_DEBUG 0

   setenv I_MPI_STATS 0

   setenv I_MPI_DEVICE sock

   setenv I_MPI_NETMASK ib0

 setenv LD_LIBRARY_PATH /opt/intel/impi/4.0.2.003/intel64/lib:$LD_LIBRARY_PATH

 setenv LD_LIBRARY_PATH /opt/intel/composerxe-2011.4.191/compiler/lib/intel64:$LD_LIBRARY_PATH

 setenv LD_LIBRARY_PATH /opt/intel/composer_xe_2013.0.079/mkl/lib/intel64:$LD_LIBRARY_PATH

 setenv MKL_SERIAL YES

 setenv MKL_NUM_THREADS 1

 setenv LD_LIBRARY_PATH /usr/local/cuda/lib64:$LD_LIBRARY_PATH

 setenv GMS_CCHEM '1'

 @ NUMGPU=1

 setenv CCHEM 'devices=0;memory=4g'

 echo " * CPU cores:   `grep cores /proc/cpuinfo | wc -l` "

 echo " * GPU devices: $NUMGPU (with settings: $CCHEM)"

 echo " "

 chdir $SCR

 set echo

 mpdboot --rsh=ssh -n $NNODES -f $HOSTFILE

  mpiexec -configfile $PROCFILE < /dev/null

 mpdallexit

 unset echo

 rm -f $PROCFILE

endif

#----------------------------------------------------------------------

echo ----- accounting info -----

if ($GDDIjob == true) cp $SCR/$JOB.F07 ~/scr/$JOB.dat

echo Files used on the master node $master were:

ls -lF $SCR/$JOB.*

rm -f  $SCR/$JOB.F*

if (-e $SCR/$JOB.V84)        mv $SCR/$JOB.V84     $USERSCR

if (-e $SCR/$JOB.V80)        rm -f $SCR/$JOB.V*

if (-e $SCR/$JOB.TEMP02)     rm -f $SCR/$JOB.TEMP*

if (-e $SCR/$JOB.orb)        mv $SCR/$JOB.orb     $USERSCR

if (-e $SCR/$JOB.vec)        mv $SCR/$JOB.vec     $USERSCR

if (-e $SCR/$JOB.mol)        mv $SCR/$JOB.mol     $USERSCR

if (-e $SCR/$JOB.molf)       mv $SCR/$JOB.molf    $USERSCR

if (-e $SCR/$JOB.mkl)        mv $SCR/$JOB.mkl     $USERSCR

if (-e $SCR/$JOB.xyz)        mv $SCR/$JOB.xyz     $USERSCR

ls $SCR/${JOB}-*.cube > $SCR/${JOB}.lis

if (! -z $SCR/${JOB}.lis) mv $SCR/${JOB}*.cube $USERSCR

rm -f $SCR/${JOB}.lis

ls $SCR/${JOB}-*.grd > $SCR/${JOB}.lis

if (! -z $SCR/${JOB}.lis) mv $SCR/${JOB}*.grd $USERSCR

rm -f $SCR/${JOB}.lis

ls $SCR/${JOB}-*.csv > $SCR/${JOB}.lis

if (! -z $SCR/${JOB}.lis) mv $SCR/${JOB}*.csv $USERSCR

rm -f $SCR/${JOB}.lis

if ($TARGET == mpi) then

   set nnodes=`wc -l $HOSTFILE`

   set nnodes=$nnodes[1]

   @ n=1

   set master=`hostname`

   set master=$master:r

   while ($n <= $nnodes)

      set host=`sed -n -e "$n p" $HOSTFILE`

      set host=$host[1]

      if ($host != $master) then

         echo Files used on node $host were:

         ssh $host -l $USER "ls -l $SCR/$JOB.*"

         ssh $host -l $USER "rm -f $SCR/$JOB.*"

      endif

      @ n++

   end

   rm -f $HOSTFILE

   if ($?I_MPI_STATS) then

      if ($I_MPI_STATS > 0) mv $SCR/stats.txt ~/$JOB.$NCPUS.stats

   endif

endif

date

time

exit