0
$\begingroup$

I am trying to get the following function, to plot with a log y axis. However, when I replace Plot with LogPlot, the plot is not computed correctly. Instead, a plot with incorrect x-axis is returned.

Constants

au = QuantityMagnitude[UnitConvert[Quantity[1, "AstronomicalUnit"], "Meters"]]; 
c = QuantityMagnitude[UnitConvert[Quantity[1, "SpeedOfLight"], "MetersPerSecond"]]; 
Qpr = 1; 
Lsun = QuantityMagnitude[UnitConvert[Quantity[1, "SolarLuminosity"], "Watts"]]; 
Rsun = QuantityMagnitude[UnitConvert[Quantity[1, "SolarRadius"], "Meters"]]; 
Msun = QuantityMagnitude[UnitConvert[Quantity[1, "SolarMass"], "Kilograms"]]; 
G = QuantityMagnitude[UnitConvert[Quantity[1, "GravitationalConstant"], ("Meters"^2*"Newtons")/"Kilograms"^2]]; 
year = QuantityMagnitude[UnitConvert[Quantity[1, "Years"], "Seconds"]]; 
Myr = year*10^6; 
Gyr = year*10^9; 
Mwd = 0.6*Msun; 
Cst = 1.27; 
U = 1*10^17; 

Functions

L[t_] := (3.26*Lsun*(Mwd/(0.6*Msun)))/(0.1 + t/Myr)^1.18; 
Roche[dens_] := (0.65*Cst*Rsun*(Mwd/(0.6*Msun))^(1/3))/(dens/3000)^3^(-1); 
Papsis[t_] := a[t]*(1 - e[t]); 

Radiative Drag

RDdadtR\[Rho]a = -((3*L[t]*Qpr*(2 + 3*e[t]^2))/(c^2*(16*Pi*\[Rho]*Rast*a[t]*(1 - e[t]^2)^(3/2)))); 
RDdedtR\[Rho]a = -((15*L[t]*e[t])/(c^2*(32*Pi*Rast*\[Rho]*a[t]^2*Sqrt[1 - e[t]^2]))); 

RDsolR\[Rho]a = ParametricNDSolveValue[{Derivative[1][a][t] == RDdadtR\[Rho]a, Derivative[1][e][t] == RDdedtR\[Rho]a, a[0] == a0, e[0] == 0.3}, {a, e}, {t, 0, 9*Gyr}, 
    {Rast, \[Rho], a0}]; 

fRDticks = {{Automatic, Automatic}, {Charting`FindTicks[{0, 1}, {0, 1/Myr}], Automatic}}; 

Manipulate[Column[{Style["Radiative Drag Working Plot", Bold], Plot[fun[func, t]/scale[func], {t, 0, 9*Gyr}, FrameTicks -> fRDticks, 
     Epilog -> {Red, Dashed, InfiniteLine[{{0, Roche[\[Rho]]}, {10, Roche[\[Rho]]}}]}, PlotStyle -> {Directive[Blue, Thickness[0.01]]}], Style["Compiled Plot", Bold], 
    If[comp === {}, Plot[fun[func, t]/scale[func], {t, 0, 9*Gyr}, FrameTicks -> fRDticks, Epilog -> {Red, Dashed, InfiniteLine[{{0, Roche[\[Rho]]}, {10, Roche[\[Rho]]}}]}, 
      PlotStyle -> {Directive[Blue, Thickness[0.01]]}], Plot[comp/scale[func], {t, 0, 9*Gyr}, FrameTicks -> fRDticks, 
      Epilog -> {Red, Dashed, InfiniteLine[{{0, Roche[\[Rho]]}, {10, Roche[\[Rho]]}}]}, PlotStyle -> {Directive[Blue, Thickness[0.01]]}]]}], 
  {{func, 1}, {1 -> "a", 2 -> "e", 3 -> "q"}}, {{Rast, 0.005}, 0.0001, 0.1, 0.001, Appearance -> "Labeled"}, {{\[Rho], 3000}, 1000, 7000, 50, Appearance -> "Labeled"}, 
  {{a0, 10, "a0 (au)"}, 1, 20, 0.2, Appearance -> "Labeled"}, Button["Append", AppendTo[comp, fun[func, t]]], Button["Reset", comp = {}], 
  TrackedSymbols -> {func, Rast, \[Rho], a0}, Initialization :> {comp = {}, fun[sel_, t_] := Switch[sel, 1, RDsolR\[Rho]a[Rast, \[Rho], a0*au][[1]][t], 2, 
      RDsolR\[Rho]a[Rast, \[Rho], a0*au][[2]][t], 3, RDsolR\[Rho]a[Rast, \[Rho], a0*au][[1]][t]*(1 - RDsolR\[Rho]a[Rast, \[Rho], a0*au][[2]][t])], 
    scale[sel_] := Switch[sel, 1 | 3, au, 2, 1]}]

The question is- how do I get this plot to have a logarithmic y axis?

Thanks in advance.

EDIT: enter image description here

Improve this question
$\endgroup$
8
  • $\begingroup$ Your function changes very little. Therefore, a log axis looks similar to a linear axis. $\endgroup$
    – Daniel Huber
    Feb 16 at 15:02
  • $\begingroup$ @DanielHuber I still need a log axis as I'm going to be changing the parameter values such that the y axis will change quickly. LogPlot completely changes the graph and so I don't think it's the correct solution. $\endgroup$
    – testing09
    Feb 16 at 15:07
  • $\begingroup$ I can not see much of a difference when I change Plot and LogPlot. Can you give an example about what you mean? $\endgroup$
    – Daniel Huber
    Feb 16 at 15:21
  • $\begingroup$ @DanielHuber I've added a picture above. $\endgroup$
    – testing09
    Feb 16 at 17:51
  • $\begingroup$ But your plot have very little y change. In this case a log plot looks like a linear plot $\endgroup$
    – Daniel Huber
    Feb 16 at 19:14
0
$\begingroup$ 
au = QuantityMagnitude[
   UnitConvert[Quantity[1, "AstronomicalUnit"], "Meters"]];
c = QuantityMagnitude[
   UnitConvert[Quantity[1, "SpeedOfLight"], "MetersPerSecond"]];
Qpr = 1;
Lsun = QuantityMagnitude[
   UnitConvert[Quantity[1, "SolarLuminosity"], "Watts"]];
Rsun = QuantityMagnitude[
   UnitConvert[Quantity[1, "SolarRadius"], "Meters"]];
Msun = QuantityMagnitude[
   UnitConvert[Quantity[1, "SolarMass"], "Kilograms"]];
G = QuantityMagnitude[
   UnitConvert[
    Quantity[1, "GravitationalConstant"], ("Meters"^2*"Newtons")/
     "Kilograms"^2]];
year = QuantityMagnitude[UnitConvert[Quantity[1, "Years"], "Seconds"]];
Myr = year*10^6;
Gyr = year*10^9;
Mwd = 0.6*Msun;
Cst = 1.27;
U = 1*10^17;


L[t_] := (3.26*Lsun*(Mwd/(0.6*Msun)))/(0.1 + t/Myr)^1.18;
Roche[dens_] := (0.65*Cst*
     Rsun*(Mwd/(0.6*Msun))^(1/3))/(dens/3000)^3^(-1);
Papsis[t_] := a[t]*(1 - e[t]);



RDdadtR\[Rho]a = -((3*L[t]*
       Qpr*(2 + 3*e[t]^2))/(c^2*(16*Pi*\[Rho]*Rast*
         a[t]*(1 - e[t]^2)^(3/2))));
RDdedtR\[Rho]a = -((15*L[t]*
       e[t])/(c^2*(32*Pi*Rast*\[Rho]*a[t]^2*Sqrt[1 - e[t]^2])));

RDsolR\[Rho]a = 
  ParametricNDSolveValue[{Derivative[1][a][t] == RDdadtR\[Rho]a, 
    Derivative[1][e][t] == RDdedtR\[Rho]a, a[0] == a0, 
    e[0] == 0.3}, {a, e}, {t, 0, 9*Gyr}, {Rast, \[Rho], a0}];

fRDticks = {{Automatic, 
    Automatic}, {Charting`FindTicks[{0, 1}, {0, 1/Myr}], Automatic}};

