QuantumATK Forum

QuantumATK => General Questions and Answers => Topic started by: postnikov on August 4, 2009, 04:29

Title: The problem about charge being zero in the central region has been solved.
Post by: postnikov on August 4, 2009, 04:29
I have a test calculation using the atk code. I find the charge in the central region become zero even if the  paramter
    integral_lower_bound has been set as 100*Rydberg.

100 Ry should be enough for the central region, however, the charge has also been zero in the two-probe calculatio.
By the way, the circle_points has also increase (100).
Title: Re: charge becomes zero in the central region
Post by: postnikov on August 4, 2009, 04:32
I have a test calculation using the atk code. I find the charge in the central region become zero even if the  paramter
    integral_lower_bound has been set as 100*Rydberg.

100 Ry should be enough for the central region, however, the charge has also been zero in the two-probe calculatio.
By the way, the circle_points has also increase (100).


The structure has been checked carefully and it seems ok.
Title: Re: charge becomes zero in the central region
Post by: zh on August 4, 2009, 06:39
Did you try to increase the size of k-point grid and to use the DZP basis set?
Title: Re: charge becomes zero in the central region
Post by: Anders Blom on August 4, 2009, 09:19
Although it's not a really critical point, I would consider shifting the coordinates in the X direction. If you visualize the system in the Nanoscope in VNL, you will find that half of the atoms are outside the unit cell. ATK will shift them inside, but it can cause some problems, and at any rate it will be a bit tricky to visualize the results later on.

Otherwise the key to this system is most likely indeed the k-point sampling, as zh pointed out. In the X direction, that is. The system is periodic in this direction, and graphene typically requires a relatively large k-point sampling, if not for convergence then at any rate to obtain accurate results.

There is no need to increase the circle to 100 Ry, you will just lose accuracy, even with 100 points. If you fail to converge with a circle of, say, 10 Ry, then the problem lies elsewhere. Since there appears to be enough electrode layers in your setup, I think the k-point sampling will be very helpful. In addition, you can experiment with increasing the temperature.

Increasing the k-point sampling increases the computation time quite a bit if you run in serial. If you have possibility, consider running the calculation in parallel.
Title: Re: charge becomes zero in the central region
Post by: postnikov on August 4, 2009, 09:32
zh and Blom, thanks!

I am performing the  KP 10x10X200  case calculations.

The obtained results are also very bad as follows:

