0

Suppose I have the following code:

ClearAll["Global`*"]
q := 1.6*10^-19; (* Electron charge in Coulomb *)
me := 9.1*10^-31; (* Free electron rest mass in kg *)
h := 6.63*10^-34;  (* Reduced Planck's constant in J.s *)
kb := 1.38*10^-23; (* Boltzmann constant in J/K *)
FD[d_, η_] := -PolyLog[
    d + 1, -E^η];(* Defining the Fermi-Dirac integrals *)
Nc[d_, t_, α_, gs_, gv_, T_] := 
  2*gs*gv*((kb*T)/α)^(d/t)* 1/(2*π^0.5)^d*1/t*Gamma[d/t]/
   Gamma[d/2];  (* Effective band-edge DOS in d dimensions *)
n[d_, t_, α_, gs_, gv_, T_, ηF_] := 
  Nc[d, t, α, gs, gv, T]*FD[(d - t)/t, ηF]*(100)^-d; 
(* Effective SI carrier density in d dimensions in cm units*)
ηS[d_, t_, α_, gs_, gv_, T_, v_, nd_] := 
 Quiet[Chop[
     FindRoot[
      1/2*(n[d, t, α, gs, gv, T, η] + 
          n[d, t, α, gs, gv, T, η - (q*v)/(kb*T)]) == 
       nd, {η, 10000}]][[1]][[
   2]]]; (* Source Fermi Level at voltage v in d dimensions*)
J0[d_, t_, α_, gs_, gv_, T_] := 
 gs*gv*q^2/h*((kb*T)/α)^((d - 1)/t)* 1/(2*π^0.5)^(d - 1)*(
  kb*T)/q*Gamma[(d - 1 + t)/t]/Gamma[(d + 1)/2];
Jcore[d_, t_, α_, gs_, gv_, T_, v_, nd_] := 
  FD[(d - 1)/t, ηS[d, t, α, gs, gv, T, v, nd]] -FD[(
    d - 1)/t, ηS[d, t, α, gs, gv, T, v, nd] - (q*v)/(
     kb*T)];
Jall[d_, t_, α_, gs_, gv_, T_, v_, nd_] := 
  J0[d, t, α, gs, gv, T]* 
   Jcore[d, t, α, gs, gv, T, v, nd];
vinj[d_, t_, α_, gs_, gv_, T_, nd_, Eop_] := 
  Jall[d, t, α, gs, gv, T, Eop/q, nd]/(q*nd);

On plotting vinj for some parameters mentioned below you should get the following plot:

v1 = LogLinearPlot[{vinj[2, 2, h^2/(4*π^2*2*0.2*me), 2, 1, 300, 
     n2d, 0.092*q]*10^-9}, {n2d, 10^10, 10^14}, Frame -> True, 
  FrameTicks -> Automatic, 
  FrameLabel -> {"Carrier Density (1/\!\(\*SuperscriptBox[\(cm\), \(2\
\)]\))", "\!\(\*SubscriptBox[\(v\), \(inj\)]\) (\!\(\*SuperscriptBox[\
\(10\), \(7\)]\)cm/s)"}, BaseStyle -> {FontSize -> 15}, 
  PlotRange -> Full, AxesOrigin -> {Automatic, 0}, 
  AspectRatio -> GoldenRatio, PlotStyle -> {Thick, Blue}]

enter image description here

As you can see there is a maxima in the value of vinj as a function of n2d. How do I find this maxima numerically? Is there some function which can directly extract it? I tried FindMaxValue and NMaximize but both of them didn't work. Any suggestions are welcome.

A little mouse on the pampas
  • 5,760
  • 2
  • 13
  • 40
Indeterminate
  • 601
  • 3
  • 9
  • https://mathematica.stackexchange.com/a/26037/4999 might be helpful in getting FindMaxValue or FindMaximum to work. – Michael E2 Aug 30 '20 at 02:59

3 Answers3

3
ClearAll["Global`*"]
q = 1.6*10^-19;(*Electron charge in Coulomb*)
me = 
 9.1*10^-31;(*Free electron rest mass in kg*)
h = 
 6.63*10^-34;(*Reduced Planck's constant in J.s*)
kb = 
 1.38*10^-23;(*Boltzmann constant in J/K*)

FD[d_, η_] := -PolyLog[
   d + 1, -E^η];(*Defining the Fermi-Dirac integrals*)

Nc[d_, t_, α_, gs_, gv_, T_] := 
 2*gs*gv*((kb*T)/α)^(d/t)*1/(2*π^(1/2))^d*1/t*
  Gamma[d/t]/Gamma[d/2];(*Effective band-edge DOS in d dimensions*)
n[d_, t_, α_, gs_, gv_, T_, ηF_] := 
 Nc[d, t, α, gs, gv, T]*FD[(d - t)/t, ηF]*(100)^-d;
(*Effective SI carrier density in d dimensions in cm units*)

Note that since ηS uses a numeric technique (FindRoot), its arguments should be restricted to numeric values using NumericQ

ηS[d_?NumericQ, t_?NumericQ, α_?NumericQ, gs_?NumericQ, 
  gv_?NumericQ, T_?NumericQ, v_?NumericQ, nd_?NumericQ] := 
 Quiet[Chop[
     FindRoot[1/
         2*(n[d, t, α, gs, gv, T, η] + 
          n[d, t, α, gs, gv, T, η - (q*v)/(kb*T)]) == nd, {η,
        10000}]][[1]][[2]]];(*Source Fermi Level at voltage v in d dimensions*)

J0[d_, t_, α, gs, gv_, T_] := 
  gsgvq^2/h((kbT)/α)^((d - 1)/t)1/(2π^0.5)^(d - 1)(kbT)/q*
   Gamma[(d - 1 + t)/t]/Gamma[(d + 1)/2];
Jcore[d_, t_, α, gs, gv_, T_, v_, nd_] := 
  FD[(d - 1)/t, ηS[d, t, α, gs, gv, T, v, nd]] - 
   FD[(d - 1)/t, ηS[d, t, α, gs, gv, T, v, nd] - (qv)/(kbT)];
Jall[d_, t_, α, gs, gv_, T_, v_, nd_] := 
  J0[d, t, α, gs, gv, T]Jcore[d, t, α, gs, gv, T, v, nd];
vinj[d_, t_, α, gs, gv_, T_, nd_, Eop_] := 
  Jall[d, t, α, gs, gv, T, Eop/q, nd]/(qnd);

To find the maximum and argument use Maximize

Maximize[{vinj[2, 2, h^2/(4*π^2*2*(2/10)*me), 2, 1, 300, 
    n2d, (92/1000)*q]*10^-9, 10^10 < n2d < 10^14}, n2d]

(* {1.37927, {n2d -> 3.3406310^12}} )

For just the maximum use MaxValue

MaxValue[{vinj[2, 2, h^2/(4*π^2*2*(2/10)*me), 2, 1, 300, 
    n2d, (92/1000)*q]*10^-9, 10^10 < n2d < 10^14}, n2d]

(* 1.37927 *)

For just the argument of the maximum use ArgMax

ArgMax[{vinj[2, 2, h^2/(4*π^2*2*(2/10)*me), 2, 1, 300, 
    n2d, (92/1000)*q]*10^-9, 10^10 < n2d < 10^14}, n2d]

(* 3.3406310^12 )
Bob Hanlon
  • 157,611
  • 7
  • 77
  • 198
  • Thanks a lot. Though I don't understand what does NumericQ doing? – Indeterminate Aug 30 '20 at 04:19
  • NumericQ keeps the function from trying to evaluate with non-numeric arguments. See the "Applications" section of the documentation for NumericQ. Presumably, one of the reasons that you included Quiet in the function definition is that you got error messages when the function tried to evaluate with symbolic arguments. – Bob Hanlon Aug 30 '20 at 04:32
  • I see. That makes sense. Thanks a lot. – Indeterminate Aug 30 '20 at 04:45
2

This is my solution:

vinjnew[x_] := 
 Chop[vinj[2, 2, h^2/(4*\[Pi]^2*2*0.2*me), 2, 1, 300, x, 
    0.092*q]*10^-9]

vinjnew[#] & /@ Table[x, {x, 10^12, 10^13, 10^11}] // Max
vinjnew[#] & /@ Table[x, {x, 10^12, 10^13, 10^10}] // Max
vinjnew[#] & /@ Table[x, {x, 10^12, 10^13, 10^9}] // Max

The results of above codes are 1.37924, 1.37927, 1.37927 respectively. Hopefully you can try different precison like 10^8 to check the result.

I hope this can help.

Hao Wang
  • 46
  • 1
0

I used this method to find the value of x

vinjnew2[x_] := 
  Chop[vinj[2, 2, h^2/(4*\[Pi]^2*2*0.2*me), 2, 1, 300, x*10^12, 
     0.092*q]*10^-9];

Finder[x_] := Module[{para1, para2}, (
   para1 = x;
   If[Abs[vinjnew2[para1] - 1.37927] <= 0.000001, para2 = para1, 
    para2 = para1 + 0.01];
   Return[para2])];


Nest[Finder, 3, 1000]

I think you may change the precison of x by changing the vinjnew function (I set x*10^12 inside the function) and the factor in the Finder (I set +0.01 per iterartion )

Hao Wang
  • 46
  • 1