[Pw_forum] how to choose bds in PWCOND?

2010-10-01 Thread Manoj Srivastava
Dear All, 
 I am confused about choosing bds parameter in the transmission
calculation of PWCOND. My question is about example on Smogunov's website
Monatomic Ni wire with a spin reversal, which can be accessed at
http://people.sissa.it/~smogunov/ 
 The manual says- right boundary of the scatt. reg. (left one is at 0 if
prefixs is used and = bdl if prefixt is used) Now in this example bds is
chosen to be at 1.0, but if I let the system relax and atom at z=1.0 moves
away from z=1, would bds still be equal to 1.0? I am guessing so, but I am
not sure. Would someone mind explaining this?

-Manoj Srivastava 
Department of Physics, 
University of Florida, 
Gainesville, FL



[Pw_forum] How to choose the energy0 value while doing a Ballistic conductance using the code pwcond.x

2010-09-07 Thread Manoj Srivastava
Hello, 
 I just wanted to stress the fact that when you use energy0=0, your
calculation is being done at Fermi energy. You can look up the
code in do_cond.f90 for more details. nenergy is number of energy steps,
denergy is energy of each step. Again from  do_cond.f90 

!   the array of energies is automatically formed
DO ien = 1, nenergy
earr(ien) = energy0 + (ien-1)*denergy
tran_tot(ien) = 0.d0
ENDDO

HTH 
Manoj Srivastava 
Department of Physics, 
University of Florida 

On Tue, 7 Sep 2010, mohsen
modaresi wrote:

> Dear Dimpy Sharma,
> pwcond.x give the Transmission probablity as a function of energy. energy0
> determine the minimum of energy in the energy axis.
> 
> Mohsen Modarresi
> Ferdowsi University
> 
> 
> 
> On Tue, Sep 7, 2010 at 10:50 PM, Dimpy Sharma  tyndall.ie>wrote:
> 
> >
> > Hi there,
> >
> > My warm greetings to all of you. I would like to know how do we choose the
> > energy0 (inital energy ) while doing a ballistic conductance using
> > pwcond.x.As <http://pwcond.x.as/> from the example and manual it mentions
> > as the inital energy,which energy does it refer to?is this the total energy
> > of the system ?I could not follow the term nenergy as well.Any advice would
> > be appreciated.
> >
> > Thanks and regards
> >
> > Dimpy
> > Dimpy Sharma
> > ETG-Group
> > Tyndall National Institute
> > Cork
> >
> >
> > ___
> > Pw_forum mailing list
> > Pw_forum at pwscf.org
> > http://www.democritos.it/mailman/listinfo/pw_forum
> >
> >
> 



[Pw_forum] Left and right going Bloch's states from PWCOND

2009-10-30 Thread Manoj Srivastava
Just to add one more question. 
what is 'lorb' switch. I get following- 
  lorb,  &  ! if .t. calculate the scattering (or Bloch) states
so, i guess for transmission calculation i will have to set it to
't'. now- 
IF (lorb.and.ikind.eq.2) call errore('do_cond','lorb not working with
ikind = 2',1)
so, it seems that calculation for r and t for right and left going
state, when left and right leads are different doesnot work. Is that
correct? A possible work around may be rotating the whole system by 180
degree, which would be like treating left lead as right lead and
scattering region rotated. Can I do that?

Regards, 
Manoj

On Fri, 30 Oct 2009, Manoj Srivastava wrote:

> Dear Alexander,
> Thank you very much! I read the code, it seems that there are
> quiet a few changes in 4.1.1 compared to 4.0.4 :) .So, just to make sure,
> if I want to calculate r and t for left going Bloch's state as well as
> right going Bloch's state, I need to run the code two times, one time with
> 'left_to_right' switch true, and another time false. I dont need to change
> scattering region or interchange left and right leads. 
> 
> Regards, 
> Manoj Srivastava
>  On Thu, 29 Oct 2009, Alexander
> Smogunov wrote:
> 
> > Dear Manoj
> > 
> > On Wed, 2009-10-28 at 01:49 -0400, Manoj Srivastava wrote:
> > > Dear Alexander, 
> > >  You are correct, left and right moving states are in general not
> > > related. I was wrong about trying to create a pair. Thanks for
> > > clarification. 
> > >  I have a few more questions.  
> > >  1. I want to find out the reflection coefficient in the code. In
> > > transmit.f90, transmission coefficient is |tmat(ig,n)|^2, where tmat is
> > > vec1(ntran-n2d-npol*nocrosr+ig,n). (line 206, of version 4.0.5).
> > > Reflection is not explicitly calculated in the code, but it should be
> > > modulus square of vec1(2*n2d+npol*norbs+ig,n) (line 234). Am I correct?
> > 
> > in the 4.1 version there is output of both t and r.
> > 
> > 
> > > 
> > >  2. In the code reflection and transmission coefficients are calculated
> > > for one direction lets say t++ and t+-, which is transmission and
> > > reflection coefficient for Bloch' state moving from - to + , or in other
> > > words, left to right. To calculate transmission and reflection for Bloch's
> > > state going from right to left, what should I do? Should i just rotate my
> > > system 180 degree, that way leads get interchanged and also scattering
> > > area has orientation 180 different from before.
> > > 
> > 
> > again in the 4.1 version you have possibility to calculate t and r for
> > both directions, in the do_cond.f90 the transmit routine is called with 
> > an extra parameter which is true (false) if you want to calculate
> > scattering states propagating from the left (from the right). 
> > 
> > >  3. In the subroutine sunitary.f90, is raux measure of conservation of
> > > flux? If I sum all the transmission and reflecion coefficient for
> > > a particular Bloch's state, I should get 1, and computationally a number
> > > close to 1. How large could it be for S matrix not to be unitary, i m
> > > getting raux more than one e.g. 1.0011  but code still considers smatrix
> > > to be unitary. 
> > 
> > in the case when |R + T - 1| > 1.d-4 for some band, it should print out
> > the value of R+T but continues to run anyway.
> > 
> > All the best,
> > Alexander
> >  
> > 
> > 
> > > 
> > > Regards, 
> > > Manoj Srivastava
> > > University of Florida, Gainesville, FL
> > > 
> > > 
> > > 
> > >  
> > > 
> > > 
> > > On Thu, 22
> > > Oct 2009, Alexander Smogunov wrote:
> > > 
> > > > On Thu, 2009-10-22 at 11:34 -0400, Manoj Srivastava wrote:
> > > > > Dear Alexander, 
> > > > >  Thanks for your answer. I just want to make sure. Imagine we have 
> > > > > total
> > > > > number of channels in the left lead 2, so total number of Bloch's 
> > > > > state
> > > > > are 4. 2 of them left going say a and b, and 2 right going say c and 
> > > > > d.
> > > > > So, are you saying that for left going state a, the corresponding 
> > > > > right
> > > > > going state is c? Are they ordered this way?
> > > > 
> > > > what do you mean by corresponding? Left and right moving Bloch
> > > > states are in general not related one to another, you can even have
> &g

[Pw_forum] Left and right going Bloch's states from PWCOND

2009-10-28 Thread Manoj Srivastava
Dear Alexander, 
 You are correct, left and right moving states are in general not
related. I was wrong about trying to create a pair. Thanks for
clarification. 
 I have a few more questions.  
 1. I want to find out the reflection coefficient in the code. In
transmit.f90, transmission coefficient is |tmat(ig,n)|^2, where tmat is
vec1(ntran-n2d-npol*nocrosr+ig,n). (line 206, of version 4.0.5).
Reflection is not explicitly calculated in the code, but it should be
modulus square of vec1(2*n2d+npol*norbs+ig,n) (line 234). Am I correct?

 2. In the code reflection and transmission coefficients are calculated
for one direction lets say t++ and t+-, which is transmission and
reflection coefficient for Bloch' state moving from - to + , or in other
words, left to right. To calculate transmission and reflection for Bloch's
state going from right to left, what should I do? Should i just rotate my
system 180 degree, that way leads get interchanged and also scattering
area has orientation 180 different from before.

 3. In the subroutine sunitary.f90, is raux measure of conservation of
flux? If I sum all the transmission and reflecion coefficient for
a particular Bloch's state, I should get 1, and computationally a number
close to 1. How large could it be for S matrix not to be unitary, i m
getting raux more than one e.g. 1.0011  but code still considers smatrix
to be unitary. 

Regards, 
Manoj Srivastava
University of Florida, Gainesville, FL



 


On Thu, 22
Oct 2009, Alexander Smogunov wrote:

> On Thu, 2009-10-22 at 11:34 -0400, Manoj Srivastava wrote:
> > Dear Alexander, 
> >  Thanks for your answer. I just want to make sure. Imagine we have total
> > number of channels in the left lead 2, so total number of Bloch's state
> > are 4. 2 of them left going say a and b, and 2 right going say c and d.
> > So, are you saying that for left going state a, the corresponding right
> > going state is c? Are they ordered this way?
> 
> what do you mean by corresponding? Left and right moving Bloch
> states are in general not related one to another, you can even have
> different number of them ... Only if you have some symmetry S which
> brings kz to -kz conserving k_parallel, then the state with 
> \psi_{-kz} will be S \psi{kz}. This is true for example at 2D G point
> when you have time reversal operation.
> 
> Now the code simply arranges the propagating states in the order of 
> increasing |k_z|...
> 
> Regards, Alexander
>  
> 
> > 
> > Regards, 
> > Manoj
> > 
> >  
> >  On Thu, 22 Oct 2009, Alexander
> > Smogunov wrote:
> > 
> > > Dear Manoj.
> > > 
> > > The output of complex k vectors is performed in
> > > summary_band.f90 routine. If you want to see all 
> > > the complex k vectors, not only propagating ones, 
> > > you can change at the end of this routine:
> > > 
> > > ---
> > >   do i = 1, nchanl
> > > WRITE( stdout,'(3f12.7)')  DBLE(kvall(i)), AIMAG(kvall(i)), eev
> > >   enddo
> > > ---
> > > 
> > > to
> > > ---
> > > 
> > >   do i = 1, 2*nstl
> > > WRITE( stdout,'(3f12.7)')  DBLE(kvall(i)), AIMAG(kvall(i)), eev
> > >   enddo
> > > ---
> > > 
> > > Altogether there are 2*nstl (or 2*nstr) Bloch states in the left 
> > > (or right) lead. First half, [1,nstl], are propagating or decaying to
> > > the right states, another half, [nstl+1,2*nstl], - propagating or
> > > decaying to the left. In each group, first nchanl states are propagating
> > > states.
> > > 
> > > The propagating states are normalized by the current and are arranged in
> > > the above order at the end of jbloch.f90 routine, after the following
> > > lines:
> > > 
> > > !
> > > ! Right ordering (+, >, -, <)
> > > !
> > > 
> > > 
> > > 
> > > Notice, that in the last versions the code gives in output
> > > both propagating to the right and to the left states.
> > > 
> > > Hope this helps,
> > > Alexander
> > > 
> > >  
> > > 
> > > 
> > > 
> > > 
> > > 
> > > 
> > > On Wed, 2009-10-21 at 14:13 -0400, Manoj Srivastava wrote:
> > > > Dear All, 
> > > >  I am trying to figure out the left and right going Bloch's states in 
> > > > the
> > > > lead from PWCOND. For a given (kx,ky)and energy we get kz. The code only
> > > > prints out Bloch's state moving in one direction. eg. in one of the
> > > > calculation- 
&g

[Pw_forum] Left and right going Bloch's states from PWCOND

2009-10-22 Thread Manoj Srivastava
Dear Alexander, 
 Thanks for your answer. I just want to make sure. Imagine we have total
number of channels in the left lead 2, so total number of Bloch's state
are 4. 2 of them left going say a and b, and 2 right going say c and d.
So, are you saying that for left going state a, the corresponding right
going state is c? Are they ordered this way?

Regards, 
Manoj

 
 On Thu, 22 Oct 2009, Alexander
Smogunov wrote:

> Dear Manoj.
> 
> The output of complex k vectors is performed in
> summary_band.f90 routine. If you want to see all 
> the complex k vectors, not only propagating ones, 
> you can change at the end of this routine:
> 
> ---
>   do i = 1, nchanl
> WRITE( stdout,'(3f12.7)')  DBLE(kvall(i)), AIMAG(kvall(i)), eev
>   enddo
> ---
> 
> to
> ---
> 
>   do i = 1, 2*nstl
> WRITE( stdout,'(3f12.7)')  DBLE(kvall(i)), AIMAG(kvall(i)), eev
>   enddo
> ---
> 
> Altogether there are 2*nstl (or 2*nstr) Bloch states in the left 
> (or right) lead. First half, [1,nstl], are propagating or decaying to
> the right states, another half, [nstl+1,2*nstl], - propagating or
> decaying to the left. In each group, first nchanl states are propagating
> states.
> 
> The propagating states are normalized by the current and are arranged in
> the above order at the end of jbloch.f90 routine, after the following
> lines:
> 
> !
> ! Right ordering (+, >, -, <)
> !
> 
> 
> 
> Notice, that in the last versions the code gives in output
> both propagating to the right and to the left states.
> 
> Hope this helps,
> Alexander
> 
>  
> 
> 
> 
> 
> 
> 
> On Wed, 2009-10-21 at 14:13 -0400, Manoj Srivastava wrote:
> > Dear All, 
> >  I am trying to figure out the left and right going Bloch's states in the
> > lead from PWCOND. For a given (kx,ky)and energy we get kz. The code only
> > prints out Bloch's state moving in one direction. eg. in one of the
> > calculation- 
> > k//=(0.375,-0.375)
> >  Nchannels of the left tip =1
> > k1(2pi/a)k2(2pi/a)   E-Ef (eV)
> > 
> >0.3157801   0.000   0.000
> > 
> > Now if I want Bloch's state moving in right as well as left direction, I
> > can go to kbloch.f90 subroutine, and print out all the eigen values of
> > AX=exp(ikd)BX, and out of those the ones with real solution would be our
> > Bloch's state, so I get for each channel two solutions- 
> >  kval (-0.275409421993275,1.823688001395235E-010)
> >  kval (0.315780119742506,-3.611201785292708E-012)
> > 
> > To figure out the direction, I can calculate current associated with these
> > Bloch's sate and if the current is +ive it is right moving , and if '-'ive
> > its left moving Bloch's state. I can print out current from jbloch.f90
> > subroutine which are - 
> >  current eigenvalue -1.86502143831863  1.59149029314457
> > 
> > So, clearly the first state with kval=-0.2754094 is left moving and the
> > other one right moving. Upto here its clear to me how to identify left and
> > right moving states.
> > 
> > I get confused when for a given (kx,ky,E), I have more than one Bloch'
> > state. In another calculation where i get multiple Bloch's state- 
> >  Nchannels of the left tip =5
> > k1(2pi/a)k2(2pi/a)   E-Ef (eV)
> > 
> >   -0.0746301   0.000   0.000
> >0.1205527   0.000   0.000
> >0.3112908   0.000   0.000
> >0.4200218   0.000   0.000
> >   -0.4935150   0.000   0.000
> > 
> > so i did the same trick i did above to first print out kz and then
> > current, which gives me - 
> >  kval (-0.420023481074359,1.979595081419732E-010)(call it a)
> >  kval (0.420023367986768,2.500979698670295E-011)  (b)
> >  kval (-0.306507431678779,-1.236804629184431E-011)(c)
> >  kval (-0.125376071175573,-6.134512510438736E-011)(d)
> >  kval (-7.945001124706894E-002,6.683546930037856E-011)(e)
> >  kval (0.106554601758169,-6.427946951285107E-011) (f)
> >  kval (8.866867725358024E-002,8.342250371574646E-011) (g)
> >  kval (0.32514672671,1.260810749228185E-011)  (h)
> >  kval (-0.488725859521576,1.769197678346003E-010) (i)
> >  kval (0.479509832763231,1.765499400037283E-010)  (j)
> > 
> > current eigenvalue  -9.31389492882581   -1.24296522993488
> >   -1.21324078359658   -1.11950286753963   -1.08166842367443
> >1.081874821648641.119731465842631.21295295042188
> >1.242800315349409.313897787790
> > 
> &

