Voltage sweep (I-V curve)

[Anders Blom]

One of the fundamentally important applications of ATK is to produce an I-V curve for a two-probe system. There is, however, no simple way to do it without quite a bit of scripting, as noted in another post on this Forum.

Here is a script/module ivcurve.py (attached) that takes care of that, at least in a quick-and-dirty way. The script is not a masterpiece of software design, and does no error checking etc, but ... it works :-)

The idea is to minimize the code the user has to write in his own script. The only required code to add is something along the lines of

Code:

import ivcurve
   
voltages=[0.0,0.1,0.2,0.3]*Volt

ivcurve.runIVcurve (
    twoprobe_configuration,
    two_probe_method,
    runtime_parameters,
    voltages,
    vnl_filename='myfile.vnl', sample_name='mysample',
    current_k_point_sampling = (1,1),
    current_number_of_points = 100
  )

iv = ivcurve.extractIVcurveFromVNLFile('myfile.vnl','mysample')

ivcurve.plotIVCurve(iv,'iv.png')

Hopefully the above is relatively self-explanatory, but we'll review some details below. Note how the variable names above are chosen to match those used by Virtual NanoLab when it produces NanoLanguage code for a two-probe system. Therefore it is very easy to just include the code above at the end of a script produced by Virtual NanoLab to obtain an I-V curve.

For a more complete example, see the attached script lih2li_iv.py. It is based on the Li-H2-Li system from the manual, just because it's relatively quick to compute, so it's easy to test the I-V curve script using this system. To run the example you should have the geometry file lih2li.vnl too; for simplicity it is also attached.

The code above should be inserted into a script, after the usual part which defines the TwoProbeConfiguration and the TwoProbeMethod, and also the run-time parameters. Then, insert the code snipped above at the bottom of the script, replacing the execution statement (executeSelfConsistentCalculation()) if present. If you used your own variable names for the configuration, method and run-time parameters, adjust accordingly.

The code snippet above

• defines the voltages for which to compute the current
• runs a sweep with one calculation for each voltage bias
• produces a plot, stored in "iv.png" of the I-V curve

Important notes:

• The bias values should come in sequence and start with 0 V, since for each bias, the previous converged calculation is used to initialize the calculation.
• An alternative way to set up the bias is to use numpy; to sweep from 0.0 to 1.0 V in steps of 0.1 V, use

Code:

import numpy
voltages = numpy.arange(0.,1.01,0.1)*Volt

• If the runtime parameters contain a checkpoint filename, the value of the bias will be added before the extension to create a unique NetCDF file for each bias. Thus, if you specified file.nc for a run with bias 0.0, 0.1 and 0.2 V, you will get three nc files file-0.0.nc, file-0.1.nc, and file-0.2.nc.
• The runtime parameters argument is required, but you can give None if you don't want any NetCDF files.
• All the computed values of the current will be stored in the VNL file specified, under the given sample name. You can inspect these values in Virtual NanoLab afterwards too.
• Remember to provide the correct integration parameters for the current (k-point sampling, in particular!).
• The plot can be customized to your liking. ATK actually has a built-in plotting engine (matplotlib), and you can produce very beautiful plots directly in NanoLanguage by using it. I'll make a separate post on that below!
• You can change the extension of the image file from "png" to any other supported file format (such as "eps").
• To print the I-V values, add this code to the end of the script:

Code:

print 'Bias (V)\tCurrent (A)\n' + '-'*40
for i in iv:
    print i[0],'\t\t',i[1]

It's not unlikely that I overlooked some detail in this script, which causes it not to work in a particular situation... But please try it out, and together we can update it to make it a powerful NanoLanguage tool!

Edit: Updated version of script, which is safe for running in parallel.

[serhan]

Thank you very much Dr. Blom. I will try as soon as possible:)

Regards
Serhan

[Anders Blom]

I mentioned in the original post that you can do nice plotting with ATK + matplotlib. A specific example of this is how the I-V script plots the results, and we can present that in a bit more detail, for those interested.

Remember from the actual I-V calculation that we have stored the values of bias and current in an array called "iv".

The simplest version of a plotting code would look like this:

Code:

import matplotlib
matplotlib.use('Agg')
import pylab as P
from numpy import array

