User Tools

Site Tools


Zdenka Chromcova 2010/03/18 11:23


Switching the classic MD on/off

The using of classic potential could be switched on/off in file :

  • the implicit settings is off: iclassMD is not defined in or iclassicMD = 0
  • classic MD on: iclassicMD = 1

Settings of parameters for classic MD

in directory Cdata are two types of files:

  1. Cdata/ (define parameters of simulation)
  2. Cdata/NAME-OF-POTENTIAL-DEFINED-IN-FILE-ABOVE.dat (defines the parameters of potential)

In this file is defined which potential will be used and other parameters used in classical MD only.

  • Defines which potential will be used

Format of file

  all Lennard-Jones

or there could be defined potential for each pair: (not useful at this moment):

  H-H Lennard-Jones
  Si Lennard-Jones
  Si-H Lennard-Jones

Note: all pair potentials have to be of the same type. Thus this second type of definition of potential is not very usefull at this moment.

In the directory Cdata have to exist file with parameters of potential with name NAME-OF-POTENTIAL.dat (in case above the program will look for file Lennard-Jones.dat).

  • Defines frequency of outputs in MD simulation.

The implicit value is each 1000th cycle. The example of the definition:

  freq_of_outputs 100

Files with potentials

The parameters of potential is defined in files with name: Cdata/NAME-OF-POTENTIAL.dat . The string NAME-OF-POTENTIAL have to match with string which define potential in file .

In this file are parameters for pair potential for different types of atoms in shape: Atom1-Atom2 params = [ params separate by comas ] . Comments starts with character “#”. The example for for Lennard-Jones potential is below:


  #                   epsilon         ro          cutoff  
  # shape of potential( 4epsilon( (ro/r)**12 - 2(ro/r)**6) cutoff = 2.5 ro/2**(1/6)
  H-H     params = [      2   ,       2.4     ,   9.35 ]
  Si-Si   params = [      2   ,       1.4     ,   9.35 ]
  Si-H    params = [      2   ,       1.9     ,   9.35 ]

In these files could be defined more pairs than will be used in simulations. If the program do not find an pair defined in basisfile it ends with error message which pair is missing.

Supported potentials

At this moment there is only Lennard-Jones potential. There going to be addaed also RGL, EAM, MEAM and Tersoff potential.

List of potentials:


The input files for a nanocluster of Pd13 is here: pd.tar.gz. The initial configuration (in.bas file) is in global minima according a program for MD of Ficcardo Ferrando using the same shape of RGL potential with the same parameters. The value of total energy in the global minima is -42.09318eV (according Riccardo's program). The shape of the cluster in GM is icosahedral.

New files and changes in Fireball:

New files

Because of classic MD there were added these files:

  • DASSEMBLERS/getforces_classic.f90
  • READFILES/readdata_classicMD.f90
  • MODULES/classicMD.f90
  • MODULES/classicMEAMforces.f90 routines used in empirical potentials
  • NEIGHBOR/find_neigh_max_class.f90
  • INITIALIZERS/init_getforces.f90 assignment of the pointer on the function which give the force (one of dft or the classic)
  • DASSEMBLERS/getforces_classic_RGL.f90
  • DASSEMBLERS/getforces_classic_Tersoff.f90 (not yet in use)

Details of files:


Here is defined the main procedure for classic_forces which call the relevant procedure for given potential. The “relevant procedures” will be defined in the individual files in future.


Definition of global variable “*Potential”, user-defined types and interface for subroutines and functions.

NEIGHBORS/neighbors.f90, MODULES/neighbor_map.f90, ALLOCATE/allocate_neigh.f90, NEIGHBOR/find_neigh_max_class.f90, INITIALIZERS/initbasics.f90

if iclassicMD = 1 there is allocate and initialized also new arrays:

  !CHROM neighbor for classical MD simulation (MODULES/neighbor_map.f90)
  integer :: neigh_max_class = 0 !size of arrays
  integer, dimension(:,:),   allocatable     ::  neigh_classic
  integer, dimension(:),     allocatable     ::  neighn_classic
  integer, dimension(:,:),   allocatable     ::  neigh_b_classic

The arrays of quantum neighbors are kept unchanged to should be used in future eventually (i.e. there should be implemented the combination of classic and quantum simulation)

Changes in existing files

The code which was added or changed because of classic MD (by Zdenka Chromcova) is between tags:

  • begin: !CHROM
  • end: !END CHROM

The list of files which were changed:

  • Makefile
  • LOOPS/main_loop_MD.f90 ( defined loop for classic MD )
  • INITIALIZERS/initbasics.f90 (if iclassicMD = 1: added filling of classical neighbors)
  • ALLOCATE/allocate_neigh.f90 (if iclassicMD = 1:allocation of classical neighbors)
  • MODULES/neighbor_map.f90 (added definition of classical neighbors)
  • MODULES/options.f90
  • NEIGHBORS/neighbors.f90 (added actualization of classicl neighbor if iclassicMD = 1. Note: if iclassicMD>1 both classic and quantum arrays for neighbors are actualized)
  • READFILES/readparam.f90
  • READFILES/readdata.f90
  • MD/move-ions.f90 ( A) added variable vmax (maximum if the velocity) and in case that vmax*dt>2 && iclassicMD>0 it break the program. B) if iclassicMD>0 the outputs are with period given in C) some do-loops was joined together (because of computer time))

i.e. this (2x3xN loops)

vatom(:,1:natoms) = vatom(:,1:natoms)*vscale 
xdot(1,:,1:natoms) = xdot(1,:,1:natoms)*vscale 

was rewritten as (Nx3 loops):

do iatom = 1, natoms 
  do k=1,3 
    vatom(k,iatom) = vatom(k,iatom)*vscale 
    xdot(1,k,iatom) = xdot(1,k,iatom)*vscale 
  • DASSEMBLERS/getforces.f90 Not more used. Now getforces is not the function, but the pointer on the function (defined in MODULES/forces.f90).
  • MODULE/forces.f90 (added row: “procedure (), pointer :: getforces” )

—- Download (a bit old verion: January 2010)

ToDo list (zdenka)

  • output on the screen and into the files is not done in each step (some intervals if iclassicMD==1, definition in
  • pointer on subroutines which calculates the forces (not the case statement)
  • if velocity*dt>2 then it write out the warrning that dt is too small (in case iclassicMD>0 only)
  • made better structure of type(POTENTIALTYPE) (One type of potential for all atoms - delete the components ptype and type in the definition of type(POTENTIALTYPE) )
  • add RGL potencial (the total energy of cluster is the same as the code of Riccardo Ferrando gives)
  • verify the biatomical clusters! (comape with Riccardo's code)
  • add Tersoff potential
  • add EAM potencial (Talad Rahman, Tapio A. Nissila)
  • add MEAM potencial
  • global optimization/Metadynamics
  • freq_of_outputs - to write out something on the end of simulation (if freq_of_outputs > num of cycles in the simulation)

Optimazlization of the code

  • controll that the do-loops are written economically (ie: ratom=a*ratom;veloc=b*veloc; ⇒ do i=1,natoms; do j=1,3; ratom(j,i)=a*ratom(j,i),veloc(j,i)=a*veloc(j,i);enddo; enddo;)
  • compare the time for 1 cycle in Fireball and Riccardo's code
classical_md.txt · Last modified: 2011/02/18 13:13 (external edit)