[Pw_forum] Left and right going Bloch's states from PWCOND

2009-10-21 Thread Manoj Srivastava
Dear All, 
 I am trying to figure out the left and right going Bloch's states in the
lead from PWCOND. For a given (kx,ky)and energy we get kz. The code only
prints out Bloch's state moving in one direction. eg. in one of the
calculation- 
k//=(0.375,-0.375)
 Nchannels of the left tip =1
k1(2pi/a)k2(2pi/a)   E-Ef (eV)

   0.3157801   0.000   0.000

Now if I want Bloch's state moving in right as well as left direction, I
can go to kbloch.f90 subroutine, and print out all the eigen values of
AX=exp(ikd)BX, and out of those the ones with real solution would be our
Bloch's state, so I get for each channel two solutions- 
 kval (-0.275409421993275,1.823688001395235E-010)
 kval (0.315780119742506,-3.611201785292708E-012)

To figure out the direction, I can calculate current associated with these
Bloch's sate and if the current is +ive it is right moving , and if '-'ive
its left moving Bloch's state. I can print out current from jbloch.f90
subroutine which are - 
 current eigenvalue -1.86502143831863  1.59149029314457

So, clearly the first state with kval=-0.2754094 is left moving and the
other one right moving. Upto here its clear to me how to identify left and
right moving states.

I get confused when for a given (kx,ky,E), I have more than one Bloch'
state. In another calculation where i get multiple Bloch's state- 
 Nchannels of the left tip =5
k1(2pi/a)k2(2pi/a)   E-Ef (eV)

  -0.0746301   0.000   0.000
   0.1205527   0.000   0.000
   0.3112908   0.000   0.000
   0.4200218   0.000   0.000
  -0.4935150   0.000   0.000

so i did the same trick i did above to first print out kz and then
current, which gives me - 
 kval (-0.420023481074359,1.979595081419732E-010)(call it a)
 kval (0.420023367986768,2.500979698670295E-011)  (b)
 kval (-0.306507431678779,-1.236804629184431E-011)(c)
 kval (-0.125376071175573,-6.134512510438736E-011)(d)
 kval (-7.945001124706894E-002,6.683546930037856E-011)(e)
 kval (0.106554601758169,-6.427946951285107E-011) (f)
 kval (8.866867725358024E-002,8.342250371574646E-011) (g)
 kval (0.32514672671,1.260810749228185E-011)  (h)
 kval (-0.488725859521576,1.769197678346003E-010) (i)
 kval (0.479509832763231,1.765499400037283E-010)  (j)

current eigenvalue  -9.31389492882581   -1.24296522993488
  -1.21324078359658   -1.11950286753963   -1.08166842367443
   1.081874821648641.119731465842631.21295295042188
   1.242800315349409.313897787790

So, the first 5 are left moving and rest are right moving. But I dont know
the pairs. for example for left moving state a, what is the corresponding
right moving state whether its f or g ... j ?

Any help would be appreciated. 

Regards, 
Manoj Srivastava
University of Florida, Gainesville. 



[Pw_forum] Relaxation of low symmetry lattices

2009-10-09 Thread Manoj Srivastava
Dear Lorenzo and Gabriele, 
 Thank you very much for your suggestions. I will try them out. 

Regards, 
Manoj

On Fri, 9 Oct 2009, Gabriele Sclauzero wrote:

> 
> 
> Manoj Srivastava wrote:
> > Dear users,
> >  Thank you all for quick reply. Through all the replies I am under the
> > impression that the forces can never be actually zero, even for perfect
> > FCC lattice. I dont know the exact reason yet, whether it is just
> > numerical or something bad on my side. So, I would like to give it one
> > more shot, and i have some questions on that.
> > 
> >  Duy Le: I will try using crystal and report to the forum. Although
> > relaxation stops after completing one scf cycle, so atoms actually dont
> > move, but still nonzero forces.
> > 
> >  S. Baroni: When you use cubic symmetry, do you get forces on all atoms
> > zero. Did you do the calculation or is this your intution? also, can one
> > use cubic symmetry for low symmetry unit cell as in my case? 
> 
> As I wrote you in my previous mail, it is not that you use cubic symmetry, 
> but you start 
> from cubic symmetry operations to find the actual symmetries (by removing 
> those which do 
> not hold anymore depending on the lattice and atomic positions).
> You could do the exercise to find out by hand what the symmetry should be in 
> your case and 
> check that the symmetry operations correspond to those found by pw.x
> 
> 
> > 
> >  which is also one order of magnitude smaller than before, so i am
> > guessing there might me some other parameters that are not set up
> > correctly and causing symmetry breaking? Would someone mind to tell me? I
> > dont know much about FFT grid, Isn't FFT dimension determined by the KE
> > cutoff? Is there anyother way of changing it?
> 
> I don't know what the cause could be, anyway, since you seem to be very 
> curious and to 
> have some spare time for testing, I can suggest the following: in forces.f90 
> uncomment the 
> section between #if defined (DEBUG) and #endif and look at the various 
> contributions to 
> the forces (this is present in the CVS version, I don't know if also in the 
> released one). 
> By looking at the single terms you might find why they do not cancel exactly.
> 
> Cheers
> 
> GS
> 
> > 
> > Regards, 
> > Manoj Srivastava
> > University of Florida
> > Gainesville, FL
> > 
> 
> -- 
> 
> 
> o  o
> | Gabriele Sclauzero, PhD Student  |
> | c/o:   SISSA & CNR-INFM Democritos,  |
> |via Beirut 2-4, 34014 Trieste (Italy) |
> | email: sclauzer at sissa.it |
> | phone: +39 040 3787 511  |
> | skype: gurlonotturno |
> o  o
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] Relaxation of low symmetry lattices

2009-10-09 Thread Manoj Srivastava
Dear users,
 Thank you all for quick reply. Through all the replies I am under the
impression that the forces can never be actually zero, even for perfect
FCC lattice. I dont know the exact reason yet, whether it is just
numerical or something bad on my side. So, I would like to give it one
more shot, and i have some questions on that.

 Duy Le: I will try using crystal and report to the forum. Although
relaxation stops after completing one scf cycle, so atoms actually dont
move, but still nonzero forces.

 S. Baroni: When you use cubic symmetry, do you get forces on all atoms
zero. Did you do the calculation or is this your intution? also, can one
use cubic symmetry for low symmetry unit cell as in my case? I increased
the kinetic energy cutoff to a really high value 150Ry !!  forces are -
 atom   1 type  1   force =-0.0074   -0.03150.0310
 atom   2 type  1   force = 0.0084   -0.03910.0354
 atom   3 type  1   force =-0.01280.0470   -0.0401
 atom   4 type  1   force = 0.01110.0044   -0.0083
 atom   5 type  1   force = 0.00080.0193   -0.0181

which is one order of magnitude smaller than before :) but the KE cutoff
is ridiculously high! I also increase k points with keeping KE cutoff 50
Ry, and that gave me -

 atom   1 type  1   force =-0.01030.0086   -0.0038
 atom   2 type  1   force = 0.0535   -0.0155   -0.0025
 atom   3 type  1   force = 0.06430.0039   -0.0244
 atom   4 type  1   force =-0.0603   -0.01410.0310
 atom   5 type  1   force =-0.04710.0171   -0.0003 

 which is also one order of magnitude smaller than before, so i am
guessing there might me some other parameters that are not set up
correctly and causing symmetry breaking? Would someone mind to tell me? I
dont know much about FFT grid, Isn't FFT dimension determined by the KE
cutoff? Is there anyother way of changing it?

Gabriele: Would you mind to tell me how can i get rid of FFT
incompatibility issue?

Andrea: I changed it from e(-8) to e(-10), still no improvement.

Once again, thank you all for your help. 

Regards, 
Manoj Srivastava
University of Florida
Gainesville, FL



On Thu, 8 Oct 2009, Lorenzo Paulatto wrote:

> In data 08 ottobre 2009 alle ore 10:46:13, Gabriele Sclauzero  
>  ha scritto:
> > if you use cubic it should find 2 more, but they are not compatible with  
> > your FFT grid (see at the beginning of the output).
> 
> The code will only use the 2 additional operations if nr3 is a multiple of  
> nr1 and nr2, which is weird as it is, by default, much smaller in your  
> case. Anyway, if you set nr3 equal to nr1 it will use all the symmetry  
> operations.
> 
> For the symmetry experts: could you please check if the code behaviour is  
> correct?
> 
> cheers
> 
> -- 
> Lorenzo Paulatto
> SISSA  &  DEMOCRITOS (Trieste)
> phone: +39 040 3787 511
> skype: paulatz
> www:   http://people.sissa.it/~paulatto/
> 
>  *** save italian brains ***
>   http://saveitalianbrains.wordpress.com/
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] Relaxation of low symmetry lattices

2009-10-07 Thread Manoj Srivastava
Dear PWSCF users, 
 I have a question regarding relaxation of low symmetry lattices. I did
relaxation of 5 atoms FCC unit cell of copper. This is pure copper, no
interface, no defect or anything. The unit cell is just not conventional,
but a different way of choosing bulk FCC. For the default value of
convergence the calculation finishes at the end of 1st scf cycle, but I
see different forces when I describe nosym=.TRUE. or cell_parameters
{cubic}. Following is the input file.


calculation='relax'
pseudo_dir = '/home/manoj/espresso-4.0.4/pseudo',
outdir='./',
prefix='lcu',
wf_collect=.TRUE.
 /
 
ibrav =0,
celldm(1)=6.81650937063832
nat= 5,
ntyp= 1,
ecutwfc =50.0,
occupations='smearing',
smearing='gaussian',
degauss=0.02,
ecutrho=500
 /

conv_thr = 1.0e-8
mixing_beta=0.7
/

ion_dynamics='bfgs'
pot_extrapolation='second_order'
wfc_extrapolation='second_order'
/
ATOMIC_SPECIES
 Cu 63.55  Cu.pz-d-rrkjus.UPF
ATOMIC_POSITIONS
  Cu  0.0   0.0   0.0
  Cu  0.5   0.5   0.0
  Cu  0.0   1.0   0.0
  Cu  1.0   1.0   0.0
  Cu  0.5   1.5   0.0

K_POINTS (automatic)
8 8 8 0 0 0
CELL_PARAMETERS {cubic}
 1.5   0.5  0.0
-0.5   1.5  0.0
 0.0   0.5  0.5
 
Result of relaxation : The covergence is achieved with follwing forces on
the atoms. 
Forces acting on atoms (Ry/au):

   Forces acting on atoms (Ry/au):

 atom   1 type  1   force = 0.0.0.
 atom   2 type  1   force = 0.05670.0185   -0.0172
 atom   3 type  1   force = 0.11550.0377   -0.0372
 atom   4 type  1   force =-0.1155   -0.03770.0372
 atom   5 type  1   force =-0.0567   -0.01850.0172


I realized that I should not have chosen cubic symmetry group, as my unit
cell is not cubic. So, as described in the input file parameters, for the
low symmetry lattices, I specifed nosym=.TRUE., and again the relaxation
stopped after 1 scf cycle, acheiving convergence. But the forces on the
atoms are -  

Forces acting on atoms (Ry/au):

 atom   1 type  1   force = 0.13550.0524   -0.0516
 atom   2 type  1   force =-0.1493   -0.03950.0395
 atom   3 type  1   force = 0.0253   -0.00430.0045
 atom   4 type  1   force = 0.08200.0230   -0.0239
 atom   5 type  1   force =-0.0935   -0.03160.0315

So, still forces are .001 ev/a.u, which i think is higher for a perfect
FCC lattice. I can specify stricter convergence of force, but that would
move atoms around and break symmetry. So, as long as the lattice parameter
is correct, I think for perfect FCC lattice, forces on atoms should be
smaller than what they are presently. I have checked the lattice constant
for bulk FCC with a different calculation, so lattice constant is not an
issue. Another possibility is symmetry group, as you can see with the same
input variables, changing symmetry group changes results.I was wondering
if there is a way i can choose correct symmetry group for this low
symmetry lattice which might give me smaller forces? Are there any other
parameters I need to specify? On a broad view I am curious to know how to
do relaxation of low symmetry lattices.

Regards, 
Manoj Srivastava
University of Florida
Gainesville, USA



[Pw_forum] relaxating only z coordinate

2009-09-24 Thread Manoj Srivastava
Thanks Sagar!

On Thu, 24 Sep 2009, ambavale sagar wrote:

> Dear Manoj

i.e now in your input file the atomic position card will look like:

ATOMIC POSITIONS
Si 000001
Si0.50.50.5001
Si0.70.70.7100

This atomic position card allow you to displace  first two atoms in z direction 
while the third Si in x direction. 


Regards
Sagar Ambavale
The M. S. University of Baroda
India
===


Message: 8
Date: Wed, 23 Sep 2009 17:50:22 -0400 (EDT)
From: Manoj Srivastava <ma...@phys.ufl.edu>
Subject: [Pw_forum] relaxating only z coordinate?
To: PWSCF Forum 
Message-ID:

Content-Type: TEXT/PLAIN; charset=US-ASCII

Dear PWSCF users, 
This might have been answered previously, but I was unable to find
it. How can one relax only z coordinate of atoms inside unit cell using
PW? In other words how can I constrain x and y coordinates of atoms during
relaxation?


Regards, 
Manoj Srivastava
University of Florida 
Gainesville, FL 



--