X = array(iv)[:,0]
Y = array(iv)[:,1]
P.plot(X,Y)
 
P.xlabel('Bias (Volt)')
P.ylabel('Current (Ampere)')
P.title('I-V curve for Li-H2-Li')

P.savefig('lih2li_iv.png')

This looks very much like Matlab, and should be trivial to understand. The plot is attached as "lih2li_iv.png", and we can see that there are some things we would like to improve one (no plot points, scaling of the y-axis values).

Now, in the actual I-V script we go a step beyond this, since we wish to give the user control over things like the output file name, plot title, axis labels, etc, but also

• scaling (currents are in microamperes)
• plot symbols and colors
• axis min/max
• ... and perhaps something more

There are many ways to solve that, but one simple way is to use a dictionary. The script contains an internal dict(ionary) the defines the default values, and we then allow the user to pass another dict to the plot function, and use it to overwrite the defaults with any user-defined parameters.

A simple example of this technique is

Code:

# Default values
parameters = { 'a' : 5, 'b' : 3 }
# User-defined parameters
my_parameters = { 'b' : 8 }
# Apply user parameters
parameters.update(my_parameters)

After that, the dict "parameters" will have a=5 (default) and b=8 (user-defined).

We use this in the plot function of the I-V script as follows. The default values are:

Code:

plot_params = {
        # Blue circles with a line is default plot style
        'plot_color' : 'b',   
        'plot_symbol' : 'o-',
        'ylabel' : 'Current (micro-Ampere)',
        'xlabel' : 'Bias (V)',
        'grid' : False,
        'title' : 'I-V curve',
        'ymin' : None,
        'ymax' : None,
        # Scaling factor for the y-axis
        'yscale' : 1.e6
    }

To result of using these default can be seen in the original post (updated!).

The user can then make his own preferences, like

Code:

    my_plot_params = {
        # I want a green line with X plot symbols
        'plot_color' : 'g',
        'plot_symbol' : 'x-',
        # No scaling, please
        'yscale' : 1
        'ylabel' : 'Current (Ampere)',
        # I want a grid
        'grid' : False,
        'title' : 'I-V curve',
    }

He then passes this to the plot script:

Code:

ivcurve.plotIVCurve(iv,'lih2li_iv_green.png',my_plot_params)

The result is also attached, as "lih2li_iv_green.png".

As mentioned above already, more details on matplotlib is available on their webpage.

Hopefully this demonstrates some of the powerful stuff you can do with plotting inside ATK!

[zdhlover]

I have compute successfully with your script, then I create a script through the Nanolanguage Scripter, and I try my best to modify the script to get its' I-V curve as mentioned above ,but I couldn't get it.Could you helpe me ?Give me a simple demonstration with my script.

from ATK.TwoProbe import *
from ATK.MPI import processIsMaster

# Generate time stamp
if processIsMaster():
    import platform, time
    print '#',time.ctime()
    print '#',platform.node(),platform.platform()+'\n'

# Opening vnlfile
if processIsMaster(): file = VNLFile('D:/fe-mgo-fe.vnl')

# Scattering elements and coordinates
scattering_elements = [Iron,      Iron,      Iron,      Iron,     
                       Oxygen,    Oxygen,    Oxygen,    Oxygen,   
                       Magnesium, Magnesium, Magnesium, Magnesium,
                       Iron,      Iron,      Iron,      Iron]
scattering_coordinates = [[  0.7165 ,   0.7165 ,   1.433  ],
                          [  2.1495 ,   2.1495 ,   2.866  ],
                          [  0.7165 ,   0.7165 ,   4.299  ],
                          [  2.1495 ,   2.1495 ,   5.732  ],
                          [  2.1495 ,   2.1495 ,   7.932  ],
                          [  0.7165 ,   0.7165 ,  10.12735],
                          [  2.1495 ,   2.1495 ,  12.3227 ],
                          [  0.7165 ,   0.7165 ,  14.51805],
                          [  0.7165 ,   0.7165 ,   7.932  ],
                          [  2.1495 ,   2.1495 ,  10.12735],
                          [  0.7165 ,   0.7165 ,  12.3227 ],
                          [  2.1495 ,   2.1495 ,  14.51805],
                          [  0.7165 ,   0.7165 ,  16.71805],
                          [  2.1495 ,   2.1495 ,  18.15105],
                          [  0.7165 ,   0.7165 ,  19.58405],
                          [  2.1495 ,   2.1495 ,  21.01705]]*Angstrom
       