#-------------------------------------------------------------------------------
# ----------------------------------------------------------------
# TwoProbe Calculation
# ----------------------------------------------------------------
# sc  0 : q =  576.00000 e
#-------------------------------------------------------------------------------
# Mulliken Population for sc 0
#-------------------------------------------------------------------------------
48  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
....
63  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
......
89  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
90  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
91  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
92  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
93  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
94  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
95  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
96  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
97  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
98  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
......
124  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
125  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
126  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
127  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
128  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
129  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
130  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
131  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
132  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
133  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
134  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
135  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
136  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
137  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
138  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
........
168  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
169  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
170  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
171  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
172  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
173  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
174  C     q =  4.00000 [ s:  2.000, py:  0.667, pz:  0.667, px:  0.667 ]
......
189  C     q =  7.98406 [ s:  1.994, py:  1.997, pz:  1.997, px:  1.995 ]
190  C     q =  7.67588 [ s:  1.920, py:  1.948, pz:  1.869, px:  1.939 ]
191  C     q =  7.98417 [ s:  1.994, py:  1.998, pz:  1.997, px:  1.995 ]
#-------------------------------------------------------------------------------
# Total Charge =  927.80288
#-------------------------------------------------------------------------------
# sc  2 : q =   -2.32246 e  dRho =  7.9168E+00
#-------------------------------------------------------------------------------
# Mulliken Population for sc 2
#-------------------------------------------------------------------------------
48  C     q = -0.00113 [ s: -0.000, py: -0.000, pz:  0.000, px: -0.001 ]
49  C     q = -0.20335 [ s: -0.056, py: -0.023, pz: -0.087, px: -0.037 ]
50  C     q = -0.00110 [ s: -0.001, py: -0.000, pz:  0.000, px: -0.001 ]
51  C     q = -0.20361 [ s: -0.057, py: -0.023, pz: -0.088, px: -0.037 ]
52  C     q = -0.00110 [ s: -0.001, py: -0.000, pz:  0.000, px: -0.001 ]
53  C     q = -0.20364 [ s: -0.057, py: -0.023, pz: -0.087, px: -0.037 ]
54  C     q = -0.00110 [ s: -0.001, py: -0.000, pz:  0.000, px: -0.001 ]
55  C     q = -0.20364 [ s: -0.057, py: -0.023, pz: -0.087, px: -0.037 ]
56  C     q = -0.00110 [ s: -0.001, py: -0.000, pz:  0.000, px: -0.001 ]
57  C     q = -0.20361 [ s: -0.057, py: -0.023, pz: -0.088, px: -0.037 ]
58  C     q = -0.00113 [ s: -0.000, py: -0.000, pz:  0.000, px: -0.001 ]
59  C     q = -0.20335 [ s: -0.056, py: -0.023, pz: -0.087, px: -0.037 ]
60  C     q = -0.00001 [ s: -0.000, py:  0.000, pz: -0.000, px: -0.000 ]
61  C     q = -0.00017 [ s: -0.000, py: -0.000, pz: -0.000, px: -0.000 ]
62  C     q = -0.00001 [ s: -0.000, py:  0.000, pz: -0.000, px: -0.000 ]
63  C     q = -0.00016 [ s: -0.000, py: -0.000, pz: -0.000, px: -0.000 ]
64  C     q = -0.00001 [ s: -0.000, py:  0.000, pz: -0.000, px: -0.000 ]
65  C     q = -0.00016 [ s: -0.000, py: -0.000, pz: -0.000, px: -0.000 ]
66  C     q = -0.00001 [ s: -0.000, py:  0.000, pz: -0.000, px: -0.000 ]
67  C     q = -0.00016 [ s: -0.000, py: -0.000, pz: -0.000, px: -0.000 ]
68  C     q = -0.00001 [ s: -0.000, py:  0.000, pz: -0.000, px: -0.000 ]
69  C     q = -0.00016 [ s: -0.000, py: -0.000, pz: -0.000, px: -0.000 ]
70  C     q = -0.00001 [ s: -0.000, py:  0.000, pz: -0.000, px: -0.000 ]
71  C     q = -0.00017 [ s: -0.000, py: -0.000, pz: -0.000, px: -0.000 ]
72  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
73  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
74  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
75  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
76  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
77  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
78  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
79  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
80  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
81  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
82  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
83  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
84  C     q =  0.00000 [ s:  0.000, py:  0.000, pz:  0.000, px:  0.000 ]
......
175  C     q = -0.00001 [ s: -0.000, py:  0.000, pz: -0.000, px: -0.000 ]
176  C     q = -0.00014 [ s: -0.000, py: -0.000, pz: -0.000, px: -0.000 ]
177  C     q = -0.00001 [ s: -0.000, py:  0.000, pz: -0.000, px: -0.000 ]
178  C     q = -0.00015 [ s: -0.000, py: -0.000, pz: -0.000, px: -0.000 ]
179  C     q = -0.00001 [ s: -0.000, py:  0.000, pz: -0.000, px: -0.000 ]
180  C     q = -0.18116 [ s: -0.044, py: -0.023, pz: -0.080, px: -0.034 ]
181  C     q = -0.00068 [ s: -0.000, py:  0.000, pz:  0.000, px: -0.001 ]
182  C     q = -0.18142 [ s: -0.044, py: -0.023, pz: -0.080, px: -0.034 ]
183  C     q = -0.00068 [ s: -0.000, py:  0.000, pz:  0.000, px: -0.001 ]
184  C     q = -0.18187 [ s: -0.044, py: -0.023, pz: -0.080, px: -0.034 ]
185  C     q = -0.00067 [ s: -0.000, py:  0.000, pz:  0.000, px: -0.001 ]
186  C     q = -0.18187 [ s: -0.044, py: -0.023, pz: -0.080, px: -0.034 ]
187  C     q = -0.00067 [ s: -0.000, py:  0.000, pz:  0.000, px: -0.001 ]
188  C     q = -0.18142 [ s: -0.044, py: -0.023, pz: -0.080, px: -0.034 ]
189  C     q = -0.00068 [ s: -0.000, py:  0.000, pz:  0.000, px: -0.001 ]
190  C     q = -0.18116 [ s: -0.044, py: -0.023, pz: -0.080, px: -0.034 ]
191  C     q = -0.00068 [ s: -0.000, py:  0.000, pz:  0.000, px: -0.001 ]
#-------------------------------------------------------------------------------
# Total Charge =   -2.32246
#-------------------------------------------------------------------------------
Title: Re: charge becomes zero in the central region
Post by: Anders Blom on August 4, 2009, 09:39
There is no need to use more than 1 k-point in the Y direction, it just makes the calculation run slower, without changing or improving the results in any way!

