\documentclass{article}

\usepackage{graphicx}
\usepackage{setspace}
\usepackage{listings}
\usepackage{color}
\usepackage{amsmath}
\usepackage{float}
\usepackage{caption}
\usepackage{subcaption}
\usepackage[margin=0.75in]{geometry}

\renewcommand{\thesubsection}{\indent(\alph{subsection})}

\definecolor{dkgreen}{rgb}{0,0.6,0}
\definecolor{gray}{rgb}{0.5,0.5,0.5}
\definecolor{mauve}{rgb}{0.58,0,0.82}

%\lstset{language=Matlab,%
%    %basicstyle=\color{red},
%    breaklines=true,%
%    morekeywords={matlab2tikz},
%    keywordstyle=\color{blue},%
%    morekeywords=[2]{1}, keywordstyle=[2]{\color{black}},
%    identifierstyle=\color{black},%
%    stringstyle=\color{mylilas},
%    commentstyle=\color{mygreen},%
%    showstringspaces=false,%without this there will be a symbol in the places where there is a space
%    numbers=left,%
%    numberstyle={\tiny \color{black}},% size of the numbers
%    numbersep=9pt, % this defines how far the numbers are from the text
%    emph=[1]{for,end,break},emphstyle=[1]\color{red}, %some words to emphasise
%    %emph=[2]{word1,word2}, emphstyle=[2]{style},    
%}
\lstset{basicstyle=\small,
        keywordstyle=\color{mauve},
        identifierstyle=\color{dkgreen},
        stringstyle=\color{gray},
        numbers=left,
        xleftmargin=5em
        }

\title{ECE 456 - Problem Set 1}
\date{2021-02-06}
\author{David Lenfesty \\ lenfesty@ualberta.ca
    \and Phillip Kirwin \\ pkirwin@ualberta.ca}

\setcounter{tocdepth}{2} % Show subsections

