2
$\begingroup$

I am trying to plot the derivative of function g[x] below where g[x] is defined as the root of another equation. However, I am getting nonsensical answers for the derivative.

Any advice would be much appreciated!

Simplified code is below:

jRa[T3_] := 10^-5/(1 + 5*10^-8/T3) 

jRb[x_, T3_] := 
 x*1.3*10^-4/(10^-4 + 
     x*(1 + 10^-6/(10^-5 - T3)))

g[x_?NumericQ] := (T3 /. 
    FindRoot[{jRa[T3] == jRb[x, T3]}, {{T3, 10^-7}}, 
     MaxIterations -> 1000])/Ttot

s1 = Plot[(g[x]), {x, 10^-7, 2*10^-5}, PlotRange -> All]

s2 = Plot[g'[x], {x, 10^-7, 2*10^-5}, PlotRange -> All]

g'[.9*10^-5]
$\endgroup$
3
  • $\begingroup$ Welcome to Mathematica.SE! I have formatted the question in the way that is customary on this site. You can click the EDIT link above to see how I did this so you can do it yourself next time. Also, can you please indicate what answers you are getting yourself, and (if necessary) explain why they are wrong? It might happen that the results are not the same on others' machines, so it is good practice to do this. Also, to make the code complete, provide the numerical value of Ttot (I know that this is just a small oversight). $\endgroup$
    – Szabolcs
    Oct 30, 2013 at 20:41
  • $\begingroup$ Now in this case, I can see the problem on my machine. The derivative Mma computes, using g'[x], is clearly incorrect and it's negative where it should be positive. It may be because of the unusual magnitude of the numbers you used, but it's clearly incorrect on my machine. The ND function from the <<NumericalCalculus` package is usually much more reliable. However, in this case it also gives a wrong result on my machine! I do not know why. Looks like a bug to me. $\endgroup$
    – Szabolcs
    Oct 30, 2013 at 20:43
  • $\begingroup$ For other SE members, here's what I get for the function and its derivatives using various methods: dropbox.com/s/7tnzgjfjt60ykpz/derivatives.pdf (Ttot = 1.0) $\endgroup$
    – Szabolcs
    Oct 30, 2013 at 20:46

2 Answers 2

2
$\begingroup$

It is a pain to see ND behave in a erratic manner for this example! However a quick fix is also within reach! We can use finite difference here. For more detail about this check here. In the code supplied by the OP I have assumed Ttot = 1.0.

n = 1000;
h = 2*10^-5/n;
grid = h Range[0, n];
(* Interpolating function of your derivative*)
fixder[fun_, degree_, grid_] := 
Interpolation@(Transpose@{grid, 
  NDSolve`FiniteDifferenceDerivative[Derivative[degree], grid,Map[fun, grid]]});
(* the first derivative  and the second one*)
{d1,d2} = fixder[g, #, grid]&/@{1,2};
GraphicsRow[{s1, Plot[d1[x], {x, 0, 2*10^-5}, PlotRange -> All], 
 Plot[d2[x], {x, 0, 2*10^-5}, PlotRange -> All]}, ImageSize -> 600]

enter image description here

$\endgroup$
2
$\begingroup$

You can use implicit differentiation to find the derivative of g:

Ttot = 1.0;

dg = (Dt[T3]/Ttot)/Dt[x] /. 
    First@Solve[Dt[jRa[T3] == jRb[x, T3]], Dt[T3]] /. T3 -> g[x] // Simplify

(* (99.6169 (5.*10^-8 + 1. g[x])^2 (1.*10^-10 - 0.00002 g[x] + 
     1. g[x]^2))/(3.83142*10^-21 + 8.42912*10^-17 x + 
   4.66092*10^-13 x^2 + (-7.66284*10^-16 - 1.6092*10^-11 x + 
      1.53257*10^-8 x^2) g[x] + (3.83142*10^-11 + 7.66284*10^-7 x + 1. x^2) g[x]^2) *)

Plot[g[x], {x, 10^-7, 2*10^-5}, PlotRange -> All]
Plot[dg, {x, 10^-7, 2*10^-5}, PlotRange -> All]

Mathematica graphics

Mathematica graphics

$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

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