# Set up left electrode
left_electrode_elements = [Iron, Iron, Iron, Iron]
left_electrode_coordinates = [[ 0.7165,  0.7165,  0.7165],
                              [ 2.1495,  2.1495,  2.1495],
                              [ 0.7165,  0.7165,  3.5825],
                              [ 2.1495,  2.1495,  5.0155]]*Angstrom

left_electrode_cell = [[ 2.866,  0.   ,  0.   ],
                       [ 0.   ,  2.866,  0.   ],
                       [ 0.   ,  0.   ,  5.732]]*Angstrom

# Set up right electrode
right_electrode_elements = [Iron, Iron, Iron, Iron]
right_electrode_coordinates = [[ 0.7165,  0.7165,  0.7165],
                               [ 2.1495,  2.1495,  2.1495],
                               [ 0.7165,  0.7165,  3.5825],
                               [ 2.1495,  2.1495,  5.0155]]*Angstrom

right_electrode_cell = [[ 2.866,  0.   ,  0.   ],
                        [ 0.   ,  2.866,  0.   ],
                        [ 0.   ,  0.   ,  5.732]]*Angstrom

# Set up electrodes
left_electrode_configuration = PeriodicAtomConfiguration(
    left_electrode_cell,
    left_electrode_elements,
    left_electrode_coordinates
    )

right_electrode_configuration = PeriodicAtomConfiguration(
    right_electrode_cell,
    right_electrode_elements,
    right_electrode_coordinates
    )

# Set up two-probe configuration
twoprobe_configuration = TwoProbeConfiguration(
    (left_electrode_configuration,right_electrode_configuration),
    scattering_elements,
    scattering_coordinates,
    electrode_repetitions=[[1,1],[1,1]],
    equivalent_atoms=([0,0],[3,15])
    )
if processIsMaster(): nlPrint(twoprobe_configuration)
if processIsMaster(): file.addToSample(twoprobe_configuration, 'twoprobe_configuration')

######################################################################
# Central region parameters
######################################################################
exchange_correlation_type = LDA.PZ

iteration_mixing_parameters = iterationMixingParameters(
    algorithm = IterationMixing.Pulay,
    diagonal_mixing_parameter = 0.1,
    quantity = IterationMixing.Hamiltonian,
    history_steps = 6
)

electron_density_parameters = electronDensityParameters(
    mesh_cutoff = 100.0*Rydberg
)

basis_set_parameters = basisSetParameters(
    type = DoubleZetaPolarized,
    radial_sampling_dr = 0.001*Bohr,
    energy_shift = 0.01*Rydberg,
    delta_rinn = 0.8,
    v0 = 40.0*Rydberg,
    charge = 0.0,
    split_norm = 0.15
)

iteration_control_parameters = iterationControlParameters(
    tolerance = 1e-005,
    criterion = IterationControl.TotalEnergy,
    max_steps = 100
)

electrode_voltages = (0.0,0.0)*Volt

two_probe_algorithm_parameters = twoProbeAlgorithmParameters(
    electrode_constraint = ElectrodeConstraints.Off,
    initial_density_type = InitialDensityType.EquivalentBulk
)

energy_contour_integral_parameters = energyContourIntegralParameters(
    circle_points = 50,
    integral_lower_bound = 7.0*Rydberg,
    fermi_line_points = 10,
    fermi_function_poles = 4,
    real_axis_infinitesimal = 0.01*electronVolt,
    real_axis_point_density = 0.02*electronVolt
)

two_center_integral_parameters = twoCenterIntegralParameters(
    cutoff = 2500.0*Rydberg,
    points = 1024
)

######################################################################
# Left electrode parameters
######################################################################
left_electrode_electron_density_parameters = electronDensityParameters(
    mesh_cutoff = 100.0*Rydberg
)

left_electrode_iteration_control_parameters = iterationControlParameters(
    tolerance = 1e-005,
    criterion = IterationControl.TotalEnergy,
    max_steps = 100
)