So, try something like 10x1x200 instead.
Title: Re: charge becomes zero in the central region
Post by: postnikov on August 4, 2009, 09:49
thanks again!

Now, I am trying the 16x1x200,
and the integral_lower_bound has been set as 10*Rydberg, circle_points = 100
Title: Re: charge becomes zero in the central region
Post by: postnikov on August 4, 2009, 12:22
thanks again!

Now, I am trying the 16x1x200,
and the integral_lower_bound has been set as 10*Rydberg, circle_points = 100

in  the case of 16x1x200, the charge variation in the calculation:

# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  576.00000
# Total Charge =  927.83657
# Total Charge =   -2.32268
# Total Charge =   38.14169
# Total Charge = 1150.31016
# Total Charge =  120.60376
# Total Charge =   97.44324
# Total Charge =   -5.98985
# Total Charge =   -5.60166
# Total Charge =   11.28651
# Total Charge =  739.07469
# Total Charge =  889.81866
# Total Charge =  765.27955
# Total Charge =  885.32093
# Total Charge =  215.58614
# Total Charge =  231.24017
# Total Charge =  375.51966
# Total Charge =  550.63754
# Total Charge =  481.29176
# Total Charge =  689.45979
# Total Charge =  195.55278
# Total Charge =  326.48935
# Total Charge =  589.78465
# Total Charge =  485.65247
# Total Charge =  428.85457
# Total Charge =  436.84078
# Total Charge =  377.56452

Title: Re: charge becomes zero in the central region
Post by: postnikov on August 4, 2009, 12:24


in  the case of 10x1x200, the charge variation in the calculation:


# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  192.00000
# Total Charge =  576.00000
# Total Charge =  927.80347
# Total Charge =   -2.32289
# Total Charge =   38.15299
# Total Charge = 1150.30976
# Total Charge =  120.64785
# Total Charge =   97.38514
# Total Charge =   -5.95919
# Total Charge =   -5.59841
# Total Charge =   38.45065
# Total Charge =  768.23992
# Total Charge =  996.36881
# Total Charge =  913.11851
# Total Charge = 1108.56319
# Total Charge =  929.81157
# Total Charge =  345.36823
# Total Charge =  360.73868
# Total Charge =  141.53380
# Total Charge =  610.70247
# Total Charge =  565.37767
# Total Charge =  552.51291
# Total Charge =  581.16821
# Total Charge =  569.63964
# Total Charge =  538.48649
# Total Charge =  548.40624
# Total Charge =  571.78094
# Total Charge =  574.82316
# Total Charge =  577.37777
# Total Charge =  580.78817
# Total Charge =  580.76599
# Total Charge =  581.45682
Title: Re: charge becomes zero in the central region
Post by: Anders Blom on August 4, 2009, 12:53
Try increasing the temperature, and lowering the diagonal mixing to 0.05. Is this SingleZeta? I'd use SingleZetaPolarized.
Title: Re: charge becomes zero in the central region
Post by: zh on August 4, 2009, 13:33
How did you specify the "electrode_parameters" in the "TwoProbeMethod()"? Could you paste it here? If the left and right electrodes are identical to each other, one had better specify the same entity to the "electrode_parameters" in the "TwoProbeMethod()", as shown in follows:
electrode_params = ...
basis_set = ...
initial_density = ...
twoprobe_method = TwoProbeMethod(
    electrode_parameters=(electrode_parms, electrode_parms),
    basis_set_parameters=basis_set,
    electron_density_parameters=initial_density,
    ......
    )

Although the definitions of "left_electrode_params" and "right_electrode_params" in the followings example are identical to each other, it is not recommended to specify the "electrode_parameters" in the "TwoProbeMethod()" in the following way:
left_electrode_params = ...
right_electrode_params=
basis_set = ...
initial_density = ...
twoprobe_method = TwoProbeMethod(
    electrode_parameters=(left_electrode_parms, right_electrode_parms),
    basis_set_parameters=basis_set,
    electron_density_parameters=initial_density,
    ......
    )

