Hello Anton,

I checked my code, I think the problem is still about the modes of the leads.
I can modify the phase of the modes in the leads. When I call the modes (using lead.info), I obtain modified wave functions but the scattering matrix is the same. I guess,the scattering matrix is still calculated with original wave functions.

I posted my code below.
Thank you

Nur

import sys as systemCommands
import kwant
import numpy as np

def modes(self, energy, args=()):
prop_modes_fix, stab_modes_fix = \
size_fix   = prop_modes_fix.wave_functions.shape[1]
fix = []
for i in range(size_fix):
fix0 = prop_modes_fix.wave_functions[0][i]
fix.append(fix0)
prop_modes, stab_modes = \
wf_size = prop_modes.wave_functions.shape[0] #number of orbital
num_modes   = prop_modes.wave_functions.shape[1]#number of modes (2*num_modes)
## change of the phase
for i in range(wf_size):
for j in range(num_modes):
f = fix[j]
prop_modes.wave_functions[i,j]=prop_modes.wave_functions[i,j]*(prop_modes.wave_functions[i,j]/f)
if prop_modes.wave_functions[i,j].real >= 0 :
prop_modes.wave_functions[i,j] = -1 * prop_modes.wave_functions[i,j]
else :
continue
return prop_modes, stab_modes
def make_system(a, t, L, f, r, rcc, kso, w):
lat = kwant.lattice.square(a)
sys = kwant.Builder()
def myShape(pos):
(x, y) = pos
if (L/2-f<= y < L/2+f) and (0<=x <2):
return True
elif (L/2-f<= y < L/2+f) and (L-1<= x< L+1):
return True
elif ( 2 <= x < L-1 ) and ( 0 <= y < L-1) and rcc ** 2 <(x ** 2 + y ** 2) and (3*r ** 2 < ((x-(L/2)) ** 2 + (y-(L/2)) ** 2)):
return True
def onsite(site1, tm):## Electric field in x direction
x1 =site1.pos[0]
y1 = site1.pos[1]
return 4*t  + np.sin(w*tm)*w*x1*kso
sys[lat.shape(myShape, (0, L/2))] = onsite
sys[kwant.builder.HoppingKind((1, 0), lat, lat)] = -t
sys[kwant.builder.HoppingKind((0, 1), lat, lat)] = -t
lead0 = kwant.Builder(kwant.TranslationalSymmetry((-a, 0)))
(x, y) = pos
return ( L/2-f <= y < L/2+f)

lead0[kwant.builder.HoppingKind((1, 0), lat, lat)] = -t
lead0[kwant.builder.HoppingKind((0, 1), lat, lat)] = -t
lead1 = kwant.Builder(kwant.TranslationalSymmetry((a, 0)))
(x, y) = pos
return ( L/2-f <= y < L/2+f)

lead1[kwant.builder.HoppingKind((1, 0), lat, lat)] = -t
lead1[kwant.builder.HoppingKind((0, 1), lat, lat)] = -t
fsys = sys.finalized()
###
###
return fsys

def main():
## w = frequency
## kso = Rashba spin orbit coupling
sys = make_system(a=1, t=1, L=50, f=2, r=6, rcc=5, kso = 0.3/50, w = 0.001)
energies=[0.00001 * i+0.40 for i in xrange(300)]
tm =1000 #time parameter
for energy in energies:
S = kwant.smatrix(sys, energy, args=[tm])
wf =mode[0].wave_functions
print wf
if __name__ == '__main__':
main()

On 24 January 2016 at 22:52, Anton Akhmerov wrote:
Hi Nur,

I don't think you need to do anything inside the scattering region to
modify the phase of the scattering matrix. Are you sure this is the
problem?

Best,
Anton

On Sun, Jan 24, 2016 at 6:04 PM, Fahriye Nur Gursoy
<fahriyetortop@sabanciuniv.edu> wrote:
> Hello all,
>
> My question is about modifying the phase of wave functions. I want to
> calculate derivative of scattering matrix with respect to energy and they
> come out to be discontinuous because of the phases.
>
> To remedy this I need to modify these phases. Following the suggestion of
> Anton Akhmerov in this mail group I understood how to modify the phase of
> modes in the  leads.
>
> But I also need to change the phase of wave functions in the scattering
> region to obtain continuous scattering matrices at different energies.
>
> I wonder if there is an another way to call wave functions in the scattering
> region (apart from kwant.solvers.default.wave_function) or a method to
> modify these phases.
>
> Thank you in advance.
>
> Nur