6.3. Thermoelectric transport
The Seebeck coefficient (\(S\)) reflects the ability of the device to produce a current induced by a temperature gradient [DDV] . At low temperature the \(S\) reads approximately [LWC]:
where:
\(\tau\) is the electronic transmission,
\(T\) is the absolute temperature in Kelvin,
\(k_{B}\) the Boltzmann constant,
\(E_{F}\) the Fermi energy.
This expression shows that \(S\) is related not only to the value of \(\tau\) but also to its slope at the Fermi energy.
In this tutorial, we will calculate the thermoelectric properties in a device structure using an armchair graphene nanoribbon (AGNR) model. Since we have already built the 8AGNR we’ll use it to investigate the thermoelectric transport properties (Fig. 6.3.1) using NanoDCAL.
Fig. 6.3.1 Schematic view of 8AGNR transport geometry.
The calculation steps are:
Build:
8AGNR unit cell;
Expand the primitive 8AGNR unit cell;
Convet to device structure;
Self-consistent calculation:
Generate the input files;
SCF of electrodes;
NEGF-SCF calculation for the scattering region.
Transmission calculation.
Seebeck coefficient calculation.
6.3.1. Build 8AGNR device structure
For practical purposes of this tutorial, we will consider a pristine 8AGNR, periodic along the z direction, as we showed in previous tutorial (Thermoelectric properties of materials). This device can be built using Device Studio modeling as shown in Fig. 6.3.2:
Open Device Studio software, go to the menu bar and click Build ➟ Gr_Nanoribbon;
In the Gr_Nanoribbon window, keep uncheck the Passivate dangling bonds with hydrogen checkbox option, set n = 0 and m = 4, click Preview and Build the structure.
In the menu bar click Build ➟ Redefine Crystal;
In the cell vector set u = (1a, 0b, 0c), v = (0a, 1b, 0c), and w = (0a, 0b, 8c), click Preview or Build to expand the cell 8 times in the c direction, click on Preview and Build the structure.
In the menu bar click Build ➟ Convert to Device, check Left electrode and Right electrode options, click Preview and Build the structure.
Fig. 6.3.2 Graphical user interface of Device Studio to build, expand the unit cell, and convert to device.
The 8AGNR device Fig. 6.3.1 can be divided into three regions: two electrodes (left and right) and a central region. The left and right electrodes are semi-infinite with periodic boundary conditions along the \(\pm\) z-direction and don’t differ from each other. The central region is terminated by buffer regions (electrode copies) and a scattering region in between. The calculation of a two-probe system typically consists of the following steps:
SCF calculation for the electrodes;
NEGF-SCF calculation for the scattering region.
6.3.2. Generate the input files
The required input files for these steps could be generated using Device Studio with NanoDCAL module simulator, as follows:
In the menu bar of Device Studio, click Simulator ➟ NanoDCAL ➟ SCF Calculation to pop up the interface shown in Fig. 2.1.3;
Users can set calculation parameters related to SCF-NEGF method;
In this example the default parameters will be used.
Click Generate files;
The Device Studio will create a new folder that contain the scf.input files of LeftElectrode, RightElectrode and Device.
Fig. 6.3.3 NanoDCAL module simulator.
Now, we are ready to start the calculations. But, before doing that we’ll briefly present the simulation parameters specified in the generated files.
Electrode SCF input file
%%What quantities should be calculated
calculation.name = scf
%Basic setting
calculation.occupationFunction.temperature = 100
calculation.realspacegrids.E_cutoff = 80 Hartree
calculation.xcFunctional.Type = LDA_PZ81
calculation.k_spacegrids.number = [ 1 1 100 ]'
system.centralCellVectors = [[32.678 0 0]' [0 20 0]' [0 0 4.26]']
system.spinType = NoSpin
%Iteration control
calculation.SCF.monitoredVariableName = {'rhoMatrix','hMatrix','totalEnergy','bandEnergy','gridCharge','orbitalCharge'}
calculation.SCF.convergenceCriteria = {1e-04,1e-04,[],[],[],[]}
calculation.SCF.maximumSteps = 200
calculation.SCF.mixMethod = Pulay
calculation.SCF.mixRate = 0.1
calculation.SCF.mixingMode = H
calculation.SCF.startingMode = H
%calculation.SCF.donatorObject = NanodcalObject.mat
%Basic set
system.neutralAtomDataDirectory = '../'
system.atomBlock = 16
AtomType OrbitalType X Y Z
C LDA-DZP 11.42000000 10.00000000 2.84000000
C LDA-DZP 11.42000000 10.00000000 1.42000000
C LDA-DZP 12.64975607 10.00000000 0.71000000
C LDA-DZP 12.64975607 10.00000000 3.55000000
C LDA-DZP 13.87951215 10.00000000 2.84000000
C LDA-DZP 13.87951215 10.00000000 1.42000000
C LDA-DZP 15.10926822 10.00000000 0.71000000
C LDA-DZP 15.10926822 10.00000000 3.55000000
C LDA-DZP 16.33902429 10.00000000 2.84000000
C LDA-DZP 16.33902429 10.00000000 1.42000000
C LDA-DZP 17.56878037 10.00000000 0.71000000
C LDA-DZP 17.56878037 10.00000000 3.55000000
C LDA-DZP 18.79853644 10.00000000 2.84000000
C LDA-DZP 18.79853644 10.00000000 1.42000000
C LDA-DZP 20.02829251 10.00000000 0.71000000
C LDA-DZP 20.02829251 10.00000000 3.55000000
end
The keywords specify the following:
calculation.nameThis input parameter specifies the task that NanoDCAL should perform.calculation.occupationFunction.temperatureElectronic temperature of the system, in the unit of K.calculation.realspacegrids.E_cutoffThe equivalent energy cut-off of the grid density.calculation.xcFunctional.TypeThe type of exchange-correlation functional to be used.calculation.k_spacegrids.numberThe small k-space grid number in each direction. It is used to divide the Brillouin zone into the input number of small grids.system.centralCellVectorsThis parameter defines the three base vectors of the central cell, either the primitive or super cell for a periodic system or the central scattering region for an open system, which is with a shape of parallelepiped.system.spinTypeDetermines the spin description.calculation.SCF.monitoredVariableNameName list of variables to be monitored.calculation.SCF.convergenceCriteriaThe SCF loop will stop when all monitored variables satisfy this criteria.calculation.SCF.maximumStepsName list of variables to be monitored.calculation.SCF.convergenceCriteriaThe maximum steps in the SCF cycle.calculation.SCF.monitoredVariableNameName list of variables to be monitored.calculation.SCF.mixMethodWith this parameter an user may select a built-in mixer.calculation.SCF.mixRateA parameter to control the mixing rate with the method selected by the parameter calculation.SCF.mixMethod.calculation.SCF.mixingModeIn the SCF iterations, the mixing value could be any of Hamiltonian matrix, rho matrix, and real space density.calculation.SCF.startingModeThe SCF cycle can start from different physical quantities.calculation.SCF.donatorObjectThe value of the starting physical quantity will be taken from the donatorObject.system.neutralAtomDataDirectoryDirectory containing neutral atom (basis) filessystem.atomBlockThis block defines the number, type, orbital, and xyz positions of each atom in the system.
Note
The left and right electrodes are the same in this pristine device. Thus, the calculations will be performed only with the left electrode.
Central region SCF input file
Again, we present a briefly description of the scf.input from the central region.
%%What quantities should be calculated
calculation.name = scf
%Basic setting
calculation.occupationFunction.temperature = 100
calculation.realspacegrids.E_cutoff = 80 Hartree
calculation.xcFunctional.Type = LDA_PZ81
calculation.k_spacegrids.number = [ 1 1 1 ]'
%Description of electrode
system.numberOfLeads = 2
system.typeOfLead1 = left
system.voltageOfLead1 = 0
system.objectOfLead1 = ../LeftElectrode/NanodcalObject.mat
system.typeOfLead2 = right
system.voltageOfLead2 = 0
system.objectOfLead2 = ../LeftElectrode/NanodcalObject.mat
%Contour integral
%calculation.complexEcontour.lowestEnergyPoint = 1.5 Hartree
calculation.complexEcontour.numberOfPoints = 40
calculation.realEcontour.interval = 0.0272114
calculation.realEcontour.eta = 0.0272114
system.centralCellVectors = [[32.678 0 0]' [0 20 0]' [0 0 25.56]']
system.spinType = NoSpin
%Iteration control
calculation.SCF.monitoredVariableName = {'rhoMatrix','hMatrix','totalEnergy','bandEnergy','gridCharge','orbitalCharge'}
calculation.SCF.convergenceCriteria = {1e-04,1e-04,[],[],[],[]}
calculation.SCF.maximumSteps = 300
calculation.SCF.mixMethod = Pulay
calculation.SCF.mixRate = 0.1
calculation.SCF.mixingMode = H
calculation.SCF.startingMode = H
calculation.SCF.donatorObject = NanodcalObject.mat
calculation.SCF.maximumTime = 10
%Basic set
system.neutralAtomDataDirectory = '../'
system.atomBlock = 96
AtomType OrbitalType X Y Z
C LDA-DZP 11.42000000 10.00000000 2.84000000
C LDA-DZP 11.42000000 10.00000000 7.10000000
C LDA-DZP 11.42000000 10.00000000 11.36000000
C LDA-DZP 11.42000000 10.00000000 15.62000000
C LDA-DZP 11.42000000 10.00000000 19.88000000
C LDA-DZP 11.42000000 10.00000000 24.14000000
C LDA-DZP 11.42000000 10.00000000 1.42000000
C LDA-DZP 11.42000000 10.00000000 5.68000000
C LDA-DZP 11.42000000 10.00000000 9.94000000
C LDA-DZP 11.42000000 10.00000000 14.20000000
C LDA-DZP 11.42000000 10.00000000 18.46000000
C LDA-DZP 11.42000000 10.00000000 22.72000000
C LDA-DZP 12.64975607 10.00000000 0.71000000
C LDA-DZP 12.64975607 10.00000000 4.97000000
C LDA-DZP 12.64975607 10.00000000 9.23000000
C LDA-DZP 12.64975607 10.00000000 13.49000000
C LDA-DZP 12.64975607 10.00000000 17.75000000
C LDA-DZP 12.64975607 10.00000000 22.01000000
C LDA-DZP 12.64975607 10.00000000 3.55000000
C LDA-DZP 12.64975607 10.00000000 7.81000000
C LDA-DZP 12.64975607 10.00000000 12.07000000
C LDA-DZP 12.64975607 10.00000000 16.33000000
C LDA-DZP 12.64975607 10.00000000 20.59000000
C LDA-DZP 12.64975607 10.00000000 24.85000000
C LDA-DZP 13.87951215 10.00000000 2.84000000
C LDA-DZP 13.87951215 10.00000000 7.10000000
C LDA-DZP 13.87951215 10.00000000 11.36000000
C LDA-DZP 13.87951215 10.00000000 15.62000000
C LDA-DZP 13.87951215 10.00000000 19.88000000
C LDA-DZP 13.87951215 10.00000000 24.14000000
C LDA-DZP 13.87951215 10.00000000 1.42000000
C LDA-DZP 13.87951215 10.00000000 5.68000000
C LDA-DZP 13.87951215 10.00000000 9.94000000
C LDA-DZP 13.87951215 10.00000000 14.20000000
C LDA-DZP 13.87951215 10.00000000 18.46000000
C LDA-DZP 13.87951215 10.00000000 22.72000000
C LDA-DZP 15.10926822 10.00000000 0.71000000
C LDA-DZP 15.10926822 10.00000000 4.97000000
C LDA-DZP 15.10926822 10.00000000 9.23000000
C LDA-DZP 15.10926822 10.00000000 13.49000000
C LDA-DZP 15.10926822 10.00000000 17.75000000
C LDA-DZP 15.10926822 10.00000000 22.01000000
C LDA-DZP 15.10926822 10.00000000 3.55000000
C LDA-DZP 15.10926822 10.00000000 7.81000000
C LDA-DZP 15.10926822 10.00000000 12.07000000
C LDA-DZP 15.10926822 10.00000000 16.33000000
C LDA-DZP 15.10926822 10.00000000 20.59000000
C LDA-DZP 15.10926822 10.00000000 24.85000000
C LDA-DZP 16.33902429 10.00000000 2.84000000
C LDA-DZP 16.33902429 10.00000000 7.10000000
C LDA-DZP 16.33902429 10.00000000 11.36000000
C LDA-DZP 16.33902429 10.00000000 15.62000000
C LDA-DZP 16.33902429 10.00000000 19.88000000
C LDA-DZP 16.33902429 10.00000000 24.14000000
C LDA-DZP 16.33902429 10.00000000 1.42000000
C LDA-DZP 16.33902429 10.00000000 5.68000000
C LDA-DZP 16.33902429 10.00000000 9.94000000
C LDA-DZP 16.33902429 10.00000000 14.20000000
C LDA-DZP 16.33902429 10.00000000 18.46000000
C LDA-DZP 16.33902429 10.00000000 22.72000000
C LDA-DZP 17.56878037 10.00000000 0.71000000
C LDA-DZP 17.56878037 10.00000000 4.97000000
C LDA-DZP 17.56878037 10.00000000 9.23000000
C LDA-DZP 17.56878037 10.00000000 13.49000000
C LDA-DZP 17.56878037 10.00000000 17.75000000
C LDA-DZP 17.56878037 10.00000000 22.01000000
C LDA-DZP 17.56878037 10.00000000 3.55000000
C LDA-DZP 17.56878037 10.00000000 7.81000000
C LDA-DZP 17.56878037 10.00000000 12.07000000
C LDA-DZP 17.56878037 10.00000000 16.33000000
C LDA-DZP 17.56878037 10.00000000 20.59000000
C LDA-DZP 17.56878037 10.00000000 24.85000000
C LDA-DZP 18.79853644 10.00000000 2.84000000
C LDA-DZP 18.79853644 10.00000000 7.10000000
C LDA-DZP 18.79853644 10.00000000 11.36000000
C LDA-DZP 18.79853644 10.00000000 15.62000000
C LDA-DZP 18.79853644 10.00000000 19.88000000
C LDA-DZP 18.79853644 10.00000000 24.14000000
C LDA-DZP 18.79853644 10.00000000 1.42000000
C LDA-DZP 18.79853644 10.00000000 5.68000000
C LDA-DZP 18.79853644 10.00000000 9.94000000
C LDA-DZP 18.79853644 10.00000000 14.20000000
C LDA-DZP 18.79853644 10.00000000 18.46000000
C LDA-DZP 18.79853644 10.00000000 22.72000000
C LDA-DZP 20.02829251 10.00000000 0.71000000
C LDA-DZP 20.02829251 10.00000000 4.97000000
C LDA-DZP 20.02829251 10.00000000 9.23000000
C LDA-DZP 20.02829251 10.00000000 13.49000000
C LDA-DZP 20.02829251 10.00000000 17.75000000
C LDA-DZP 20.02829251 10.00000000 22.01000000
C LDA-DZP 20.02829251 10.00000000 3.55000000
C LDA-DZP 20.02829251 10.00000000 7.81000000
C LDA-DZP 20.02829251 10.00000000 12.07000000
C LDA-DZP 20.02829251 10.00000000 16.33000000
C LDA-DZP 20.02829251 10.00000000 20.59000000
C LDA-DZP 20.02829251 10.00000000 24.85000000
end
The keywords are specifies the following:
system.numberOfLeadsThe number of leads (probes) in the system.system.typeOfLead1The value must be given for each of the leads in the system. It describes the relative position of the lead to the central scattering region.system.voltageOfLead1The bias voltage applied to the lead, in the unit of Volt. The change of chemical potential of the lead due to the bias voltage is -eV.system.objectOfLead1The calculation of a system with leads can be performed in two steps. The first step is to perform calculation for each lead which is considered equivalently as an infinite periodic system. In the second step, the physical properties (potential, density etc.) of the leads are used as boundary conditions applied to the central scattering region of the device. In this scenario, this input value is normally a filename of the calculated and saved results of the leads.calculation.complexEcontour.lowestEnergyPointThe lowest energy point on the complex energy contour.calculation.complexEcontour.numberOfPointsNumber of energy points used on the complex energy contour for integrating the equilibrium part of the lesser Green’s function.calculation.realEcontour.intervalEnergy interval used to determine the number of energy points for integrating the non-equilibrium part of lesser Green’s function on the real energy axis.calculation.realEcontour.etathe small eta used in the calculation of self-energy and/or Green’s function along the real energy contour.
6.3.3. SCF calculations
To perform the calculation download the input files:
scf.input file from the Left electrode
8AGNR.L-elec.scf.input;scf.input file from the Right electrode
8AGNR.R-elec.scf.input;scf.input file from the Central region
8AGNR.CentralRegion.scf.input.
The users can also use Device Studio to peform the calculations, as following:
On the Device Studio navigate to the Project panel:
Go to the LeftElectrode folder ➟ right click over scf.input file ➟ run to pop up the Run dialog window, as shown in Fig. 2.1.4;
In the Run window, choose the Gateway location, number of cores and press run;
Once the calculation has completed, the results will be returned to the LeftElectrode folder.
Attention
The electrode calculation will take few minute. After the job finished continue the steps.
Go to the Device folder ➟ right click over scf.input file ➟ click on Open with to pop up the input file dialog.
Replace the path of Right electrode in system.objectOfLead2 ➟ system.objectOfLead2 = ../LeftElectrode/NanodcalObject.mat ➟ press save and close the window.
Back to the Device folder ➟ right click over scf.input ➟ run to pop up the Run dialog window ➟ choose the Gateway location, number of cores and press run.
Fig. 6.3.4 Run dialog window of Device Studio.
Important
Several ground state properties of the system could also be calculated from the obtained SCF results. In fact, the electronic properties of the electrode could be relevant to further analys of transport properties. The users can visit our previous tutorials calculate and plot
Band sctructure properties;
Charge density;
Electrostatic potencial.
Once the self-consistent two-probe calculation finish. To obtained the thermoelectric transport properties, We must have to peforme the transmission calculation, which is the most fundamental quantum transport property.
6.3.4. Transmission calculation
The transmission.input file can be created following the steps below:
In menu bar of Device Studio, click on Simulator ➟ NanoDCAL ➟ Analysis to set physical quantity, as shown in Fig. 6.3.5;
In the Analysis window, on the top left choose the Device option;
Navigate to the Analysis panel (left side) ➟ select Transmission ➟ click on the Right arrow button to transfer Transmission calculator to the Calculation Selected panel.
The TransmissionSpectrum options are shown in the right side of the Analysis window. Edit the k-point sampling set to n1 = 1, n2 = 1, n3 = 1 and Energy range set to -1.5, 1.5.
click on Generate files.
Fig. 6.3.5 Physical properties interface of the Device Studio.
The Device Studio will create a new file in the Device directory that will be the Transmission input file.
The generated file could be downloaded or users can use Device Studio to peforme the calculation, as we shown above (run SCF).
Before starting the calculation, we’ll briefly present the simulation parameters specified in the Transmission file.
Transmission input file
%%What quantities should be calculated
system.object = NanodcalObject.mat
calculation.name = transmission
calculation.transmission.kSpaceGridNumber = [ 1 1 100 ]'
calculation.transmission.energyPoints = [-2.00:0.005:2.00]
%calculation.transmission.plot = true
calculation.control.xml = true
The keywords specifies the following:
calculation.transmission.kSpaceGridNumberthe small k-space grid number in each direction which, together with kSpaceGridShift, are used to produce the parameter kSpacePoints.calculation.transmission.energyPointsThe energy points at which the transmission will be calculated. Note that the energy values are measured from chemical potential of a lead having zero applied voltage.calculation.bandStructure.coordinatesOfTheSymmetryKPointsThe k-space fractional coordinates of those high symmetry k-points listed in symmetryKPoints, each column corresponds to one of them.
Transmission spectrum
Next, you should analyses the results from the produced transmission dataset (.xml and .mat) that can then be plotted using the MATLAB or Device Studio by:
In menu bar of Device Studio, click on Simulator ➟ NanoDCAL ➟ Analysis Plot ➟ import Transmission.xml file.
As expected for a perfect 1D system, the transmission spectrum in Fig. 6.3.6 exhibits a sequence of steps with several integer steps (conducting channels) and narrow jump points, which is consistent with previous published results [DZD].
Fig. 6.3.6 The calculated transmission coefficient for the 8AGNR device.
The transmission spectrum is clearly non-zero in the vicinity of the Fermi level, which is related with the unsatisfied valence at the edge of the 8AGNR that produces states around the Fermi energy owing to unpaired electrons.
The key features and interpretations of the transmission spectrum include:
The transmission function of the material describes how electrons or charge carriers propagate through the material at different energies.
Identify peaks or features in the transmission spectrum that correspond to specific electronic states or energy levels in the material.
Determine the overall transmission properties of the material, such as the transmission coefficient at different energy levels.
6.3.5. Seebeck coefficient calculation
The Seebeck.input file can be created following the steps below:
In menu bar of Device Studio, click on Simulator ➟ NanoDCAL ➟ Analysis to set physical quantity, as shown in Fig. 6.3.7;
In the Analysis window, on the top left choose the Device option;
Navigate to the Analysis panel (left side) ➟ select Seebeck ➟ click on the Right arrow button to transfer Seebeck calculator to the Calculation Selected panel.
click on Generate files.
Fig. 6.3.7 Physical properties interface of the Device Studio.
The Device Studio will create new files in the Device directory that will be the Seebeck input file.
The generated input file could be downloaded or users can use Device Studio to peforme the calculations, as we shown above (run SCF).
Before starting the calculation, we’ll briefly present the simulation parameters specified in the Seebeck and Seebeck files.
Seebeck input file
system.object = NanodcalObject.mat
calculation.name = SeebeckCoefficient
calculation.Seebeck.temperature = [50:50:300]
calculation.Seebeck.chemicalPotentialEnergies = [-1:0.01:1]
%calculation.Seebeck.plot = true
calculation.control.xml = true
The keywords specifies the following:
calculation.Seebeck.temperatureThe Seebeck coefficients will be calculated respectively when the temperature is set to those input values. It is in the unit of Kelvin. The value of Boltzmann constant \(k = 8.61733326 \times 10^{-5}\) eV/K, \(k = 27.2113834\) eV/Hartree, \(k = 3.1668151 \times 10^{-5}\) Hartree/Kelvin.calculation.Seebeck.chemicalPotentialEnergiesThe Seebeck coefficients will be calculated respectively when the chemical potential energy is set to those input values. It is in the unit of eV.
Important
First, run the simulation with the transmission input file
Transmission.input. This calculation will produced transmission dataset (.xml and .mat).Second, run the simulation with the Seebeck input file
Seebeck.input. This calculation will produced Seebeck coefficient dataset (.xml and .mat).
Seebeck coefficient
Next, you should analyses the results from the produced transmission dataset (.xml and .mat) that can then be plotted using the MATLAB or Device Studio by:
In menu bar of Device Studio, click on Simulator ➟ NanoDCAL ➟ Analysis Plot ➟ import SeebeckCoefficient.xml file.
As expected for a perfect armchair graphene nanoribbon, the Seebeck spin in Fig. 6.3.8 is zero due to the absence of edge states.
Fig. 6.3.8 The calculated transmission coefficient for the 8AGNR device.
The Seebeck coefficient increases as the temperature increases from 50K to 300K. The coefficient \(S\) is quite large, in the energy range of the conducting channels. For instance, the \(S\) reach to \(65\mu V/K\) around the Fermi energy, where transmission coefficients about \(T(E)=3\).
The key features and interpretations of the Seebeck coefficient include:
Compare the features of the transmission spectrum with the Seebeck coefficient to understand how the material’s electronic properties influence its thermoelectric behavior.
Investigate the behavior of the Seebeck coefficient at different temperature ranges and its impact on the overall thermoelectric performance of the material.
Consider the anisotropic nature of the material and how it influences the Seebeck coefficient.
Analyze how different types and concentrations of dopants impact the Seebeck coefficient and overall thermoelectric performance.
Compare the Seebeck coefficient with other thermoelectric properties, such as electrical conductivity and thermal conductivity. Understand how these properties interact and contribute to the overall efficiency of the thermoelectric material.
Tip
As an additional task the users should be able to:
Calculate the Seebeck coefficient for 8AGNR passivated with hydrogen atom.
Explore the Seebeck properties for different ribbon width \(N=3p\), \(N=3p+1\) and \(N=3p+2\).
Dubi and M. Di Ventra. Thermospin effects in a quantum dot connected to ferromagnetic leads. Phys. Rev. B 80 (2009), p. 119902.
Y.-S. Liu, X.-F. Wang, and F. Chi. Non-magnetic doping induced a high spin-filter efficiency and large spin Seebeck effect in zigzag graphene nanoribbons. J. Mater. Chem. C, 1 (2013), p. 8046.
M.-M. Dubois, Z. Zanolli, X. Declerck, and J.-C. Charlier. Electronic properties and quantum transport in Graphene-based nanostructures. Eur. Phys. J. B 72 (2009), p. 1.