left_electrode_brillouin_zone_integration_parameters = brillouinZoneIntegrationParameters(
    monkhorst_pack_parameters = (5, 5, 100)
)

left_electrode_iteration_mixing_parameters = iterationMixingParameters(
    algorithm = IterationMixing.Pulay,
    diagonal_mixing_parameter = 0.1,
    quantity = IterationMixing.Hamiltonian,
    history_steps = 6
)

left_electrode_eigenstate_occupation_parameters = eigenstateOccupationParameters(
    temperature = 1300.0*Kelvin
)

######################################################################
# Collect left electrode parameters
######################################################################
left_electrode_parameters = ElectrodeParameters(
    brillouin_zone_integration_parameters = left_electrode_brillouin_zone_integration_parameters,
    electron_density_parameters = left_electrode_electron_density_parameters,
    eigenstate_occupation_parameters = left_electrode_eigenstate_occupation_parameters,
    iteration_mixing_parameters = left_electrode_iteration_mixing_parameters,
    iteration_control_parameters = left_electrode_iteration_control_parameters
)

######################################################################
# Right electrode parameters
######################################################################
right_electrode_electron_density_parameters = electronDensityParameters(
    mesh_cutoff = 100.0*Rydberg
)

right_electrode_iteration_control_parameters = iterationControlParameters(
    tolerance = 1e-005,
    criterion = IterationControl.TotalEnergy,
    max_steps = 100
)

right_electrode_brillouin_zone_integration_parameters = brillouinZoneIntegrationParameters(
    monkhorst_pack_parameters = (5, 5, 100)
)

right_electrode_iteration_mixing_parameters = iterationMixingParameters(
    algorithm = IterationMixing.Pulay,
    diagonal_mixing_parameter = 0.1,
    quantity = IterationMixing.Hamiltonian,
    history_steps = 6
)

right_electrode_eigenstate_occupation_parameters = eigenstateOccupationParameters(
    temperature = 1300.0*Kelvin
)

######################################################################
# Collect right electrode parameters
######################################################################
right_electrode_parameters = ElectrodeParameters(
    brillouin_zone_integration_parameters = right_electrode_brillouin_zone_integration_parameters,
    electron_density_parameters = right_electrode_electron_density_parameters,
    eigenstate_occupation_parameters = right_electrode_eigenstate_occupation_parameters,
    iteration_mixing_parameters = right_electrode_iteration_mixing_parameters,
    iteration_control_parameters = right_electrode_iteration_control_parameters
)

######################################################################
# Initialize self-consistent field calculation
######################################################################
two_probe_method = TwoProbeMethod(
    electrode_parameters = (left_electrode_parameters,right_electrode_parameters),
    exchange_correlation_type = exchange_correlation_type,
    iteration_mixing_parameters = iteration_mixing_parameters,
    electron_density_parameters = electron_density_parameters,
    basis_set_parameters = basis_set_parameters,
    iteration_control_parameters = iteration_control_parameters,
    energy_contour_integral_parameters = energy_contour_integral_parameters,
    two_center_integral_parameters = two_center_integral_parameters,
    electrode_voltages = electrode_voltages,
    algorithm_parameters = two_probe_algorithm_parameters
)
if processIsMaster(): nlPrint(two_probe_method)

runtime_parameters = runtimeParameters(
    verbosity_level = 10,
    checkpoint_filename = 'D:/fe-mgo-fe.nc'
)

# Perform self-consistent field calculation
scf = executeSelfConsistentCalculation(
    twoprobe_configuration,
    two_probe_method,
    runtime_parameters = runtime_parameters
)

######################################################################
# Calculate physical properties
######################################################################
current = calculateCurrent(
    self_consistent_calculation = scf,
    brillouin_zone_integration_parameters = brillouinZoneIntegrationParameters((1, 1)),
    green_function_infinitesimal = 1.0e-5*electronVolt,
    number_of_points = 100
)
if processIsMaster(): nlPrint(current)
if processIsMaster(): file.addToSample(current, 'twoprobe_configuration', 'Current')

[Anders Blom]

It should be quite easy to do. Let's see, step by step.

1. Insert the line

Code:

import ivcurve