.
Title: Re: charge becomes zero in the central region
Post by: ipsecog on August 4, 2009, 13:40
Well, the second way is how VNL produces the scripts, whether the setup is homogeneous or not...!

I doubt changing the setup this way will change the convergence, but if the system is homogeneous (left electrode == right electrode), then the calculation will run faster if you set it up using the same electrode parameters, since it uses FFT rather than the multigrid method for the electrostatics.

This should be the default for homogeneous systems, but it isn't. So, you have to edit the script by hand, after VNL has provided the "template".
Title: Re: charge becomes zero in the central region
Post by: postnikov on August 4, 2009, 14:42
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('sample_tp.vnl')

# Scattering elements and coordinates
scattering_elements = [Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
  .....
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon,
                       Carbon, Carbon, Carbon, Carbon]
scattering_coordinates = [[  7.755,   0.   ,  10.99 ],
                          [  8.46 ,   0.   ,   9.769],
                          [ 10.575,   0.   ,  10.99 ],
                          [  9.87 ,   0.   ,   9.769],
....
                          [ 14.805,   0.   ,  18.316],
                          [ 14.1  ,   0.   ,  17.095],
                          [ 16.215,   0.   ,  18.316],
                          [ 16.92 ,   0.   ,  17.095],
                          [ 19.035,   0.   ,  18.316],
                          [ 18.33 ,   0.   ,  17.095],
                          [  7.755,   0.   ,  20.759],
                          [  8.46 ,   0.   ,  19.538],
                          [ 10.575,   0.   ,  20.759],
                          [  9.87 ,   0.   ,  19.538],
                          [ 11.985,   0.   ,  20.759],
                          [ 12.69 ,   0.   ,  19.538],
                          [ 14.805,   0.   ,  20.759],
                          [ 14.1  ,   0.   ,  19.538],
                          [ 16.215,   0.   ,  20.759],
                          [ 16.92 ,   0.   ,  19.538],
                          [ 19.035,   0.   ,  20.759],
                          [ 18.33 ,   0.   ,  19.538],
                          [  7.755,   0.   ,  23.201],
                          [  8.46 ,   0.   ,  21.98 ],
                          [ 10.575,   0.   ,  23.201],
                          [  9.87 ,   0.   ,  21.98 ],
                          [ 11.985,   0.   ,  23.201],
                          [ 12.69 ,   0.   ,  21.98 ],
                          [ 14.805,   0.   ,  23.201],
                          [ 14.1  ,   0.   ,  21.98 ],
                          [ 16.215,   0.   ,  23.201],
                          [ 16.92 ,   0.   ,  21.98 ],
                          [ 19.035,   0.   ,  23.201],
 ...
                          [ 16.215,   0.   ,  28.085],
                          [ 16.92 ,   0.   ,  26.864],
                          [ 19.035,   0.   ,  28.085],
                          [ 18.33 ,   0.   ,  26.864],
                          [  7.755,   0.   ,  30.527],
                          [  8.46 ,   0.   ,  29.306],
                          [ 10.575,   0.   ,  30.527],
                          [  9.87 ,   0.   ,  29.306],
                          [ 11.985,   0.   ,  30.527],
                          [ 12.69 ,   0.   ,  29.306],
                          [ 14.805,   0.   ,  30.527],
                          [ 14.1  ,   0.   ,  29.306],
                          [ 16.215,   0.   ,  30.527],
                          [ 16.92 ,   0.   ,  29.306],
                          [ 19.035,   0.   ,  30.527],
                          [ 18.33 ,   0.   ,  29.306],
                          [  7.755,   0.   ,  32.97 ],
                          [  8.46 ,   0.   ,  31.748],
                          [ 10.575,   0.   ,  32.97 ],
                          [  9.87 ,   0.   ,  31.748],
                          [ 11.985,   0.   ,  32.97 ],
                          [ 12.69 ,   0.   ,  31.748],
                          [ 14.805,   0.   ,  32.97 ],
                          [ 14.1  ,   0.   ,  31.748],
                          [ 16.215,   0.   ,  32.97 ],
                          [ 16.92 ,   0.   ,  31.748],
                          [ 19.035,   0.   ,  32.97 ],
                          [ 18.33 ,   0.   ,  31.748],
                          [  7.755,   0.   ,  35.412],
                          [  8.46 ,   0.   ,  34.191],
                          [ 10.575,   0.   ,  35.412],
.....
                          [ 16.92 ,   0.   ,  36.633],
                          [ 19.035,   0.   ,  37.854],
                          [ 18.33 ,   0.   ,  36.633]]*Angstrom
        