\begin{document}

    \doublespacing
    \pagenumbering{gobble}
    \maketitle
    \newpage

    \singlespacing
    \pagenumbering{arabic}

    \section*{Question 1}
        \subsection*{(a)}
        Beginning with the following two equations:
        \begin{equation}
            \label{eq:N_old}
            N = \int_{-\infty}^{\infty}\frac{\gamma_1 f_1(E) + \gamma_2 f_2(E)}{\gamma_1 + \gamma_2} D(E-U) dE,
        \end{equation}

        \begin{equation}
            \label{eq:I_old}
            I = \frac{q}{\hbar}\frac{\gamma_1 \gamma_2}{\gamma_1 + \gamma_2} \int_{-\infty}^{\infty}[f_1(E) - f_2(E)] D(E-U) dE,
        \end{equation}
        
        and changing the variable of integration to $E' = E - U$:
        
        \begin{equation*}
            N = \int_{-\infty}^{\infty}\frac{\gamma_1 f_1(E' + U) + \gamma_2 f_2(E' + U)}{\gamma_1 + \gamma_2} D(E') dE',
        \end{equation*}

        \begin{equation*}
            I = \frac{q}{\hbar}\frac{\gamma_1 \gamma_2}{\gamma_1 + \gamma_2} \int_{-\infty}^{\infty}[f_1(E' + U) - f_2(E' + U)] D(E') dE'.
        \end{equation*}
        
        Replacing $E' \rightarrow E$, we obtain equations \ref{eq:N_new} and \ref{eq:I_new}.

        \begin{equation}
            \label{eq:N_new}
            N = \int_{-\infty}^{\infty}\frac{\gamma_1 f_1(E + U) + \gamma_2 f_2(E + U)}{\gamma_1 + \gamma_2} DE dE,
        \end{equation}

        \begin{equation}
            \label{eq:I_new}
            I = \frac{q}{\hbar}\frac{\gamma_1 \gamma_2}{\gamma_1 + \gamma_2} \int_{-\infty}^{\infty}[f_1(E + U) - f_2(E + U)] DE dE,
        \end{equation}

        \subsection*{(b)}

        For the provided constants, the plots of the number of channel electrons and the channel current follow:


        \begin{figure}[H]
            \centering
            \begin{subfigure}{0.5\textwidth}
                % Mess around with widths later
                \centering
                \includegraphics[width=\textwidth]{q1b_electrons.png}
                \caption{Plot of channel electrons vs. drain voltage.}
                \label{fig:q1b_electrons}
                % Note that the comment after \end{subfigure} is required for side by side figures
            \end{subfigure}%
            \begin{subfigure}{0.5\textwidth}
                % Mess around with widths later
                \centering
                \includegraphics[width=\textwidth]{q1b_current.png}
                \caption{Plot of channel current vs. drain voltage.}
                \label{fig:q1b_current}
            \end{subfigure}

            \caption{Number of electrons and current versus drain voltage.}
        \end{figure}

            Below is our code. Note that some variable names are different from those in the example code.
            \begin{lstlisting}[language=Matlab]
            
clear all;

%% Constants

% Physical constants
hbar = 1.052e-34;

% Single-charge coupling energy (eV)
U_0 = 0.25;
% (eV)
kBT = 0.025;
% Contact coupling coefficients (eV)
gamma_1 = 0.005;
gamma_2 = gamma_1;
gamma_sum = gamma_1 + gamma_2;
% Capacitive gate coefficient
a_G = 0.5;
% Capacitive drain coefficient
a_D = 0.5;
a_S = 1 - a_G - a_D;

% Central energy level
mu = 0;

% Energy grid, from -1eV to 1eV
NE = 501;
E = linspace(-1, 1, NE);
dE = E(2) - E(1);
% TODO name this better
cal_E = 0.2;

% Lorentzian density of states, normalized so the integral is 1
D = (gamma_sum / (2*pi)) ./ ( (E-cal_E).^2 + (gamma_sum/2).^2 );
D = D ./ (dE*sum(D));

% Reference no. of electrons in channel
N_0 = 0;

voltages = linspace(0, 1, 101);

% Terminal Voltages
V_G = 0;
V_S = 0;

for n = 1:length(voltages)
    % Set varying drain voltage
    V_D = voltages(n);

    % Shifted energy levels of the contacts
    mu_1 = mu - V_S;
    mu_2 = mu - V_D;

    % Laplace potential, does not change as solution is found (eV)
    % q is factored out here, we are working in eV
    U_L = - (a_G*V_G) - (a_D*V_D) - (a_S*V_S);

    % Poisson potential must change, assume 0 initially (eV)
    U_P = 0;

    % Assume large rate of change
    dU_P = 1;

    % Run until we get close enough to the answer
    while dU_P > 1e-6
        % source Fermi function
        f_1 = 1 ./ (1 + exp((E + U_L + U_P - mu_1) ./ kBT));
        % drain Fermi function
        f_2 = 1 ./ (1 + exp((E + U_L + U_P - mu_2) ./ kBT));

        % Update channel electrons against potential
        N(n) = dE * sum( ((gamma_1/gamma_sum) .* f_1 + (gamma_2/gamma_sum) .* f_2) .* D);

        % Re-update Poisson portion of potential
        tmpU_P = U_0 * ( N(n) - N_0);
        dU_P = abs(U_P - tmpU_P);

        % Unsure why U_P is updated incrementally, perhaps to avoid oscillations?
        %U_P = tmpU_P;
        %U_P = U_P + 0.1 * (tmpU_P - U_P);
    end

    % Calculate current based on solved potential.
    % Note: f1 is dependent on changes in U but has been updated prior in the loop
    I(n) = q * (q/hbar) * (gamma_1 * gamma_1 / gamma_sum) * dE * sum((f_1-f_2).*D);

end


%%Plotting commands

figure(1); 
h = plot(voltages, N,'k'); 
grid on;
set(h,'linewidth',[2.0]); 
set(gca,'Fontsize',[18]);
xlabel('Drain voltage [V]'); 
ylabel('Number of electrons');

figure(2); 
h = plot(voltages, I,'k'); 
grid on;
set(h,'linewidth',[2.0]); 
set(gca,'Fontsize',[18]);
xlabel('Drain voltage [V]'); 
ylabel('Current [A]');
            \end{lstlisting}

        \subsection*{(c)}

        \begin{figure}[H]
            \centering
            \begin{subfigure}{0.5\textwidth}
                % Mess around with widths later
                \centering
                \includegraphics[width=\textwidth]{q1c_1.png}
                \caption{Plot of channel electrons vs. drain voltage.}
                \label{fig:q1c_1}
                % Note that the comment after \end{subfigure} is required for side by side figures
            \end{subfigure}%
            \begin{subfigure}{0.5\textwidth}
                % Mess around with widths later
                \centering
                \includegraphics[width=\textwidth]{q1c_2.png}
                \caption{Plot of channel electrons vs. drain voltage.}
                \label{fig:q1c_1}
                % Note that the comment after \end{subfigure} is required for side by side figures
            \end{subfigure}
            \begin{subfigure}{0.5\textwidth}
                % Mess around with widths later
                \centering
                \includegraphics[width=\textwidth]{q1c_3.png}
                \caption{Plot of channel electrons vs. drain voltage.}
                \label{fig:q1c_1}
                % Note that the comment after \end{subfigure} is required for side by side figures
            \end{subfigure}%
            \begin{subfigure}{0.5\textwidth}
                % Mess around with widths later
                \centering
                \includegraphics[width=\textwidth]{q1c_4.png}
                \caption{Plot of channel electrons vs. drain voltage.}
                \label{fig:q1c_1}
                % Note that the comment after \end{subfigure} is required for side by side figures
            \end{subfigure}
            \begin{subfigure}{0.5\textwidth}
                % Mess around with widths later
                \centering
                \includegraphics[width=\textwidth]{q1c_5.png}
                \caption{Plot of channel electrons vs. drain voltage.}
                \label{fig:q1c_1}
                % Note that the comment after \end{subfigure} is required for side by side figures
            \end{subfigure}%

            \caption{Visual representation of the fermi functions of the contacts and channel.}
        \end{figure}



% TODO appendix for Part C code

\end{document}