3
$\begingroup$

I have a function ($a(x)=0.00005x^4-0.01x^3+0.68x^2-23.087x+322$) which defines a profile desired.

I need to get a new function that is shifted relative to the first function using a list (offset=${10,9.2,8.06,7.09,6.12,5.15,4.18,3.21,2.24,1.27,0.3}$) that determines points least for this new function. This means that my new function should go in addition to these points and may not be below these values. Also should not allow the existence of any inflection point.

Through my code is possible to have na idea of what I need, but how can see, I couldn't avoid the inflection points and my function is not well as all least points.

Another issue: I this code inters = Solve[a[x] == b, x], how so I use the value rational and positive? The code inters = x /. inters[[4]], the value $4$ insert manually.

Clear["Global`*"];
a[x_]:=0.00005x^4-0.01x^3+0.68x^2-23.087x+322;
ent=Range[0,20,2];
xinit=a[#]&/@ent//MatrixForm;
point1={ent[[#]],xinit[[1,#]]}&/@Range[Length[ent]];
offset={10,9.2,8.06,7.09,6.12,5.15,4.18,3.21,2.24,1.27,0.3};
newdata=ent[[#]]+offset[[#]]&/@Range[Length[ent]];
point2={newdata[[#]],xinit[[1,#]]}&/@Range[Length[ent]];
b=Fit[point2,{1,x,x^2,x^3,x^4},x];
inters=Solve[a[x]==b,x]
inters=x/. inters[[4]];
Plot[{a[x],b},{x,0,60},PlotRange->{{0,inters},{0,point1[[1,2]]}},ImageSize->1100,AspectRatio->1,Epilog->{Red,PointSize[0.01],Point[{point1}],Point[{point2}],Black,Point[{x,b}/.Solve[D[b,{x,2}]==0]],Purple,Dashed,Line[{point1[[#]],point2[[#]]}]&/@Range[11]}]

enter image description here

$\endgroup$
2
  • $\begingroup$ In short: your offsets are horizontal offsets, as opposed to vertical or orthogonal offsets? $\endgroup$ Commented Jul 4, 2016 at 18:17
  • $\begingroup$ Yes. These are horizontal. $\endgroup$
    – LCarvalho
    Commented Jul 4, 2016 at 20:25

1 Answer 1

2
$\begingroup$

Updated based on comments:

Clear["Global`*"];

curve[x_] = 0.00005 x^4 - 0.01 x^3 + 0.68 x^2 - 23.087 x + 322;
offset = {10, 9.2, 8.06, 7.09, 6.12, 5.15, 4.18, 3.21, 2.24, 1.27, 0.3};

ent = Range[0, 20, 2];

curveList = Partition[Riffle[ent, curve[ent]], 2];
curveOffset = Partition[Riffle[ent + offset, curve[ent]], 2];

penalty[residuals_?VectorQ, scale_: 10] := 
  scale*Length@residuals*(1 - Sign@Min[residuals]);

scale = Max[curveOffset] - Min[curveOffset];
model = a x^4 + b x^3 + c x^2 + d x + e;
fit = FindFit[curveOffset, model, {a, b, c, d, e}, x, 
   NormFunction -> (Norm[#1] + penalty[#1, scale] &), 
   Method -> "NMinimize", MaxIterations -> 10000];
modelf = Function[{x}, Evaluate[model /. fit]]

{tl, yl} = Transpose[curveOffset];
residuals = yl - Map[modelf, tl];
ListPlot[residuals, Filling -> Axis]

Show[
 ListPlot[{curveList, curveOffset}],
 Plot[modelf[x], {x, Min[ent + offset], Max[ent + offset]}]
 ]

Residuals plotted to show horizontal offset points are the minimum surface.

enter image description here

142.235 + 149.992 x - 22.4808 x^2 + 1.12361 x^3 - 0.0192096 x^4

Original and offset data overlaid with the above fit equation:

enter image description here

Reference: How to fit a function to data so that the fit is always greater than or equal to the data?

$\endgroup$
5
  • $\begingroup$ It is necessary that the curve is to the right of the minimum points. This list that I received define the lowest possible values ​​X may shift these values ​​to be larger so that there is inflection points $\endgroup$
    – LCarvalho
    Commented Jul 4, 2016 at 20:46
  • $\begingroup$ Sorry if the grammar is bad. English is not my native language. $\endgroup$
    – LCarvalho
    Commented Jul 4, 2016 at 20:50
  • $\begingroup$ @LeandroMacieldeCarvalho Let me know if the newly edited answer is what you wanted. $\endgroup$
    – Young
    Commented Jul 4, 2016 at 22:43
  • $\begingroup$ When I executed your code occured the following error: The model ax^4+bx^3+cx^2+dx+e is linear in the parameters \ {a,b,c,d,e}, but a nonlinear method or non-Euclidean norm was \ specified, so nonlinear methods will be used. Would be that´s right? $\endgroup$
    – LCarvalho
    Commented Jul 5, 2016 at 21:32
  • 1
    $\begingroup$ It's just a warning, but it's ok $\endgroup$
    – Young
    Commented Jul 5, 2016 at 21:39

Your Answer

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

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