0

NDsolve gives me a solution to a certain NL PDE. I want to evaluate the rate of exponential decay of the maxima of the solutions in time. So, I plot

Max1[t_?NumericQ] := First[Maximize[sol1[[1]][x, t], 0 <= x <= L, x]]
data1 = Table[{y, LMax1[y]}, {y, T0, Tfit}];
efit1 = FindFit[data1, Exp[a*t], {a}, t];
Plot[{Max1[y], Exp[efit1[1]*y]}, {y, T0, Tfit}, PlotRange -> All, Frame -> True]

and it looks like this:

Now, since I would like to fit this exponentially, or take the logarithm to evaluate the slope, how can I smooth that little jump which creates me a lot of trouble?

Here is what I get if I take the Log:

LMax1[t_?NumericQ] := Log[Abs[First[Maximize[sol1[[1]][x, t], 0 <= x <= L, x]]] - 1]
data1 = Table[{y, LMax1[y]}, {y, T0, Tfit}];
fit1 = Fit[data1, {1, y}, y];
Plot[{LMax1[y], fit1}, {y, T0, Tfit}]

enter image description here

usumdelphini
  • 996
  • 4
  • 14
  • You could try the HPFilter package, see this thread: http://mathematica.stackexchange.com/q/31629/131 (esp. the link to the updated package). – Yves Klett Nov 18 '13 at 12:18
  • 1
    ... and please add all code to reproduce your problem, like e.g. sol1. – Yves Klett Nov 18 '13 at 12:20
  • Do you expect the little jump at all ? Can it be that you need to make the solution "smoother" ? – b.gates.you.know.what Nov 18 '13 at 12:20
  • There was a topic about taking the function without sudden drops/peaks but I can't find it :/ – Kuba Nov 18 '13 at 12:23
  • @b.gatessucks I wouldn't expect that jump, but apart from the jump everything else suits my expections! – usumdelphini Nov 18 '13 at 12:23
  • @YvesKlett Just not to make the post too heavy, here is my nb, to reproduce the problem: https://www.dropbox.com/s/ek544yaceiwi9o8/mod_3.nb – usumdelphini Nov 18 '13 at 12:25
  • Could you perhaps post a small MWE to get a similar plot? Would be much more convenient and reproducible. – Yves Klett Nov 18 '13 at 12:27
  • @YvesKlett I hardly see how to write that, because sol1 comes from a complicated nonlinear PDE... sorry about that – usumdelphini Nov 18 '13 at 12:32
  • 1
    @YvesKlett i = Import["http://i.stack.imgur.com/u3y75.jpg"]; i1 = ImageTake[ColorReplace[i, Black], {20, 300}, {50, 500}]; l1 = PixelValuePositions[i1, Blue, .9]; ListPlot@l1 – Dr. belisarius Nov 18 '13 at 12:44
  • @YvesKlett I tried to apply the HPFilter, but when I call it with the code: testdata = Accumulate[RandomVariate[NormalDistribution[0, 0.5], 200]]; Get["HPFilter", Path ->"/Users/dc13661/Dropbox/Phd/1_High_density/Codes"] hptest = HPFilter[testdata, 1600]; ListLinePlot[{testdata, hptest}]

    It gives me an empty plot :(

     – usumdelphini Nov 18 '13 at 12:58
  • This works for me: testdata = Accumulate[RandomVariate[NormalDistribution[0, 0.5], 200]]; Get["http://www.verbeia.com/mathematica/mma/HPFilter.m"] hptest = HPFilter[testdata, 1600]; ListLinePlot[{testdata, hptest}] – Yves Klett Nov 18 '13 at 14:06
  • As it has been already proposed above you might think on improving the solving operator. If this does not work, I would say that the answer to your question is not that much depends on Mathematica, but on the nature of the problem. Why do not you take the smooth part of the curve (that is at x>0.05)? If the curve behavior at x>0.05 meets your expectations (which should be based on the problem nature), that will do. – Alexei Boulbitch Nov 18 '13 at 14:50
  • @AlexeiBoulbitch I cannot take the last part of the curve as it is, because when I take the Log, since everything is of order $10^{-7}$, Mathematica goes crazy... and I sincerely do not know what it does with the Log of a very small number :( – usumdelphini Nov 18 '13 at 15:58
  • @YvesKlett I get Get::noopen: Cannot open http://www.verbeia.com/mathematica/mma/HPFilter.m. >> – usumdelphini Nov 18 '13 at 16:00

0 Answers0