How to get polynom formula from drawing?

Question

I would like to have an analytic formula for certain functions that I can draw with a pen. For example, if I want a function like that then I can model it with the function $$ V(x) = (x-a)^2 (x-b)^2.$$ In this way I can model multiple minima on the same height. But how can I model a potential with multiple minima with different height and some high barriers, like this:

EDIT

I tried the method from mathreadler trying to construct a polynomial that looks like this: which is specified on nine points $x_1 = -2.7 ,\dots ,x_9 = 2$. I am looking for a polynom that satisfies $p(x_i) = a_i$ where $a_i$ can be found in the picture and $p'(x_i) = 0$. Also, according to the wish that $x_i$ is minimum or maximum I have conditions on $p''(x_i)$. These are 27 conditions, so that I have to find a polynom of degree 26 that fulfills these conditions. I wrote a short matlab code to construct the polynom as follows:

function findNicePoly()

h1=@(x)([x^26, x^25 , x^24 , x^23, x^22 , x^21 , x^20 ,x^19 ,x^18 ,x^17 ,x^16 ,x^15 , x^14 , x^13, x^12 ,x^11, x^10, x^9 ,x^8 , x^7 ,x^6 ,x^5 , x^4 ,x^3, x^2  ,x^1 , 1 ]);
h2=@(x)([26*x^25, 25*x^24 , 24*x^23 , 23*x^22, 22*x^21 , 21*x^20 , 20*x^19 ,19*x^18 ,18*x^17 ,17*x^16 ,16*x^15 ,15*x^14 , 14*x^13 , 13*x^12, 12*x^11 ,11*x^10, 10*x^9, 9*x^8 ,8*x^7 , 7*x^6 ,6*x^5 ,5*x^4 , 4*x^3 ,3*x^2, 2*x^1  ,1 , 0 ]);
h3=@(x)([25*26*x^24, 24*25*x^23 , 23*24*x^22 , 22*23*x^21, 21*22*x^20 , 20*21*x^19 , 19*20*x^18 ,18*19*x^17 ,17*18*x^16 ,16*17*x^15 ,15*16*x^14 ,14*15*x^13 , 13*14*x^12 , 12*13*x^11, 11*12*x^10 ,10*11*x^9, 9*10*x^8, 8*9*x^7 ,7*8*x^6 , 6*7*x^5 ,5*6*x^4 ,4*5*x^3 ,3*4*x^2 ,2*3*x^1, 2  ,0 , 0 ]);


M = zeros(27,27);

M(1,:)=h1(-2.7);
M(2,:)=h1(-2.5);
M(3,:)=h1(-2);
M(4,:)=h1(-1.5);
M(5,:)=h1(-1);
M(6,:)=h1(-0.5);
M(7,:)=h1(-0.3);
M(8,:)=h1(1);
M(9,:)=h1(2);

M(10,:)=h2(-2.7);
M(11,:)=h2(-2.5);
M(12,:)=h2(-2);
M(13,:)=h2(-1.5);
M(14,:)=h2(-1);
M(15,:)=h2(-0.5);
M(16,:)=h2(-0.3);
M(17,:)=h2(1);
M(18,:)=h2(2);

M(19,:)=h3(-2.7);
M(20,:)=h3(-2.5);
M(21,:)=h3(-2);
M(22,:)=h3(-1.5);
M(23,:)=h3(-1);
M(24,:)=h3(-0.5);
M(25,:)=h3(-0.3);
M(26,:)=h3(1);
M(27,:)=h3(2);

b=[1, 4, 2,3,2,3,1,8,1, 0,0,0,0,0,0,0,0,0, 1,-1,2,-2,2,-2,20,-50,10]'

a=M\b;
a

x=[-3:0.01:3];
y=polyval(a,x);
plot(x,y)

but the ploted polynom does not look as I wished for.

Is there a better method, or am I doing something completely wrong?

Answer 1

You can set that the first derivative should be $0$ at the points you are interested in. This gives one linear equation in the coefficients of the polynomial for each point.

For example third order polynomial $p(x)$ should have extremum at point $x_0$: $$p(x) = ax^3+bx^2+cx+d$$

$$p'(x_0) = 3a(x_0)^2+2bx_0 + c = 0$$

