This is the equation that I want to plot in MATLAB -

I want to sweep the variable V_{bi} from `-5`

to `+5`

and then plot the variable `T`

for each individual value of V_{bi} (essentially `T`

vs. `Vbi`

) for `Nd = Na =`

10^{15}, 10^{16} and 10^{17}.

I know V_{bi} can be swept by creating a vector:

```
Vbi = -5 : 0.001 : 5;
```

But I am not sure how to go about solving this problem, as I have never faced something like this before. Can someone please advice on how to code this?

Simplest is to use a loop through all values you want, and use `fzero`

to find all values for the non-linear implicit function:

```
% Values for Na/Nd you want
Na = [1e15 1e16 1e17];
Nd = Na;
% Vbi values you want to sweep
VbiValues = -5 : 0.001 : 5;
% Initialize outputs
Tout = zeros(numel(VbiValues), numel(Na));
% Solve equation for each value of Vbi
for ii = 1:numel(VbiValues)
for jj = 1:numel(Na)
% Re-define equation for new value of Vbi, Na, Nd
eq = @(T) -T + 11603.24*VbiValues(ii) ./ log(Na(jj)*Nd(jj)/2.798e39 * 300./T .* exp(13452./T));
% Solve it
Tout(ii,jj) = fzero(eq, 10);
end
end
% plot results
plot(...
VbiValues, Tout(:,1), 'r', ...
VbiValues, Tout(:,2), 'g', ...
VbiValues, Tout(:,3), 'b')
L = legend(...
'$N_A = N_D = 10^{15}$',...
'$N_A = N_D = 10^{16}$',...
'$N_A = N_D = 10^{17}$');
set(L, 'Interpreter', 'LaTeX');
```

Alternatively, if you have the optimization toolbox, you can use `fsolve`

:

```
function topLevelFunction
% Values for Na/Nd you want
Na = [1e15 1e16 1e17];
Nd = Na;
% Vbi values you want to sweep
VbiValues = -5 : 0.01 : 5;
% Use fsolve() to solve the system in one go
Tout = fsolve(@(T)F(T, VbiValues, Nd), ones(numel(Nd),numel(VbiValues)));
% plot results
plot(...
VbiValues, Tout(1,:), 'r', ...
VbiValues, Tout(2,:), 'g', ...
VbiValues, Tout(3,:), 'b')
L = legend(...
'$N_A = N_D = 10^{15}$',...
'$N_A = N_D = 10^{16}$',...
'$N_A = N_D = 10^{17}$');
set(L, 'Interpreter', 'LaTeX');
end
function Fvals = F(Ttrials, Vbis, Nad)
% Output function values for
% - all Vbi (each column is a different Vbi)
% - all Na/Nd (each row is a different Na/Nd)
Fvals = -Ttrials + bsxfun(@rdivide, ...
11603.24*Vbis(:).', ...
log( bsxfun(@times, Nad(:).^2/2.798e39, 300./Ttrials .* exp(13452./Ttrials)) ));
end
```

Note that one of the two solutions (or both :) still contains an error, as the two plots are not the same, but I'm sure this will help you get started.