right after the two import statements at the top of your script. The reason we do this early in the script is because if it fails (for some reason), we will know immediately. The file needs to be in the same directory as your script.

2. Remove all code in your script after and including the line

Code:

# Perform self-consistent field calculation

3. Insert this code at the end of the script instead:

Code:

voltages=[0.0,0.1,0.2,0.3]*Volt

ivcurve.runIVcurve (
    twoprobe_configuration,
    two_probe_method,
    runtime_parameters,
    voltages,
    vnl_filename='FeMgOFe_iv.vnl', sample_name='FeMgOFe',
    current_k_point_sampling = (10,10),
    current_number_of_points = 100
  )

4. Change the voltages to whatever range you prefer, as well as the file and sample name.

Notes:

In your original script, the k-point sampling for the current calculation was 1x1. A characteristic feature of the FeMgOFe system is that you need very, very many k-points to get an accurate current for the minority transmission in the parallel electrode confiuration. I have used 10x10 in the example code above, but you should really check this number, and find out how many you need. Just do it for the zero bias case, that should be fine.

However, and this is crucial: the way you have set up the system, it is not spin-polarized. For this you must add an initial spin to the electron density parameters. For more information, see the tutorial on FeMgO.

Another point is that this calculation will be extremely time-consuming. However, by running it on a cluster you can probably save a factor of 10 or more in time. In that case, you may have to take some special care about the path and import statements etc. However, I note that you are using Windows, and currently there is not parallel Windows version of ATK. So, I really hope you have opportunity to access a parallel Linux machine, or you may have to wait a week or so for the results of this calculation!

[zdhlover]

First,thank you for  Anders Blom's timely and enthusiastic reply. I have a little trick can enjoy with others:you can put the script-----  ivcurve.py  in the Atomistix\atk\lib\site-packages files, then you can site and use it in each script.But if you want to ues the ivcurve.py in the liux system ,you hve to copy the script with your script in the same file.

[artingchen]

I edit the script according to the what  you describe. once the job run,  an import error which is NO module named ivcurve will come out. what's the problem? could you help me?
my script is edited as follows:
from ATK.TwoProbe import *
from ATK.MPI import processIsMaster
import numpy

# Generate time stamp
if processIsMaster():
    import platform, time
    print '#',time.ctime()
    print '#',platform.node(),platform.platform()+'\n'

# Opening vnlfile
if processIsMaster(): file = VNLFile('C:/Users/617/lih2li.vnl')

# Scattering elements and coordinates
scattering_elements = [lithium,  Lithium,  Lithium,  Hydrogen,
                       Hydrogen, Lithium,  Lithium,  Lithium]
scattering_coordinates = [[  0.   ,   0.   ,   0.   ],
                          [  0.   ,   0.   ,   2.9  ],
                          [  0.   ,   0.   ,   5.8  ],
                          [  0.   ,   0.   ,   8.166],
                          [  0.   ,   0.   ,   8.97 ],
                          [  0.   ,   0.   ,  11.336],
                          [  0.   ,   0.   ,  14.236],
                          [  0.   ,   0.   ,  17.136]]*Angstrom
       

electrode_elements = [lithium, Lithium, Lithium, Lithium]
electrode_coordinates = [[ 4.35,  4.35,  0.  ],
                         [ 4.35,  4.35,  2.9 ],
                         [ 4.35,  4.35,  5.8 ],
                         [ 4.35,  4.35,  8.7 ]]*Angstrom

electrode_cell = [[  8.7,   0. ,   0. ],
                  [  0. ,   8.7,   0. ],
                  [  0. ,   0. ,  11.6]]*Angstrom

# Set up electrodes
electrode_configuration = PeriodicAtomConfiguration(
    electrode_cell,
    electrode_elements,
    electrode_coordinates
    )

# Set up two-probe configuration
twoprobe_configuration = TwoProbeConfiguration(
    (electrode_configuration,electrode_configuration),
    scattering_elements,
    scattering_coordinates,
    electrode_repetitions=[[1,1],[1,1]],
    equivalent_atoms=([0,0],[3,7])
    )
if processIsMaster(): nlPrint(twoprobe_configuration)
if processIsMaster(): file.addToSample(twoprobe_configuration, 'lih2li')

