21.08.2026 Version 4.9.2
^^^^^^^^^^^^^^^^^^^^^^^^
1. Change of the format of [currents.dat] in option EA
   [file: density.f90]
2. Improved the reading of [current.dat] file for the EA option so that
   the order of grid points does not matter anymore.
   [file: density.f90]

06.11.2025 Version 4.9
^^^^^^^^^^^^^^^^^^^^^^^
1. A bug in interpolate_EA found and fixed.
   [file: device.f90]
2. Changed x11 to q1 in plotting since gnuplot does not longer support x11.
   [file: plotting.f90]
3. Used abs(grid1(i1,i2,i3)) in density.f90 to read the current in case it is negative.
   [file: density.f90]

14.09.2023 Version 4.81
^^^^^^^^^^^^^^^^^^^^^^^
1. VASP changed their output , so some changes needed to be implemented to read OUTCAR.
   [file: do_param.f90]

07.04.2023 Version 4.80
^^^^^^^^^^^^^^^^^^^^^^^
1. A bug corrected in Ci and Cs options of the EA routine related to an incorrect value of the current at the boundaries
   of the region within the plane. Might be useful in the non-periodic case (iperiodic=0).
   [file: stm_EA.f90, device.f90]
2. A possibility implemented that enanbles one, when the lattice is square, to extend the image across the axes,
   i.e. the image is treated as a upper-right quarter of a rectangular and extended to the whole of it (option EI
   in the EA model).
   [file: stm_EA.f90]

07.04.2023 Version 4.70
^^^^^^^^^^^^^^^^^^^^^^^
1. A very significant rewrite of the EA routine: 1D, 2D and 3D cases; a direct interpolation along Z.
   [file: stm_EA.f90, density.f90, bastr.f90, device.f90]

05.04.2023 Version 4.60
^^^^^^^^^^^^^^^^^^^^^^^
1. NEW: current within the EA model can be calcualted either within the WBA for the tip (current) or with an explicit
   gamma(w) of the tip (tip_current).
   [files: directory EA_model]

01.04.2023 Version 4.54
^^^^^^^^^^^^^^^^^^^^^^^
1. The directory EA_model modified, specifically some of the scripts changed and the README file updated.

29.03.2023 Version 4.53
^^^^^^^^^^^^^^^^^^^^^^^
1. NEW: constant height mode is implemented in STM_EA. Also, a bug is fixed related to a shift when going along Z.
   Exponential interpolation is arranged as a separate routine here.
   [file: stm_EA.f90]

27.03.2023 Version 4.52
^^^^^^^^^^^^^^^^^^^^^^^
1. Optioon EA debugged and seemded to be fuly working, many bugs fixed.
   [files: device.f90, stm_EA.f90]

The suggested procedure for calculating the STM image (constant current mode) within the EA model is given below:

>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> STM CALCULATION INSTUCTIONS  (NOT CHIMERA) <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<

This is how the files for the STM calculation within the Efective Atom model need to be prepared:

(A) prepare input files with different positions of the Effective Atom (local platform)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

1. start from the <input.inp> file that contains He1 species above the surface somewhere

2. create a copy of that file <1.inp> in which the poisions x,y,z of the He1 are replaced by XXX YYY ZZZ, e.g.

       He1   0.0 0.0 3.5

is replaced by

       He1   XXX YYY ZZZ

3. Run tetr and read <input.inp>, go to EA menu option and run it. There:

     - define the box within which the grid is going to be defined; this is done by defining:
         Cp - the corner point of the box
         Up - the upper Z height of the box (Zupper)
         Vp - the vectors along the two in-plane sides of the box, both starting from Cp
              The 3rd vector of the box is defined perpendicular to the plane of the surface and between Cp(3) and Zupper

     - define the grid within the box
     - using S will generates all directories with files <xyz> contanining coordinates of the He1 species. It will also place
       there an input file as well with the correct geometry, however, in the case of CP2K this file needs editing and hence
       cannot eb directly used.
     - tetr will also create the file <for_lev00.dat> with the information on the grid chosen and the box to be used by <lev00>
       when previewing the image.

4. Edit the script <run_prep_input> to include all grid points starting from 0.

5. Run the script. It will place in each directory N1.N2.N3 of the grid the correct CP2K input file <input.inp> that has the
   required position of the He1 species. 

(B) Run the CP2K simulations and prepare the necessary file for the current simulation (HPC platform)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

6. tar/gzip the directories and scp the archive to an HPC.

7. Run all simulations using <job_run_cp2k> script, e.g.

#!/bin/bash -l
#$ -P Gold
#$ -A KCL_Kantorovitch
#$ -l h_rt=12:00:00
#$ -N cp2k_test
#$ -pe mpi 20
#$ -cwd

module purge
module load gerun
module load gcc-libs
module load compilers/gnu/4.9.2
module load mpi/openmpi/3.1.4/gnu-4.9.2
module load openblas/0.3.7-openmp/gnu-4.9.2
module load cp2k/7.1/ompi/gnu-4.9.2

export OMP_NUM_THREADS=1


for d1 in 0 1 2 3 4
do
  for d2 in 0 1 2 3 4
  do
    for d3 in 0 1 2 3 4
    do
        cd $d1.$d2.$d3

#..... run CP2K
    
        gerun cp2k.popt -inp input.inp > output.out  

#..... collect all file names of the KS Hamiltonian

        ls Gold-KS_MAT-1_0_*Log > names.dat

#...... determine the KS Hamiltonian file from the last iteration

        ~/bin/lastname > 0
        echo 'last name found is ' `cat 0`

#...... move it to a fixed name; remove the rest of the files

        mv `cat 0` Gold-KS_MAT_last.Log
        rm Gold-KS_MAT-1*.Log
        rm *bak*
        rm cp2k*

        cd ..

    done
  done
done


This script uses a F90 code <lastname> (the source - <names.f90>) that picks up the KS Hamiltonian file from the last iteration
- in the case if the print out is made every iteration. If only the file from the last iteration is printed by CP2K, the script
requires a simple modification (to move that file into the fixed name <Gold-KS_MAT_last.Log>).

8. Extract the KS Hamiltonian file and the eigenvectors/eigenvalues by running the script <run_python_script>. It relies on a few
   python scripts:    <ks_ham.cp2k.py>,  <wf_coeff_cp2k.py> and <wfn_gs.py>.
   This script creates in each directory two files: <tot_ks.dat> with the full Hamiltonian matrix (we if fact only need its 1st row
   or column, so this step needs to be modified) and <eigenvalues_eigenvectors.dat> that contains the LCAO coefficients and KS
   energies.

9. Run <>run_extract_clean> scripts that creates a temporary directory <tmp> with all information needed to calculate the current
   (hopefully, this directory is smaller in size than the directory with all CP2K runs).

10. Go into <tmp> and tar/gzip everything. S-copy it to the local computer platform to the same directory where all the input files
    were prepared (see A). Remove <tmp> entirey. Clean the directory with all the CP2K files if necessary.