electrode_elements = [Carbon, Carbon, Carbon, Carbon,
                      Carbon, Carbon, Carbon, Carbon,
...
                      Carbon, Carbon, Carbon, Carbon,
                      Carbon, Carbon, Carbon, Carbon,
                      Carbon, Carbon, Carbon, Carbon]
electrode_coordinates = [[  7.755,   0.   ,   1.221],
                         [  8.46 ,   0.   ,   0.   ],
                         [ 10.575,   0.   ,   1.221],
                         [  9.87 ,   0.   ,   0.   ],
                         [ 11.985,   0.   ,   1.221],
                         [ 12.69 ,   0.   ,   0.   ],
                         [ 14.805,   0.   ,   1.221],
 ....
                         [  7.755,   0.   ,   3.663],
                         [  8.46 ,   0.   ,   2.442],
                         [ 10.575,   0.   ,   3.663],
                         [  9.87 ,   0.   ,   2.442],
                         [ 11.985,   0.   ,   3.663],
                         [ 12.69 ,   0.   ,   2.442],
                         [ 14.805,   0.   ,   3.663],
                         [ 14.1  ,   0.   ,   2.442],
                         [ 16.215,   0.   ,   3.663],
                         [ 16.92 ,   0.   ,   2.442],
                         [ 19.035,   0.   ,   3.663],
                         [ 18.33 ,   0.   ,   2.442],
                         [  7.755,   0.   ,   6.105],
                         [  8.46 ,   0.   ,   4.884],
                         [ 10.575,   0.   ,   6.105],
                         [  9.87 ,   0.   ,   4.884],
                         [ 11.985,   0.   ,   6.105],
                         [ 12.69 ,   0.   ,   4.884],
                         [ 14.805,   0.   ,   6.105],
 ...
                         [  9.87 ,   0.   ,   7.327],
                         [ 11.985,   0.   ,   8.548],
                         [ 12.69 ,   0.   ,   7.327],
                         [ 14.805,   0.   ,   8.548],
                         [ 14.1  ,   0.   ,   7.327],
                         [ 16.215,   0.   ,   8.548],
                         [ 16.92 ,   0.   ,   7.327],
                         [ 19.035,   0.   ,   8.548],
                         [ 18.33 ,   0.   ,   7.327]]*Angstrom

electrode_cell = [[ 12.78774,   0.     ,   0.     ],
                  [  0.     ,  30.     ,   0.     ],
                  [  0.     ,   0.     ,   9.844  ]]*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],[47,143])
    )
if processIsMaster(): nlPrint(twoprobe_configuration)
if processIsMaster(): file.addToSample(twoprobe_configuration, 'sample_tp')

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

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

electron_density_parameters = electronDensityParameters(
    mesh_cutoff = 150.0*Rydberg
)

basis_set_parameters = basisSetParameters(
    type = SingleZetaPolarizedPolarized,
    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.NeutralAtom
)

energy_contour_integral_parameters = energyContourIntegralParameters(
    circle_points = 100,
    integral_lower_bound = 10*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 = (16, 1, 200)
)

left_electrode_iteration_mixing_parameters = iterationMixingParameters(
    algorithm = IterationMixing.Pulay,
    diagonal_mixing_parameter = 0.05,
    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 = (16, 1, 200)
)

right_electrode_iteration_mixing_parameters = iterationMixingParameters(
    algorithm = IterationMixing.Pulay,
    diagonal_mixing_parameter = 0.05,
    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 = None
)

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






This is the used script file in my ongoing calculation, in which the diagonal mixing has been changed into 0.05 and the SingleZetaPolarized orbital has been adopted.

In the previous calculation, the diagonal mixing is 0.1 and the SingleZeta orbital is used.
Title: Re: charge becomes zero in the central region
Post by: postnikov on August 4, 2009, 14:52
It is just a simple test calculation using ATK code.
However, I think it may be a good choice for atk users and developers to check electron loss problem.
I have found there are many posts about such a problem in this forum.

I have also test the case, in which electrode_constraint has been set as ElectrodeConstraints.RealSpaceDensity.
In such a case, everything is Ok.
Title: Re: charge becomes zero in the central region
Post by: Anders Blom on August 4, 2009, 15:26
You are right about the charge loss being a problem, but most of the time the solution is quite simple (increase the k-points sampling or the number of electrode layers).