Message: 9
Date: Thu, 24 Sep 2009 06:29:11 +0800
From: lan haiping <lanhaip...@gmail.com>
Subject: Re: [Pw_forum] relaxating only z coordinate?
To: PWSCF Forum 
Message-ID:

Content-Type: text/plain; charset="iso-8859-1"

Hi, Please check INPUT_PW document,  you will find below description:
if_pos(1), if_pos(2), if_pos(3) INTEGER  *Default:* 1

component i of the force for this atom is multiplied by if_pos(i),
which must be either 0 or 1.  Used to keep selected atoms and/or
selected components fixed in meta-dynamics, neb, smd, MD dynamics or
structural optimization run.

If you want to fix some atom or some direction , just set if_pos(i) =0 .

Regards,

Hai-Ping


  From cricket scores to your friends. Try the Yahoo! India Homepage! 
http://in.yahoo.com/trynew



[Pw_forum] relaxating only z coordinate?

2009-09-24 Thread Manoj Srivastava
Thanks Hai-Ping!

-Manoj

On Thu, 24 Sep 2009, lan haiping wrote:

> Hi, Please check INPUT_PW document,  you will find below description:
> if_pos(1), if_pos(2), if_pos(3) INTEGER  *Default:* 1
> 
> component i of the force for this atom is multiplied by if_pos(i),
> which must be either 0 or 1.  Used to keep selected atoms and/or
> selected components fixed in meta-dynamics, neb, smd, MD dynamics or
> structural optimization run.
> 
> If you want to fix some atom or some direction , just set if_pos(i) =0 .
> 
> Regards,
> 
> Hai-Ping
> 
> On Thu, Sep 24, 2009 at 5:50 AM, Manoj Srivastava  phys.ufl.edu>wrote:
> 
> > Dear PWSCF users,
> >  This might have been answered previously, but I was unable to find
> > it. How can one relax only z coordinate of atoms inside unit cell using
> > PW? In other words how can I constrain x and y coordinates of atoms during
> > relaxation?
> >
> >
> > Regards,
> > Manoj Srivastava
> > University of Florida
> > Gainesville, FL
> >
> > ___
> > Pw_forum mailing list
> > Pw_forum at pwscf.org
> > http://www.democritos.it/mailman/listinfo/pw_forum
> >
> 
> 
> 
> -- 
> Hai-Ping Lan
> Department of Electronics ,
> Peking University , Bejing, 100871
> lanhaiping at gmail.com, hplan at pku.edu.cn
> 



[Pw_forum] relaxating only z coordinate?

2009-09-23 Thread Manoj Srivastava
Dear PWSCF users, 
 This might have been answered previously, but I was unable to find
it. How can one relax only z coordinate of atoms inside unit cell using
PW? In other words how can I constrain x and y coordinates of atoms during
relaxation?
 

Regards, 
Manoj Srivastava
University of Florida 
Gainesville, FL 



[Pw_forum] How to control number of slabs in pwcond calculation?

2009-09-11 Thread Manoj Srivastava
Dear Hui Wang, 
 I dont think there is a straight forward way of doing this, at least not
that I know of. (Community members are most welcome to correct me on
this).changing the number of slabs along z direction is determined using
the unit cell length along z direction. Here is what I found upon some
grep. I tried to show here only when we have only one region present, lets
say left lead.
 from cond_out.f90
 do k=1, nrzl
   write(stdout,'(2x,i3,2x,3f7.4,3x,80i1)')   &
 k,zl(k),zl(k+1),zl(k+1)-zl(k),(crosl(iorb,k),iorb=1,norbl)
 enddo
so if you wish to change the number of slabs along z direction, you need
to change nrzl and so on 
 from init_cond.f90
line 211 
 nrzl = nrzreg(1)
line 86
nrz1=0.0
do iz = 2, 5
 naux=(dwid(iz)+dz1*0.5d0)/dz1-nrz1
 nrzreg(iz-1) = naux
  ...
   
 we have two options, dz1 and dwid(2)
line 44
 nrztot = nr3s
  if(nrztot/2*2.eq.nrztot) nrztot = nrztot+1
  zlen = at(3,3)
  dz1 = zlen/nrztot 
 zlen=z component of a_3 lattice vector. nr3s=
smooth fft dimension in z direction, so if you could change nr3s, you
might be able to change nrzl. 
 now, dwid(2)- line 80
 dwid(2) = bd2
line 60 
 if(flag.eq.'l') then
bd2 = bdl
 from the input file description - 
bdl: right boundary of the left lead (left one is supposed to be at 0)

I dont know how one can change smooth fft dimesion in z direction. I
believe its done in SCF part not in PWCOND. 

Regards, 
Manoj

On Fri, 11 Sep 2009, xirainbow wrote:

> Dear all:  I find the following content at
> "/espresso-4.0.4/examples/example12/results/al.cond.out" file.
> *" nrz  =   21*
> *.*
> *k slabz(k)  z(k+1) crossing(iorb=1,norb)*
> *1   0. 0.0673 0.0673   *
> *.*
> *   20   1.2793 1.3467 0.0673   0111*
> *   21   1.3467 1.4140 0.0673     "*
> 
>  I want to know, at *.cond.in input file, how to control the number of
> slabs along z direction in left, right lead and scattering region,
> respectively.
>  By the way, what is the meaning of "crossing(iorb=1,norb)"?
> 
>  Thanks in advance ??
> 
> 
> 
> Hui Wang
> School of physics, Nankai University, Tianjin, China
> 



[Pw_forum] How to figure out machine accuracy

2009-09-04 Thread Manoj Srivastava
Dear Lorenzo, 
 Thank you for your reply. I figured out my problem. You are right, it is
not coming from lapack routine. 

Regards, 
Manoj

On Fri, 4 Sep 2009, Lorenzo Paulatto wrote:

> In data 04 settembre 2009 alle ore 15:54:38, Manoj Srivastava  
>  ha scritto:
> >  The input is CMACH, I can specify various input variables for 'CMACH' to
> > get output, but what variable is returned as output, as there is only one
> > in this subroutine?  Would you mind having a look at this- its in
> > espresso4.0.4/flib/dlamch.f
> 
> 
> Dear Manoj,
> QE does not use simple LAPACK routines for most of its algebra.  
> Hamiltonian diagonalization is computed using a custom implementation of  
> the iterative Davidson (or conjugate-gradient) algorithm that only uses  
> LAPACK to diagonalize a reduced hamiltonian during each iteration.
> 
> Said that, you can control the accuracy to which eigenvalues are computed  
> cannot be controlled directly by the user (without modifying the code) and  
> is decreased automatically as the self-consistency approaches convergence.
> 
> The formula used is the following:
>   ethr = MIN( ethr, 0.1D0*dr2 / MAX( 1.D0, nelec ) )
> where dr2 is the estimated scf accuracy.
> As you can see you can increase the accuracy (at convergence) of  
> eigenvalues and eigenvectors by increasing the convergence threshold.
> 
> 
> In a more general frame, I think you could disclose some more information  
> on your "symmetry problem"; personally I doubt it is caused by  
> insufficient diagonalization accuracy.
> 
> best regards
> 
> -- 
> Lorenzo Paulatto
> SISSA  &  DEMOCRITOS (Trieste)
> phone: +39 040 3787 511
> skype: paulatz
> www:   http://people.sissa.it/~paulatto/
> 
>  *** save italian brains ***
>   http://saveitalianbrains.wordpress.com/
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] How to figure out machine accuracy

2009-09-04 Thread Manoj Srivastava
Dear PWSCF users and developers, 
 I am wondering how can one find machine accuracy? I am having
some symmetry problem, and for that I want to solve eigen value
problem with the best accuracy. In the lapack library it says-
ABSTOL  (input) DOUBLE PRECISION
*  The absolute error tolerance for the eigenvalues.
*  An approximate eigenvalue is accepted as converged
*  when it is determined to lie in an interval [a,b]
*  of width less than or equal to
*
*  ABSTOL + EPS *   max( |a|,|b| ) ,
*
*  where EPS is the machine precision.  If ABSTOL is less than
*  or equal to zero, then  EPS*|T|  will be used in its place,
*  where |T| is the 1-norm of the tridiagonal matrix obtained
*  by reducing A to tridiagonal form.
*
*  Eigenvalues will be computed most accurately when ABSTOL is
*  set to twice the underflow threshold 2*DLAMCH('S'), not zero.
*  If this routine returns with INFO>0, indicating that some
*  eigenvectors did not converge, try setting ABSTOL to
*  2*DLAMCH('S').

the last paragraph is useful. so i looked up DLAMCH-
  DOUBLE PRECISION FUNCTION DLAMCH( CMACH ) 
 The input is CMACH, I can specify various input variables for 'CMACH' to
get output, but what variable is returned as output, as there is only one
in this subroutine?  Would you mind having a look at this- its in
espresso4.0.4/flib/dlamch.f

espresso4.0.4 can be replaced by whatever version you are using. 

Regards, 
Manoj Srivastava
University of Florida
Gainesville, USA



[Pw_forum] question about PWCOND

2009-08-16 Thread Manoj Srivastava
Dear Alexander and PWSCF users,
 I have a question about an expression in the scatter_forw.f90 subroutine,
specifically about the following expression - 
 IF (ABS(AIMAG(zk(lam, kp))).LT.eps) THEN
   f1(lam,iorb)=-ezk(lam)*CONJG(di(iorb,lam))*zk2(lam)
   f2(lam,iorb)=-ezk(lam)*CONJG(ci(iorb,lam))*zk2(lam)
  ELSE
   f1(lam,iorb)=-CONJG(ci(iorb,lam))*zk2(lam)
   f2(lam,iorb)=-CONJG(di(iorb,lam))*zk2(lam)
  ENDIF 
 Why does this expression have if/else statement depending on imaginary
part of zk? I did some algebra and realized that f1 for (ABS(AIMAG(zk(lam,
kp))).LT.eps) is related with f1 for (ABS(AIMAG(zk(lam, kp))).GT.eps) by
changing zk to -zk. Same is true for f2. It is defined in previous
subroutine that
  zk=sqrt(E-E_lambda)
where E is given energy and E_lambda is eigenvalue of 2D problem. 
Square root gives us two solution, one would be zk and
another -zk, and I dont understand how do you choose which solution to use
in the expression for f1 and f2, depending on imaginary part of zk. 
Would anyone mind to explain? 

Regards, 
Manoj



[Pw_forum] Request for Pseudopotentials Si and H

2009-08-14 Thread Manoj Srivastava
Hey SaptariBshi, 
 You can find them in /pseudo directory under espresso, e.g. 
espresso-4.0.4/pseudo


-Manoj Srivastava
University of Florida 
Gainesville, FL

On Fri, 14 Aug 2009, Shaptrishi Sharma wrote:

> Hi Quantum espresso users,
> 
> I need the  following pseudopotential file for my simulation  Si.vbc.UPF and
> H.vbc.UPF . Howeevr I could not find them in PWScf, so can anybody suggest
> me from where can I download these pseudopotentials. Please help me!!
> Thanks
> 
> Shaptarishi.
> 



[Pw_forum] generating k point weights

2009-08-10 Thread Manoj Srivastava
Dear All, 
 Thank you very much for all the help :)

Regards,
Manoj Srivastava
University of Florida
Gainesville, FL
On Mon, 10 Aug 2009, Gabriele Sclauzero wrote:

> Just one additional note of caution: kpoints.x reduces the number of k-points 
> (and compute 
> weights) according to the symmetry of the bravais lattice only, while the 
> subroutine 
> kpoint_grid in PW/ used by pw.x takes into account the crystal symmetry 
> (which can be 
> lower than the lattice symmetry if you have more than one atom per cell or 
> non-collinear 
> magnetism)
> 
> HTH
> 
> GS
> 
> Manoj Srivastava wrote:
> > Dear Lorenzo, 
> >  Thanks a lot. I will look into the code. In the mean time I found a paper
> > which talks about generation of special k points and evaluation of their
> > weight. Its by Chadi & Cohen, Phys. Rev. B, 8,5747-5753 (1973). They have
> > some examples in it too.
> > 
> > Regards, 
> > Manoj Srivastava
> > University of Florida, 
> > Gainesville, FL
> > 
> >  On Mon, 3 Aug 2009, Lorenzo Paulatto wrote:
> > 
> >> In data 31 luglio 2009 alle ore 16:25:39, Manoj Srivastava  
> >>  ha scritto:
> >>
> >>> Dear All,
> >>>  I was wondering if someone can tell me, how to generate k point weights
> >>> in the BZ. Reference of a paper would be great!
> >> Dear Manoj,
> >> you can find the algorithm tha tgenerates points and weights in  
> >> pwtools/kpoints.f
> >> The obvious reference is the original article from Monkhorst and Pack:  
> >> Phys. Rev. B 13, 5188 - 5192 (1976), although I don't know if they give an 
> >>  
> >> explicit formulation for the weights there (weights depend on how many  
> >> symmetry-equivalent copies of a point are present in the set, not on the  
> >> generation algorithm).
> >>
> >> cheers
> >>
> >>
> >> -- 
> >> Lorenzo Paulatto
> >> SISSA  &  DEMOCRITOS (Trieste)
> >> phone: +39 040 3787 511
> >> skype: paulatz
> >> www:   http://people.sissa.it/~paulatto/
> >>
> >>  *** save italian brains ***
> >>   http://saveitalianbrains.wordpress.com/
> >> ___
> >> Pw_forum mailing list
> >> Pw_forum at pwscf.org
> >> http://www.democritos.it/mailman/listinfo/pw_forum
> >>
> > 
> > ___
> > Pw_forum mailing list
> > Pw_forum at pwscf.org
> > http://www.democritos.it/mailman/listinfo/pw_forum
> > 
> 
> -- 
> 
> 
> o  o
> | Gabriele Sclauzero, PhD Student  |
> | c/o:   SISSA & CNR-INFM Democritos,  |
> |via Beirut 2-4, 34014 Trieste (Italy) |
> | email: sclauzer at sissa.it |
> | phone: +39 040 3787 511  |
> | skype: gurlonotturno |
> o  o
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] Fermi enery

2009-08-03 Thread Manoj Srivastava
Dear Mahasin, 
 'grep' is very useful to find stuff around, for example if you have a
SCF output file, just grep Fermi.

-manoj

On Mon, 3 Aug 2009, MAHASIN ALAM wrote:

> Dear Friends,
> ?
> Can any one tell me how to get the Fermi Energy from PWSCF (QE) result files. 
> Thank you.
> ?
> Warm regards
> ?
> Mahasin





[Pw_forum] generating k point weights

2009-08-03 Thread Manoj Srivastava
Dear Lorenzo, 
 Thanks a lot. I will look into the code. In the mean time I found a paper