11. Run a separate calculation without the He atom, just the surface, and extract the value of the chemical potential under: 'Fermi energy',
    e.g.

    grep Fermi output.out 
    Fermi energy:                                                -0.12405316266654

(C) Claculation of the current (local platform)
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^

12. Un-tar the archive file. You should have in each subdirectory N1.N2.N3 of the grid: <input.inp> of CP2K, as well as the KS Hamiltonian
    matrix and the eigenvalues/vectors in the files  <tot_ks.dat> and <eigenvalues_eigenvectors.dat>, respectively.

13. Run the script <run_current> that would create the file <current.dat> with all the values of the current from each of the directory,
    something like this:

# plot of the current vs z of the He atom 
0 0 0 6.318714964813027E-002
0 0 1 0.413616560810557
0 0 2 4.233095158921797E-003
0 0 3 3.688743327850073E-006
0 0 4 2.338825239116674E-010
0 1 0 6.306286248128982E-002
.................

    where in the first three columns the frid is shown, the current - in the 4-th column. The input parameters
    for the calculations (temperature, Gamma, EA on-site energy and the Fermi energy of the surface) are in <input.dat> file,
    while the bias voltage value is given explcitily in the script <run_current> in the command:

 curr tot_ks.dat eigenvalues_eignevectors.dat 1.0 > ../tmp

    (1.0 Volt in this example).

14. Run lev00:
     - choose C (CP2K),
     - give the name of the CP2K input file,
     - skip reading the CP2K output file
     - skip reading the CP2K cube file (the electronic density)
     - choose D option (twice)
     - choose  EA option of the density menu, where: use I to set up the value of the current and Ci to calculate the image, Vi to preview it.
       Spectroscopy: using XY, Cs, Vs. All other options are self-explnatory (e.g., Ex for the exponential fit).

>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> E N D <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<

23.03.2023 Version 4.51
^^^^^^^^^^^^^^^^^^^^^^
1. The option EA almost fuly rewritten. Two options are suypported: (i) the region where the current was calcualted
   can be periodically continued within the XY plane (e.g. the calcualtion is done of a perfect surface in which case the
   CP2K unit cell is much larger as the surface primitive cell to avoid interact5ion between He atoms), and (ii) cannot. 
   [files: density.f90, stm_TH.f90, device.f90, code.f90, read_density_cp2k.f90,stm_EA.f90]

03.03.2023 Version 4.5
^^^^^^^^^^^^^^^^^^^^^^
1. NEW: an infrastructure for the Effective Atom (EA) model has been implemented (probably still buggy).
   A current file is read in by choosing CP2K. Now there are two options: to read the (spin-)density or the current file.
   The current file has a format containing: a new (less dense) grid (NGX_,NGY_,NGZ_), the lattice vectros (must be identical
   to the system lattice vectors), and the values of the current on that grid. Then, the same routine as stm_TH is run
   to calculate the constant current image using the new grid.
   [files: density.f90, stm_TH.f90, device.f90, code.f90, read_density_cp2k.f90]

16.11.2022 Version 4.41
^^^^^^^^^^^^^^^^^^^^^^
1. Units in which the electrostatic potential from VASP is viewed along the line / plane (using LOCPOT) is changed to eV
   (which is the original LOCPOT unit). Previously it was divided by the unit cell volume VOLC.
   [files: density.f90, plot_add.f90]

10.10.2022 Version 4.4
^^^^^^^^^^^^^^^^^^^^^^
1. A very significant rewrite of the stm_TH routine to improve its stability; many bugs corrected and performce improved.
   In particular:
      a. Interpolation routine substantially improved.
      b. Option Ez removed.
      c. BT option helped by showing the actual z window
      d. Spectrosopy option: interpolatoin (line) is shown together with the raw data (points)
      e. Exponential interpolation can be done from any value of the DOS (option Ds)
      f. Items in the menu reshuffled
      g. Number of points in spectroscopy is set to ngrid(3), cannot be changed (option Np revmoed)
      h. Calculation of zz1 correspoding to the given DOS value is forced withitn the Bot-Top window
   [file: stm_TH.f90]
2. A bug in density.f90 corrected
   [file: density.f90]

>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>> STM CALCULATION INSTUCTIONS  (NOT CHIMERA) <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<

I. Instructions on using lev00 for the calcualtion of the STM image 
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ 
The steps below are advised for the calculation of the constant LDOS STM image:

1. Crereate a grid file corresponding to the desired voltage (PARCHG for VASP)
2. Run lev00 (V -> OUTCAR -> POSCAR -> D -> P)
3. Go to TH menu option.
4. In there:

       1 - specify the X,Y centre of your image (need to give Z as well, not used)
       2 - specify the X,Y sizes of the image (in A)
       4 - specify the grid (e.g., 100 100 100)
      BT - choose the Z window for the image; this should be in the vacuum gap between 
           the top of the atoms and up to the middle of the gap; this info is provided
       I - choose the value of the LDOS for which the image (x,y,z) is yo be calculated
           within the specified interval of LDOS values
      Ci - calculate the image
      Vi - preview the image; use 3D to select the desired mode ("3D; in colour" is a good one)
      Fz - udjust the z-window for the image and rerun Vi
      Ps - create a postscript file of the image

5. If larger z values are desired than available (up to the middle of the gap) or in order to
   alleviate the noise at small LDOS values (larger z), use exponential interpolation:
   
      Repeat the first 3 steps from item 4
      Ex - enable the fitting
      Ds - choose the value of the LDOS below which the interpolation is to be used;
           for LDOS values larger than this value the VASP data will be used as without
	   interpolation
      Run Ci, Vi, Ps as before.

6. Spectroscopy at up to 8 points is possible; this is useful also to check the character of the LDOS(z)
   for particular x,y points and also for checking the fitting procedure (choosing the values for BT)

      Repeat the first 3 steps from item 4 (the z grid value givess the number of points along z)
      Nl - choose the number of lateral points (up to 8)
      XY - specify positions of all points
      Sl - preview the points if desired
      zW - z-window for showing atoms
      Cs - perform the calculation
      Vs - preview curves, both the actual data (dots) and the interplation with the exponents (lines)

7. There are also additional options:

      Co - preview all atoms
      An - specify positions in A, fractional coordinates, atomic numbers
      
>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<

II. Instructions on using lev00 for the calcualtion of the STM dI/dV image 
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^ 
The steps below are advised for the calculation of the dI/dV image following the tip trajectory at the same
driving voltage at which the derivative is calculated (constant LDOS) [VASP only]

1. Create two PARCHG and PARCH1 files, for the given voltage and for a slightly different one
2. Run lev00 (V -> OUTCAR -> POSCAR -> D -> PP); note that PP is to be used
3. Go to T1 menu option

Operation there is practically identical to that of the TH men. There is one additional option:
   dV - voltage difference

