[gpaw-users] Problem simulating STM images
Susi Toma
toma.susi at aalto.fi
Tue Sep 24 12:31:03 CEST 2013
Hello again,
I have a few further questions.
First of all, my graphene unit cell is not a square, but orthorhombic, and furthermore, it is 8 unit cells in one direction and 6 in the other (made with the attached script). The stm image simulation, however, yields square images, resulting in a significant and non-trivial distortion compared to the actual real-space structure. This makes interpretation of the images rather hard.
Can the stm functionality somehow create non-square images? Or if not, is there a way to automatically overlay the structure on the images?
Many thanks,
Toma
On 5.9.2013, at 15.28, Toma Susi <toma.susi at aalto.fi<mailto:toma.susi at aalto.fi>> wrote:
Hi,
After Jussi kindly updated the gpaw-env module on the CSC Taito cluster, the simulations now run correctly following Jens' instructions. I'm now getting reasonable looking STM image simulations for my system, with no problems with either of the previous errors regardless of the values I use for z, bias or z0.
With any luck we'll have a chance to compare the simulations to experimental measurements that we are hoping to do soon.
Thanks a lot for all your help!
Best,
Toma
On 30.8.2013, at 10.34, Jens Jørgen Mortensen <jensj at fysik.dtu.dk<mailto:jensj at fysik.dtu.dk>> wrote:
Den 29-08-2013 08:58, Susi Toma skrev:
On 29.8.2013, at 8:54, Jens Jørgen Mortensen wrote:
Den 28-08-2013 14:46, Susi Toma skrev:
Hi again,
I managed to get the stm tutorial to run on one of the older CSC clusters, and am in contact with them to sort out the issues on the newer ones. However, I am still having trouble getting the stm scan of my particular system.
I'm closely following the tutorial, but since in my case I have a planar system, I have k-points only in the x and y directions. Thus if I try to run the tutorial's code "stm = STM(atoms, symmetries=[0, 1, 2])" on my .gpw file, I get the following error:
>>> c = stm.get_averaged_current(bias, z)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/home/u22/tsusi/ase/ase/dft/stm.py", line 95, in get_averaged_current
self.calculate_ldos(bias)
File "/home/u22/tsusi/ase/ase/dft/stm.py", line 78, in calculate_ldos
ldos += ldos.transpose((1, 0, 2)).copy()
ValueError: operands could not be broadcast together with shapes (96,72,56) (72,96,56) (96,72,56)
However, changing the symmetries to what I presume to be the correct "stm = STM(atoms, symmetries=[0, 1])" and then continuing gets me one step further with a successful run of the stm.get_averaged_current() command. However, the next step crashes again:
>>> c = stm.get_averaged_current(bias, z)
>>> h = stm.scan(bias, c)
Traceback (most recent call last):
File "<stdin>", line 1, in <module>
File "/home/u22/tsusi/ase/ase/dft/stm.py", line 119, in scan
heights[i] = find_height(a, current, h)
File "/home/u22/tsusi/ase/ase/dft/stm.py", line 170, in find_height
c2, c1 = ldos[n:n + 2]
ValueError: need more than 1 value to unpack
I have defined the variables bias = 1 and z = 3, so the error doesn't seem to be due to an undefined variable. What am I doing wrong?
>From looking at the code, it would seem that your current is too low. Is z=3 Ang the correct place in your unit cell? Let's say your unit cell is 20 Ang high and your atoms are distributed from 7 to 13 Ang and you want to scan at an approximate height of 4 Ang above the surface, then you should use z=13+4 in your call to stm.get_averaged_current() and then play with the value you get.
As you can see you are the first one to use the code for real systems :-) so your feedback is very important. If you give me a script that can calculate the wave functions for your system, then I could give it a shot also and fix some of the strange error messages you got.
Jens Jørgen
Sure, I'd love to help benchmark the code with a real system.
My graphene plane is at z = 0, so I'd assumed z = 3 is reasonable (btw, is it advisable to have planes/slabs at the center of unit cells? I haven't had issues with this in other calculations).
I always do that and the STM code assumes that you do, but with the latest version of ASE (that I have just checked in) you can specify the hight where you want to scan from - default is to start from the top of the cell and search downwards until you find the correct current.
In your case, do this:
h = stm.scan(bias, c, z0=6.0)
where z0 is somewhere in the vacuum where the current is below c.
If you get a "Tip crash" you should lower the current.
In any case, I tried centering the atoms, but it didn't make a difference.
Strange. That should have worked.
One other thing: From your GPAW text output:
Symmetries present: 1
25 k-points: 5 x 5 x 1 Monkhorst-Pack grid
13 k-points in the Irreducible Part of the Brillouin Zone
you can see that there are no symmetries found (only the identity) so you should use:
stm = STM(atoms, symmetries=[])
or just "stm=STM(atoms)".
On last thing: The wave function file is quite big, so you may want to not write it. Instead you can calculate the current directly in the script where you calculate the wave functions and then write that to a file - takes much less space:
from ase.dft.stm import STM
from ase.parallel import rank
stm = STM(atoms)
stm.calculate_ldos(1.0)
if rank == 0:
stm.write('ldos.1.0.pckl')
Then in your plotting script you can use:
stm = STM('ldos.1.0.pckl')
and play with z, c and z0.
Jens Jørgen
The current does perhaps seem a bit low, but playing with z and bias I got the averaged current even up to 0.38 – no difference, same ValueError. Incidentally, at one point I accidentally played with the bias without doing the "c = stm.get_averaged_current(bias, z)" step in between, and for some values managed to reproduce the tip crash error. No successful calculation, though.
I've attached the relaxed POSCAR file and a script for getting the wavefunctions; I've been running relaxations with 128 processors and 2048 MB memory each, but you might get away with a bit less. The script takes the POSCAR structure file as argument. I also included the stm script, same deal there.
-Toma
Many thanks,
Toma
On 27.8.2013, at 15:40, Marcin Dulak wrote:
On 08/27/2013 01:44 PM, Susi Toma wrote:
On 27.8.2013, at 14:34, Marcin Dulak wrote:
On 08/27/2013 12:14 PM, Susi Toma wrote:
On 27.8.2013, at 13:03, Marcin Dulak wrote:
Hi,
On 08/27/2013 11:06 AM, Susi Toma wrote:
Hi,
I've done as you suggested, and am now running the latest snv version. However, now I get an error on the numpy version number: "RuntimeError: module compiled against API version 7 but this version of numpy is 6". This does not
i guess you using a different GPAW installation now, compiled against a different version of numpy (1.7) than the available (or incorrectly loaded) one (1.6).
Check the loaded version of numpy with: python -c "import numpy; print numpy.__version__"
Revert to the old GPAW installation, and try changing to the latest ASE only:
export PYTHONPATH=~/myase:$PYTHONPATH
Hi,
So, the way GPAW works on the Taito cluster is to load the gpaw-env module, which results in a number modules being loaded:
----------------------------------------------------------------------------
Setting up new environment, removing all currently loaded modules
----------------------------------------------------------------------------
Loading new modules:
gcc/4.7.2 intelmpi/4.1.0 ase/3.6.0 openblas/0.2.6 python/2.7.3 hdf5-par/1.8.10 gpaw/0.9.0
Due to MODULEPATH changes the following have been reloaded:
1) intelmpi/4.1.0
maybe it's the ase module that loads numpy-1.7?
At that point do:
python -c "import numpy; print numpy.__version__, numpy"
module show ase
module show python
module show gpaw
module show gpaw-env
Here's the output after loading gpaw-env – it seems it loads numpy 1.6.2:
[tsusi at taito-login3 ~]$ python -c "import numpy; print numpy.__version__, numpy"
1.6.2 <module 'numpy' from '/appl/nano/gpaw/numpy/lib/python/numpy/__init__.pyc'>
[tsusi at taito-login3 ~]$ module show ase
------------------------------------------------------------
/appl/modulefiles/Linux/ase/3.6.0.lua:
------------------------------------------------------------
help([[ This loads the Atomic Simulation Environment.
Version 3.6.0
Modifies: PATH, PYTHONPATH
]])
prepend_path("PATH", "/appl/nano/ase/3.6.0/bin")
prepend_path("PYTHONPATH", "/appl/nano/ase/3.6.0/lib/python")
prepend_path("PATH", "/appl/nano/ase/3.6.0/bin")
[tsusi at taito-login3 stm]$ module show python
------------------------------------------------------------
/appl/modulefiles/MPI/gcc/4.7.2/intelmpi/4.1.0/python/2.7.3.lua:
------------------------------------------------------------
prereq("openblas/0.2.6")
family("python")
prepend_path("MANPATH", "/appl/opt/python/2.7.3-gcc-shared/share/man")
prepend_path("PATH", "/appl/opt/python/2.7.3-gcc-shared/bin")
prepend_path("LD_LIBRARY_PATH", "/appl/opt/python/2.7.3-gcc-shared/lib")
prepend_path("LD_LIBRARY_PATH", "/appl/opt/mysql/mysql-5.6.10-linux-glibc2.5-x86_64/lib")
[tsusi at taito-login3 stm]$ module show gpaw
------------------------------------------------------------
/appl/modulefiles/MPI/gcc/4.7.2/intelmpi/4.1.0/gpaw/0.9.0.lua:
------------------------------------------------------------
help([[ This loads GPAW real-space DFT program
Version 0.9.0
Modifies: PYTHONPATH, PATH, GPAW_PYTHON
]])
prereq("ase/3.6.0", "openblas/0.2.6", "python/2.7.3", "hdf5-par/1.8.10")
prepend_path("PATH", "/appl/nano/gpaw/0.9.0/bin")
prepend_path("PYTHONPATH", "/appl/nano/gpaw/0.9.0/lib/python")
prepend_path("PYTHONPATH", "/appl/nano/gpaw/numpy/lib/python")
setenv("GPAW_SETUP_PATH", "/appl/nano/gpaw/setups/0.8.7929")
[tsusi at taito-login3 stm]$ module show gpaw-env
------------------------------------------------------------
/appl/modulefiles/Linux/gpaw-env/0.9.0.lua:
------------------------------------------------------------
help([[ This module loads the gpaw environment version 0.9.0
Note that this module removes all other modules during load
]])
If I don't unload the loaded ase/3.6.0, the system uses it:
Python 2.7.3 (default, Mar 18 2013, 15:03:18)
[GCC 4.7.2] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> import ase
>>> print ase
<module 'ase' from '/appl/nano/ase/3.6.0/lib/python/ase/__init__.pyc'>
However, having made the suggested myase installation, when I unload the gpaw-env's ase module, the system resorts to using the latest svn version I installed:
[tsusi at taito-login4 stm]$ module unload ase
print numpy information again after unloading ase:
python -c "import numpy; print numpy.__version__, numpy"
After unloading ase:
[tsusi at taito-login3 stm]$ python -c "import numpy; print numpy.__version__, numpy"
1.6.2 <module 'numpy' from '/appl/nano/gpaw/numpy/lib/python/numpy/__init__.pyc'>
So still the old version.
[tsusi at taito-login4 stm]$ module list
Currently Loaded Modules:
1) gpaw-env/0.9.0 3) intelmpi/4.1.0 5) python/2.7.3 7) gpaw/0.9.0
2) gcc/4.7.2 4) openblas/0.2.6 6) hdf5-par/1.8.10
[tsusi at taito-login4 stm]$ python
Python 2.7.3 (default, Mar 18 2013, 15:03:18)
[GCC 4.7.2] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> import ase
>>> print ase
<module 'ase' from '/homeappl/home/tsusi/myase/ase/__init__.pyc'>
you don't need to unload ase, just prepend the PYTHONPATH, after loading of all modules:
export PYTHONPATH=~/myase:$PYTHONPATH
Normally you can have all the module loading and exporting of variables in the batch system script
that you submit, and if the only code you are using is GPAW you can even add all these to your login script (~/.bashrc),
so you have all the settings ready at login and when a batch job starts.
Ah, right, doing the export manually after loading of modules overrides the path, of course. However, this still doesn't help: the loaded numpy version doesn't work with ase/svn.
I think I have to contact the administrators and ask them to update numpy, or...?
the question if where numpy-1.7 comes from, and I think the "RuntimeError: module compiled against API version 7 but this version of numpy is 6"
will only show with GPAW, not ASE. At this point:
ls -al `which gpaw-python`
gpaw-python
python -c "import ase; print ase" # we know that works and points to myase
ls -ald $(dirname `python -c "import ase; print ase.__file__"`)
python -c "import gpaw; print gpaw"
Does the original gpaw-python work still without unloading ase or exporting myase PYTHONPATH?
Marcin
________________________________
"Love is the only emotion that enhances our intelligence."
-Humberto Maturana
________________________________
http://physics.aalto.fi/personnel/?id=322
http://mostlyphysics.wordpress.com/
"Love is the only emotion that enhances our intelligence."
-Humberto Maturana
________________________________
http://physics.aalto.fi/personnel/?id=322
http://mostlyphysics.wordpress.com<http://mostlyphysics.wordpress.com/>
"Love is the only emotion that enhances our intelligence."
-Humberto Maturana
________________________________
http://mostlyphysics.wordpress.com
"Love is the only emotion that enhances our intelligence."
-Humberto Maturana
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://listserv.fysik.dtu.dk/pipermail/gpaw-users/attachments/20130924/dbca9c8e/attachment-0001.html>
-------------- next part --------------
A non-text attachment was scrubbed...
Name: make_gra.py
Type: text/x-python-script
Size: 605 bytes
Desc: make_gra.py
URL: <http://listserv.fysik.dtu.dk/pipermail/gpaw-users/attachments/20130924/dbca9c8e/attachment-0001.bin>
More information about the gpaw-users
mailing list