Run time of ldos

[esp]

OK, I am posting this with full details in hopes that you can help:   I have been having a hard time getting ldos to finish ... (it runs forever):   I am probably doing something wrong, but please take a look.  I have no errors, and my logs are updated so I know it is running, but no matter what I do it does not finish! .. I tried parallelization with a 24 core machine, and my latest attempt is a simple single core, running just 1 energy point.  It has been running for 8 hours already ... is it almost done?   Last time I ran a parallelized version with -n 12 it ran for 2 days and never completed.  No errors.   Here is the script I am running:  Attached also is the log file

Code

# file:  run.sh
# run single instance
./QuantumWise/atk-11.8.2/atkpython/bin/atkpython ./makeCfg.py
./QuantumWise/atk-11.8.2/atkpython/bin/atkpython ./gnr_12_0_AGNR_z10.py > log1.txt

makeCfg is a script that makes a simple GNR FET type device, attached. parms.py is also attached which has some parameters.   The second file is where the trouble is (code below).  Note that ldos_energies = numpy.linspace(0,0,1)*eV, (I am running only 1 energy point only because I was hoping it would complete within a shorter amount of time, so I could estimate how long an actual calc would take)

Code

def my_DeviceLCAOCalculator():

	# -------------------------------------------------------------
	# Calculator
	# -------------------------------------------------------------

	#----------------------------------------
	# Numerical Accuracy Settings
	#----------------------------------------
	left_electrode_numerical_accuracy_parameters = NumericalAccuracyParameters(
		k_point_sampling=(kPntsX, kPntsY, kPntsZ),
	)

	right_electrode_numerical_accuracy_parameters = NumericalAccuracyParameters(
		k_point_sampling=(kPntsX, kPntsY, kPntsZ),
	)

	#----------------------------------------
	# Iteration Control Settings
	#----------------------------------------
	left_electrode_iteration_control_parameters = IterationControlParameters()
	right_electrode_iteration_control_parameters = IterationControlParameters()
	device_iteration_control_parameters = IterationControlParameters()

	if SCF == 1:

		#----------------------------------------
		# Electrode Calculators
		#----------------------------------------
		left_electrode_calculator = LCAOCalculator(
			numerical_accuracy_parameters=left_electrode_numerical_accuracy_parameters,
			# set these for SCF
			iteration_control_parameters=left_electrode_iteration_control_parameters,
		)

		right_electrode_calculator = LCAOCalculator(
			numerical_accuracy_parameters=right_electrode_numerical_accuracy_parameters,
			# set these for SCF
			iteration_control_parameters=right_electrode_iteration_control_parameters,
		)

		#----------------------------------------
		# Device Calculator
		#----------------------------------------
		calculator = DeviceLCAOCalculator(
			# set these for SCF
			iteration_control_parameters=device_iteration_control_parameters,
			electrode_calculators=
				[left_electrode_calculator, right_electrode_calculator],
		)

	else:
		#----------------------------------------
		# Electrode Calculators
		#----------------------------------------
		left_electrode_calculator = LCAOCalculator(
			numerical_accuracy_parameters=left_electrode_numerical_accuracy_parameters,
		)

		right_electrode_calculator = LCAOCalculator(
			numerical_accuracy_parameters=right_electrode_numerical_accuracy_parameters,
		)

		#----------------------------------------
		# Device Calculator
		#----------------------------------------
		calculator = DeviceLCAOCalculator(
			electrode_calculators=
				[left_electrode_calculator, right_electrode_calculator],
		)

	return calculator


def doLDOS():
	configuration = nlread(outFileName, object_id="cfg")[0]	
	calculator = my_DeviceLCAOCalculator()
	configuration.setCalculator(calculator)			
	configuration.update()

	# Energies at which to compute the LDOS; 
	# arguments to linspace are (min, max, number of points)
	energies = ldos_energies

	# calculate ldos for each energy in the list
	for energy in energies:
		LDOS = LocalDeviceDensityOfStates(
			configuration=configuration,
			energy=energy,
			kpoints=MonkhorstPackGrid(3,3),
			contributions=All,
			energy_zero_parameter=AverageFermiLevel,
			infinitesimal=1e-06*eV,
			self_energy_calculator=KrylovSelfEnergy(),
			spin=Spin.Sum,
		)
		nlsave(outFileName + '_ldos.nc', LDOS, object_id="LDOS at %s" % energy)


if __name__ == "__main__":

	from parms import *

	try:

		if processIsMaster():
			print "Master node"
		else:
			print "Slave node"
							
		doLDOS()
		print "doLDOS done"

       except:
              print "error", sys.exc_info()[0]

[kstokbro]

I tried running the script on my laptop with 11.8.2, it took 7 minutes, I have an i7 dual core.

I changed the script slightly I attach input and output.
Let me know if you can run my script

[esp]

Yes it did complete now.  Thank you. 

[esp]

