Ask Your Question
0

Error in polynomial interpolation

asked 2020-04-30 12:13:48 +0100

Cyrille gravatar image

updated 2020-05-01 07:29:36 +0100

I have edited this question to insist on what goes wrong

The following code ends with an indetermined parameter.

# Polynomial interpolation
# Passing the matrix
A = Matrix([[0,0],[0.5,0.25],[0.6,0.36],[2,4]])
var('a0, a1, a2, a3, a4')
B =[[A[i][1], a3*A[i][0]^3+ a2*A[i][0]^2+ a1*A[i][0]+ a0] for i in range(0,4)]
C=[[(B[i][0]-B[i][1])^2] for i in range(0,4)]
# sum of the square of errors
ssr=sum(C[i][0] for i in range(4))
# First order condition of optimality
ssr_a0=ssr.diff(a0)
ssr_a1=ssr.diff(a1)
ssr_a2=ssr.diff(a2)
ssr_a3=ssr.diff(a3)
ssr_a4=ssr.diff(a4)
# solution
sol=solve([ssr_a0==0,ssr_a1==0,ssr_a2==0,ssr_a3==0,ssr_a4==0], a0, a1, a2, a3, a4, solution_dict=True)
sol

The proposed solution being [{a0: 0, a1: 0, a2: 1, a3: 0, a4: r3}].

But the ssr_ai $i = 0,\ldots, 3$ are all linear. So sol should gives the same result that the following code

C=matrix([[8, 6.2, 9.22, 16.682], [6.2, 9.22,16.68,32.3842],[9.22,16.682,32.3842,64.21802],[16.68,32.3842,64.21802,128.124562]])
b=vector([9.22,16.68,32.3842,64.21902])
C.solve_right(b)

So why it is not the case ? (There is one and only one solution)

This lead to a second question : for the ssr_ai how can I select the parameters to construct the C matrix and the b vector not by hand as I have been obliged to do ?

edit retag flag offensive close merge delete

Comments

1

Please add what are your expectations. Deciphering the code allows to understand the expected result, but not what is in your mind.

Sébastien gravatar imageSébastien ( 2020-04-30 13:02:04 +0100 )edit

I should obtain a value for a0, a1, a3, and a4. So what is this result a4 : r3 ? What means r3 ? I already encounter such a result but I have forgoten how to resolve it.

Cyrille gravatar imageCyrille ( 2020-04-30 17:21:47 +0100 )edit

r3 is a parameter in the space of solution. In the documentation of solve, you can see this example: "If there is a parameter in the answer, that will show up as a new variable. In the following example, "r1" is an arbitrary constant (because of the "r")":

sage: forget()
sage: x, y = var('x,y')
sage: solve([x+y == 3, 2*x+2*y == 6],x,y)
[[x == -r1 + 3, y == r1]]
Sébastien gravatar imageSébastien ( 2020-04-30 18:03:52 +0100 )edit

Thanks Sebastien for your comment. But this is weird. My system is linear so either there is no solution either there is one. I have tried to substitute find_root() to solve but I suspect that it is only programmed to solve one variable equation. The solution would be to write the system in matrix notation but I do not know how to select the parameters to wtrite a matrix.

Cyrille gravatar imageCyrille ( 2020-05-01 07:00:28 +0100 )edit

"But this is weird. My system is linear so either there is no solution either there is one." This is not true as a system of linear equations can have 0, 1 or an infinity of solutions. See the example in my previous comment.

Sébastien gravatar imageSébastien ( 2020-05-01 09:34:10 +0100 )edit

1 Answer

Sort by » oldest newest most voted
1

answered 2020-05-01 09:53:45 +0100

Sébastien gravatar image

updated 2020-05-01 13:48:09 +0100

If you look at the system of linear equations, you can see that it really has 4 equations with 4 unknowns (the 5th equation is 0=0 and the fifth unknown a4 just doesn't appear)

sage: [ssr_a0==0,ssr_a1==0,ssr_a2==0,ssr_a3==0,ssr_a4==0]
[8*a0 + 6.20*a1 + 9.220000*a2 + 16.68200*a3 - 9.2200000 == 0,
 6.20000*a0 + 9.2200000*a1 + 16.6820000*a2 + 32.3842000*a3 - 16.682000 == 0,
 9.2200*a0 + 16.6820000*a1 + 32.384200*a2 + 64.21802000*a3 - 32.38420000 == 0,
 16.68200*a0 + 32.384200*a1 + 64.2180200*a2 + 128.12456200*a3 - 64.2180200 == 0,
 0 == 0]

That system has a solution involving only the first four unknown, which explains what you obtain:

[{a0: 0, a1: 0, a2: 1, a3: 0, a4: r3}]

Indeed, we can see that the coefficient of a2 matches the constant terms. Alternatively, you can just remove a4 in the call to solve since that variable does not appear in the system:

sage: solve([ssr_a0==0,ssr_a1==0,ssr_a2==0,ssr_a3==0,ssr_a4==0], a0, a1, a2, a3, solution_dict=True)
[{a0: 0, a1: 0, a2: 1, a3: 0}]

EDIT (answering a question in the comments): To transform the system of linear equations to matrices and solve the system, you may do:

sage: SSR = [ssr_a0,ssr_a1,ssr_a2,ssr_a3,ssr_a4]
sage: A = matrix([[ssr_a.coefficient(a) for a in [a0,a1,a2,a3]] for ssr_a in SSR])
sage: b = vector([-ssr_a.subs(a0=0,a1=0,a2=0,a3=0) for ssr_a in SSR])
sage: A \ b
(-4.81096644000283e-16, 8.08282642120849e-14, 0.999999999999817, 7.15728811560799e-14)
sage: A.change_ring(QQ) \ b.change_ring(QQ)  # or in the ring of rationals
(0, 0, 1, 0)
edit flag offensive delete link more

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 2020-04-30 12:13:48 +0100

Seen: 630 times

Last updated: May 01 '20