Let us know when you get a parameter setting that works, it will, as you say, be very useful as an instructive system for other users. :)

I should also say, that when calculating a similar system, we found that you need a really high k-point sampling for the transmission spectrum to get accurate results. Therefore I'm not surprised if many k-points are also needed for the SCF loop.

In general, SingleZeta should be reserved for cases when the system uses so much memory that it's the only option, or if you have already tested a small system with this basis set and obtained accurate results. I wouldn't expect many systems to pass this test, however, the SingleZeta set is simply too limited for any serious work.
Title: Re: charge becomes zero in the central region
Post by: postnikov on August 4, 2009, 15:43
You are right about the charge loss being a problem, but most of the time the solution is quite simple (increase the k-points sampling or the number of electrode layers).

Let us know when you get a parameter setting that works, it will, as you say, be very useful as an instructive system for other users. :)

I should also say, that when calculating a similar system, we found that you need a really high k-point sampling for the transmission spectrum to get accurate results. Therefore I'm not surprised if many k-points are also needed for the SCF loop.

In general, SingleZeta should be reserved for cases when the system uses so much memory that it's the only option, or if you have already tested a small system with this basis set and obtained accurate results. I wouldn't expect many systems to pass this test, however, the SingleZeta set is simply too limited for any serious work.


I am expecting you can give me some advice on this problem, and I am confused on this problem as other atk users.
Therefore, could you summarize the general advice for our atk users to solve this problem?
Title: Re: charge becomes zero in the central region
Post by: Anders Blom on August 4, 2009, 15:55
There are 4 things to start with (in this order), if the calculation converges to zero charge:


If the problems still persist, and you have checked that the system is periodic in the Z direction (if not, turn off the equivalent bulk calculation as use the NeutralAtom "initial density" setting), it's time for some more advanced attempts:


A lot of this is documented in various Forum posts, but it's a good idea to collect it. Also, most of it is described in the tutorials (see http://www.quantumwise.com/publications/tutorials), most notably the one on "Convergence Tricks" ;)
Title: Re: charge becomes zero in the central region
Post by: zh on August 5, 2009, 03:14
Well, the second way is how VNL produces the scripts, whether the setup is homogeneous or not...!

I doubt changing the setup this way will change the convergence, but if the system is homogeneous (left electrode == right electrode), then the calculation will run faster if you set it up using the same electrode parameters, since it uses FFT rather than the multigrid method for the electrostatics.

This should be the default for homogeneous systems, but it isn't. So, you have to edit the script by hand, after VNL has provided the "template".

The second way will lead to the system being treated as a heterogeneous system and then the multigrid method is used to solve the Poisson's equation.  The experience and several failure cases suggest that it's better to do in the first way where the FFT method will be used to solve the Poisson's equation.  In some cases where the system is indeed homogeneous (i.e., the left and right electrodes are identical to each other), these two different ways for the calculation setup can not give identical results, e.g., the detailed shape and characteristics of the transmission spectrum.


Title: Re: charge becomes zero in the central region
Post by: zh on August 5, 2009, 03:32
...
# Set up electrodes
electrode_configuration = PeriodicAtomConfiguration(
    electrode_cell,
    electrode_elements,
    electrode_coordinates
    )

# Set up two-probe configuration
twoprobe_configuration = TwoProbeConfiguration(
    (electrode_configuration,electrode_configuration),
  ....
    )
...
two_probe_method = TwoProbeMethod(
    electrode_parameters = (left_electrode_parameters,right_electrode_parameters),
.......
)
.....

From the definition of "twoprobe_configuration", the left and right electrodes are exactly identical to each other, so it is safe to define the "two_probe_method" in the following way:
two_probe_method = TwoProbeMethod(
    electrode_parameters = (left_electrode_parameters,left_electrode_parameters),
.......
)
or
two_probe_method = TwoProbeMethod(
    electrode_parameters = (right_electrode_parameters,right_electrode_parameters),
.......
)
.
In this way, the FFT method will be used to solve the Poisson equation and it will be more safe to make the self-consistent calculations be quickly converged.

Title: Re: charge becomes zero in the central region
Post by: postnikov on August 5, 2009, 03:38
...
# Set up electrodes
electrode_configuration = PeriodicAtomConfiguration(
    electrode_cell,
    electrode_elements,
    electrode_coordinates
    )

