Show Posts

This section allows you to view all posts made by this member. Note that you can only see posts made in areas you currently have access to.


Messages - ipsecog

Pages: 1 2 [3]
31
Well, at least in NanoLanguage you can use the small script provided in the examples directory. It is called import_xyz.py and it's quite straight-forward to use; in fact, and example (import_xyz-test.py) is also included in the same directory.

For convenience, I copy the example here:

Code
from ATK.KohnSham import *

# Import XYZ file
from import_xyz import *
(elements,positions) = importXYZFile ('3-octanol.xyz')

# Create NanoLanguage molecule
molecule = MoleculeConfiguration(elements, positions*Angstrom)

# Save in VNL file   
vnlfile = VNLFile('3-octanol.vnl')
vnlfile.addToSample(new_vnl_molecule,'3-octanol')

Just remember to copy the file import_xyz.py to the same directory as your script, or add these lines at the top:

Code
import sys
sys.path.append('c:/Program Files/Atomistix/ATK 2008.02/examples')

(Yeah, I run ATK on Windows! the performance is actually not much worse, except there is no MPI parallelization of course...)

32
Putting data into VNL files is simple in NanoLanguage, using the addToSample() function, and appears to be a convenient way of storing data calculated by ATK. However, getting it out again can be a challenge (which kind of makes the idea of using it for storing data a bit less useful).

Of course, Virtual NanoLab can read the data and plot it, but what if you want the actual data, the raw numbers? That is, the data you had access to before you put it into the VNL file (typically via the toArray() method on the data object, e.g. the transmission spectrum).

Well, it turns out that the vnl files have a rather simple but undocumented structure; they are actually just Python zipfiles, and the data is stored in nested dicts with quite obvious labels! Therefore, it's relatively simple to extract the data, as long as you know what you are doing.

Disclaimer: The scripts below use undocumented features, and may not work in future (and/or some old) releases of ATK.

It's best to do some of these operations in interactive mode, but there is an initial part, which is in the attached script. Download it and run it as

Code
atk -i extract_from_vnl_file.py

Note the "-i" option; this will cause ATK to execute the script and then remain in interactive mode. The script will ask for the filename (because it cannot be given as a script parameter when you want to run in interactive mode); be sure to specify it with full path if it's not in the current directory. (On most platforms you can use the Tab key to complete the file name!)

At this point the script has created a variable (a dictionary) called samples which contains all the data. In this top-level dict, the keys correspond to the sample names. You can list the sample names:

Code
print samples.keys()

Let us assume that the sample we are interested in is called "mySample". To drill down to the data, give the command

Code
data = samples['mySample']['NanoLanguage']

Now, you can see which data has been stored on this sample:

Code
print data.keys()

Let us assume that we are looking to extract a transmission spectrum. In that case, you will hopefully see (perhaps among many other things) a key called "Transmission Spectrum Transmission Spectrum" in the printed list. However, it depends a bit on the version of ATK you used to create the VNL file. In some older versions it's just "Transmission Spectrum" and in some really old ones it might even be "TransmissionSpectrum". So, it's best to inspect the data to be sure!

Now we can extract the data:

Code
trans = data["Transmission Spectrum Transmission Spectrum"]

