Solve system of two variables involving 3rd order polynomials

[BasicPoke]

Hello, I am trying to understand how to solve a system of two variables (let's say s and p representing two physical quantities), where there is a third-order polynomial representing each. I'm not sure I am describing this correctly in words, but here is the system I need to solve.

P = a0 + a1*p + a2*p^2 + a3*p^3
a0 = b00 + b01*s + b02*s^2 + b03*s^3
a1 = b10 + b11*s + b12*s^2 + b13*s^3
a2 = b20 + b21*s + b22*s^2 + b23*s^3
a0 = b30 + b31*s + b32*s^2 + b33*s^3

I have seen some solutions to this mentioning matrix inversion, for example here:
automotive pressure sensor calibration

I would like to see an explanation of how to solve for the bxx constants. Any help is appreciated.
Thanks
Ron

 

[Ackbach]

Just to be clear here: you're trying to solve for $s$ and $p$? Everything else is known, correct? By the way, in $\LaTeX$, you could write your system as follows:
\begin{align*}
P&=a_{0}+a_{1}p+a_{2}p^{2}+a_{3}p^{3} \\
a_{0}&=b_{00}+b_{01}s+b_{02}s^{2}+b_{03}s^{3} \\
a_{1}&=b_{10}+b_{11}s+b_{12}s^{2}+b_{13}s^{3} \\
a_{2}&=b_{20}+b_{21}s+b_{22}s^{2}+b_{23}s^{3} \\
a_{3}&=b_{30}+b_{31}s+b_{32}s^{2}+b_{33}s^{3}.
\end{align*}

 

[I like Serena]

Hello, I am trying to understand how to solve a system of two variables (let's say s and p representing two physical quantities), where there is a third-order polynomial representing each. I'm not sure I am describing this correctly in words, but here is the system I need to solve.

P = a0 + a1*p + a2*p^2 + a3*p^3
a0 = b00 + b01*s + b02*s^2 + b03*s^3
a1 = b10 + b11*s + b12*s^2 + b13*s^3
a2 = b20 + b21*s + b22*s^2 + b23*s^3
a0 = b30 + b31*s + b32*s^2 + b33*s^3

I have seen some solutions to this mentioning matrix inversion, for example here:
automotive pressure sensor calibration

I would like to see an explanation of how to solve for the bxx constants. Any help is appreciated.
Thanks
Ron

Welcome to MHB, BasicPoke! :)

I believe you're talking about a cubic surface?

A point $\mathbf{\vec P}$ can then be written in matrix notation as:
$$\mathbf{\vec P} = \mathbf{\vec P}(p,s) = \begin{pmatrix}1 &p&p^2&p^3\end{pmatrix}
\begin{bmatrix}\mathbf{\vec b_{00}} & \dots & \mathbf{\vec b_{03}} \\
\vdots & & \vdots \\
\mathbf{\vec b_{30}} & \dots & \mathbf{\vec b_{33}} \end{bmatrix}
\begin{pmatrix}1 \\s\\s^2\\s^3\end{pmatrix}$$
where each of the coefficients in the matrix is a 3-dimensional vector.

Are you familiar with this notation?

 

[BasicPoke]

Just to be clear here: you're trying to solve for $s$ and $p$? Everything else is known, correct? By the way, in $\LaTeX$, you could write your system as follows:
\begin{align*}
P&=a_{0}+a_{1}p+a_{2}p^{2}+a_{3}p^{3} \\
a_{0}&=b_{00}+b_{01}s+b_{02}s^{2}+b_{03}s^{3} \\
a_{1}&=b_{10}+b_{11}s+b_{12}s^{2}+b_{13}s^{3} \\
a_{2}&=b_{20}+b_{21}s+b_{22}s^{2}+b_{23}s^{3} \\
a_{3}&=b_{30}+b_{31}s+b_{32}s^{2}+b_{33}s^{3}.
\end{align*}

That looks much better Ackbach. I don't know what LATEX is. It looks like I cannot see all of the buttons on this editor so can't really format anything. No, I need to find the b constants b00, b01, b02, etc. All of this is used to calibrate an instrument that measures p and s from non-ideal sensors, and must calculate the actual P. P depends on p and s. For calibration, we have a system that can control p and s. We measure 4 values of p at 4 different values of s. Sorry I am being mysterious. The result is a 3D surface for P, that is dependent on the measurement of p and s.

Physical calibration process:

The 3rd order polynomial for P is: P = a0 + a1*p + a2*p^2 + a3*p^3
To correct for s, a 3rd order polynomial in s should be run for each a0.

Control s of the system at the minimum design range s, record s0
Control p at minimum design range p, record p00
Control p at 1/3 design range p, record p10
Control p at 2/3 design range p, record p20
Control p maximum design range p, record p30

Control s of the system at 1/3 max design range s, record s1
Control p at minimum design range p, record p01
Control p at 1/3 design range p, record p11
Control p at 2/3 design range p, record p21
Control p maximum design range p, record p31

...two more sets to capture in all: s0-s3, all 16 values of pxx

Using these measurements, there is a way to solve for the bxx constants, involving matrix inversion, but I don't quite understand it. I hope this is making sense.
Ron

 

[BasicPoke]

Welcome to MHB, BasicPoke! :)