Manipulate[
 Column[{Style["Radiative Drag Working Plot", Bold], 
   LogPlot[fun[func, t]/scale[func], {t, 0, 9*Gyr}, 
    FrameTicks -> fRDticks, 
    Epilog -> {Red, Dashed, 
      InfiniteLine[{{0, Roche[\[Rho]]}, {10, Roche[\[Rho]]}}]}, 
    PlotStyle -> {Directive[Blue, Thickness[0.01]]}], 
   Style["Compiled Plot", Bold], 
   If[comp === {}, 
    LogPlot[fun[func, t]/scale[func], {t, 0, 9*Gyr}, 
     FrameTicks -> fRDticks, 
     Epilog -> {Red, Dashed, 
       InfiniteLine[{{0, Roche[\[Rho]]}, {10, Roche[\[Rho]]}}]}, 
     PlotStyle -> {Directive[Blue, Thickness[0.01]]}], 
    LogPlot[comp/scale[func], {t, 0, 9*Gyr}, FrameTicks -> fRDticks, 
     Epilog -> {Red, Dashed, 
       InfiniteLine[{{0, Roche[\[Rho]]}, {10, Roche[\[Rho]]}}]}, 
     PlotStyle -> {Directive[Blue, Thickness[0.01]]}]]}], {{func, 
   1}, {1 -> "a", 2 -> "e", 3 -> "q"}}, {{Rast, 0.005}, 0.0001, 0.1, 
  0.001, Appearance -> "Labeled"}, {{\[Rho], 3000}, 1000, 7000, 50, 
  Appearance -> "Labeled"}, {{a0, 10, "a0 (au)"}, 1, 20, 0.2, 
  Appearance -> "Labeled"}, 
 Button["Append", AppendTo[comp, fun[func, t]]], 
 Button["Reset", comp = {}], 
 TrackedSymbols -> {func, Rast, \[Rho], a0}, 
 Initialization :> {comp = {}, 
   fun[sel_, t_] := 
    Switch[sel, 1, RDsolR\[Rho]a[Rast, \[Rho], a0*au][[1]][t], 2, 
     RDsolR\[Rho]a[Rast, \[Rho], a0*au][[2]][t], 3, 
     RDsolR\[Rho]a[Rast, \[Rho], a0*au][[1]][
       t]*(1 - RDsolR\[Rho]a[Rast, \[Rho], a0*au][[2]][t])], 
   scale[sel_] := Switch[sel, 1 | 3, au, 2, 1]}]
Improve this answer
$\endgroup$
2
  • $\begingroup$ the x-axis has changed from Myr to seconds. $\endgroup$
    – testing09
    Feb 16 at 19:44
  • $\begingroup$ Where is the difference? $\endgroup$
    – Daniel Huber
    Feb 16 at 20:00

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.