>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<

III. Instructions on using lev00 for the calcualtion of the STM dI/dV image 
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^  
The steps below are advised for the calculation of the dI/dV image following the tip trajectory at one
driving voltage (constant LDOS), while the derivative is calculated at another [VASP only].

1. Create PARCHG style file with VASP and link it to the name PARCHG  - to record the tip trajectory. 
2. Follow all steps from I; files ldos.dat and ldos_settings.dat are created containing the trajectory.
   Rename PARCHG to keep it.
3. Then, you create two other PARCHG style files corresponding to two different close voltages, call one again
   PARCHG and another PARCH1 (use logical links).
4. Rerun lev00 (V -> OUTCAR -> POSCAR -> D -> PP); note that PP is to be used
   as well as PARCHG and PARCH1 files, for the given voltage and for a slightly different one
5. Run lev00 (V -> OUTCAR -> POSCAR -> D -> PP)
6. Go to T2 menu option

Operation there is practically identical to that of the TH menu although not options are available.
There is one additional option:
   dV - voltage difference

>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>><<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<


07.10.2022 Version 4.33
^^^^^^^^^^^^^^^^^^^^^^
1. A new option Ez in the TH option (TH STM) that enables calculation of the STM image using only an exponential fit,
   for all z values.
   [file: stm_TH.f90]

22.09.2021 Version 4.32
^^^^^^^^^^^^^^^^^^^^^^
1. A bug corrected.
   [file: mainmenu.f90]
2. When entering the Ge menu, all atoms are chosen as a default.
   [file: mainmenu.f90]   
3. Ge menu: added option Aa to choose whether in option Cf to show all atoms with all their forces
   or only atom with foerces larger than the minimum one. Also, the name of the VASP output file is not asked
   anymore, it is taken from the previously entered.
   [file: mainmenu.f90]   

22.09.2021 Version 4.3
^^^^^^^^^^^^^^^^^^^^^^
1. NEW: dI/dV image can be executed now within the constant current mode (VASP)
   with the driving voltage being different from the one at which the differentiation
   wrt to the voltage is made. 
   [files: stm_TH.f90]

The procedure:

   (a) run VASP to create PARCHG corresponding to the voltage used to drvie the tip (tip trajectory voltage)
   (b) run V -> D -> P -> TH. Files ldos.dat and  and ldso_settings.dat are written.
   (c) run VASP again to create two PARCHG files for two values of the voltage (close enough, for differentiation);
   (d) set the first one as PARCHG and the second as PARCH1.
   (e) run V -> D -> PP -> T2. The image (via Ci and Vi) will be plotted corresponding to dI/dV.


21.09.2021 Version 4.2
^^^^^^^^^^^^^^^^^^^^^^
1. NEW: dI/dV image can be executed now within the constant current mode (VASP)
   [files: stm_TH.f90, density.f90, read_density_vasp.f90]

The procedure:

  (a) starting from the relaxed geometry, perform a VASP calcualtion to have WAVECAR at an
      extended set of k-points:

ISTART=1
ICHARGE=11

      Keep OUTCAR and of course the geometry file.

  (b) generate the partial charge dnesity between two values of the voltage V (between the Fermi energy
      and the energy either above [if VV>0] or below [if V<0] it):

ISTART=1
ICHARGE=11
LPARD = .TRUE.
MBMOD = -2               !  MBMOD=-3
EINT = -0.6312 0.6688    !  EINT = -0.3 0.0 (wrt Fermi energy)

      A PARCHG file is to be produced.

This voltage is assumed to be the one that was used to determine the tip height.

  (c) Repeat (b) for another value of the voltage close to the first one (for numerical differentiation).
      The produced PARCHG file rename as PARCH1. It is assumed here that THE SAME tip position
      was used in the experiment for calculating dI/dV. In other words, PARCHG from step (b) is used
      for determining the tip height and for the main LDOS value, while at step (c) the height is the same,
      but the voltage is different.

  (d) In lev00: V -> D -> PP -> T1. You can use Ex (exponential fit) with this option.

9.09.2021 Version 4.1
^^^^^^^^^^^^^^^^^^^^^^
1. TH option (STM): a new exponential fit routine
   [files: stm_TH.f90]

9.09.2021 Version 4.02
^^^^^^^^^^^^^^^^^^^^^^
1. A bug in a format statement corrected when printing a selceted geometry within the Ge menu
   [file: mainmenu.f90]
2. A format of reading VASP eigenvalues from OUTCAR adjusted to that used currently in VASP.
   [file: do_param.f90]

5.06.2017 Version 4.01
^^^^^^^^^^^^^^^^^^^^^^^
1. NEW: atomic numbers can be speciefied across the bred system (e.g. in choosing direction
   of plotting).
   [files: plot_add.f90, device.f90, density.f90, stm_TH.f90]

10.05.2017 Version 4.00
^^^^^^^^^^^^^^^^^^^^^^^
1. NEW: LEV00 is bound to work seemlessly with CHIMERA. Initially, a geom.xyz file
   is written, so that CHIMERA can show the system on its screen. Then, in the density options
   pL, pP, CS, cA, gO and TH the menu items B0 and Bb are added so that it is possible to show
   more than a single unit cell. Every time a change of the cell box is made, geom.xyz is updated
   automatically (no W item is needed as compared to TETR).  
   [files: lev00.f90, param.f90, device.f90, plot_add.f90, density.f90, stm_TH.f90]