# Set up two-probe configuration
twoprobe_configuration = TwoProbeConfiguration(
    (electrode_configuration,electrode_configuration),
  ....
    )
...
two_probe_method = TwoProbeMethod(
    electrode_parameters = (left_electrode_parameters,right_electrode_parameters),
.......
)
.....

From the definition of "twoprobe_configuration", the left and right electrodes are exactly identical to each other, so it is safe to define the "two_probe_method" in the following way:
two_probe_method = TwoProbeMethod(
    electrode_parameters = (left_electrode_parameters,left_electrode_parameters),
.......
)
or
two_probe_method = TwoProbeMethod(
    electrode_parameters = (right_electrode_parameters,right_electrode_parameters),
.......
)
.
In this way, the FFT method will be used to solve the Poisson equation and it will be more safe to make the self-consistent calculations be quickly converged.



Thanks! The input file has been changed
Title: Re: charge becomes zero in the central region
Post by: zh on August 5, 2009, 03:56
Wish this hint to solve your problem. Actually, I have also encountered this kind of problem before. Sometimes, tuning the k-point grid, basis set or other parameters may not work, although these parameters are indeed important for the safe convergence of self-consistent calculation. 
Title: Re: charge becomes zero in the central region
Post by: postnikov on August 5, 2009, 09:29
Anders Blom, the problem has still not been solved.

1) About the Electrode depth, in this case, the unit cell is four zigzag periods and it should be deep enough;
2) As for k-point sampling, the K-Point (x, y, z) has been varied in the range (1,1,100),(10,1,100),(10,1,200),(16,1,200).
So, I think the K-P should be enough, especially in case of (16,1,200)
3) As to Temperature, I have tried two cases. One is 300K, the other is 500K.
4) For Mixing parameters, the diagonal_mixing_parameter has changed into 0.05, and the history_steps is 6.
5) About contour parameters, the contour radius is set as 6, or 7, or 10 Ry, and the number of points is varied from

50 to 70.
6) About the mesh cut-off, 150 Ry is always adopted. I think it should be enough.
7) As to the FFT method, zh's advice is adopted.
8) The SZP basis has also used.

However, the problem is still not solved.
The charge is always changed into about zero in the 2nd or 3rd step of two-probe part caculation.

It is very strange!!!!

Any more advice ??

Moderator edit: Cleaned up stray bb codes
Title: Re: The problem about charge being zero in the central region has not been solved.
Post by: zh on August 5, 2009, 10:42
If possible, could you pass your script file to me  and let me check it? Did you tune the "electrode_constraint" in the "twoProbeAlgorithmParameters()"? Please see the explanation on the "twoProbeAlgorithmParameters()" in manual here:
http://www.quantumwise.com/documents/manuals/ATK-2008.10/ref.twoprobealgorithmparameters.html (http://www.quantumwise.com/documents/manuals/ATK-2008.10/ref.twoprobealgorithmparameters.html)
Title: Re: The problem about charge being zero in the central region has not been solved.
Post by: Anders Blom on August 5, 2009, 10:49
I agree most parameters appear to be well chosen.

However, the coordinates are still outside the unit cell, which is never a good idea.

As for increasing the temperature, 300 or 500 K is just about the same; normally, it takes 1200 to 2000 to get an effect from this parameter.
Title: Re: The problem about charge being zero in the central region has not been solved.
Post by: postnikov on August 5, 2009, 10:57
I agree most parameters appear to be well chosen.

However, the coordinates are still outside the unit cell, which is never a good idea.

As for increasing the temperature, 300 or 500 K is just about the same; normally, it takes 1200 to 2000 to get an effect from this parameter.


OK, I will move the coordinate origin and test again.
Title: Re: The problem about charge being zero in the central region has not been solved.
Post by: Anders Blom on August 5, 2009, 10:58
One other thing, not really related to the convergence itself, but still worth considering.

In a sense you have a too large system. Since it's periodic in the X direction, you only need to include one period in this direction, that is 4 atom layers. This reduces the number of atoms but a factor of 3, and since most numerical things in ATK scale as N^3 with the number of atoms N, this could mean a reduction in calculation time of an order of magnitude or more.

Perhaps you are preparing for running ribbons etc later, but as long as you're trying to figure the parameters out, using a smaller system will make the calculations go that much faster, thus saving yourself a lot of time in this test phase.