By typing "trans" and hitting the Tab key, we now find that the so-created "trans" data object has several methods and properties. (You can also give the command

Code
dir(trans)

to see the list.

What we are interested in here are the energies E and the values of the transmission at these energies, T(E). Here's how to access it:

Code
E = t.energy()._dataarray_
T = t.average_transmission()

And now we can just print the transmission spectrum (or do whatever you want with it):

Code
for i in range(len(E)):
    print E[i],T[i]

Note that if the calculation was spin-polarized, the transmission has a bit different shape.

Using the general method outlined above it should be possible to extract any kind of data stored in a VNL file. Do note, however, that each specific data type is stored in different ways, and you have to figure out the (undocumented) methods and properties for each one by using the "dir()" command quite a lot...

To summarize, here is a complete script that extracts and prints the transmission spectrum from a VNL file, by asking the user for the file and sample name (this script is also attached to this post as extract_transmission_from_vnl_file.py):

Code
import zipfile, cPickle, sys

filename = raw_input('Please enter VNL file name: ')
try:
    f = zipfile.ZipFile(filename,'r')
except:
    print 'Unable to locate VNL file "%s"' % filename
    sys.exit(1)

# Restore the samples from the VNL file into a dict
samples = {}
for zinfo in f.infolist():
    s = f.read(zinfo.filename)
    try:
        obj = cPickle.loads(s)
        samples[obj.name()] = {}
        for h in obj.history():
            samples[obj.name()][h.name] = h.resultsample.properties()
    except:
        pass

print 'The following samples are present in this VNL file:'
print samples.keys()

sname = raw_input('Please enter sample name to extract (default=first sample): ')
if sname.strip()=="":
    data = samples[samples.keys()[0]]['NanoLanguage']
else:
    data = samples[sname]['NanoLanguage']

# Some older versions used "Transmission Spectrum" or "TransmissionSpectrum"
trans = data['Transmission Spectrum Transmission Spectrum']
print 'Energy (eV)\tTotal transmission'
print '----------------------------------'
energies = trans.energy()._dataarray_
# Check for energy unit, if you want to be sure
# e_unit = t.energy().units()
T = trans.average_transmission()
for i in range(len(energies)):
    print '%g\t%g' % (energies[i],T[i])

33
Scripts, Tutorials and Applications / Re: Density difference
« on: December 3, 2008, 15:11 »
Another example: How to calculate the spin polarization.

This script performs a (simplified) self-consistent calculation of an oxygen (O2) molecule and afterwards computes and stores in the VNL file the spin polarization density (spin up density minus spin down density).

Code
from ATK.KohnSham import *

# Define the geometry for O2 molecule:
elements = [Oxygen, Oxygen]
coordinates =  [ (0.0,0.0,0.6) , (0.0,0.0,-0.6) ]*Ang
o2 = MoleculeConfiguration(elements, coordinates)

# Set the initial density with a maximum scaled spin:
initial_electron_density = electronDensityParameters(
    initial_scaled_spin=(1,1)
    )

# Define the calculation method:
dft_method=KohnShamMethod(
        electron_density_parameters=initial_electron_density,
        basis_set_parameters = basisSetParameters (type = SingleZeta)
        )

import ATK
ATK.setCheckpointFilename("o2.nc")
ATK.setVerbosityLevel (level = 1)

# And execute the calculation:
scf=dft_method.apply(o2)

def calculateSpinDensityDifference (spin_density):

    import copy
   
    n_up = spin_density.toArray()[0]
    n_dn = spin_density.toArray()[1]
    new_density = copy.copy(spin_density)
    new_density._ElectronDensity__electron_density_data = (n_up-n_dn)
    return new_density

density = calculateElectronDensity(scf)
dens_diff = calculateSpinDensityDifference(density)

f = VNLFile('o2.vnl')
f.addToSample(o2,'Oxygen molecule')
f.addToSample(dens_diff,'Oxygen molecule','Spin polarization density')
f.addToSample(density,'Oxygen molecule')

The result, plotted in VNL is attached!

34
Scripts, Tutorials and Applications / Re: Density difference
« on: December 3, 2008, 15:08 »
Here is an example, using the function above! This is based on the famous Li-H2-Li example in the manual, and intends to calculate and plot the change in the self-consistent electron density at 0.7 V bias compared to the zero-bias case.

Code
from ATK.TwoProbe import *

def calculateDensityDifference (d1,d2):

    import copy
   
    n1 = d1.toArray()
    n2 = d2.toArray()
    density = copy.copy(d1)
    density._ElectronDensity__electron_density_data = (n1-n2)
    return density

scf1 = restoreSelfConsistentCalculation ('lih2li-bias-0.0.nc')
scf2 = restoreSelfConsistentCalculation ('lih2li-bias-0.7.nc')
density1 = calculateElectronDensity(scf1)   # 0.0 V bias
density2 = calculateElectronDensity(scf2)   # 0.7 V bias
diff_density = calculateDensityDifference(density1,density2)

f = VNLFile('lih2li.vnl')
f.addToSample(density1,'lih2li','Density 0 V')
f.addToSample(density2,'lih2li','Density 1 V')
f.addToSample(diff_density,'lih2li','Difference density')

The resulting plot from VNL is shown in the attached figure! (Coming as soon as I've calculated it!)

35
Scripts, Tutorials and Applications / Density difference
« on: December 3, 2008, 15:06 »
NanoLanguage doesn't allow you to subtract two ElectronDensity objects from each other (it does allow you to add them, however), which means it's not immediately obvious how to calculate e.g. the spin polarization density.

By using some undocumented tricks, it is however possible to do it relatively easily anyway.

Disclaimer: The scripts below use undocumented features, and may not work in future (and/or some old) releases of ATK.

The steps involved are:

1) Extract the actual array data from the ElectronDensity objects. This is easy, just use the toArray() method on the object.
2) Calculate the difference. This is also easy, because the arrays (not the objects!) support the minus operation.
3) Store the difference. This is the tricky bit, but all we need to know is that the density data is stored in a private property called _ElectronDensity__electron_density_data.

The second part of the trick is to use the Python module copy to create a new ElectronDensity object (since this class has no constructor exposed in NanoLanguage).

Putting this together in a utility function, we obtain

Code
def calculateDensityDifference (d1,d2):

    import copy
   
    n1 = d1.toArray()
    n2 = d2.toArray()
    density = copy.copy(d1)
    density._ElectronDensity__electron_density_data = (n1-n2)
    return density

Calling this with two densities (d1 and d2), the function returns their difference in an ElectronDensity object, which you can put into a VNL file just as normal, and plot in Virtual NanoLab.

Pages: 1 2 [3]