which talks about generation of special k points and evaluation of their
weight. Its by Chadi & Cohen, Phys. Rev. B, 8,5747-5753 (1973). They have
some examples in it too.

Regards, 
Manoj Srivastava
University of Florida, 
Gainesville, FL

 On Mon, 3 Aug 2009, Lorenzo Paulatto wrote:

> In data 31 luglio 2009 alle ore 16:25:39, Manoj Srivastava  
>  ha scritto:
> 
> > Dear All,
> >  I was wondering if someone can tell me, how to generate k point weights
> > in the BZ. Reference of a paper would be great!
> 
> Dear Manoj,
> you can find the algorithm tha tgenerates points and weights in  
> pwtools/kpoints.f
> The obvious reference is the original article from Monkhorst and Pack:  
> Phys. Rev. B 13, 5188 - 5192 (1976), although I don't know if they give an  
> explicit formulation for the weights there (weights depend on how many  
> symmetry-equivalent copies of a point are present in the set, not on the  
> generation algorithm).
> 
> cheers
> 
> 
> -- 
> Lorenzo Paulatto
> SISSA  &  DEMOCRITOS (Trieste)
> phone: +39 040 3787 511
> skype: paulatz
> www:   http://people.sissa.it/~paulatto/
> 
>  *** save italian brains ***
>   http://saveitalianbrains.wordpress.com/
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] generating k point weights

2009-07-31 Thread Manoj Srivastava
Dear All, 
 I was wondering if someone can tell me, how to generate k point weights
in the BZ. Reference of a paper would be great! 

Regards, 
Manoj



[Pw_forum] pwcond-problem with linear system

2009-07-22 Thread Manoj Srivastava
Dear Sagar, 
 The error 'Problems with the linear system' generally gets
removed if you change 'ecut2d' parameter, which is the 2D energy cutoff
and determines how many plane waves would you have in your calculation. I
am guessing that the other problem is also coming from energy convergence.
I suggest to check that.

-Manoj 

On Wed, 22 Jul 2009, ambavale sagar wrote:

> Dear all,
I am running a pwcond program for system made of Al electrodes and vaccum (> 10 
angstrom) in between. Even for such a large vaccum I get nonzero values (10e-1 
to 10e-4) of transmission coefficient (T) at higher energies (> 2 eV +Ef). for 
energies < 2 eV, T become almost zero. After 0.8 eV +Ef energy the program got 
terminated giving error: error # 5411
problems with the linear system. 
What might be a source of error in my case? Can I substract this contribution 
of leads from the T values I got for electrode-molecule-electrode system to get 
accurate T in case of molecule attached to electrodes?
The other question is little sidelined. Why Al is not used as electrode in most 
of the studies? In other words, what are the drawbacks of using Al as an 
electrode?

regards
Sagar Ambavale
PhD student
M.S. University of Baroda
India



  See the Webs breaking stories, chosen by people like you. Check out 
Yahoo! Buzz. http://in.buzz.yahoo.com/



[Pw_forum] transmission calculation

2009-06-04 Thread Manoj Srivastava
Dear PWSCF users and developers,
 I have been trying to do the conductance calculation for the twin
boundary of Cu. I compared my results with previously done calculations
and for most K points the transmission coefficients match up very well.
But there are some k points that I get transmission coefficients more than
1! This is not reasonable, as the maximum value of transmission
coefficient could be 1.I am getting 1.000345, and similar for other k
points. Do you think transmission coefficients being more than 1 is just
numerical error of the code or there is something wrong with my input
file? Attached is the input file for lead and scatterign region -

SCF for the LEFT lead-

calculation='scf'
pseudo_dir = '/home/manoj/espresso-4.0.1/pseudo',
outdir='./',
prefix='lcu',
 /
 
ibrav =0,
celldm(1)=4.82,
nat= 3,
ntyp= 1,
ecutwfc =150.0,
occupations='smearing',
smearing='gaussian',
degauss=0.02,
ecutrho=500
 /
 
conv_thr = 1.0e-8
mixing_beta=0.7
 /
ATOMIC_SPECIES
 Cu 63.55  Cu.pz-d-rrkjus.UPF
ATOMIC_POSITIONS
  Cu  0.0   0.00.0
  Cu -0.5   0.2886751340.81649658
  Cu -1.0   0.5773502691.632993162

K_POINTS (automatic)
6 6 8 1 1 1
CELL_PARAMETERS {hexagonal}
 1.0   0.00.0
-0.5   0.8660254037844386 0.0
 0.0   0.0 2.449489743
>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
SCF for the RIGHT lead-

calculation='scf'
pseudo_dir = '/home/manoj/espresso-4.0.1/pseudo',
outdir='./',
prefix='rcu',
 /
 
ibrav =0,
celldm(1)=4.82,
nat= 3,
ntyp= 1,
ecutwfc =150.0,
occupations='smearing',
smearing='gaussian',
degauss=0.02,
ecutrho=500
 /
 
conv_thr = 1.0e-8
mixing_beta=0.7
 /
ATOMIC_SPECIES
 Cu 63.55  Cu.pz-d-rrkjus.UPF
ATOMIC_POSITIONS
  Cu 0.0   0.00.0
  Cu 0.5  -0.2886751340.81649658
  Cu 1.0  -0.5773502691.632993162

K_POINTS (automatic)
6 6 8 1 1 1
CELL_PARAMETERS {hexagonal}
 1.0  0.00.0
-0.5  0.8660254037844386  0.0
 0.0   0.0 2.449489743

>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
SCF for SCATTERING region-

calculation='scf'
pseudo_dir = '/home/manoj/espresso-4.0.1/pseudo',
outdir='./',
prefix='scu',
 /
 
ibrav =0,
celldm(1)=4.82,
nat= 12,
ntyp= 1,
ecutwfc =150.0,
occupations='smearing',
smearing='gaussian',
degauss=0.02,
ecutrho=500
 /
 
conv_thr = 1.0e-8
mixing_beta=0.7
/
ATOMIC_SPECIES
 Cu 63.55  Cu.pz-d-rrkjus.UPF
ATOMIC_POSITIONS
  Cu  0.0   0.0 0.0
  Cu -0.5   0.288675134 0.816496581
  Cu -1.0   0.577350269 1.632993162
  Cu -0.5   0.288675134 2.449489743
  Cu  0.0   0.0 3.265986324
  Cu  0.5  -0.288675134 4.082482905
  Cu  1.0  -0.577350269 4.898979486
  Cu  1.5  -0.866025404 5.715476066
  Cu  2.0  -1.154700538 6.531972647
  Cu  1.5  -0.866025404 7.348469228
  Cu  1.0  -0.577350269 8.164965809
  Cu  0.5  -0.288675134 8.981462390

K_POINTS (automatic)
6 6 4 1 1 1
CELL_PARAMETERS {hexagonal}
 1.0   0.0 0.0
-0.5   0.866025403 0.0
 0.0   0.0 9.797958971

>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>
input file for conductance part-


outdir='./'
prefixs='scu'
prefixl='lcu'
prefixr='rcu'
tran_file ='trans.twin'
ikind=2
ecut2d=200
energy0=0.0
denergy=0.0
ewind=101.d0
epsproj=1.d-6
nz1=25
bds=3.265986324

 /
  6
  0.500 -0.500   1
  4.2584095E-17  0.500   1
  0.250000.000   1
  0.167  0.167   1
 -1.2775231E-16 -0.250   1
  0.250 -0.251   1
1
0.0


Regards, 
Manoj Srivastava 
University of Florida, 
Gainesville, FL



[Pw_forum] transmission calculation

2009-04-30 Thread Manoj Srivastava
Dear Gabriele,
Thanks for answering. 
On Thu, 30 Apr 2009, Gabriele Sclauzero wrote:

> Manoj Srivastava wrote:
> > Dear PWSCF users and developers,
> >  I have tried to do a transmission calcuation and not getting right
> > results. I have tried to figure out the reason behind this and realized
> > that in my set up for some reason the Bloch waves in left lead get changed
> > resulting wrong transmission. By this I meant that, if I just try to
> > calculate complex band of the left lead, I have a kz for each kx,ky and E
> > but when I am doing transmission calculation, the Bloch's wave in the same
> > lead has different kz for the same kx,ky and E, which does not make sense
> > as why should the Bloch's wave in lead get affected by the presence of
> > scattering region? For the right lead kz remains the same in both cases
> > though, which I believe how it should be. 
> 
> Does it happen also when you do not use G_perp basis reduction? Can you try 
> again the CBS 
> and the transmission calculation without specifying ecut2d, nz1 and putting 
> ewind=1.d3?

Yes, I have tried this. No difference, the kz value for left lead gets
changed, when i switch from ikind=0 to ikind=2.

> 
> Better, first check the size of the reduced basis using
> 
> grep n2d  
> 
> and see if they are very different (or if the second is too low with respect 
> to the size 
> of the whole basis set).
> Also try to understand which of the two results (kz vectors) is correct, by 
> comparing the 
> real part of kz with the band structure from pw.x

The one with ikind=0 is correct.

> 
> 
> GS
> 
> 
> >  I did this for several cases, I interchanged right and left leads,
> > increased number of atoms in the scattering region, chose scattering
> > region in different way, and everytime kz of left lead changes. Following
> > is my input file for the scattering region - 
> > 
> > scf part -
> > 
> > 
> > calculation='scf'
> > pseudo_dir = '/home/manoj/espresso-4.0.1/pseudo',
> > outdir='./',
> > prefix='cus',
> >  /
> >  
> > ibrav =0,
> > celldm(1)=4.82,
> > nat= 10,
> > ntyp= 1,
> > ecutwfc =50.0,
> > occupations='smearing',
> > smearing='gaussian',
> > degauss=0.02,
> > ecutrho=500
> >  /
> >  
> > conv_thr = 1.0e-8
> > mixing_beta=0.7
> > /
> > ATOMIC_SPECIES
> >  Cu 63.55  Cu.pz-d-rrkjus.UPF
> > ATOMIC_POSITIONS
> >   Cu  0.0   0.0 0.0
> >   Cu -0.5   0.288675134 0.816496581
> >   Cu -1.0   0.577350269 1.632993162
> >   Cu -1.5   0.866025404 2.449489743
> >   Cu  0.0   0.0 3.265986324
> >   Cu  0.5  -0.288675134 4.082482905
> >   Cu  1.0  -0.577350269 4.898979486
> >   Cu  1.5  -0.866025404 5.715476066
> >   Cu  1.0  -0.577350269 6.531972647
> >   Cu  0.5  -0.288675134 7.348469228
> >  K_POINTS (automatic)
> > 6 6 4 1 1 1
> > CELL_PARAMETERS {hexagonal}
> >  1.0   0.0 0.0
> > -0.5   0.866025403 0.0
> >  0.0   0.0 8.16496581
> > 
> > 
> > 
> > conductance part -
> > 
> > 
> > outdir='./'
> > prefixl='cul'
> > prefixr='cur'
> > prefixs='cus'
> > tran_file ='trans.twin'
> > ikind=2
> > ecut2d=50
> > energy0=0.0
> > denergy=0.0
> > ewind=5.d0
> > epsproj=1.d-6
> > nz1=25
> >  /
> >   2
> >   0.00   0.00   1
> >   0.500 -0.500   1
> > 
> > 1
> > 0.0
> > 
> > Thank you for your attention. 
> > 
> > Regards, 
> > Manoj Srivastava
> > Graduate Student
> > Department of Physics
> > University of Florida, Gainesville, FL
> > 
> > ___
> > Pw_forum mailing list
> > Pw_forum at pwscf.org
> > http://www.democritos.it/mailman/listinfo/pw_forum
> > 
> 
> -- 
> 
> 
> o  o
> | Gabriele Sclauzero, PhD Student  |
> | c/o:   SISSA & CNR-INFM Democritos,  |
> |via Beirut 2-4, 34014 Trieste (Italy) |
> | email: sclauzer at sissa.it |
> | phone: +39 040 3787 511  |
> | skype: gurlonotturno |
> o  o
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] transmission calculation

2009-04-30 Thread Manoj Srivastava
Dear Lorenzo, 
 Thank you for answering. Layers 4 and 5 are supposed be like that, it is
done on purpose. The file I provided was not for the lead but for the
scattering region. I have visualized my structure, it looks ok.

Regards, 
Manoj

On Thu, 30 Apr 2009, Lorenzo Paulatto wrote:

> In data 30 aprile 2009 alle ore 13:23:25, Manoj Srivastava  
>  ha scritto:
> >  I did this for several cases, I interchanged right and left leads,
> > increased number of atoms in the scattering region, chose scattering
> > region in different way, and everytime kz of left lead changes. Following
> > is my input file for the scattering region -
> 
> I'm no expert in transimission calculation, but I have had a look at your  
> scf input file for the lead region, using XCrysDen, and I've noticed tha  
> tlayers 4 and 5 are aligned in a strange way, reducing the system  
> symmetry. Are you sure you have done it on purpose? Have you tried  
> checking your input files with XCD?
> 
> best regards
> 
> 
> -- 
> Lorenzo Paulatto
> SISSA  &  DEMOCRITOS (Trieste)
> phone: +39 040 3787 511
> skype: paulatz
> www:   http://people.sissa.it/~paulatto/
> 
>  *** save italian brains ***
>   http://saveitalianbrains.wordpress.com/
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] transmission calculation

2009-04-30 Thread Manoj Srivastava

Dear PWSCF users and developers,
 I have tried to do a transmission calcuation and not getting right
results. I have tried to figure out the reason behind this and realized
that in my set up for some reason the Bloch waves in left lead get changed
resulting wrong transmission. By this I meant that, if I just try to
calculate complex band of the left lead, I have a kz for each kx,ky and E
but when I am doing transmission calculation, the Bloch's wave in the same
lead has different kz for the same kx,ky and E, which does not make sense
as why should the Bloch's wave in lead get affected by the presence of
scattering region? For the right lead kz remains the same in both cases
though, which I believe how it should be. 
 I did this for several cases, I interchanged right and left leads,
increased number of atoms in the scattering region, chose scattering
region in different way, and everytime kz of left lead changes. Following
is my input file for the scattering region - 

scf part -


calculation='scf'
pseudo_dir = '/home/manoj/espresso-4.0.1/pseudo',
outdir='./',
prefix='cus',
 /
 
ibrav =0,
celldm(1)=4.82,
nat= 10,
ntyp= 1,
ecutwfc =50.0,
occupations='smearing',
smearing='gaussian',
degauss=0.02,
ecutrho=500
 /
 
conv_thr = 1.0e-8
mixing_beta=0.7
/
ATOMIC_SPECIES
 Cu 63.55  Cu.pz-d-rrkjus.UPF