Note that decreasing the cell size requires a large k-point sampling, so perhaps now 30-40 k-points in X will be needed. However,
a) the calculation time only scales linearly in k-points, so you still win a lot, and
b) the calculation is parallelized over the k-point but not over geometry, so replacing atoms by k-points mean a large advantage when running in parallel
Title: Re: The problem about charge being zero in the central region has not been solved.
Post by: postnikov on August 5, 2009, 11:18
One other thing, not really related to the convergence itself, but still worth considering.

In a sense you have a too large system. Since it's periodic in the X direction, you only need to include one period in this direction, that is 4 atom layers. This reduces the number of atoms but a factor of 3, and since most numerical things in ATK scale as N^3 with the number of atoms N, this could mean a reduction in calculation time of an order of magnitude or more.

Perhaps you are preparing for running ribbons etc later, but as long as you're trying to figure the parameters out, using a smaller system will make the calculations go that much faster, thus saving yourself a lot of time in this test phase.

Note that decreasing the cell size requires a large k-point sampling, so perhaps now 30-40 k-points in X will be needed. However,
a) the calculation time only scales linearly in k-points, so you still win a lot, and
b) the calculation is parallelized over the k-point but not over geometry, so replacing atoms by k-points mean a large advantage when running in parallel


I have two questions about this post.

1) As you know, the approx. in the atk code requires the unit cell should be large enough to avoid the next neighboring interaction. Then, in my opinion, the size of unit cell in the x, y, and z direction should be enough. Are my understanding right?

2)Even, in such a test case, the unit cell size in the x direction can be reduced corresponding, however if I want to deal with such a big (or more big ) system as this case, how can I solve such a charge missing problem?

By the way, I will move the coordinate origin as you adivced.

Many thanks for your and zh's help!
Title: Re: The problem about charge being zero in the central region has not been solved.
Post by: Anders Blom on August 5, 2009, 11:22
1) The overlap condition you mention only applies in the Z direction! In XY ATK uses normal periodic boundary conditions, but in Z the special open boundary conditions lead to the requirement of non-overlap beyond the next cell.

2) I doubt the origin of the non-convergence problem is related to the X/Y directions, it's something to do with the transport in Z. So reducing the cell you just save time while figuring out the proper parameters. Moreover, if you later move to ribbons, these are anyway terminated in the X direction too, so you will then only need 1 k-point there too.

Title: Re: The problem about charge being zero in the central region has not been solved.
Post by: postnikov on August 6, 2009, 04:17
If possible, could you pass your script file to me  and let me check it? Did you tune the "electrode_constraint" in the "twoProbeAlgorithmParameters()"? Please see the explanation on the "twoProbeAlgorithmParameters()" in manual here:
http://www.quantumwise.com/documents/manuals/ATK-2008.10/ref.twoprobealgorithmparameters.html (http://www.quantumwise.com/documents/manuals/ATK-2008.10/ref.twoprobealgorithmparameters.html)


Hi, ZH, do you check it? What is your checking result?
Title: Re: The problem about charge being zero in the central region has not been solved.
Post by: zh on August 6, 2009, 04:48
Sorry, I did not get your input script file.
Title: Re: The problem about charge being zero in the central region has not been solved.
Post by: postnikov on August 6, 2009, 08:10
Sorry, I did not get your input script file.
Hi, Zhu feng, the script has been sent to you.
Many thanks for your help!
Title: Re: The problem about charge being zero in the central region has not been solved.
Post by: postnikov on August 7, 2009, 05:27

  After the coordinates have been shifted, all of atoms are placed in unit cell.
Everything is OK. The error should be come from the atk coordinate-shifting.

 Zhufeng and Blom, many thanks for your help!
Title: Re: The problem about charge being zero in the central region has been solved.
Post by: jjhskang on September 5, 2009, 02:28
Hi!

I tried many possible combinations of (N_k points, mesh cut_off, LDA/GGA, Temperature, Integral lower bound), al of which fail to avoid the disappearance of charge. As I describe before, only constrained SCF calculation for the real space density gives nonzero chage! You'd better tryo to fix this problem.
Title: Re: The problem about charge being zero in the central region has been solved.
Post by: Anders Blom on September 5, 2009, 11:40
The post above refers to http://quantumwise.com/forum/index.php?topic=354.msg2037#msg2037 I suppose. We'll have a look at the Au 111 wire system and see if we can find a working parameter set.