######################################################################
# Central region parameters
######################################################################
exchange_correlation_type = LDA.PZ

iteration_mixing_parameters = iterationMixingParameters(
    algorithm = IterationMixing.Pulay,
    diagonal_mixing_parameter = 0.1,
    quantity = IterationMixing.Hamiltonian,
    history_steps = 6
)

electron_density_parameters = electronDensityParameters(
    mesh_cutoff = 150.0*Rydberg
)

basis_set_parameters = basisSetParameters(
    type = DoubleZetaPolarized,
    radial_sampling_dr = 0.001*Bohr,
    energy_shift = 0.01*Rydberg,
    delta_rinn = 0.8,
    v0 = 40.0*Rydberg,
    charge = 0.0,
    split_norm = 0.15
)

iteration_control_parameters = iterationControlParameters(
    tolerance = 1e-005,
    criterion = IterationControl.TotalEnergy,
    max_steps = 100
)

electrode_voltages = (0.0,0.0)*Volt

two_probe_algorithm_parameters = twoProbeAlgorithmParameters(
    electrode_constraint = ElectrodeConstraints.Off,
    initial_density_type = InitialDensityType.EquivalentBulk
)

energy_contour_integral_parameters = energyContourIntegralParameters(
    circle_points = 30,
    integral_lower_bound = 3*Rydberg,
    fermi_line_points = 10,
    fermi_function_poles = 4,
    real_axis_infinitesimal = 0.01*electronVolt,
    real_axis_point_density = 0.02*electronVolt
)

two_center_integral_parameters = twoCenterIntegralParameters(
    cutoff = 2500.0*Rydberg,
    points = 1024
)

######################################################################
# Left electrode parameters
######################################################################
left_electrode_electron_density_parameters = electronDensityParameters(
    mesh_cutoff = 150.0*Rydberg
)

left_electrode_iteration_control_parameters = iterationControlParameters(
    tolerance = 1e-005,
    criterion = IterationControl.TotalEnergy,
    max_steps = 100
)

left_electrode_brillouin_zone_integration_parameters = brillouinZoneIntegrationParameters(
    monkhorst_pack_parameters = (1, 1, 500)
)

left_electrode_iteration_mixing_parameters = iterationMixingParameters(
    algorithm = IterationMixing.Pulay,
    diagonal_mixing_parameter = 0.1,
    quantity = IterationMixing.Hamiltonian,
    history_steps = 6
)

left_electrode_eigenstate_occupation_parameters = eigenstateOccupationParameters(
    temperature = 300.0*Kelvin
)

######################################################################
# Collect left electrode parameters
######################################################################
left_electrode_parameters = ElectrodeParameters(
    brillouin_zone_integration_parameters = left_electrode_brillouin_zone_integration_parameters,
    electron_density_parameters = left_electrode_electron_density_parameters,
    eigenstate_occupation_parameters = left_electrode_eigenstate_occupation_parameters,
    iteration_mixing_parameters = left_electrode_iteration_mixing_parameters,
    iteration_control_parameters = left_electrode_iteration_control_parameters
)

######################################################################
# Right electrode parameters
######################################################################
right_electrode_electron_density_parameters = electronDensityParameters(
    mesh_cutoff = 150.0*Rydberg
)

right_electrode_iteration_control_parameters = iterationControlParameters(
    tolerance = 1e-005,
    criterion = IterationControl.TotalEnergy,
    max_steps = 100
)

right_electrode_brillouin_zone_integration_parameters = brillouinZoneIntegrationParameters(
    monkhorst_pack_parameters = (1, 1, 500)
)

right_electrode_iteration_mixing_parameters = iterationMixingParameters(
    algorithm = IterationMixing.Pulay,
    diagonal_mixing_parameter = 0.1,
    quantity = IterationMixing.Hamiltonian,
    history_steps = 6
)

right_electrode_eigenstate_occupation_parameters = eigenstateOccupationParameters(
    temperature = 300.0*Kelvin
)