ATOMIC_POSITIONS
  Cu  0.0   0.0 0.0
  Cu -0.5   0.288675134 0.816496581
  Cu -1.0   0.577350269 1.632993162
  Cu -1.5   0.866025404 2.449489743
  Cu  0.0   0.0 3.265986324
  Cu  0.5  -0.288675134 4.082482905
  Cu  1.0  -0.577350269 4.898979486
  Cu  1.5  -0.866025404 5.715476066
  Cu  1.0  -0.577350269 6.531972647
  Cu  0.5  -0.288675134 7.348469228
 K_POINTS (automatic)
6 6 4 1 1 1
CELL_PARAMETERS {hexagonal}
 1.0   0.0 0.0
-0.5   0.866025403 0.0
 0.0   0.0 8.16496581



conductance part -


outdir='./'
prefixl='cul'
prefixr='cur'
prefixs='cus'
tran_file ='trans.twin'
ikind=2
ecut2d=50
energy0=0.0
denergy=0.0
ewind=5.d0
epsproj=1.d-6
nz1=25
 /
  2
  0.00   0.00   1
  0.500 -0.500   1

1
0.0

Thank you for your attention. 

Regards, 
Manoj Srivastava
Graduate Student
Department of Physics
University of Florida, Gainesville, FL



[Pw_forum] PWCOND

2009-04-22 Thread Manoj Srivastava
Dear Alexander, 
 Thank you for your quick response. One follow up question on this. So,
how do you describe your left and right leads, if they are different? If
the leads are different then I need to do two calculations for the lead,
one for each leads. They are semi-infinite, so periodicity in one is in +z
direction and while in the other it is in -z direction.

Regards, 
Manoj


 On
Wed, 22 Apr 2009, Alexander Smogunov wrote:

> Dear Manoj.
The pwcond code was written assuming that the unit cell
starts at z=0 and goes in positive z direction, the direction
of transport. Quite probable that when you invert the
cell it does not work properly... 
Alexander


?? ??, 21/04/2009 ?? 23:44 -0400, Manoj Srivastava ??:
> Dear PWSCF users, 
>  I have a question regarding real band structure calculated from
> PWCOND. Complex band code requires kx,ky and E as input and gives kz
> as output. I have noticed that the answer for kz depends whether you have
> unit cell vectors as (a1,a2,a3) or (a1,a2,-a3). I did a calculation for
> Cu(001) with 2 atoms tetragonal unit cell, with the following input file - 
> 
> & control
> calculation='scf'
> pseudo_dir = '/home/manoj/espresso-4.0.1/pseudo',
> outdir='./',
> prefix='cu',
>  /
>  
> ibrav = 0,
> celldm(1)=4.7588286373854648
> nat= 2,
> ntyp= 1,
> ecutwfc =150.0,
> occupations='smearing',
> smearing='gaussian',
> degauss=0.02,
> ecutrho=400
>  /
>  
> conv_thr = 1.0e-8
> mixing_beta=0.7
>  /
> ATOMIC_SPECIES
>  Cu 63.55  Cu.pz-d-rrkjus.UPF
> ATOMIC_POSITIONS
>   Cu 0.0 0.0 0.0
>   Cu 0.5 0.5 0.707106781
> K_POINTS (automatic)
> 8 8 8  1 1 1
> CELL_PARAMETERS
> 1.0 0.0 0.0
> 0.0 1.0 0.0
> 0.0 0.0 1.414213562
> 
> 
> outdir='./'
> prefixl='cu'
> band_file ='bands.cu'
> ikind=0
> energy0=2.d0
> denergy=-0.5
> ewind=104.d0
> epsproj=1.d-6
>  /
> 1
> 0.0 0.0 1.0
> 24
> 
> This gives me reasonable answer and it matches up with existing results. 
> Now when I changed a3 to -a3, so that the unit cell now has become 
> (a1,a2,-a3), the answer does not match up with the previous one. 
> Just to give you one concrete example- 
> at Fermi energy for kx,ky (0,0) for unit cell (a1,a2,a3) we get kz=-0.1571
> while for (a1,a2,-a3) we get kz=0.07323. 
> I understand that for two different set of unit cell vectors we should not
> in general expect k_z to match, but this is a very special example, as it
> is just mirror symmetry about xy plane. so shouldn't  kz in one set have
> value -kz in other?
> Thanks for your attention. 
> 
> Regards, 
> Manoj Srivastava
> Department of Physics 
> University of Florida
> Gainesville, FL  
> 
> 
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
-- 
Alexander Smogunov 
Email: smogunov at sissa.it
Home page: http://people.sissa.it/~smogunov


___
Pw_forum mailing list
Pw_forum at pwscf.org
http://www.democritos.it/mailman/listinfo/pw_forum



[Pw_forum] PWCOND

2009-04-22 Thread Manoj Srivastava
Dear PWSCF users, 
 I have a question regarding real band structure calculated from
PWCOND. Complex band code requires kx,ky and E as input and gives kz
as output. I have noticed that the answer for kz depends whether you have
unit cell vectors as (a1,a2,a3) or (a1,a2,-a3). I did a calculation for
Cu(001) with 2 atoms tetragonal unit cell, with the following input file - 

& control
calculation='scf'
pseudo_dir = '/home/manoj/espresso-4.0.1/pseudo',
outdir='./',
prefix='cu',
 /
 
ibrav = 0,
celldm(1)=4.7588286373854648
nat= 2,
ntyp= 1,
ecutwfc =150.0,
occupations='smearing',
smearing='gaussian',
degauss=0.02,
ecutrho=400
 /
 
conv_thr = 1.0e-8
mixing_beta=0.7
 /
ATOMIC_SPECIES
 Cu 63.55  Cu.pz-d-rrkjus.UPF
ATOMIC_POSITIONS
  Cu 0.0 0.0 0.0
  Cu 0.5 0.5 0.707106781
K_POINTS (automatic)
8 8 8  1 1 1
CELL_PARAMETERS
1.0 0.0 0.0
0.0 1.0 0.0
0.0 0.0 1.414213562


outdir='./'
prefixl='cu'
band_file ='bands.cu'
ikind=0
energy0=2.d0
denergy=-0.5
ewind=104.d0
epsproj=1.d-6
 /
1
0.0 0.0 1.0
24

This gives me reasonable answer and it matches up with existing results. 
Now when I changed a3 to -a3, so that the unit cell now has become 
(a1,a2,-a3), the answer does not match up with the previous one. 
Just to give you one concrete example- 
at Fermi energy for kx,ky (0,0) for unit cell (a1,a2,a3) we get kz=-0.1571
while for (a1,a2,-a3) we get kz=0.07323. 
I understand that for two different set of unit cell vectors we should not
in general expect k_z to match, but this is a very special example, as it
is just mirror symmetry about xy plane. so shouldn't  kz in one set have
value -kz in other?
Thanks for your attention. 

Regards, 
Manoj Srivastava
Department of Physics 
University of Florida
Gainesville, FL  




[Pw_forum] input file for transmission

2009-04-14 Thread Manoj Srivastava
Dear Gabriele, 

Thank you very much quick reply. I have  some follow up questions on this. 

On Tue, 14 Apr 2009, Gabriele Sclauzero wrote:

> 
> Manoj Srivastava wrote:
> > Dear PWSCF users and developers, 
> >  I wish to do a transmission calculation and confused about the input
> > file. I have a question on example 12 of the package, where transmission
> > of monoatomic Al wire with a H atom adsorbed on the side is done. The SCF
> > run in the device region is done with (some part of input file is given
> > below)
> 
> it is usually called the "scattering region"
> 
> 
> > 
> > 
> > ibrav = 6,
> > celldm(1) =12.0,
> > celldm(3) =1.875,
> > 
> > where the atomic postions of different atoms is 
> > Al 0. 0. 0.
> > Al 0. 0. 0.375
> > Al-0.02779870 0. .75537515
> > H  0.19269012 0. .9375
> > Al-0.02779870 0. 1.11962485
> > Al 0. 0. 1.5
> > 
> > So, looking at the z coordinate of above system, we notice that device
> > region is periodic with period 1.875. 
> 
> You are right, the scattering region is a periodic system, since pwscf always 
> uses PBCs. 
> The lenght of the scattering region id 12.0*1.875 a.u. and contains all 5 Al 
> and the H 
> impurity.

But that is exactly the trouble I am having in this. The scattering
region should not be periodic like leads, as now we have infinite
scattering region! But PWSCF always has PBC, so we should have a large 3rd
lattice vector to make the scattering region practically finite. I dont
see any super cell here. a_3 is just 1.875*a_0, where a_0 is lattice
constant. The atomic postions are all in a_0 unit, which makes me believe
that it is a practically infinite system with a_3=1.875*a_0. 
 
> 
> > So, physically we are solving for an
> > infinite device region, but in the physical setting of a transmission
> > problem leads are semi-infinite and device is finite. Shouldn't we use
> > some kind of vacum, i.e. taking 3rd lattice vector large, which
> > effectively would represent the finite device region? 
> 
> 
> > Also how much part
> > of the leads should be taken as part of device region, 
> 
> I don't understand this point. The leads are conceptually different thing 
> than the 
> scattering region. The lead is a periodic unit of the "bulk" region (in this 
> case an 
> infinitely long monatomic wire) and it is used to compute the generalized 
> Bloch states, 
> which in turn are propagated in the scattering region. 

 In the above example, in principle we can have one atom H as scattering
region, and Al wire as left and righ leads, but we have taken few Al atoms
with H and treated it as scattering region. Thats what I
meant by how much part of leads should be taken as scattering region.

> 
> > Is there some kind
> > of convergence criterion? Is it like keep increasing part of lead in the
> > device reion till further increase does not substantial change device
> > behavior, e.g. Bloch's state?
> 
> There is a main convergence criterion (though I don't understand if you are 
> actually 
> refering to this). You have to increase the scattering region, adding more Al 
> atoms in the 
> wire, such that the complex band structure with real wave-vectors computed 
> using the 
> leftmost periodic unit of the wire included in the supercell (the H impurity 
> being in the 
> middle of the s.c.) converges to the band structure of an impurity-free wire, 
> obtained for 
> instance from a pwscf calculation (or from a pwcond calculation with a 1 atom 
> cell 
> containing an Al atom).
  
> 
> To do this you can use pwcond with
> ...
> prefixt='prefix of the scattering region'
> bdl=ratio between the lenght of the periodic unit and celldm(1)
> ikind=0
> band_file='name of file containing the CBS'
> ...
> 
> 
> then compare the real bands (contained in .re) with those from 
> pwscf (obtained 
> using the 1 atom cell).
> 
> Also convergence of the transmission with the lenght of the scattering region 
> can be used, 
> but it is quite more cheap to check convergence of CBS (which can also help 
> to understand 
> if everything is going fine), and when the CBS of your "bulk" region (leads) 
> is correctly 
> reproduced the transmission should be converged as well.
> 
> HTH
> 
> GS
> 
> 
> 
> > 
> > Regards, 
> > Manoj Srivastava
> > Ph.D. student
> > Department of Physics
> > University of Florida, Gainesville, FL
> > 
> > 
> > 

[Pw_forum] input file for transmission

2009-04-14 Thread Manoj Srivastava
Dear PWSCF users and developers, 
 I wish to do a transmission calculation and confused about the input
file. I have a question on example 12 of the package, where transmission
of monoatomic Al wire with a H atom adsorbed on the side is done. The SCF
run in the device region is done with (some part of input file is given
below)


ibrav = 6,
celldm(1) =12.0,
celldm(3) =1.875,

where the atomic postions of different atoms is 
Al 0. 0. 0.
Al 0. 0. 0.375
Al-0.02779870 0. .75537515
H  0.19269012 0. .9375
Al-0.02779870 0. 1.11962485
Al 0. 0. 1.5

So, looking at the z coordinate of above system, we notice that device
region is periodic with period 1.875. So, physically we are solving for an
infinite device region, but in the physical setting of a transmission
problem leads are semi-infinite and device is finite. Shouldn't we use
some kind of vacum, i.e. taking 3rd lattice vector large, which
effectively would represent the finite device region? Also how much part
of the leads should be taken as part of device region, Is there some kind
of convergence criterion? Is it like keep increasing part of lead in the
device reion till further increase does not substantial change device
behavior, e.g. Bloch's state?

Regards, 
Manoj Srivastava
Ph.D. student
Department of Physics
University of Florida, Gainesville, FL




[Pw_forum] Transmission calculation using PWCOND

2009-04-04 Thread Manoj Srivastava
Thank you Alexander!

-Manoj
On Fri, 3 Apr 2009, Alexander Smogunov wrote:

> Dear Manoj
> 
> On Wednesday 01 April 2009 22:57, Manoj Srivastava wrote:
> > Dear PWSCF users and developers,
> >  I have been trying to understand the PWCOND code for transmission
> > calculation, basically 'transmit.f90' file. The wavefunction in the left
> > lead, right lead and scattering region is given by -
> >
> >psi_k+sum_{k'}r_{kk'}psi_k'
> >   psi= sum_{n} a_n psi_n+ sum_{\alpha l m} C_{\alpha l m}psi_{\alpha l m}
> >sum_{k'}t_{kk'}psi_k'
> >
> > The above expression is from Choi & Ihm's paper, on which I believe PWCOND
> > is based on. Now undetermined coefficients in above expression are
> > r_{kk'}, a_n, C_{\alpha,l,m}, t_kk'. I am confused on code evaluates these
> > coefficients.  My understanding is you first get the wavefunction in the
> > leads and then solve scattering region the same way to get the
> > wavefunction in scattering region, so at this point you have already
> > determined a_n and C_{\alpha,l,m}. At last do boundary matching condition
> > between leads and scattering region to get reflection and transmission
> > coefficient. Am I right? 
> No, the scattering state (SS)  is completely determined by unknown 
> coefficients which are: {r_kk', a_n, a_alpha, t_kk'}, the first defines the 
> SS in the left lead, the two next in the scattering region, and the last one 
> in the right lead. What you have found before are just basis solutions phi_n 
> and phi_alpha of the Sroedinger eq. in the scattering region over which you 
> expand the SS with coefficients (unknown yet) a_n and a_alpha.
> 
> Hope this helps,
> Alexander
> 
> > My thinking is based on simple quantum mechanics 
> > scattering problem where we solve schrodinger's equation in each region
> > and then do boundary matching to get reflection and transmission
> > coefficients.
> >
> > Now looking at the subroutine 'transmit.f90', this does not appear to be
> > the case. It seems as if we are solving the Schrodinger's equation in the
> > scattering region and doing boundary matching at the same time. Using a
> > big matrix and solving something like AX=B. Would anybody mind explaing
> > what's going on in subroutine transmit.f90?
> >
> > Regards,
> > Manoj Srivastava
> > Physics Graduate Student
> > University Of Florida, Gainesville, FL
> >
> > ___
> > Pw_forum mailing list
> > Pw_forum at pwscf.org
> > http://www.democritos.it/mailman/listinfo/pw_forum
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] Transmission calculation using PWCOND