I cannot get this to converge self-consistently .. it runs for days on end .. i have tried to increase k points, tweak other parms .... cannot converge

[kstokbro]

try to use the following,
IterationControlParameters(number_of_history_steps=20, tolerance=1e-4)

This will be the new default in 12.2, let me know if you system still does not converge, then it will our first example of a non-converging system with these parameters!

[esp]

I am pretty sure I did try this and it ran for like 12 hours without ending .. i will try again tonight (let it run while i sleep)

I will post tomorrow

[esp]

Ok I think I got it to converge now . it is still running a transmission calc, but the self-consistent part is done i think ... it took like 4 hours ... my new parms file is below: Thank you ... i think it was a matter of k-points and other parms but at least for now it is ok ....

Code

# -------------------------------------------------------------
# setup parms
# -------------------------------------------------------------

from NanoLanguage import *
import numpy
import ase

nc_files = "../nc_files/"
plots = "../plots/"

# the number of carbon atoms to replace, every <dopingLvl> atoms
dopingLvl = 8

outFileName = "gnr_12_0_AGNR_z10.nc"
outFileName2 = outFileName + ".ivout.nc"
gate_voltage_list = [-0.2 + 0.05*n for n in range(0,9,1)]*Volt
bias_voltage_list = [-0.2 + 0.05*n for n in range(0,9,1)]*Volt
transSpecEnergies = numpy.linspace(-1,1,100)*eV
kPntsX = 1
kPntsY = 3
kPntsZ = 80
ldos_energies = numpy.linspace(-1.0,1.0,11)*eV
iteration_control_parameters = IterationControlParameters(number_of_history_steps=20, tolerance=1e-4)
numerical_accuracy_parameters = NumericalAccuracyParameters(k_point_sampling=(kPntsX, kPntsY, kPntsZ))

# bandstructure parms
pps=20 	# points_per_segment
bafl=All	# bands_above_fermi_level

SCF = 1 # self consistent or not

[Nordland]

A few points.

• Always have a odd number of k-points in the z-axis, thus you are then guaranteed have included the Gamma-point which is important.
• I haven't seen your structure, but unless it is periodic in the y-axis as well, it is only a waste of CPU cycles to have 3 points in that direction. So if the system is a graphene ribbon, you can reduce the calculation time to 1/3

[esp]

Thank you again, I do need to review that theory so I understand fully but for now I will try what you said. .. thank you ... Question:  wouldn't a GNR ribbon width of say 6 hex rings in the y direction be considered periodic?

[Anders Blom]

It's periodic in Y if there is no vacuum around it in this direction. I think you're trying nanoribbons, right? They are only periodic in one direction, the transport direction (Z).

[esp]

Well I suppose it depends on how much length you define as periodic, which is why I specified something like 6 rings in the y direction .. surely the physical structure is repetitive .. but how long must it be for you to consider it periodic ... am I understanding the issue correctly?  ... yes I understand transport would be in the z direction, but there is still a periodic potential in z and y isnt there?    What I am saying is that, for example, suppose the size of the "ribbon" was not a ribbon really, but equal in z and y directions, say 12 and 12 square .... would you still say it is only periodic in z because that is the transport direction?    Or is it that because we are assumming a two probe system with Vd and Vs potentials, that in a direction y the potential is constant, as far as Vds is concerned?  ..

[Anders Blom]

I guess what you are saying is, how long should the period length be (in y) in order for 1 k-point to be good enough? Probably longer than 6 rings.

We assume no periodicity in Z for transport calculations. If you do, you cannot apply a finite bias.

[esp]

hmm ... not sure if there are typos or i am misunderstanding: 
You are saying 1 k point is not enough then? .. i think you meant to say "shorter than" 6 rings? and then, I think you meant to say "we assume no periodicity in Y", not "Z" above in the last sentence?

...

I think what you are saying is that when considering a nanoribbon, and when using k points like (1,1,100), it should in fact be a ribbon, where z >> y ... is that the right idea?

thank you again

[Anders Blom]

No typos

Here's the thing: we always assume periodic boundary conditions in X and Y - and never in Z.

In X and Y, you can basically have two situations: the structure can be finite (like a nanotube, or a ribbon). Since the model is anyway periodic, you add some vacuum so the periodic copies don't interact. In this case, 1x1x100 would be the appropriate k-point sampling.

If you have periodicity in X and Y, like if the electrodes are formed by a cleaved crystal surface. Then you need NAxNBx100 k-points, and how many NB and NB depends on the length of the period - the shorter the period length, the more k-points. So 1 might be still ok if the period is really long, which usually means you make several repetition of the system.

You can also have a mixture in X and Y, like finite in X and periodic in Y. This would be the same of a graphene sheet in the YZ plane. Then 1xNBx100 would be an appropriate sampling, with the same rules applying for NB as just above.

So in short: for a graphene sheet - 1x9x100 (or 12, or 15, etc). For a ribbon, 1x1x100.