I believe you're talking about a cubic surface?

A point $\mathbf{\vec P}$ can then be written in matrix notation as:
$$\mathbf{\vec P} = \mathbf{\vec P}(p,s) = \begin{pmatrix}1 &p&p^2&p^3\end{pmatrix}
\begin{bmatrix}\mathbf{\vec b_{00}} & \dots & \mathbf{\vec b_{03}} \\
\vdots & & \vdots \\
\mathbf{\vec b_{30}} & \dots & \mathbf{\vec b_{33}} \end{bmatrix}
\begin{pmatrix}1 \\s\\s^2\\s^3\end{pmatrix}$$
where each of the coefficients in the matrix is a 3-dimensional vector.

Are you familiar with this notation?

I think you are onto something there. This looks similar to some of the documentation I have, except I am not familiar with vectors being elements of a matrix. If this were multiplied out, it is looking familiar. I was hoping someone has this process written out and doesn't have to go thru it from scratch for me. Anything you can show me would be great.
Ron

 

[I like Serena]

I think you are onto something there. This looks similar to some of the documentation I have, except I am not familiar with vectors being elements of a matrix. If this were multiplied out, it is looking familiar. I was hoping someone has this process written out and doesn't have to go thru it from scratch for me. Anything you can show me would be great.
Ron

Yeah. I know. ;)

Moving on, your 16 calibration points look in matrix notation like:
$$
\begin{bmatrix}\mathbf{\vec P_{00}} & \dots & \mathbf{\vec P_{03}} \\
\vdots & & \vdots \\
\mathbf{\vec P_{30}} & \dots & \mathbf{\vec P_{33}} \end{bmatrix}
=
\begin{pmatrix}1 &p_0&p_0^2&p_0^3 \\
\vdots &\vdots & \vdots & \vdots \\
1 &p_3&p_3^2&p_3^3 \\
\end{pmatrix}

\begin{bmatrix}\mathbf{\vec b_{00}} & \dots & \mathbf{\vec b_{03}} \\
\vdots & & \vdots \\
\mathbf{\vec b_{30}} & \dots & \mathbf{\vec b_{33}} \end{bmatrix}

\begin{pmatrix}1 & \dots &1 \\
s_0 & \dots & s_3 \\
s_0^2 & \dots & s_3^2 \\
s_0^3 & \dots & s_3^3 \\
\end{pmatrix}
$$
Note that your p values would have only 1 index instead of 2.

The solution is:
$$
\begin{bmatrix}\mathbf{\vec b_{00}} & \dots & \mathbf{\vec b_{03}} \\
\vdots & & \vdots \\
\mathbf{\vec b_{30}} & \dots & \mathbf{\vec b_{33}} \end{bmatrix}
=
\begin{pmatrix}1 &p_0&p_0^2&p_0^3 \\
\vdots & \vdots & \vdots & \vdots \\
1 &p_3&p_3^2&p_3^3 \\
\end{pmatrix}^{-1}

\begin{bmatrix}\mathbf{\vec P_{00}} & \dots & \mathbf{\vec P_{03}} \\
\vdots & & \vdots \\
\mathbf{\vec P_{30}} & \dots & \mathbf{\vec P_{33}} \end{bmatrix}