This equation is linear in a,b, and c. So if you gather enough constraints you can just build a large enough equation system and it will work.

---

Edit: you may want to ensure that you actually get a minima, you can do that by encouraging the values of the second derivative to be lower or higher than 0. But this can also be done with the same kind of approach: set $p''(x_0) = \pm \epsilon_0$ with the sign depending on maximum or minimum. Since derivative is linear with respect to the coefficients this will fit nicely into the linear equation system.

---

Another edit Let us define C to be the matrix $\bf C = \left[x^0,x^1,\cdots,x^n\right]$, $\bf x$ being the column vector of x:s for extrema. This way if we have a polynomial with coefficients into a vector $\bf v$ the values of $p(x_k) = ({\bf Cv})_k$

Now if we minimize $$O({\bf v}) = \|{\bf Wv}\|_2^2+ \|{\bf CDv}\|_2^2 + \|{\bf CD^2v} - {\bf e}\|_2^2$$ the first term is a regularization term with diagonal weight matrix $\bf W$ which pushes coefficients towards 0. Second term pushes the derivatives towards 0. Third term pushes second derivatives towards a vector $\bf e$ where each element is $\pm \epsilon$ with the sign depending on maxima or minima. We can also add a fourth term $\|{\bf Cv} - {\bf d} \|_2^2$, where $\bf d$ is a vector of values such that $p(x_k) = {\bf d}_k$ is encouraged (which you asked in the comment).

Just to show that it works. Here we want minima at 1,3,5 and maxima at 2,4. We don't put any constraints on the function values ( although we could do that if we want to ).

Answer 2

The V(x) you gave touches the x-axis but your graph is above x-axis.

I have not "designed " roots like this before, but believe it is plausible. I outline a numerical procedure to start with.

The way you have sketched, all the roots must be complex.

They are type 1 and type 2 complex roots ordered alternately in a waveform well above the x-axis.

If the roots are like

$$ a \pm i A, b\pm i B, c \pm i C, d \pm i D, e \pm i E $$ etc.,let the sums and products be $ s_i, p_i$ for

$ f(x)= (x^2- s_1 x +p_1)(x^2- s_2 x +p_2)(x^2- s_3 x +p_3)(x^2- s_4 x +p_4).. $

The usual check for two type of roots is:

Vanishing first and positive value for second derivative for type 1 roots, and

Vanishing first and negative value for second derivative for type 2 roots. This is not a root and in fact can be discarded.

Order them in ascending value of real parts and slightly perturb the real/imaginary part in the neighborhood changing the coefficients as needed using a graphing calculator.

This could give the polynomial graph approximately as needed.

However an analytical approach can be developed with some modification in line with the above (I 've never done this before and there could be a minor pitfall).

If 4 minima and 3 maxima are seen in the graph above x-axis, then start with an $8^{th}$ order polynomial. Use the condition that all the roots be complex and so on as above.

As a simple example sketching non-polynomial $ y(x)= 2 + \cos \pi x , $ upto a sufficiently large degree, complex roots with real parts near odd integers can be seen.

EDIT1:

f(x)= (2 x^4 - x^3 - 20 x^2 - 10 x + 100)/10

has complex roots very approximately at real part $x=\pm$ 2 expected by graphical

inspection. The computed roots are ~ $ ( -2.335 + 1.33 i), (2.585,.491 i). $

EDIT 2:

Eighth order may be tough to start with.The fourth order polynomial can be be design/tailored using Mathematica code as above described using DE :

$$ y'''(x) = A\, 4! $$

You can choose where you seed your complex roots.

f[x_] := (2 x^4 - x^3 - 20 x^2 - 10 x + 100)/10 g1 = Plot[f[x], {x,-4, 4}, GridLines -> Automatic, PlotStyle -> {Thick, Magenta}] NSolve[f[x] == 0, {x}] NDSolve[{D[Y[x], {x, 4}] == 4.8, Y[-2] == 8, Y'[-2] == 0, Y[2.5] == 2, Y'[2.5] == 0}, Y, {x, -4, 4}] y[x_] = Y[x] /. First[%]; g2 = Plot[y[x], {x, -4, 4}, GridLines -> Automatic, PlotStyle -> {Thick, Purple}] Show[g1, g2]

You can automate it for higher orders next following such a procedure.There is reasonable good match between them I hope.