2. NEW:`when CHIMERA=True, then in the gO option, when a box is chosen, all atoms in the box are
   highlighted in the CHIMEERA window. Note that only those atoms are shown which 
   are within the Bb extension box. If not all atoms of the chosen box are shown, choose a bigger 
   Bb extension. Next, when W option of gO is hit, the electron density is shown in the CHIMERA
   window together with the atoms in the chosen box. In this case Bb box is ignored. Atoms within
   the Bb option can be shown again by toggling the corresponding option on the lev00 CHIMERA window.
   [files: density.f90, device.f90]
3. NEW: an add-on optoin in Ge for VASP to check forces on atoms. One can set up a minimum value of the
   force, and then all atoms with the forces on them larger than this value will be printed on the screen.
   Atoms which are frozen (at least one of the X,Y,Z components) are ignored.
   [file: mainmenu.f90]
4. NEW: an add-on option to Ge: with Cm one can compare two geometries along the relaxation path.
   [file: mainmenu.f90]
5. NEW: an add-on for Ge: plot the max force vs. iteration number.
   [file: mainmenu.f90]
6. NEW: add-on options in Ge: one can indicate geometries numbers between which all plots are made (energy, force,
   atomic positions).
   [file: mainmenu.f90, plotting.f90]

30.04.17 Version 3.51
^^^^^^^^^^^^^^^^^^^^^^
1. NEW: add-on to option Ge of VASP: option PE allows previewing all energies along te relaxation/MS process, while S
   enables one to save VASP_5.2 styl;e inpnut file from any iteration presented in OUTCAR
   [file: tools_strings.f90, mainmenu.f90]

09.04.17 Version 3.50
^^^^^^^^^^^^^^^^^^^^^^
1. A bug in stm_TH fixed: the number of simultaneous 1D plots cannot be more than 8 (it was 9).
   [file: stm_TH.f90]

02.03.16 Version 3.49
^^^^^^^^^^^^^^^^^^^^^^
1. NEW: a new option for VASP enabling a simple analysis of MD runs and/or geometry relaxations:
   - one can create a single XYZ file geometries.xyz containing all geometries one after another:
   - one can choose atoms of interest (by tagging), and then display (using gnuplot) their trajectories, separately
     for X, Y and Z components.
   [files: mainmenu.f90]

27.05.15 Version 3.48
^^^^^^^^^^^^^^^^^^^^^^
1. A bug [wrong length of a character in the call of plot_bunch()] fixed.
   [file: plotting.f90]

07.10.14 Version 3.47
^^^^^^^^^^^^^^^^^^^^^^
1. A bug in closeing unit 15 was fixed, related to reading vasp input file: in case some information is missing, the unit
   was not clsed.
   [file: read_vasp_geom.f90]

07.10.14 Version 3.46
^^^^^^^^^^^^^^^^^^^^^^
1. NEW: some improvements to the STM routine, namely, a new option Ex is introduced which enables an exponential
   fit for the LDOS at each lateral grid point. This is used for the value of the current (=LDOS) which is smaller
   then the available lower bound, i.e. some image can still be constructed in this case giving larger Z on the 3rd axis.
   [file: stm_TH.f90] 

25.08.14 Version 3.45
^^^^^^^^^^^^^^^^^^^^^^
1. For VASP: before, the option PD was going ahead even if the file PROCAR was not found. A warning was
   given but the code erroneously proceeded to PDOS regardless. Now fixed.
   [file: mainmenu.f90]
2. For VASP: in its new version the PROCAR file for the PDOS contains an extra F component
   as well. Although it is presently ignored (so that the total PDOS, TOT, is larger than S+P+D), 
   the PROCAR file is now read in correctly, so that PDOS can be displayed (without the F 
   component though - to be rectified in the future).
   [file: read_vasp+psi2.f90]
3. English.
   [file" choose_tasks.f90] 
4. For plotting of the PDOS: format of band numbers enlarged from i3 to i5. Plot() and Plot1() 
   are more general and accept any length of title_pl(), may be useful for the future.
   [files: dos_add.f90, prep_dos.f90,plotting.f90]

28.01.14 Version 3.44
^^^^^^^^^^^^^^^^^^^^^^
1. NEW: merger of density files is possible using the manipulation menu. In this case
   the coefficients c1=c2=1.0 and atomic positions are merged together from both density files.
   It is assumed (i.e. not checked) that in both density files atoms are ordered by species.
   Also, it is not checked if there are equivalent atoms in both sets. So far, this has only been
   done for CP2K.
   [files: code.f90, manip_dens.f90, read_density_cp2k.f90]
2. Bugs related to allocation/deallocation fixed.
   [file: read_density_cp2k.f90]
3. Atoms are shifted together with the density when option Sf is applied
   [file: density.f90]

14.01.14 Version 3.43
^^^^^^^^^^^^^^^^^^^^^^
1, NEW: writing general (cube) density files in the density menu has been added for CP2K.
   [files: density.f90, write_dens.f90]
2. Density for a shifted system: a help prompt changed to that for the 'shifted system'
   [file: density.f90]
3. When choosing the input cube file during initialisation (CP2K), all files with 
extension .cube will be shown.
   [file: get_param_cp2k.f90]

03.12.13 Version 3.42
^^^^^^^^^^^^^^^^^^^^^^
1. Some improvements in reading CP2K input files
   [file: read_cp2k_geometry.f90]

15.10.13 Version 3.41
^^^^^^^^^^^^^^^^^^^^^^
1. When reading CP2K density, only a warning is given if the lattice vectors mismatch;
there is also a possiblity to continue regadless.
   [file: read_density_cp2k.f90]
2. A stupid (!) check on the number of atoms is removed from the routine related to
   manipulation of CP2K densities.
   [file: read_density_cp2k.f90]
3. Recognizes both @include and @INCLUDE commands in the CP2K input files
   [file: get_parama_cp2k.f90] 
4. A bug fixed (segmentation fault: VKPT, WTKPT - used before allocation)
   [files: lev00.f90, get_param_cp2k.f90]
5. Temporary file tmp.inp created while reading in the CP2K input file is deleted after reading.
   [[file: read_density_cp2k.f90]

1.09.13 Version 3.4
^^^^^^^^^^^^^^^^^^^^^^
1. NEW: added charge and spin densities support for CP2K
   [files: read_cp2k_geometry.f90, lev00.f90, mainmenu.f90, get_param_cp2k.f90, density.f90, code.f90,
           manip_dens.f90, read_density_cp2k.f90]
2. Reogranisation of the files related to reading the density.
   [files: read_density_siesta.f90, read_density_vasp.f90]
3. Line skipping modified to register only matching at the beginning of lines.
   [file: tools_strings.f90]

14.12.12 Version 3.35
^^^^^^^^^^^^^^^^^^^^^^
1. An unimportant bug in a TH STM menu fixed.
   [file: stm_TH.f90]
2. A warning message added under BT option in TH STM when the LDOS window cannot be found.
   [file: stm_TH.f90]

18.11.12 Version 3.34
^^^^^^^^^^^^^^^^^^^^^^
1. an issue related to reading species number in a new version of vasp fixed.
   [file: do_param.f90]
2. an issue related to reading vasp v. > 5.2 density for the manipulation options 
   resolved.
   [file: manip_dens.f90, read_density.f90]
3. a bug is plot_bunch() related to plotting several curves fixed
   [file: plotting.f90]
4. README: the history revised' the latest changes are shown at the front of the file

27.02.12 Version 3.33
^^^^^^^^^^^^^^^^^^^^^^
1. wrong titles for X and Y axes in the STM spectroscopy; now swapped, so that they are correct.
   [file: stm_TH.f90]
2. wrong title for Z axis in 3D plots in STM routine; fixed.
   [file: stm_TH.f90]
3. error trap in inputing LDOS value in the STM menu (option I) has been put in
   [file: stm_TH.f90]
4. difiver was missing in the format in the PLot3D routine, now fixed
   [file: plotting.f90]

09.02.12 Version 3.32
^^^^^^^^^^^^^^^^^^^^^^
1. Added an oprion El to the TH menu to enable/disable empty lines in the ldos.dat file for the LDOS
   [file: stm_TH.f90]
2. 'set size ratio 1' command added to test.gnu for 3D plots
   [file: plotting.f90]

07.02.12 Version 3.31
^^^^^^^^^^^^^^^^^^^^^^
1. NEW: 8 options in total are provided for the 3D plots, some of them are really nice!
   [files: menu.f90, stm_TH.f90, plotting.f90]
2. It has been specifically checked that constant LDOS image can also be created from PARCHG file
   using gOpenmol menu and the value of the current chosen beforehand with the TH menu. This way it is possible
   also to have atoms displayed as well under the image.

04.02.12 Version 3.30
^^^^^^^^^^^^^^^^^^^^^^
1. A missing command for 3D gnuplot plots was added (this one is needed for the Gnuplot 4.4)
   [file: plotting.f90]
2. Some cleaning of variables, etc.
   [files:lev_coulomb.f90, plot_add.f90, device.f90, prep_dos.f90 ]
3. NEW: Tersoff-Hamann routine (option TH in the D menu) added for working with STM images and spectroscopy. 
   It was so far only tested for VASP. It may work with SIESTA as well, but I need to check the SIESTA LDOS
   file first.
   This is how the STM image (spectroscopy) to be calculated:
      (i) run VASP normally to relax the structure and get CHGCAR; 
     (ii) run VASP again with more k-points and the following options in INCAR:

ISTART =1        [restart (0,1,2)]
PREC=HIGH
NSW = 0
ICHARG=11        [keeps the density fixed (DOS or BAND options)]
LPARD=.TRUE.     [for the local DOS]
MBMOD = -3       [in EINT the energy window is given wrt E_Fermi]
EINT = -1.0 0.0  [energy window: in this case all occ states 1 eV below the Fermi]

    (iii) The run above should create PARCHG file which lev00 can read
     (iv) run lev00 (you need OUTCAR, CONTCAR, PARCHG): 
          V -> OUTCAR -> CONTCAR -> D -> P -> TH
      (v) in the TH menu one can do either spectroscopy in up to 9 lateral points
          or a 2D/3D plot of the constant LDOS topography
    
     [files: stm_TH.f90, plotting.f90, density.f90]
4. In a number of places numerical menu items were changed into letters.
   [file: density.f90]

12.11.11 Version 3.29
^^^^^^^^^^^^^^^^^^^^^^
1. Gnuplot 4.4 (the latest at this time) has a different syntax, correspondingly the gnuplot feature did not work
   with the new version. This has been fixed in this version of lev00. Hopefully, backward compatibility is also
   ensured.
   [file: plotting.f90]

26.10.11 Version 3.28
^^^^^^^^^^^^^^^^^^^^^^
1. Vasp 5.2 (or even 4.0 ?) may print out in OUTCAR eigenvalues at every iteration which before would result in
   read in only the very first set (not yet converged). This is now fixed: it searches for the very last iteration
   and only then it read the KS eigenvalues for each k-point.
   [file: do_param.f90]
2. Option 1 in the PDOS menu ("Choose/Add a single PDOS task") was not actually working (likely, it was not finished
   when new features ST and CT were added). Now fixed: using this option one can choose a single task or append to the
   existing list (e.g. after 1 or CT). 
   [file: prep_dos.f90]
3. A bug fixed related to the tasks and going to a new PICTURE.
   [file: prep_dos.f90]
4. px,py, pz contributions for PDOS from the vasp PROCAR file were read in a wrong order, now fixed; no effect this
   had on the P projection as only the sum of them was used anyway.
   [file: read_vasp_pis2.f90]

10.09.11 Version 3.27
^^^^^^^^^^^^^^^^^^^^^^
1. Bug fixed in DOS/PDOS routines (a little mess with allocations/deallocations)
   [files: prep_dos.f90, read_vasp_psi2.f90, do_tar_improved ]

20.07.11 Version 3.26
^^^^^^^^^^^^^^^^^^^^^^
1. An inconsistency with spins fixed (jspin was made to be either 1 or 2, but never 0)
   [file: mainmenu.f90]
2. A bug in checking siesta forces (option Ck): in some cases siesta prints forces differently, and that was now taken
   into account.
   [file: mainmenu.f90]

11.07.11 Version 3.25
^^^^^^^^^^^^^^^^^^^^^^
1. Relaxation flags are read in for any siesta input and printed initially with the atomic coordinates.
   [files: atoms.f90, read_siesta_input.f90, lev00.f90, read_vasp_geom.f90]
2. NEW: an option  Ck is added in the main menu for siesta to check the last atomic forces and whether the relaxation
   successfuly completed.
   [files: mainmenu.f90]

10.07.11 Version 3.24
^^^^^^^^^^^^^^^^^^^^^^
1. Atomic number and species are displayed in the Tasks menu of PD in the case of vasp,
   i.e. whne showing all projected DOS tasks, one can see the actual atoms numbers and their
   species.
   [files: dos_inc.f90, prep_dos.f90, read_vasp_psi2.f90]
2. NEW: Instead of calculating PDOS in the PD menu for a single task (e.g. s-projection on a particular atom),
   it is now possible to choose various sets of atoms and calculate either t,s,p or d -type PDOS
   on all these atoms, i.e. the sum of the elementary PDOS. The whole menu is revamped completely. 
   This is achieved with option CT. By repeatedly invoking this option, a new subset of PDOS tasks can be added,
   i.e. a complex list of tasks can be created if needed. A single task can be chosen/added as well using the old
   option 1.
   [files: prep_dos.f90, dos_add.f90, read_vasp_psi2.f90, dos_inc.f90, prep_disp.f90, choose_tasks.f90, Makefile ]

06.07.11 Version 3.23
^^^^^^^^^^^^^^^^^^^^^^
1. Bugs in calculating total and projected DOS for siesta from the EIG and PDOS files were fixed.
   [files: siesta_eig.f90, read_siesta_pdos.f90]
2. NEW: in calculating total and projected DOS from siesta EIG and PDOS files it is now
   possible to choose spin in the case of a spin-polarised calculation
   [files: mainmenu.f90, siesta_eig.f90, read_siesta_pdos.f90]
3. In addition, all the menus related to DOS are better presented (revamped)
   [files: get_param_siesta.f90, mainmenu.f90]
4. It is indicated in the case of the tetrahedron method for the DOS if there is no spread
   of eigenvalues (energies) wrt to the k-points
   [file: dos.add.f90]

06.11.10 Vewrsion 3.22
^^^^^^^^^^^^^^^^^^^^^^
1. A bug in a format statement fixed.
   [file: density.f90]

29.05.10 Vewrsion 3.21
^^^^^^^^^^^^^^^^^^^^^^
1. differences in the format of OUTCAR, CHGCAR and PARCHG in vasp version 5.2 with the previous
   versions were addressed. 
   [files: do_param.f90, density.f90, read_density.f90]
2. Format of on-screen statements when reporting the reading from a file changed ("... done 10%", etc.)
   [files: read_density.f90, manip_dens.f90, simulate.f90, write_dens.f90]
3. New species added - the whole Mendellev's Table, all 112 of them!
   [files: read_vasp_geom.f90, mendeleev.f90]

15.03.10 Vewrsion 3.20
^^^^^^^^^^^^^^^^^^^^^^
1. A bug in calling routines find_string fixed (rewind at the beginning was missing)
   [files: get_param_siesta.f90]

10.01.10 Vewrsion 3.19
^^^^^^^^^^^^^^^^^^^^^^
1. Allocation of Z_atom was moved at the beginning of the routine density() since some compilers
   give an error after hitting 4 in the type of density menu for VASP.
   [file: density.f90]

21.12.09 Version 3.18
^^^^^^^^^^^^^^^^^^^^^
1. Reading of RWIGS from the output of newer versions of vasp
2. Working out the number of species is done differently due to a different
   pseudopotential information vasp may produce
3. Errors in reading CHGCAR0-like files of VASP are made specific to help to identify them.
   [files: do_param.f90, read_vasp_geom.f90, read_density.f90]

18.06.09 Version 3.17
^^^^^^^^^^^^^^^^^^^^^
1. A bug in reading lattice vectors from siesta density file fixed.
   [files: manip_dens.f90, read_density.f90]

28.05.09 Version 3.16
^^^^^^^^^^^^^^^^^^^^^
1. A bug (since version 3.12) fixed (atomic positions read in from SIESTA were wrongly multiplied by
   a conversion factor au_A).
   [file: read_siesta_input.f90]

22.05.09 Version 3.15
^^^^^^^^^^^^^^^^^^^^^
1. NEW: a possibility added to the manipulation menu to check density files even without
   manipulating them,.
   [files: manip_dens.f90, read_density.f90]

20.05.09 Version 3.14
^^^^^^^^^^^^^^^^^^^^^
1. NEW: manipulation option has been taken out of the density menu and moved into the very
   first menu. This allows one to manipulate densities without reading in all the geometry, etc.
   files. This is especially useful for generating a density difference file.
   [files: mainmenu.f90, code.f90, read_density.f90, manip_dens.f90]
2. seed string is initialised before it is read in from a siesta input file; strangely, not doing so
   caused a problem on some platforms.
   [file: get_param_siesta.f90]

08.05.09 Version 3.13
^^^^^^^^^^^^^^^^^^^^^
1. NEW: total DOS can be plotted from siesta [seed].EIG file as well.
   [files: mainmenu.f90, read_siesta_pdos.f90, plotting.f90, 
           read_siesta_input.f90, siesta_eig.f90, Makefile]
2. A bug fixed.
   [file: read_siesta_pdos.f90]
3. Rearrangement of siesta PDOS menus.
   [file: read_siesta_pdos.f90]

06.05.09 Version 3.12
^^^^^^^^^^^^^^^^^^^^^
1. NEW: information about the Fermi energy and HOMO/LUMO energies is displayed from the
   siesta output and EIG files in the PDOS option. This addition is not yet complete,
   as spin-polarised calculation is not yet handled.
   [files: read_siesta_input.f90, get_param_siesta.f90]
2. An option TD added to the siesta PD menu to show or not the total DOS.
   [file: read_siesta_pdos.f90, plotting.f90]
3. NEW: in the case of SIESTA, XV file is not longer needed, only fdf file is now used
   for all necessary information concerning arrays dimensions, geometry, etc.
   [files: get_param_siesta.f90, read_siesta_input.f90]

30.04.09 Version 3.11
^^^^^^^^^^^^^^^^^^^^^
1. NEW: projected DOS for SIESTA. One can read in siesta PDOS files and plot them.
   Up to 20 PDOS curves can be plotted alongside the total DOS, for each PDOS curve
   a sum over a subset of atoms is done. One can choose different subsets for each PDOS.
   TTag routine from tetr (updated and suplied with an extra feature of choosing atoms
   by species, option [Sp]) has been implemented for selecting atoms.
  [files: ttag.f90, read_siesta_pdos.f90, device.f90, mainmenu.f90, tools_strings.f90, 
          plotting.f90,Makefile]
2. some minor bugs (mainly in writing on the screen)
  [files: read_siesta_input.f90, lev00.f90]

28.11.08 Version 3.10
^^^^^^^^^^^^^^^^^^^^^
1. spin components of PARCHG and LOCPOT files (if exist) can also be read in the same manner
   as of CHGCAR (in the case of VASP)
  [files: mainmenu.f90, density.f90, read_density.f90]
  
27.11.08 Version 3.09
^^^^^^^^^^^^^^^^^^^^^
1. bug fixes in the routines working on the DOS
  [files: prep_dos.f90, prep_disp.f90, dos_add.f90]
2. To avoid a crap produced sometimes by gnuplot when it sees the double precision numbers,
   we write to the file for gnuplot when working on the DOS in a simpler format.
  [files: prep_dos.f90, prep_disp.f90, tools.f90]

01.07.08 Version 3.08
^^^^^^^^^^^^^^^^^^^^^
1. bug fixes 
  [files: menu.f90, read_density.f90]

09.06.08 Version 3.07
^^^^^^^^^^^^^^^^^^^^^
1. density mainpulation routine bug: the density is saved much too late; fixed
  [file: manip_dens.f90]
2. a check on the species when reading vasp density files for manipulation is
   handled with more freedom - only warnings are given, but one may proceed in
   the calculation even if the files have different species or their numbers;
   files are considered inconsistent if only lattice vectors and the grid are
   different; no more check is performed against the geometry file read in 
   initially when staring lev00
  [file: read_density.f90]

10.04.08 Version 3.06
^^^^^^^^^^^^^^^^^^^^^
1. a bug in the routine reading siesta density for manipulation fixed
  [file:  read_density.f90]
2. a bug in a write statement
  [file: get_param_siesta.f90]
3. a bug in writing siesta density fixed
  [file: write_dens.f90]
4. a bug in defining A_au, au_A fixed
  [file:  read_density.f90]

28.02.08 Version 3.05
^^^^^^^^^^^^^^^^^^^^^
1. when reading species names from the OUTCAR file, only the first two characters
   before the underscore character are kept; the rest is ignored - no error message
   is given anymore
   [file: do_param.f90]
2. Species La is added.
   [file: mendeleev.f90]
3. format of the charge density is changed when saved in the vasp format to improve
   transferability
   [files: read_density.f90, write_dens.f90, manip_dens.f90] 
4. format improved when writing number of grid points in calculating charges in spheres
   [file: plot_add.f90]
4.12.07 Version 3.04
^^^^^^^^^^^^^^^^^^^^^

1. a check is performed on consistency of the information in te brill.dat file with the
   number of k-points from the initialisation (NKPTS); DOS will not be calculated if 
   inconsistency is found.
   [files: prep_dos.f90, prep_disp.f90]
2. a bug in reading unformatted siesta RHO files fixed
   [file: read_density.f90]

10.11.07 Version 3.03
^^^^^^^^^^^^^^^^^^^^^
1. a less stringent procedure for reading vasp input files implemented to allow some 
   variations in the vasp output in OUTCAR; in particular, after a non-scf run, band
   energies are not given, so that lev00 proceeds without the DOS option; logical 
   variable (yesen) is intorduced to distinguish between the existence of k-points 
   (yeskp) and the KS energies.
   [files: mainmenu.f90, lev00.f90, read_vasp_geom.f90, do_param.f90, get_param_siesta.f90,
           read_siesta_input.f90, kpoints.f90]

31.10.07 Version 3.02
^^^^^^^^^^^^^^^^^^^^^
1. negative VASP densities can also be read in, density is not set up to zero if
   negative (to support density differences)
   [file: read_density.f90]

17.10.07 Version 3.01
^^^^^^^^^^^^^^^^^^^^^
1. bug fixes
2. detailed account of reading SIESTA input is given on the screen for
convenience
   [files: read_siesta_input.f90, tools_strings.f90, hat.f90, hould.f90,
           lev00.f90, mainmenu.f90, read_vasp_geom.f90, seed.f90, 
           lev_coulmb.f90, get_param_siesta.f90]

16.10.07 Version 3.0
^^^^^^^^^^^^^^^^^^^^
This is MAJOR DEVELOPMENT: it is now f90 code with memory allocation, so that there is
no need to run do_param and get_param_siesta to recompile the code for the specific
system. It is now only required a single compilation during installation.
The codes  do_param and get_param_siesta are now completely redundant together
with the lev00.comp script!

NEW features: it

1. combines do_param and get_param_siesta routines with lev00 into a single package
2. all array boundaries are worked out automatically from the VASP or SIESTA output 
   and/or input files, where appropriate, and
3. all arrays depending on the system size paramaters are allocated and deallocated
   as appropriate.
4. other f90 features are used including modules instead of common blocks.
5. CASTEP and CETEP capabilities have been removed as redundant.
6. some cleaning up of the routines.
7. some little changes of the menu options around the code.
8. all (except that in the bastr.f90 routine) 'stop' commands have been
   removed so the the run could be safely contiued even if something goes wrong
   (e.g. a file does not exist or incompatible).

At this stage I tried to keep the code unchanged whenever possible to ease its
debugging and testing. Therefore, no new major features in this version. 

   [files: all; new files added ]

2.04.2007 Version 2.46.3
^^^^^^^^^^^^^^^^^^^^^^^^^
1. a bug in bastr.f corrected (interestingly, this one was fixed a long time
   ago in tetr!)
   [file: bastr.f]

12.03.2007 Version 2.46.2
^^^^^^^^^^^^^^^^^^^^^^^^^
1. Format statements in the DOS routines were changed to allow for up to 99999
   bands to be displayed
   [files: lev00.f, prep_dos.f]
2. get_param_siesta modified to read the main output file for the grid and
   spin if the RHO file is not available; also, one small bug fixed
   [file: get_param_siesta.f]

10.02.2007 Version 2.46.1
^^^^^^^^^^^^^^^^^^^^^^^^^
1. Bug fixes in reading unformatted siesta density
   [files: manip_dens.f, read_density.f]
2. Bug fixed in get_param_siesta: stopped if bands not found; now fixed
   [file: SIESTA_Tools/get_param_siesta.f]
3. New design for the hat (when lev00 starts), similar to that of TETR.
   tmake script creates build.hat file which is attached to hat.f
   [files: lev00.f, hat.f, make00.lev, tmake]

7.01.2007 Version 2.46
^^^^^^^^^^^^^^^^^^^^^^^
1. NEW: get_pram_siesta has been modified to read in NBANDS and NKPTS from siesta DOS run 
   (if the new option BandPoints in on), as well as band energies into band.out file,  
   to allow for DOS support in lev00
   [file: get_param_siesta.f]
2. NEW: corresponding changes have been done in lev00 routines as well to support total DOS
   calculation with siesta: presently supported menu options 0 and 3 of DOS (smearing 
   tests and total DOS)
   [files: input00.f, lev00.f]
3. character*(*) construction introduced in all routines to handle variable string lengths
   [files: tools_strings.f, density.f, dipole.f, input00.f, lev00.f, manip_dens.f, 
           plot_add.f, simulate.f] 
4. NEW: it is possible to set the upper limit to the energy (x axis) when previewing and
   of the DOS graphs with gnuplot
   [files: prep_dos.f, prep_disp.f, plot_add.f]
5. an old bug in working out physical bands, when the smearing is on, fixed
   [file: dos_add.f]

24.7.2005 Version 2.45
^^^^^^^^^^^^^^^^^^^^^^
1. the main menu corrected
   [file: lev00.f]]
2. English in several menus and some unimportant bugs (in calling plot() routine)
   [file: menu.f, prep_dos.f, prep_disp.f, plot_add.f, lev_coulomb,f, dipole.f]
3. several little bugs in the plane menu (affected only the menu look) fixed
   [file: plot_add.f]
4. some cosmetic changes in charge, dipole and quadrupole routines menus
   [files: dipole.f, plot_add.f]
5. a tiny uninportant bug in calling Plot_Explr() routine fixed
   [file: plot_add.f]
6. option 7 of the density menu is linked to CASTEP only
   [file: density.f]
7. removed 'island' message (given if option 4 of the main menu) for VASP and
   SIESTA
   [file: input00.f]
8. option 4 of the main menu is directly linked to CASTEP/CETEP only; in the
   case of VASP, the file name is asked (CHGCAR, PARCHG or LOCPOT) before
   reading in
   [file: lev00.f, density.f, input00.f]
9. manipulation option rewritten from scratch; includes SIESTA as well this time
   [files: lev00.f, manip_dens.f, density.f, read_density.f, make00.lev, tools_strings.f]
10. option 4 -> 4  removed (writing partial density)
   [file: prep_pl0.f]
11. option 8 of the main menu (shift density) was moved to the density menu as
   option 12, which allowed to generalise it for any code
   [files: lev00.f, density.f, make00.lev]
   [file removed: shift_charge.f]

2.2.2005 Version 2.44
^^^^^^^^^^^^^^^^^^^^^
1. Makefile make00.lev updated.
2. XY option in gOpenMol menu of the dnesity menu: default changed to keeping
   atoms inside the box.
   [file: density.f]
3. a bug in for_gOpenMol() routine (wrong calculation of center()) corrected.
   [file: density.f]
4. a very subtle bug (thanks to Chris!) in for_gOpenMol() routine, that
   resulted in a shift of the countour surfaces with respect to the geometry, 
   corrected.
   [file: density.f]

19.1.2005 Version 2.43
^^^^^^^^^^^^^^^^^^^^^^
1. can read unformatted siesta RHO (electron/spin density) files as well; detects 
   automatically if the file is formatted or unformatted. [Note: may not work
   if the unformatted RHO file was produced on another platform.]
   [files: read_density.f, SIESTA_Tools/get_param_siesta.f]
2. a debugging printout removed from vasp_dens() routine.
   [file: read_density.f]

10.1.2005 Version 2.42
^^^^^^^^^^^^^^^^^^^^^^
1. another format for reading k-points from the OUTCAR file implemented as the
one existed before did not seem to appear in the OUTCAR file  in some cases
   [file: input00.f]

24.11.2004  Version 2.41
^^^^^^^^^^^^^^^^^^^^^^^^
1. some minor changes to comply with Intel fortran-90 compiler
   [files: almost all]
2. make00.lev changed
3. new file tools_strings.f created that contains all the strings routines. A trick 
   in those to port an variable lenths string has been removed to comply with fortran
   rules.
   [file: inpu00.f]
4. PROCAR file with projected DOS is now working with option LORBIT=11 of VASP. Still, 
   only total p and d contributions are worked out by lev00, although all components of
   those are available in PROCAR. Can easily be taken into account as well, but later.
   [file: input00.f]

8.10.2004  Version 2.4
^^^^^^^^^^^^^^^^^^^^^^
1. cosmetic changes in get_all_from_SIESTA  (au_A,A_au)
   (file: input00.f)
2. similar cosmatic changes in siesta_dens 
   (file: read_density.f)
3. multiplication factor for the density used in line/plane options (multcon) 
   has been changed from 100 to 1, since the density is given in actual units
   (A^-3) 
   (files: menu.f, plot_add.f)
4. an option [Ap] is added to the plane menu to show (X,Y) positions of all atoms
   with respect to the local coordinate system fixed with the plane including
   the out-of-plane distance
   (files: density.f, plot_add.f)
5. the title line in the CONTCAR file is also printed on the screen for
   convenience to have the order of species 
   (file: input00.f)

30.3.2004 - 15.04.2004  Version 2.3
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
(1) a bug in density.f related to reading density of siesta (thanks to Ross)

(2) a bug in input00.f in setting up the first k-point while reading
siesta input (actually, is not used)

(3) menus in plot_add.f, simulate.f and density.f changed: 
 - to allow character inputs instead of only numbers; 
 - 'An' to choose the way coordinates are given; 
 - 'Co' is added in all menus options related to density to view atomic
   positions in fractional/Cartesian format
 - 'Q' is used to quit 
This all is in accord with the general implementation of menus in tetr.
Still, other menus will have to be modified accordingly - later!

(4) atomic species are now asked every time lev00 is started and the
coordinates are read in
(changes in: lev00.f, lev00.inc, mendeleev.inc, input00.f)

(5) density reading routines (vasp, siesta) have been moved from density.f
to a separate file read_density.f

(6) a bug in a routine reading siesta density has been corrected; it
affected only the units of the density array grid()
(the file affected: read_density.f)

(7) insignificant changes in density.f 

(8) A MAJOR NEW DEVELOPMENT: a new option added in the density menu (#10) 
allowing to write input files to gOpenMol to plot "orbitals" out of the 
density; two files are produced (both are needed): 

 - gOpenMol.cube containing the density and atomic positions in the
specific format; all atomic positions are printed (they are not actually used 
by gOpenMol which is a bit of a puzzle); the density is printed on a grid
within a box (the density box) to be specified; the box may run across 
several cells or overlap with several adjacent cells, etc.

 - gOpenMol.xyz which contains atomic positions (again!); two options
   are available (the XY option in the menu):
      - all cell atoms
      - only atoms within the density box (a larger box can also be chosen
        using the option 3 in the menu; the latter appears if the appropriate
         XY menu option has been chosen)

Files affected: density.f

What you need to do with gOpenMol:
 - read the gOpenMol.xyz file by going to: 
        File -> Import -> Coord -> Browse -> Apply
 - read in the gOpenMol.cube file by going to:
        Run -> gCube2plt/g94cub2pl (cube) -> Browse -> Apply
   which should produce an output in the window to be openned finished by 
   "Job done ..."; it produces a file gOpenMol.plt
 - then plot by reading the gOpenMol.plt:
        Plot -> Contour -> Browse -> Import File 
You will find the min and max values of the density; choose the density value
in the boxes underneath the min/max info and click Apply: the density will be
shown alongside the atomic positions in the same window; you will be able to 
rotate the image.

To make the density transparent, click on Details opposite to the chosen 
density level and change the Opacity level, then click Apply. It is possible 
to plot several density levels in the same way by filling in several density
levels and then changing their opacities.

To see atoms better, use: 
        View -> Atom Type 
and then the option Licorice with e.g. Cyl.rad = 0.1 and Sph. rad. = 0.3.
To have different radii for different species, use CPK as well.
See the gOpenMol manual for more instructions.
gOpenMol is freely available at: http://www.csc.fi/gopenmol/

1.02.2004 Version 2.2
^^^^^^^^^^^^^^^^^^^^^^
(1) a bug in lev00.f: a call to prep_plot() in potion MAP was missing
(2) added option 44 in the plane menu: after the plane is chosen
(i.e. the normal, central point and the dimensions), one can choose
another plane parallel to the existing one. Note that if the plane was
chosen using 3 points methos, you will have to choose the same plane again by
specifying the normal only (it is diplayed in the menu), so that the
3-point  method becomes off.

22.01.2004 Version 2.1
^^^^^^^^^^^^^^^^^^^^^^^
(1) do_param.f: do_param writes ISPIN variable to the param.inc file.

(2) lev00 (several files): spin 1,2 for spin-polarised calculations 
    is recognised by lev00 when doing various DOS; no need now to
    copy band.out.1 (or .2) into band.out for each spin; when going
    into any of the DOS options, lev00 asks which spin if ISPIN=2.
    In the case of a spin-polarised calculation any DOS is without
    the spin prefactor of 2.

(3) lev00.f: structure of the main routine has been simplified considerably
    and is now more transparent; some CASTEP/CETEP options are not visible
    if VASP or SIESTA are on and vice versa.

(4) when searching for a string consisting of 2 or 3 words, the string
    is broken down into words; this means that number of spaces
    between the words does not matter anymore. This is useful if
    e.g. the number of spaces was changed in the output file of the PW code
    which lev00 needs (as e.g. happenned between versions 4.4. and 4.6
    of vasp)  

19.09.2003  Version 2.0 
^^^^^^^^^^^^^^^^^^^^^^^

(1) some small bugs corrected.

(2) Added support for SIESTA (with some help from Adam Foster):

- lev00 can now fully support charge and spin density options
- due to a problem in setting users k-points to SIESTA, there
  is still no support for DOS and PDOS.

To compile lev00, you need param.inc file. The get one, run

get_param_siesta

(which is in TOOLS/SIESTA_Tools and should be compiled only once). It
will read the final geometry file *.XV and the density file *.RHO to
construct param.inc file. Then, as usual, run lev00.comp to get lev00
in your working directory.