\begin{pmatrix}1 & \dots &1 \\
s_0 & \dots & s_3 \\
s_0^2 & \dots & s_3^2 \\
s_0^3 & \dots & s_3^3 \\
\end{pmatrix}^{-1}
$$

How much do you understand of it?
And what puzzles you?

 

[BasicPoke]

I like Serena, I don't know what it means to show a vector as an element in a matrix. Is there some other way you can write this (or part of it) for me, regarding your last post.
Ron

 

[BasicPoke]

Also, I don't know what you mean here by the notations Pxx and px. There are actually 16 different values of px (small p) read from the p sensor. There are 4 actual (target) values of P.
Ron

 

[I like Serena]

I like Serena, I don't know what it means to show a vector as an element in a matrix. Is there some other way you can write this (or part of it) for me, regarding your last post.
Ron

It means the same thing as when it would simply be a number.
The notation is only a short hand for your original set of equations.

If we simplify it a bit, the following equation
$$\begin{pmatrix}x \\ y\end{pmatrix} = \begin{bmatrix}a_{00} & a_{01} \\ a_{10} & a_{11} \end{bmatrix} \begin{pmatrix}1 \\ s\end{pmatrix}$$
is the matrix notation for
$$\left\{\begin{aligned}x&=a_{00} + a_{01}s \\ y&=a_{10} + a_{11}s\end{aligned}\right.$$

Are you still with me?

Now if I write
$$\begin{pmatrix}\mathbf{\vec x} \\ \mathbf{\vec y}\end{pmatrix} = \begin{bmatrix}\mathbf{\vec a_{00}} & \mathbf{\vec a_{01}} \\ \mathbf{\vec a_{10}} & \mathbf{\vec a_{11}}\end{bmatrix} \begin{pmatrix}1 \\ s\end{pmatrix}$$
in the same way, it is the matrix notation for
$$\left\{\begin{aligned}\mathbf{\vec x}&=\mathbf{\vec a_{00}} + \mathbf{\vec a_{01}}s \\ \mathbf{\vec y}&=\mathbf{\vec a_{10}} + \mathbf{\vec a_{11}}s\end{aligned}\right.$$

Does that make sense?

 

[I like Serena]

Also, I don't know what you mean here by the notations Pxx and px. There are actually 16 different values of px (small p) read from the p sensor. There are 4 actual (target) values of P.
Ron

Well, we'll need to agree on the symbols here.
We seem to have a mix-up between lowercase and uppercase letters.

What type of data does your p sensor yield?
Are those 3-dimensional points (16 of them)?

And if you have 4 actual values of P, do you mean:
- the minimum design range for the control parameter p,
- the 1/3 design range,
- the 2/3 design range,
- and the maximum design range?

 

[I like Serena]

The confusion seems to stem from using the letters p and P for the calibrations points as well as for the control parameter.
Now I've become confused which is which.
Perhaps we can use different letters altogether?
You can decide which letters, as long as we get rid of the confusion which is which.

 

[BasicPoke]

Well, we'll need to agree on the symbols here.
We seem to have a mix-up between lowercase and uppercase letters.

What type of data does your p sensor yield?
Are those 3-dimensional points (16 of them)?

And if you have 4 actual values of P, do you mean:
- the minimum design range for the control parameter p,
- the 1/3 design range,
- the 2/3 design range,
- and the maximum design range?

The output of the p sensor is a single number, at 16 different sets of conditions. Yes, you are correct on the 4 actual values of P. I suppose we are also using P to mean a calculated value in the final system (of which there is only one).
Ron

 

[I like Serena]

The output of the p sensor is a single number, at 16 different sets of conditions. Yes, you are correct on the 4 actual values of P. I suppose we are also using P to mean a calculated value in the final system (of which there is only one).
Ron

Okay.

So we have:

1. 2 control parameters, p and s, although you also seem to refer to them as P and s.
2. 16 calibration values that are given by the p sensor, that you seem to refer to as pxx.
3. whenever the 2 control parameters are set to specific values, you can calculate the value in the final system, which you refer to as P. This would be the value the the p sensor would give for those 2 control parameter values.

Perhaps we can rename the first control parameter to, say, $r$?
Then we can refer to the actual values of the control parameters as $r_0, r_1, r_2, r_3$ respectively $s_0, s_1, s_2, s_3$.
These would then be the minimum, 1/3, 2/3, and maximum design range values.

We can use pxx (or Pxx) for the 16 calibration values of the p sensor.
And P for the calculated value in the final system for specific values of u and s.

 

[BasicPoke]

Okay.

So we have:

1. 2 control parameters, p and s, although you also seem to refer to them as P and s.
2. 16 calibration values that are given by the p sensor, that you seem to refer to as pxx.
3. whenever the 2 control parameters are set to specific values, you can calculate the value in the final system, which you refer to as P. This would be the value the the p sensor would give for those 2 control parameter values.

Perhaps we can rename the first control parameter to, say, $r$?
Then we can refer to the actual values of the control parameters as $r_0, r_1, r_2, r_3$ respectively $s_0, s_1, s_2, s_3$.
These would then be the minimum, 1/3, 2/3, and maximum design range values.

We can use pxx (or Pxx) for the 16 calibration values of the p sensor.
And P for the calculated value in the final system for specific values of u and s.

We're not quite on the same page here on some of this.

2 control parameters, p and s, although you also seem to refer to them as P and s.

My intent was to refer to them as p and s.

16 calibration values that are given by the p sensor, that you seem to refer to as pxx.

Exactly.

whenever the 2 control parameters are set to specific values, you can calculate the value in the final system, which you refer to as P.

Exactly. I said there was only one of these, which is true in the final system, but in calibration & the equations there are 4 actual values: P0, P1, P2, P3.

This would be the value the the p sensor would give for those 2 control parameter values.

Not really. The p sensor yields one of the control parameters (p).

Let's try to stick with p and s as the control parameters, since I think these are consistent in my documentation, and I am going to get really confused otherwise. Here are my definitions:

p00, p01, ..., p33 = 16 readings from p sensor during calibration. Four values of P times four values of S.
s0, s1, s2, s3 = 4 readings from s sensor during calibration.
P0, P1, P2, P3 = Actual values of the quantity P used during calibration (use all 4 at each value of S).
S0, S1, S2, S3 = Actual values of the quantity S used during calibration.

Following is one set of equations I am working with. I don't understand how we get from the verbal description of the calibration to these equations. I think I can see how to solve these. k00, k01, k02, and k03 are what we are solving for. (Sorry, I don't know how to make these display horizontally)

at 1st control point for s sensor (S0),
\begin{bmatrix}
1\: p00\: p00^{2}\: p00^{3} \\
1\: p10\: p10^{2}\: p10^{3} \\
1\: p20\: p20^{2}\: p20^{3} \\
1\: p30\: p30^{2}\: p30^{3}
\end{bmatrix}
\begin{bmatrix}
K00\\
K10\\
K20\\
K30
\end{bmatrix}
=
\begin{bmatrix}
P0\\
P1\\
P2\\
P3
\end{bmatrix}

\begin{bmatrix}
1\: s0\: s0^{2}\: s0^{3} \\
1\: s1\: s1^{2}\: s1^{3} \\
1\: s2\: s2^{2}\: s2^{3} \\
1\: s3\: s3^{2}\: s3^{3}
\end{bmatrix}
\begin{bmatrix}
k00\\
k10\\
k20\\
k30
\end{bmatrix}
=
\begin{bmatrix}
K00\\
K10\\
K20\\
K30
\end{bmatrix}

 

[I like Serena]

Here are my definitions:

p00, p01, ..., p33 = 16 readings from p sensor during calibration. Four values of P times four values of S.
s0, s1, s2, s3 = 4 readings from s sensor during calibration.
P0, P1, P2, P3 = Actual values of the quantity P used during calibration (use all 4 at each value of S).
S0, S1, S2, S3 = Actual values of the quantity S used during calibration.

To be honest, I can still make little sense of this.
To me all the letters seem to be mixed up and have multiple meanings.

Let's start with S versus s.
What do you mean by actual values of the quantity S as opposed to s sensor readings?
How are they related to each other?
What kind of units are they measured in?

 

[BasicPoke]

To be honest, I can still make little sense of this.
To me all the letters seem to be mixed up and have multiple meanings.

Let's start with S versus s.
What do you mean by actual values of the quantity S as opposed to s sensor readings?
How are they related to each other?
What kind of units are they measured in?

I am understanding better the more I stare at this problem. During calibration of sensors 's' and 'p', there is an instrument that can accurately set the physical quantities S and P. Since these are accurate, I am calling them 'actual' values. When the instrument is set to a actual values, the s and p sensors are read (which are just analog-to-digital converter readings and are not accurate because of non-idealities in the sensors). The purpose of the calibration process is to find the coefficients required to convert s and p sensor readings to accurate values of P, given any actual S and P condition. [Sorry, I am trying to not go anywhere close to anything that might be proprietary.] 16 kxy coefficients are produced by the calibration process.

One question I have for you is this: If p, P, s, and k are matrices, and we have the equation p^-1 * P = s * k

Solving for k, is this correct:

k = s^-1 * p^-1 * P

I am not sure about the order.
Thanks

 

[I like Serena]

I am understanding better the more I stare at this problem. During calibration of sensors 's' and 'p', there is an instrument that can accurately set the physical quantities S and P. Since these are accurate, I am calling them 'actual' values. When the instrument is set to a actual values, the s and p sensors are read (which are just analog-to-digital converter readings and are not accurate because of non-idealities in the sensors). The purpose of the calibration process is to find the coefficients required to convert s and p sensor readings to accurate values of P, given any actual S and P condition. [Sorry, I am trying to not go anywhere close to anything that might be proprietary.] 16 kxy coefficients are produced by the calibration process.

Right.
I'm finally starting to make sense of this as well.

In total you have 16 p measurements.
But shouldn't you then have 16 s measurements as well?
How come you only have 4?

One question I have for you is this: If p, P, s, and k are matrices, and we have the equation p^-1 * P = s * k

Solving for k, is this correct:

k = s^-1 * p^-1 * P

I am not sure about the order.
Thanks

Yes. That is the correct order.
Let me show you why:
\begin{array}{lcl}
p^{-1} \cdot P &=& s \cdot k & \qquad \text{multiply both on the left by }s^{-1} \\
s^{-1} \cdot p^{-1} \cdot P &=& s^{-1} \cdot s \cdot k & \qquad \text{replace the multiplication of $s$ with its inverse by identity} \\
s^{-1} \cdot p^{-1}\cdot P &=& k
\end{array}

 

[BasicPoke]

Right.
I'm finally starting to make sense of this as well.

In total you have 16 p measurements.
But shouldn't you then have 16 s measurements as well?
How come you only have 4?

Yes. That is the correct order.
Let me show you why:
\begin{array}{lcl}
p^{-1} \cdot P &=& s \cdot k & \qquad \text{multiply both on the left by }s^{-1} \\
s^{-1} \cdot p^{-1} \cdot P &=& s^{-1} \cdot s \cdot k & \qquad \text{replace the multiplication of $s$ with its inverse by identity} \\
s^{-1} \cdot p^{-1}\cdot P &=& k
\end{array}

Thanks for the clarification...There only 4 s measurements because of the following process. s is assumed to be constant over each 4 recordings of p, but in reality it is not exactly the same. I think we could record all 16 values of s and not increase the complexity of the calculations, since the s values are going into a 4x4 matrix anyway.

Control the device at the minimum s value of the design range, record s0
Control the device at the minimum design range p, record p00
Control the device at 1/3 of the max design range p, record p10
Control the device at 2/3 of the max design range p, record p20
Control the device at maximum design range p, record p30

Control the device at 1/3 of the max design range s, record s1
Control the device at the minimum design range p, record p01
Control the device at 1/3 of the max design range p, record p11
Control the device at 2/3 of the max design range p, record p21
Control the device at maximum design range p, record p31

Control the device at 2/3 of the max design range s, record s2
Control the device at the minimum design range p, record p02
Control the device at 1/3 of the max design range p, record p12
Control the device at 2/3 of the max design range p, record p22
Control the device at maximum design range p, record p32

Control the device at maximum design range s, record s3
Control the device at the minimum design range p, record p03
Control the device at 1/3 of the max design range p, record p13
Control the device at 2/3 of the max design range p, record p23
Control the device at maximum design range p, record p33