2009-04-01 Thread Manoj Srivastava

Dear PWSCF users and developers,
 I have been trying to understand the PWCOND code for transmission
calculation, basically 'transmit.f90' file. The wavefunction in the left
lead, right lead and scattering region is given by - 
  
   psi_k+sum_{k'}r_{kk'}psi_k'
  psi= sum_{n} a_n psi_n+ sum_{\alpha l m} C_{\alpha l m}psi_{\alpha l m}
   sum_{k'}t_{kk'}psi_k'

The above expression is from Choi & Ihm's paper, on which I believe PWCOND
is based on. Now undetermined coefficients in above expression are
r_{kk'}, a_n, C_{\alpha,l,m}, t_kk'. I am confused on code evaluates these
coefficients.  My understanding is you first get the wavefunction in the
leads and then solve scattering region the same way to get the
wavefunction in scattering region, so at this point you have already
determined a_n and C_{\alpha,l,m}. At last do boundary matching condition
between leads and scattering region to get reflection and transmission
coefficient. Am I right? My thinking is based on simple quantum mechanics
scattering problem where we solve schrodinger's equation in each region
and then do boundary matching to get reflection and transmission
coefficients.

Now looking at the subroutine 'transmit.f90', this does not appear to be
the case. It seems as if we are solving the Schrodinger's equation in the
scattering region and doing boundary matching at the same time. Using a
big matrix and solving something like AX=B. Would anybody mind explaing
what's going on in subroutine transmit.f90? 

Regards, 
Manoj Srivastava 
Physics Graduate Student 
University Of Florida, Gainesville, FL



[Pw_forum] symmetry operations in PW

2009-02-10 Thread Manoj Srivastava
Dear Alexandra and Gabriele, 
 Thank you for the help.

-manoj 
On Tue, 10 Feb 2009, Gabriele Sclauzero wrote:

> Dear Alexandra,
> 
> alexandra.carvalho at epfl.ch wrote:
> > Dear Gabriele and Manoj,
> > 
> > In fact I should clarify that I meant that rotations about directed axes 
> > belong
> > to different classes if the axes are not inverted by symmetry.
> 
> OK, I agree. Anyway when pwscf lists symmetry operations includes all 
> elements of the 
> symmetry group (which may or may not belong to the same class, this is 
> irrelevant at this 
> stage). The important point here is that those two rotations must always 
> appear together, 
> otherwise you will not obtain a group.
> 
> GS
> 
> 
> > (The word 'equivalent' was a poor choice...)
> > 
> > Alexandra
> > 
> > Quoting Gabriele Sclauzero :
> > 
> >> Dear Manoj and Alexandra,
> >>
> >>in my limited and maybe naive knowledge of group theory I think that
> >>
> >> alexandra.carvalho at epfl.ch wrote:
> >>> Dear Manoj,
> >>>
> >>> I do not know about the definitions intrinsic to pwscf,
> >>> but in group theory directed opperations are not allways equivalent.
> >>
> >> Rotation of 90 deg about z and about -z are _different_ operations (in 
> >> every
> >> group to
> >> which they belong...). Maybe pwscf uses a naive way to call symmetry
> >> operations, but you
> >> may call rotation of 90 deg about -z as "rotation of -90 deg about z" so
> >> maybe it's more
> >> clear (do you know the right hand rule?).
> >> If the rotation of 90 deg about z is present, also the rotation of -90 deg
> >> about z must be
> >> present (from the definition of group), as well as the rotation of 180 deg
> >> about z (which
> >> you can obtain by applying the rotation of 90 twice). This gives you a C_4
> >> symmetry about
> >> the z axis (then if you will add more symmetry operations you will have
> >> higher symmetry).
> >>
> >>  > In this case, the rotations about -z and z will be inequivalent
> >>> if there is no opperation of the group which reverses the sense of the z
> >> axis.
> >>
> >> I don't understand what do you mean by "equivalent/inequivalent" in this
> >> case...
> >>
> >> Hope it helps,
> >>
> >> GS
> >>
> >>> For example, in the T group there are 8 C3 opperations although there
> >>> are only four cube diagonals- see Tinkham "Group Theory and Quantum
> >>> mechanics" 2003 pp 58.
> >>>
> >>> Alexandra Carvalho
> >>> EPFL Lausanne
> >>>
> >>> Quoting Manoj Srivastava :
> >>>
> >>>> Dear PWSCF users and developers,
> >>>>  I was looking into the code where it calculates symmerty operations for
> >>>> the group. Now the symmetry operations are defined as ' 90 deg rotation -
> >>>> cart. axis [0,0,-1]', &  ' 90 deg rotation - cart. axis [0,0,1]
> >>>> I did not understand the point of defining rotation about z and -z
> >>>> axis. Is this because of some numerical constraint? Also, there are
> >>>> elements with 'theta' and '-theta' rotation in the same group. I was just
> >>>> wondering the reason behind this. Any help will be appreciated.
> >>>>
> >>>> Regards,
> >>>> Manoj Srivastava
> >>>> Department of Physics
> >>>> University of Florida
> >>>> Gainesville, FL, USA.
> >>>>
> >>>> ___
> >>>> Pw_forum mailing list
> >>>> Pw_forum at pwscf.org
> >>>> http://www.democritos.it/mailman/listinfo/pw_forum
> >>>>
> >>>
> >>> ___
> >>> Pw_forum mailing list
> >>> Pw_forum at pwscf.org
> >>> http://www.democritos.it/mailman/listinfo/pw_forum
> >>>
> >> --
> >>
> >>
> >> o  o
> >> | Gabriele Sclauzero, PhD Student  |
> >> | c/o:   SISSA & CNR-INFM Democritos,  |
> >> |via Beirut 2-4, 34014 Trieste (Italy) |
> >> | email: sclauzer at sissa.it |
> >> | phone: +39 040 3787 511  |
> >> | skype: gurlonotturno |
> >> o  o
> >> ___
> >> Pw_forum mailing list
> >> Pw_forum at pwscf.org
> >> http://www.democritos.it/mailman/listinfo/pw_forum
> >>
> > 
> > 
> > ___
> > Pw_forum mailing list
> > Pw_forum at pwscf.org
> > http://www.democritos.it/mailman/listinfo/pw_forum
> > 
> 
> -- 
> 
> 
> o  o
> | Gabriele Sclauzero, PhD Student  |
> | c/o:   SISSA & CNR-INFM Democritos,  |
> |via Beirut 2-4, 34014 Trieste (Italy) |
> | email: sclauzer at sissa.it |
> | phone: +39 040 3787 511  |
> | skype: gurlonotturno |
> o  o
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] symmetry operations in PW

2009-02-09 Thread Manoj Srivastava
Dear PWSCF users and developers, 
 I was looking into the code where it calculates symmerty operations for
the group. Now the symmetry operations are defined as ' 90 deg rotation -
cart. axis [0,0,-1]', &  ' 90 deg rotation - cart. axis [0,0,1] 
I did not understand the point of defining rotation about z and -z
axis. Is this because of some numerical constraint? Also, there are
elements with 'theta' and '-theta' rotation in the same group. I was just
wondering the reason behind this. Any help will be appreciated.

Regards, 
Manoj Srivastava
Department of Physics
University of Florida 
Gainesville, FL, USA.



[Pw_forum] error in pwcond.x calculations

2009-02-09 Thread Manoj Srivastava

Hey Erjun, Alexander, 
 I have noticed one thing that is odd. For ewind =3, or any value of it
less than 100, the code does basis reduction, which is not happening in
this case. have a look at this-  
ngper, shell number =   211   62
ngper, n2d =   211  211

So, i think something is wrong with ewind. I would suggest trying some
other values of ewind and see what happenes. 

Regards, 
Manoj Srivastava
Department of Physics 
University of Florida

 On Tue, 10 Feb 2009, ErJun Kan wrote:

> Hi, Alexander,
> 
> Thank you!
> 
> I want to calculate the conductivity of GaN:Mn film.
> In my case, I has used the leads same as the scattering region. 
> Unfortunately, it does not work. 
> 
> Bests,
> Erjun Kan
> 
> =
>  Erjun Kan
>  Department of Chemistry, North Carolina State University
>  Raleigh, North Carolina, USA
>  E-mail: erjunkan at mail.ustc.edu.cn  or
>  ekan at ncsu.edu
> = 
> 
> 
> 
> In your mail:
> >From: Alexander Smogunov 
> >Reply-To: PWSCF Forum 
> >To: pw_forum at pwscf.org
> >Subject: Re: [Pw_forum] error in pwcond.x calculations
> >Date:Mon, 9 Feb 2009 09:21:30 +0100
> >
> >Dear ErJun Kan
> > What kind of system do you study?
> > For conductance calculation you should have two semi-infinite
> > left and right metallic leads connected by some scattering region
> > which should match PERFECTLY to the leads. What are the leads and the
> > scattering region in your case?
> > Regards,
> > Alexander
> 
> 
> 
> 
> >   
> > 
> > On Saturday 07 February 2009 20:45, ErJun Kan wrote:
> > > Dear all,
> > >
> > > I have reduced the size of cell, it still does not work.
> > > what is the problem? thank you.
> > > +++
> > > +++
> > > ! SCF calculation part:
> > > # self-consistent calculation
> > >  
> > > calculation='scf'
> > > restart_mode='from_scratch',
> > > pseudo_dir = './',
> > > outdir='./'
> > > prefix='gan'
> > >  /
> > >  
> > > ibrav = 0,
> > >   nat = 36
> > > ntyp= 3,
> > > ecutwfc = 15.0,
> > > nspin = 2
> > > starting_magnetization(1)=0.7
> > > occupations='smearing',
> > > smearing='gaussian',
> > > degauss=0.01
> > >  /
> > >  
> > >  diagonalization='$diago'
> > >  mixing_mode = 'plain'
> > >  electron_maxstep = 200
> > > conv_thr = 1.0e-8
> > > mixing_beta = 0.7
> > >  /
> > > ATOMIC_SPECIES
> > >  Mn 54.938 Mn.pbe-sp-van.UPF
> > >  Ga 69.72  Ga.pbe-nsp-van.UPF
> > >  N  14.00   N.pbe-van_ak.UPF
> > > ATOMIC_POSITIONS {angstrom}
> > > Mn 1.590000.0   -6.425908496
> > > Ga 4.770000.0   -6.425908496
> > > Ga 7.950000.0   -6.425908496
> > > Ga-0.00.0   -3.671947712
> > > Ga 3.180000.0   -3.671947712
> > > Ga 6.360000.0   -3.671947712
> > > Ga 1.590000.0   -0.917986928
> > > Ga 4.770000.0   -0.917986928
> > > Ga 7.950000.0   -0.917986928
> > > Ga-0.02.59150   -5.507921568
> > > Ga 3.180002.59150   -5.507921568
> > > Ga 6.360002.59150   -5.507921568
> > > Ga 1.590002.59150   -2.753960784
> > > Ga 4.770002.59150   -2.753960784
> > > Ga 7.950002.59150   -2.753960784
> > > Ga-0.02.591500.0
> > > Ga 3.180002.591500.0
> > > Ga 6.360002.591500.0
> > > N -0.04.542899554   -5.507921568
> > > N  3.180004.542899554   -5.507921568
> > > N  6.360004.542899554   -5.507921568
> > > N  1.590004.542899554   -2.753960784
> > > N  4.770004.542899554   -2.753960784
> > > N  7.950004.542899554   -2.753960784
> > > N -0.04.5428995540.0
> > > N  3.180004.5428995540.0
> > > N  6.360004.5428995540.0
> > > N  

[Pw_forum] complex band

2008-12-26 Thread Manoj Srivastava
Dear Alexander and PWSCF users, 
 I have two questions about complex band code based on the Choi and Ihm
paper(PRB 59,2167,(1999)) 
 1. At the expressions where we have to use non local potential
W((k+G)_{\perp},tau-z-d), i can see the code only uses
W((k+G)_{\perp},tau-z). I guess this is because of periodicity of nonlocal
potential in the z direction. Am i right?
 
 2. This question is about determining the coefficients of the
wavefunctions. As we can see from equations (17) for the local part of the
potential it has N*2*N_2d number of unknown cofficients. Now using the
boundary conditions, that the wavefunction and its derivative should be
continuous between the slabs, we can determine (N-1)*2*N_2d unknownns,
which leaves 2*N_2d no. of coefficients to be determined later. The same
thing applies for the solutions of nonlocal part of the potential in
equation (24). Thus the number of linearly independent solution for the
nonlocal part is 2*N_2d. we use Bloch's condition to determine these
unknown coeffiecints. This is where my question is. When we are trying to
construct the general solution as in  equation (39), the local part has
number of unknown coeffieints 2*N_2d, which is fine. But the nonlocal part
has C_{\alpha,l,m} number of unknown coefficients, which is not equal to
2*N_2d. My question is why? Are we somehow combining the coefficients of
local part with nonlocal part? 
Any help will be appreciated. 
Happy Holidays!

Manoj Srivastava
Graduate Student 
Department of Physics 
University of Florida, Gainesville, Fl 32601 



[Pw_forum] calculation error about example12

2008-12-08 Thread Manoj Srivastava
Dear Chenweiguang, 
 The input file of example 12 is correct, at least in the version i am
using. I could run the example perfectly in version 4.0.1. Which version
are you using? Can you run any other examples?

Cheers, 
Manoj

 On Mon, 8 Dec 2008, Weiguang Chen wrote:

> Hi,
> 
> Another puzzle about example12 when I calculate the transmission of
> Al-wire, the last several lines of out-file is just as follows:
> 
>  E-Ef(ev), T(x2 spins) =0.250   4.000
>   
> %%
> 
>  from gep_x : error #38
> 
>  error on zggev
> 
>  
> %%
> 
> 
> 
>  stopping ...
> 
> p0_24174:  p4_error: : 0
> 
>  the input files of example12 is right or not  , and is there any
> tricks about transmission calculation?
> 
>  Thanks very much
> 
> -- 
> Best Wishes
> ChenWeiguang
> 
> 
> #   Chen, Weiguang
> #
> #Postgraduate,  Ph. D
> #  75 University Road, Physics Buliding  #  218
> #  School of Physics & Engineering
> #  Zhengzhou University
> #  Zhengzhou, Henan 450052  CHINA
> #
> #  Tel: 86-13203730117
> #  E-mail:chenweiguang82 at gmail.com;
> #chenweiguang82 at qq.com
> #**
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] complex band

2008-12-01 Thread Manoj Srivastava
Dear Alexander, 
 Dont worry about the previous mail, I figured it out :)

-manoj 
On Sun, 30 Nov 2008, Manoj Srivastava wrote:

> Dear Alexander, 
>  Thank you very much for your answer. It has started making sense, but
> unfortunately not a whole lot. I think I need to have complete
> understanding of this subroutine in order to have the bigger picture. In
> the other e-mail, copied and pasted as below - 
> > the expression -
> > norm2=(((il+xyk(1,ik))*bg(1,1)+(jl+xyk(2,ik))*bg(1,2))**2+ &
> > ((il+xyk(1,ik))*bg(2,1)+(jl+xyk(2,ik))*bg(2,2))**2)* &
> >
> > i believe is  |G+k|^2. How do you get above expression in the code for
> > this? 
> norm2 is indeed |G_perp+k_perp|^2, where G_perp=il*b_1+jl*b_2 and 2d
> vector 
> k_perp=xyk(1,ik)*b_1+xyk(2,ik)*b_2 and you indeed select out those 
> G_perp+k_perp which satisfy (G_perp+k_perp)^2 
> Why is K_perp in terms of b_1 and b_2? You can ofcourse take it
> like that but is there any specific reason to do that? If we look at the
> Fourier transform of wavefunction psi(r)=sum_{G}\psi(G)exp(i(K+G).r).
> here in this expression do we always take K in terms of G. I dont think
> so, K can be a label, which might be in terms of G vectors but not
> necessarily. What do you think? 
> In the input file where we enter k_x and k_y, do we enter  these  in
> terms of b1 and b2? or somewhere later in the code it gets multiplied with
> b1 and b2, if yes where? I am having tough time figuring it out :( 
> Thanks for help. 
> 
> Regards, 
> Manoj Srivastava
> Department of Physics,
> University of Florida, Gainesville
> 
> 
>  On Fri, 28 Nov 2008,
> Alexander wrote:
> 
> > Dear Manoj
> > Ngper is the total number of g_perp=G_perp+k_perp with |g_perp|^2 > which are arranged into the shells. Each shell contains g_perp having the 
> > same 
> > norm which is very useful because many calculations with g_perp depend only 
> > on their norms. The number of shells is ngpsh. gnsh(ngpsh) is the norm for 
> > each shell, ninsh(ngpsh)  is the number of g_perp in the shell.
> > The array nshell is the auxiliary array used further to arrange the g_perp 
> > vectors in the shells.
> > 
> > Hope this helps you,
> > regards, Alexander 
> > 
> > and ngpsh is the number of shells 
> > 
> > On Thursday 27 November 2008 00:26, Manoj Srivastava wrote:
> > > Dear Alexander, PWSCF users and developers,
> > >  I have a technical question in one of the subroutines used to calculate
> > > complex band. The subroutine 'init_gper.f90' under PWCOND directory is
> > > used to calculate number of G perpendicular vectors within the energy
> > > cutoff by using [k_{\perp}+G_{\perp}]^2.le.E_{cut}. I understand the
> > > counting of Gper, but i dont understand how do they construct the Gper
> > > vectors. For example what is 'ngpsh'? It is defined as 'no. of shells for
> > > G', but i really dont understand what this is. I dont understand after
> > > line 50 of the code in init_gper.f90.  Mainly what does the following
> > > piece of code do
> > >
> > >   do i=1,nrx
> > > do j=1,nry
> > >icount=0
> > >   do iw=1, ngpsh
> > >if (abs(norm2-gnorm2(iw)).gt.eps) then
> > >  icount=icount+1
> > >else
> > >  nshell(i,j)=iw
> > >endif
> > >  enddo
> > >  if (icount.eq.ngpsh) then
> > >ngpsh=ngpsh+1
> > >gnorm2(ngpsh)=norm2
> > >nshell(i,j)=ngpsh
> > >  endif
> > >   endif
> > >   enddo
> > >   enddo
> > > I understand it is not a good way of asking question, but I am completely
> > > lost. Any help will be appreciated.
> > >
> > > Regards,
> > > Manoj
> > >
> > > ___
> > > Pw_forum mailing list
> > > Pw_forum at pwscf.org
> > > http://www.democritos.it/mailman/listinfo/pw_forum
> > ___
> > Pw_forum mailing list
> > Pw_forum at pwscf.org
> > http://www.democritos.it/mailman/listinfo/pw_forum
> > 
> 
> 



[Pw_forum] possible bug in scatter_forw.f90

2008-10-20 Thread Manoj Srivastava

Dear Alexander, 
 Thank you very much for your answer. I have one more question, it seems
my list of questions is never gonna end :( 
 this is about init_gper.f90 subroutine. where you calculate
|G+k|^2 Dear Manoj
> 
> On Thursday 16 October 2008 23:15, Manoj Srivastava wrote:
> > Dear Alexandar,
> >  Thank you very much. Your answers make sense and now i have better
> > understanding of the subroutine scatter_forw.f90.  I have one question
> > about poten.f90 subroutine, where you calculate the V(p,gper) in each
> > slices of the slab(equation 15). I have a question about Gz, which is
> > 2*pi*q/L, with [-N/2]+1<=q>=[N/2]. Now in the code at the line at the line
> > 60, you have Gz multiplied with bg(3,3). I kind of understand this is
> > because 3rd direction is perpendicular to xy plane, 
> Exactly. As it is implemented now pwcond deals with systems 
> having monoclinic unit cell, i.e. such that a3 is perpendicular to both a1 
> and 
> a2 (a1 and a2 are not necessary orthogonal one to another) and the direction 
> of transport is a3. This seems to us the most general setup for transport and 
> I do not think we should consider at the moment the situation where the a3 is 
> not perpendicular to the xy plane.
> 
> Best regards, 
> Alexander
> 
> > but not totally clear 
> > about it. Would you mind explaining? my other question is when reciprocal
> > space vectors b1, b2 and b3 are not orhogonal to each other, which is
> > more general case, then should not we take bg(1,3) and bg(2,3) in there
> > too?
> >  Thank you once again for listening to my questions and i hope i am
> > not bothering you too much.
> >
> > Regards,
> > Manoj Srivastava
> > Physics Graduate Student
> > Department of Physics
> > Gainesville, FL, USA.
> >
> > have  On Wed, 15 Oct 2008, Alexander wrote:
> > > Dear Manoj
> > >
> > > On Tuesday 14 October 2008 17:49, Manoj Srivastava wrote:
> > > > Dear Alexader,
> > > >
> > > >  Thank you very much for your reply. I was out of town and could not
> > > > see the message till yesterday. I have some questions on it, it would
> > > > be great if you could answer them.
> > > >
> > > > On Fri, 10 Oct 2008, Alexander wrote:
> > > > > Dear Manoj
> > > > >
> > > > > On Friday 10 October 2008 05:20, Manoj Srivastava wrote:
> > > > > > Any follow ups on the message below? Alexander?
> > > > > >
> > > > > > Regards,
> > > > > > Manoj
> > > > > >
> > > > > > On Wed, 8 Oct 2008, Manoj Srivastava wrote:
> > > > > > > Dear All,
> > > > > > >  Is any body familiar with the scatter_forw.f90 subroutine of
> > > > > > > PWCOND? I think there is a bug in this subroutine at the place
> > > > > > > where you calculate intw2, which is z integration of nonlocal
> > > > > > > wavefunction with beta function. I have looked up the necessary
> > > > > > > formulae in the Choi & Ihm's paper (PRB 59, 2267, Jan 1999). In
> > > > > > > the paper, nonlocal wavefunction has 3 terms, one of which
> > > > > > > contains beta function, say W. The other two parts dont have any
> > > > > > > W in them, rather they are just plane wave solution. (Please have
> > > > > > > a look at equation 24 and 26 of the paper ). So, when we are
> > > > > > > doing z integration of nonlocal wavefunction with beta function
> > > > > > > W, we should have three terms, one of which should contain two W,
> > > > > > > but rest should just have one W. On the other hand in the code
> > > > > > > all three terms contain two beta functions!! (please have a look
> > > > > > > at line 228 of scatter_forw.f90). I am wondering if my
> > > > > > > understanding is right?
> > > > >
> > > > > Your understanding is wrong. On the line 228 you add to the integral
> > > > > ONLY contribution with two beta functions (from the 1st term in the
> > > > > nonlocal function, see Eq. 24) but on the line 393 the contribution
> > > > > from the 2nd term is added too. There is no 3rd term since every step
> > > > > you rotate your local solutions \psi_n in such a way that b_{\lambda
> > > > > \alpha lm} vanish (see the paragraph after Eq. 37).
> > > >
> > > >   I bel

[Pw_forum] possible bug in scatter_forw.f90

2008-10-16 Thread Manoj Srivastava
Dear Alexandar, 
 Thank you very much. Your answers make sense and now i have better
understanding of the subroutine scatter_forw.f90.  I have one question
about poten.f90 subroutine, where you calculate the V(p,gper) in each
slices of the slab(equation 15). I have a question about Gz, which is
2*pi*q/L, with [-N/2]+1<=q>=[N/2]. Now in the code at the line at the line
60, you have Gz multiplied with bg(3,3). I kind of understand this is
because 3rd direction is perpendicular to xy plane, but not totally clear
about it. Would you mind explaining? my other question is when reciprocal
space vectors b1, b2 and b3 are not orhogonal to each other, which is
more general case, then should not we take bg(1,3) and bg(2,3) in there
too? 
 Thank you once again for listening to my questions and i hope i am
not bothering you too much.

Regards, 
Manoj Srivastava 
Physics Graduate Student 
Department of Physics 
Gainesville, FL, USA.


have  On Wed, 15 Oct 2008, Alexander wrote:

> Dear Manoj
> 
> On Tuesday 14 October 2008 17:49, Manoj Srivastava wrote:
> > Dear Alexader,
> >
> >  Thank you very much for your reply. I was out of town and could not see
> > the message till yesterday. I have some questions on it, it would be great
> > if you could answer them.
> >
> > On Fri, 10 Oct 2008, Alexander wrote:
> > > Dear Manoj
> > >
> > > On Friday 10 October 2008 05:20, Manoj Srivastava wrote:
> > > > Any follow ups on the message below? Alexander?
> > > >
> > > > Regards,
> > > > Manoj
> > > >
> > > > On Wed, 8 Oct 2008, Manoj Srivastava wrote:
> > > > > Dear All,
> > > > >  Is any body familiar with the scatter_forw.f90 subroutine of PWCOND?
> > > > > I think there is a bug in this subroutine at the place where you
> > > > > calculate intw2, which is z integration of nonlocal wavefunction with
> > > > > beta function. I have looked up the necessary formulae in the Choi &
> > > > > Ihm's paper (PRB 59, 2267, Jan 1999). In the paper, nonlocal
> > > > > wavefunction has 3 terms, one of which contains beta function, say W.
> > > > > The other two parts dont have any W in them, rather they are just
> > > > > plane wave solution. (Please have a look at equation 24 and 26 of the
> > > > > paper ). So, when we are doing z integration of nonlocal wavefunction
> > > > > with beta function W, we should have three terms, one of which should
> > > > > contain two W, but rest should just have one W. On the other hand in
> > > > > the code all three terms contain two beta functions!! (please have a
> > > > > look at line 228 of scatter_forw.f90). I am wondering if my
> > > > > understanding is right?
> > >
> > > Your understanding is wrong. On the line 228 you add to the integral ONLY
> > > contribution with two beta functions (from the 1st term in the nonlocal
> > > function, see Eq. 24) but on the line 393 the contribution from the 2nd
> > > term is added too. There is no 3rd term since every step you rotate your
> > > local solutions \psi_n in such a way that b_{\lambda \alpha lm} vanish
> > > (see the paragraph after Eq. 37).
> >
> >   I believe that the equations 40, 41 and 42 implemented in the code are
> > in the momentum space. In the code, intw2 has 3 terms. We can tell this by
> > looking at the subroutine integrals.f90, where int2d is defined. Now, by
> > looking at the choi& Ihm 's paper, the very first term of psi{alpha,l,m}
> > has f{lambda,alpha,l,m}(equation26).And, if we substitute this in equation
> > 40, we should just get one term. I am just confuesd on the fact that how
> > come this one term in paper gets transformed into 3 terms of the code!
> 
> intw2(alpha,beta) is the integral of nonlocal solution (24) with w functions,
> intw2(alpha,beta) = \int_{0,d} w_{alpha}(r)*psi_{beta}(r)dr, where
> alpha,beta={alpha,l,m} in the paper. In each slab it will have TWO 
> contributions from the first and second terms of (24), the one from the third 
> term vanishes as I explained before. These two contributions are calculated
> on the lines  228  and 393 of the code. The integrals intw1
> (alpha,n)=\int_{0,d} w_{alpha}(r)*psi_{n}(r)dr are the integrals of local 
> solutions (17) with the same projector functions w. 
> The integrals intw2 and intw1 have nothing to do with boundary conditions 
> (40). They enter in (11) and are needed to construct the final solutions 
> according to Eq.(10). Substituting the form of this solution (10) in the 
> boundary conditions (40) you will be led to the ge

[Pw_forum] possible bug in scatter_forw.f90

2008-10-14 Thread Manoj Srivastava
Dear Alexader, 

 Thank you very much for your reply. I was out of town and could not see
the message till yesterday. I have some questions on it, it would be great
if you could answer them.
On Fri, 10 Oct 2008, Alexander wrote:

> Dear Manoj
> 
> On Friday 10 October 2008 05:20, Manoj Srivastava wrote:
> > Any follow ups on the message below? Alexander?
> >
> > Regards,
> > Manoj
> >
> > On Wed, 8 Oct 2008, Manoj Srivastava wrote:
> > > Dear All,
> > >  Is any body familiar with the scatter_forw.f90 subroutine of PWCOND? I
> > > think there is a bug in this subroutine at the place where you calculate
> > > intw2, which is z integration of nonlocal wavefunction with beta
> > > function. I have looked up the necessary formulae in the Choi & Ihm's
> > > paper (PRB 59, 2267, Jan 1999). In the paper, nonlocal wavefunction has 3
> > > terms, one of which contains beta function, say W. The other two parts
> > > dont have any W in them, rather they are just plane wave solution.
> > > (Please have a look at equation 24 and 26 of the paper ). So, when we are
> > > doing z integration of nonlocal wavefunction with beta function W, we
> > > should have three terms, one of which should contain two W, but rest
> > > should just have one W. On the other hand in the code all three terms
> > > contain two beta functions!! (please have a look at line 228 of
> > > scatter_forw.f90). I am wondering if my understanding is right?
> 
> Your understanding is wrong. On the line 228 you add to the integral ONLY
> contribution with two beta functions (from the 1st term in the nonlocal 
> function, see Eq. 24) but on the line 393 the contribution from the 2nd term 
> is added too. There is no 3rd term since every step you rotate your local 
> solutions \psi_n in such a way that b_{\lambda \alpha lm} vanish (see the 
> paragraph after Eq. 37).
  I believe that the equations 40, 41 and 42 implemented in the code are
in the momentum space. In the code, intw2 has 3 terms. We can tell this by
looking at the subroutine integrals.f90, where int2d is defined. Now, by
looking at the choi& Ihm 's paper, the very first term of psi{alpha,l,m}
has f{lambda,alpha,l,m}(equation26).And, if we substitute this in equation
40, we should just get one term. I am just confuesd on the fact that how
come this one term in paper gets transformed into 3 terms of the code!
Would you mind explaining? Also, the integrals that are defined in the
code are only in the slabs (nz1), so how are you taking the integration in
the region (0,d), which we need according to equations 40, 41 and 42 ?

> 
>  > > I have one more question in the later part of the same subroutine. What
> > > does the lapack subroutine
> > > ZGESV(2*n2d,2*n2d+norb*npol,amat,2*n2d,ipiv,xmat,2*n2d,info)
> > > do?
> > > I tried to look it up and i got the idea that it is trying to solve
> > > amat*x=xmat, with amat and xmat known and x unknown, and at the end of
> > > calculation it stores x in xmat.
> > > so, basically it does--  x=[(amat)^(-1)]*xmat. Am i right? So, it changes
> > > the structure of xmat from what is defined in 'constructs matrices'
> > > part. Is it correct?
> 
> Here, you are right.
> 
> > > Also, afterwards in this code where  it 'rotates integrals' is not very
> > > clear to me.
> > > Could somebody please tell me in little detail, what is going on
> > > here? 
> As I explained before rotation means that every step (at each slab)
> you make a linear combination of local solutions \psi_n with the 
> transformation matrix  h_{nn'} (see the last paragraph on the p. 2270).
> so that the new solutions have the property: 
> b_{\lambda n} = \delta_{\lambda n} and the new nonlocal solutions have
> b_{\lambda \alpha lm}=0. Doing this transformation of local and nonlocal
> solutions you should also perform the corresponding transformations
> of the integrals.
> 
> > > Also, is this subroutine written just on the basis of Choi and Ihm 
> > > paper, or are there more reference to it? If yes, would someone mind
> > > mentioning them?
> Unfortunately, the paper of Choi and Ihm is very detail, and in the 
> PWCOND code we followed it very closely, so no more references
> are needed. In our papers we just extended those ideas to ultrasoft
> pseudopotentials, magnetism, spin-orbit coupling and so on.
>   
> Hope this helped you,
> regards, Alexander 
> 
> > >
> > > Regards,
> > > Manoj Srivastava
> > > Physics Graduate Student
> > > University of Florida,
> > > Gainesville,FL, USA
> >
> > ___
> > Pw_forum mailing list
> > Pw_forum at pwscf.org
> > http://www.democritos.it/mailman/listinfo/pw_forum
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 

Once again thank you very much for your help. I appreciate it. 

Regards, 
Manoj Srivastava



[Pw_forum] possible bug in scatter_forw.f90

2008-10-10 Thread Manoj Srivastava

Any follow ups on the message below? Alexander? 

Regards, 
Manoj
On Wed, 8 Oct 2008, Manoj Srivastava wrote:

> 
> Dear All, 
>  Is any body familiar with the scatter_forw.f90 subroutine of PWCOND? I
> think there is a bug in this subroutine at the place where you calculate
> intw2, which is z integration of nonlocal wavefunction with beta function.
> I have looked up the necessary formulae in the Choi & Ihm's paper (PRB 59,
> 2267, Jan 1999). In the paper, nonlocal wavefunction has 3 terms, one of
> which contains beta function, say W. The other two parts dont have any W
> in them, rather they are just plane wave solution. (Please have a look at
> equation 24 and 26 of the paper ). So, when we are doing z integration of
> nonlocal wavefunction with beta function W, we should have three terms,
> one of which should contain two W, but rest should just have one W. On the
> other hand in the code all three terms contain two beta functions!!
> (please have a look at line 228 of scatter_forw.f90). I am wondering if my
> understanding is right? 
> I have one more question in the later part of the same subroutine. What
> does the lapack subroutine 
> ZGESV(2*n2d,2*n2d+norb*npol,amat,2*n2d,ipiv,xmat,2*n2d,info)
> do? 
> I tried to look it up and i got the idea that it is trying to solve 
> amat*x=xmat, with amat and xmat known and x unknown, and at the end of
> calculation it stores x in xmat. 
> so, basically it does--  x=[(amat)^(-1)]*xmat. Am i right? So, it changes
> the structure of xmat from what is defined in 'constructs matrices'
> part. Is it correct? 
> Also, afterwards in this code where  it 'rotates integrals' is not very
> clear to me. 
> Could somebody please tell me in little detail, what is going on
> here? Also, is this subroutine written just on the basis of Choi and Ihm
> paper, or are there more reference to it? If yes, would someone mind
> mentioning them? 
> 
> Regards, 
> Manoj Srivastava 
> Physics Graduate Student 
> University of Florida, 
> Gainesville,FL, USA
> 
> 



[Pw_forum] possible bug in scatter_forw.f90

2008-10-08 Thread Manoj Srivastava

Dear All, 
 Is any body familiar with the scatter_forw.f90 subroutine of PWCOND? I
think there is a bug in this subroutine at the place where you calculate
intw2, which is z integration of nonlocal wavefunction with beta function.
I have looked up the necessary formulae in the Choi & Ihm's paper (PRB 59,
2267, Jan 1999). In the paper, nonlocal wavefunction has 3 terms, one of
which contains beta function, say W. The other two parts dont have any W
in them, rather they are just plane wave solution. (Please have a look at
equation 24 and 26 of the paper ). So, when we are doing z integration of
nonlocal wavefunction with beta function W, we should have three terms,
one of which should contain two W, but rest should just have one W. On the
other hand in the code all three terms contain two beta functions!!
(please have a look at line 228 of scatter_forw.f90). I am wondering if my
understanding is right? 
I have one more question in the later part of the same subroutine. What
does the lapack subroutine 
ZGESV(2*n2d,2*n2d+norb*npol,amat,2*n2d,ipiv,xmat,2*n2d,info)
do? 
I tried to look it up and i got the idea that it is trying to solve 
amat*x=xmat, with amat and xmat known and x unknown, and at the end of
calculation it stores x in xmat. 
so, basically it does--  x=[(amat)^(-1)]*xmat. Am i right? So, it changes
the structure of xmat from what is defined in 'constructs matrices'
part. Is it correct? 
Also, afterwards in this code where  it 'rotates integrals' is not very
clear to me. 
Could somebody please tell me in little detail, what is going on
here? Also, is this subroutine written just on the basis of Choi and Ihm
paper, or are there more reference to it? If yes, would someone mind
mentioning them? 

Regards, 
Manoj Srivastava 
Physics Graduate Student 
University of Florida, 
Gainesville,FL, USA



[Pw_forum] some querries about pwcond

2008-10-04 Thread Manoj Srivastava
Dear Sagar, 
 i am attempting to answer one of the questions raised by you. I am not
entirely sure about the answer, I am just telling you what i have
understood and i am hoping somebody in the forum could correct me, if I am
wrong. 

On Sat, 4 Oct 2008, ambavale sagar wrote:

> Dear all,
> I have some queries about pwcond input file.
> 
> 1.how k-perpendicular points are defined. (can you give any example except 
> 'number of k-perpendicular points = 1' used in example 12) 
> 2. How weight is defined here for k-points? maximum weight?
  
 If the number of K point is one. You can choose any weight you like, it
normalizes to one. if you look at the file do_cond.f90 under PWCOND, you
would get the idea. I have never done any calculation with more than 1
k point, so i have not thought about that yet. 
 > 3.how ecut2D (2-D cut off) it effects the calculation. >
 
  ecut2D determines how many 2D reciprocal space vectors (G) (ngper)you
will be using in the calculation,|G+k|^2 Best regards
> Sagar Ambavale
> Research student
> The M.S. University of Baroda
> 
> 
> 


  Add more friends to your messenger and enjoy! Go to 
http://messenger.yahoo.com/invite/



[Pw_forum] question about complex band code

2008-08-29 Thread Manoj Srivastava

Dear Alexander and PWSCF users, 
 I have a question about the subroutine PWCOND of PWSCF, which calculates
complex bands. We are solving the generalized eigen value problem, A
X=exp(ikd) B X, where A and B and are required matrices. After the
problem is solved the information about k, which is complex eigen value is
stored in band.re , band.im and band.co files, and as the name suggests
they are purely real, purely imaginary and complex eigenvalues,
respectively. Generally, for each fixed energy we get few k's either real
or imaginary or both. 
 My question is about calculation of this complex eigen value. Total no.
of k we get from solving above GEP should be equal to the number of basis
used to describe A or B. I believe the information of k is contained in
kval(:) array defined in the compbs_2 subroutine of PWCOND. When I print
it out, I do get total number of k's which are same as the total number of
basis vectors.  This is not the case with bands.re, bands.im or bands.co
files. For each energy value they have arbitrary number of k's, why is it
so? What makes it throw some k values that it gets from 'kval'? I think, I
am missing something when the code goes from 'kval(:)' to bands.re,
bands.im and bands.co files. Does someone have any idea how does the code
go from kval(:) to write bands.re, bands.im or bands.co files? 

Regards, 
Manoj Srivastava,
Physics Graduate Student, 
University of Florida,
Gainesville, FL 32601 
USA




[Pw_forum] complex band of Cu (fwd)

2008-07-03 Thread Manoj Srivastava
Dear PWSCF users, 
my advisor is interested in the following question. Any suggestions?

-manoj 
-- Forwarded message --
Date: Thu, 03 Jul 2008 11:38:14 -0400
From: Xiaoguang Zhang <x...@ornl.gov>
Reply-To: zhangx at ornl.gov
To: Manoj Srivastava , pw_forum at pwscf.org
Cc: zhangx at ornl.gov
Subject: Re: [Pw_forum] complex band of Cu (fwd)

Thanks for the reply. Even though we're experimenting with the (100) 
direction which does allow a 2 atom/cell tetragonal unit which yields 
the correct band structure, we are really interested in the (111) 
complex bands. Why is it necessary to use a3 to define the transport 
direction? Wouldn't it be simpler and conceptually clearer to define the 
transport direction to be the norm of the base plane?
Thanks,
Xiaoguang

Manoj Srivastava wrote:
> 
> -- Forwarded message --
> Date: Wed, 2 Jul 2008 08:01:19 +0200
> From: Alexander 
> Reply-To: PWSCF Forum 
> To: pw_forum at pwscf.org
> Subject: Re: [Pw_forum] complex band of Cu
> 
> Dear Manoj 
> The PWCOND code works only for monoclinic cells
> where the a3 vector of the unit cell (defining the transport direction) is 
> ORTHOGONAL to both a1 and a2. You SHOULD therefore use the tetragonal unit 
> cell with two atoms per cell.
> Regards,
> Alexander
>  
> 
> On Tuesday 01 July 2008 21:31, Manoj Srivastava wrote:
> 
>>Dear PWSCF users,
>> I have been trying to calculate the band structure of bulk Cu. I tried
>>two different kind of unit cells, 2 atoms tetragonal and 1 atom
>>rhombohedron. In case of tetragonal cell, i get the right band structure,
>>but in the case of rhombohedron, d bands of Cu are not at the right place.
>>I believe that these calculations should not depend on the choice of unit
>>cell. I would appreiate if someone could tell me what am I doing wrong?
>>Attached is the band-structure file I get for rhombohedron unit cell and
>>my input file is as follows-
>>
>>
>>ibrav=0
>>celldm(1)=6.73
>>nat= 1,
>>ntyp= 1,
>>ecutwfc =50.0,
>>occupations= 'smearing',
>>smearing='gaussian',
>>degauss=0.01
>> /
>> 
>>conv_thr = 1.0e-8
>>mixing_beta=0.7
>> /
>> ATOMIC_SPECIES
>> Cu   63.55  Cu.pbe-paw_kj.UPF
>> ATOMIC_POSITIONS
>>  Cu 0.0 0.0 0.0
>>K_POINTS (automatic)
>>4 4 4 0 0 0
>>CELL_PARAMETERS
>>0.5 -0.5  0.0
>>0.5  0.5  0.0
>>0.5  0.0  0.5
>>
>># complex bands of cu along the 001 direction K_perp=0
>>cat > cu.cond.in << EOF
>> 
>>outdir='$TMP_DIR/'
>>prefixl='cu'
>>band_file ='bands.cu'
>>ikind=0
>>energy0=10.d0
>>denergy=-0.04d0
>>ewind=3.d0
>>epsproj=1.d-6
>> /
>>1
>>0.0 0.0 1.0
>>500
>>EOF
>>
>>Regards,
>>Manoj Srivastava
>>Graduate Student
>>University of Florida
>>Gainesville, USA.
>>
>>___
>>Pw_forum mailing list
>>Pw_forum at pwscf.org
>>http://www.democritos.it/mailman/listinfo/pw_forum
> 
> ___
> Pw_forum mailing list
> Pw_forum at pwscf.org
> http://www.democritos.it/mailman/listinfo/pw_forum
> 



[Pw_forum] complex band of Cu

2008-07-01 Thread Manoj Srivastava

Dear PWSCF users, 
 I have been trying to calculate the band structure of bulk Cu. I tried
two different kind of unit cells, 2 atoms tetragonal and 1 atom
rhombohedron. In case of tetragonal cell, i get the right band structure,
but in the case of rhombohedron, d bands of Cu are not at the right place.  
I believe that these calculations should not depend on the choice of unit
cell. I would appreiate if someone could tell me what am I doing wrong?
Attached is the band-structure file I get for rhombohedron unit cell and 
my input file is as follows-


ibrav=0
celldm(1)=6.73
nat= 1,
ntyp= 1,
ecutwfc =50.0,
occupations= 'smearing',
smearing='gaussian',
degauss=0.01
 /
 
conv_thr = 1.0e-8
mixing_beta=0.7
 /
 ATOMIC_SPECIES
 Cu   63.55  Cu.pbe-paw_kj.UPF
 ATOMIC_POSITIONS
  Cu 0.0 0.0 0.0
K_POINTS (automatic)
4 4 4 0 0 0
CELL_PARAMETERS
0.5 -0.5  0.0
0.5  0.5  0.0
0.5  0.0  0.5

# complex bands of cu along the 001 direction K_perp=0
cat > cu.cond.in << EOF
 
outdir='$TMP_DIR/'
prefixl='cu'
band_file ='bands.cu'
ikind=0
energy0=10.d0
denergy=-0.04d0
ewind=3.d0
epsproj=1.d-6
 /
1
0.0 0.0 1.0
500
EOF

Regards,
Manoj Srivastava
Graduate Student
University of Florida
Gainesville, USA.