######################################################################
# Collect right electrode parameters
######################################################################
right_electrode_parameters = ElectrodeParameters(
    brillouin_zone_integration_parameters = right_electrode_brillouin_zone_integration_parameters,
    electron_density_parameters = right_electrode_electron_density_parameters,
    eigenstate_occupation_parameters = right_electrode_eigenstate_occupation_parameters,
    iteration_mixing_parameters = right_electrode_iteration_mixing_parameters,
    iteration_control_parameters = right_electrode_iteration_control_parameters
)

######################################################################
# Initialize self-consistent field calculation
######################################################################
two_probe_method = TwoProbeMethod(
    electrode_parameters = (left_electrode_parameters,right_electrode_parameters),
    exchange_correlation_type = exchange_correlation_type,
    iteration_mixing_parameters = iteration_mixing_parameters,
    electron_density_parameters = electron_density_parameters,
    basis_set_parameters = basis_set_parameters,
    iteration_control_parameters = iteration_control_parameters,
    energy_contour_integral_parameters = energy_contour_integral_parameters,
    two_center_integral_parameters = two_center_integral_parameters,
    electrode_voltages = electrode_voltages,
    algorithm_parameters = two_probe_algorithm_parameters
)
if processIsMaster(): nlPrint(two_probe_method)

runtime_parameters = runtimeParameters(
    verbosity_level = 10,
    checkpoint_filename = 'C:/Users/617/lih2li_iv.nc'
)

voltages=[0.0,0.1,0.2,0.3]*Volt

ivcurve.runIVcurve (
    twoprobe_configuration,
    two_probe_method,
    runtime_parameters,
    voltages,
    vnl_filename='lih2li.vnl', sample_name='lih2li',
    current_k_point_sampling = (1,1),
    current_number_of_points = 100
  )

[Anders Blom]

You should bring back the statement import ivcurve.

But, before you run the calculation, copy the file ivcurve.py to C:\Program Files\QuantumWise\Virtual NanoLab 2008.10.0\atk\lib\site-packages (you probably have to create the directory first).

Note that running an I-V sweep like this is very time-consuming. The recommended way to do it, if you have possibility, is to use a 64-bit Linux computer (can be as much as twice faster) or - even better - in parallel.

[artingchen]

thank you very much !

[artingchen]

I  copy the file ivcurve.py to /ms/vnl-2008.10.0/atk/lib/site-packages, when the job run, an error of No module named ivcurve. can you help me ? the VNL is running on linux system.

[Anders Blom]

Ah, I thought you were using Windows, based on the posted script (checkpoint file C:/users...; that filename will not work on Linux, by the way).

On Linux, put the file in /ms/vnl-2008.10.0/atk/lib/python2.4/site-packages instead.

[artingchen]

I want to ask that if the  transmission spectrum T(E) itself depends strongly on the applied bias, can I use this method to calculate the I-V characteristics.

[Nordland]

In general the transmission T(E) depends strongly on the applied bias, however if the applied voltages is small, then the variation is small.
In some systems - most typical molecular junctions - the applied bias can be as much as 1 V - 2V before the error becomes deviation critical.

However linear response has some limitations:
--> The linear response current can never decrease with the applied bias.
--> The linear response current has a tendency to smear out small features.

However with that said, when I calculate an IV curve I do it with a couple of current calculation where it is self-consistent calcualted, and
then I use the linear response current to calculate alot of points in between the self-consistent sampling.

[Anders Blom]

Just to make sure we all understand each other:

If you use the main method in ATK, that is you run a self-consistent calculation for each bias, then the full dependence of the transmission on the bias will be properly accounted for. The I-V curve sweep used in the ivcurve.py module does this.

[artingchen]

hello,
I have eidted the script with this method in VNL, but  the job run with ATK using the code: [ms@localhost ~]$ mpiexec -n 4 /home/ms/vnl-2008.10.0/atk/bin/atk    /home/ms/b-pcur/lih.py
, during the calculation,  an error comes  out as follows:

Traceback (most recent call last):
  File "<string>", line 218, in ?
  File "/home/ms/vnl-2008.10.0/atk/lib/python2.4/site-packages/ivcurve.py", line 79, in runIVcurve
    f.addToSample(current,sample_name,'Current at %s V bias' % voltage.inUnitsOf(Volt))
OSError: [Errno 2] No such file or directory: 'lih.vnl'

what is the problem?
please help me to solve it. thank you !