# 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.

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.

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**.

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.name`

This input parameter specifies the task that NanoDCAL should perform.`calculation.occupationFunction.temperature`

Electronic temperature of the system, in the unit of K.`calculation.realspacegrids.E_cutoff`

The equivalent energy cut-off of the grid density.`calculation.xcFunctional.Type`

The type of exchange-correlation functional to be used.`calculation.k_spacegrids.number`

The small k-space grid number in each direction. It is used to divide the Brillouin zone into the input number of small grids.`system.centralCellVectors`

This 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.spinType`

Determines the spin description.`calculation.SCF.monitoredVariableName`

Name list of variables to be monitored.`calculation.SCF.convergenceCriteria`

The SCF loop will stop when all monitored variables satisfy this criteria.`calculation.SCF.maximumSteps`

Name list of variables to be monitored.`calculation.SCF.convergenceCriteria`

The maximum steps in the SCF cycle.`calculation.SCF.monitoredVariableName`

Name list of variables to be monitored.`calculation.SCF.mixMethod`

With this parameter an user may select a built-in mixer.`calculation.SCF.mixRate`

A parameter to control the mixing rate with the method selected by the parameter calculation.SCF.mixMethod.`calculation.SCF.mixingMode`

In the SCF iterations, the mixing value could be any of Hamiltonian matrix, rho matrix, and real space density.`calculation.SCF.startingMode`

The SCF cycle can start from different physical quantities.`calculation.SCF.donatorObject`

The value of the starting physical quantity will be taken from the donatorObject.`system.neutralAtomDataDirectory`

Directory containing neutral atom (basis) files`system.atomBlock`

This 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.numberOfLeads`

The number of leads (probes) in the system.`system.typeOfLead1`

The 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.voltageOfLead1`

The 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.objectOfLead1`

The 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.lowestEnergyPoint`

The lowest energy point on the complex energy contour.`calculation.complexEcontour.numberOfPoints`

Number of energy points used on the complex energy contour for integrating the equilibrium part of the lesser Green’s function.`calculation.realEcontour.interval`

Energy 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.eta`

the 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**.

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**.

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.kSpaceGridNumber`

the small k-space grid number in each direction which, together with kSpaceGridShift, are used to produce the parameter kSpacePoints.`calculation.transmission.energyPoints`

The 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.coordinatesOfTheSymmetryKPoints`

The 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].

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**.

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.temperature`

The 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.chemicalPotentialEnergies`

The 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.

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.