Got Questions? Get Answers.
Discover MakerZone

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi

Learn more

Discover what MATLAB® can do for your career.

Opportunities for recent engineering grads.

Apply Today

Thread Subject:
How to solve these equations?

Subject: How to solve these equations?

From: Liqing

Date: 18 Mar, 2008 16:23:03

Message: 1 of 13

a =13.9497;
b = 1.8350e+003;
c = 140.4557;
d =5.8077e-009;
k1_2 =9.3378e-007;
k2_2 =6.6098e-010;

syms x y
f1 = (b*d-k1_2*a) + (d+b+k1_2)*x + (b-d)*y + x^2 -y^2;
f2 = (c*d-k2_2*b) + (c-k2_2)*x + (c+d+k2_2)*y + x*y + y^2;
%--------------------------------------------------------------

[x,y]=solve(f1,f2);

What I got are not real numbers at all. But I'm sure these
equations should have real number solutions based on their
physical meaning.

What have I missed?

Thank you.

Subject: How to solve these equations?

From: Roger Stafford

Date: 18 Mar, 2008 16:55:03

Message: 2 of 13

"Liqing " <liqingj@gmail.com> wrote in message <froq97$n0c
$1@fred.mathworks.com>...
> a =13.9497;
> b = 1.8350e+003;
> c = 140.4557;
> d =5.8077e-009;
> k1_2 =9.3378e-007;
> k2_2 =6.6098e-010;
>
> syms x y
> f1 = (b*d-k1_2*a) + (d+b+k1_2)*x + (b-d)*y + x^2 -y^2;
> f2 = (c*d-k2_2*b) + (c-k2_2)*x + (c+d+k2_2)*y + x*y + y^2;
> %--------------------------------------------------------------
>
> [x,y]=solve(f1,f2);
>
> What I got are not real numbers at all. But I'm sure these
> equations should have real number solutions based on their
> physical meaning.
>
> What have I missed?
>
> Thank you.
------
  Unless I am greatly mistaken, you can obtain a cubic equation in y from
those two equations and there should therefore always be at least one real
root. Apparently 'solve' is not giving you all of them. Use the second
equation to solve for x in terms of y, and then substitute this in for x in the
first equation. Then multiply through by the denominator to obtain a
polynomial equation. The y^4 terms should cancel, leaving you with a cubic
equation. You can do all this with the Symbolic Toolbox, if you have it, to
obtain this cubic and then use 'roots' to solve it. That should give you all the
roots for y including at least one real one. Then use the above expression for
x in terms of y to obtain the corresponding value(s) of x.

Roger Stafford

Subject: How to solve these equations?

From: Roger Stafford

Date: 18 Mar, 2008 18:22:02

Message: 3 of 13

  If your Symbolic Toolbox is anything like mine, the result of

 [x,y]=solve(f1,f2);

should have y expressed as the solution to a cubic equation which is given
there, and x expressed in terms of y. You don't need to do the things I stated
earlier to get this cubic. You already have it from 'solve'. All you need to do is
solve it with 'roots'. I got three real solutions to it. I have no idea which of them
corresponds to physically realistic values, but none of them is complex-valued.

Roger Stafford

Subject: How to solve these equations?

From: Liqing

Date: 18 Mar, 2008 18:42:02

Message: 4 of 13

Hi Roger,

Thank you very much for the reply.

How did you get the 'roots' and the real number results?

The solution shoud be between -100 to +100 approximately.

Thank you.

LQ

"All you need to do is solve it with 'roots'. I got three
real solutions to it. I have no idea which of them
corresponds to physically realistic values, but none of them
is complex-valued."

Subject: How to solve these equations?

From: Roger Stafford

Date: 18 Mar, 2008 19:23:02

Message: 5 of 13

"Liqing " <liqingj@gmail.com> wrote in message <frp2dq$qmj
$1@fred.mathworks.com>...
> Hi Roger,
> Thank you very much for the reply.
> How did you get the 'roots' and the real number results?
> The solution shoud be between -100 to +100 approximately.
> Thank you.
> LQ
--------
  For simplicity I used

 f1 = A + B*x + C*y + x^2 - y^2;
 f2 = D + F*x + G*y + x*y + y^2;

at first for the equations. I obtained the cubic equation from the expression
for y from 'solve' in

 [x,y]=solve(f1,f2);

in terms of A, B, C, etc. On my version of the Symbolic Toolbox this is
expressed as "RootOf". You can see the cubic clearly if you use "pretty(y)". I
then substituted the numerical values you gave for a, b, c, etc. to evaluate A,
B, C, etc., then substituted these into the cubic coefficients, and finally
numerically solved the cubic using 'roots'. The three solutions were:

 y = -162.669, -140.456, and +2.615

You can get the corresponding x values from the 'solve' expression for x. Use
'pretty(x)' to see clearly how it depends on y.

Roger Stafford

Subject: How to solve these equations?

From: Liqing

Date: 18 Mar, 2008 20:55:04

Message: 6 of 13

Hi Roger,

Thank you so much! You are a genius! You saved my life by
helping me solve a great oceanographic question.

Sorry but could you send me the code? I still did not figure
out how to get the real numbers at the end. Please do not be
offended by my poor Matlab knowledge.

My email is liqingj@gmail.com.

Thank you very much.


LQ


"Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid>
wrote in message <frp4ql$qnu$1@fred.mathworks.com>...
> "Liqing " <liqingj@gmail.com> wrote in message <frp2dq$qmj
> $1@fred.mathworks.com>...
> > Hi Roger,
> > Thank you very much for the reply.
> > How did you get the 'roots' and the real number results?
> > The solution shoud be between -100 to +100 approximately.
> > Thank you.
> > LQ
> --------
> For simplicity I used
>
> f1 = A + B*x + C*y + x^2 - y^2;
> f2 = D + F*x + G*y + x*y + y^2;
>
> at first for the equations. I obtained the cubic equation
from the expression
> for y from 'solve' in
>
> [x,y]=solve(f1,f2);
>
> in terms of A, B, C, etc. On my version of the Symbolic
Toolbox this is
> expressed as "RootOf". You can see the cubic clearly if
you use "pretty(y)". I
> then substituted the numerical values you gave for a, b,
c, etc. to evaluate A,
> B, C, etc., then substituted these into the cubic
coefficients, and finally
> numerically solved the cubic using 'roots'. The three
solutions were:
>
> y = -162.669, -140.456, and +2.615
>
> You can get the corresponding x values from the 'solve'
expression for x. Use
> 'pretty(x)' to see clearly how it depends on y.
>
> Roger Stafford
>
>

Subject: How to solve these equations?

From: Liqing

Date: 18 Mar, 2008 20:55:04

Message: 7 of 13

Hi Roger,

Thank you so much! You are a genius! You saved my life by
helping me solve a great oceanographic question.

Sorry but could you send me the code? I still did not figure
out how to get the real numbers at the end. Please do not be
offended by my poor Matlab knowledge.

My email is liqingj@gmail.com.

Thank you very much.


LQ


"Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid>
wrote in message <frp4ql$qnu$1@fred.mathworks.com>...
> "Liqing " <liqingj@gmail.com> wrote in message <frp2dq$qmj
> $1@fred.mathworks.com>...
> > Hi Roger,
> > Thank you very much for the reply.
> > How did you get the 'roots' and the real number results?
> > The solution shoud be between -100 to +100 approximately.
> > Thank you.
> > LQ
> --------
> For simplicity I used
>
> f1 = A + B*x + C*y + x^2 - y^2;
> f2 = D + F*x + G*y + x*y + y^2;
>
> at first for the equations. I obtained the cubic equation
from the expression
> for y from 'solve' in
>
> [x,y]=solve(f1,f2);
>
> in terms of A, B, C, etc. On my version of the Symbolic
Toolbox this is
> expressed as "RootOf". You can see the cubic clearly if
you use "pretty(y)". I
> then substituted the numerical values you gave for a, b,
c, etc. to evaluate A,
> B, C, etc., then substituted these into the cubic
coefficients, and finally
> numerically solved the cubic using 'roots'. The three
solutions were:
>
> y = -162.669, -140.456, and +2.615
>
> You can get the corresponding x values from the 'solve'
expression for x. Use
> 'pretty(x)' to see clearly how it depends on y.
>
> Roger Stafford
>
>

Subject: How to solve these equations?

From: Roger Stafford

Date: 18 Mar, 2008 23:19:02

Message: 8 of 13

"Liqing " <liqingj@gmail.com> wrote in message <frpa78$c8s
$1@fred.mathworks.com>...
> Hi Roger,
> ......
> Sorry but could you send me the code? I still did not figure
> out how to get the real numbers at the end.
> ......
> LQ
-------
  Here are my matlab computations. Bear in mind that my matlab version is
very old, so the notation used for the Symbolic Toolbox expressions will be
somewhat different from those in newer versions.

  First I write the equations to be solved in this form for simplicity:

 f1 = 'A+B*x+C*y+x^2-y^2';
 f2 = 'D+F*x+G*y+x*y+y^2';
 [x,y]=solve(f1,f2);

The solution to y is this:

 y = RootOf((-B+C+2*G-2*F)*_Z^3 ...
           +(A-G*B-B*F+2*C*F+G^2+2*D-F^2)*_Z^2 ...
           +(2*A*F-B*G*F-B*D+C*F^2+2*G*D)*_Z ...
           +A*F^2-B*D*F+D^2)

which indicates the cubic coefficients involved in finding y. From pretty(x) (or
scratchings on the back of an envelope) can be deduced this equality:

 x = -(D+G*y+y^2)/(F+y)

  Now set about computing these expressions numerically:

 % The given constants
 a =13.9497;
 b = 1.8350e+003;
 c = 140.4557;
 d =5.8077e-009;
 k1_2 =9.3378e-007;
 k2_2 =6.6098e-010;

 % Evaluate the six above coefficients in the two equations
 A = b*d-k1_2*a;
 B = d+b+k1_2;
 C = b-d;
 D = c*d-k2_2*b;
 F = c-k2_2;
 G = c+d+k2_2;

 % The cubic coefficients from the "RootOf" expression above
 c1 = -B+C+2*G-2*F;
 c2 = A-G*B-B*F+2*C*F+G^2+2*D-F^2;
 c3 = 2*A*F-B*G*F-B*D+C*F^2+2*G*D;
 c4 = A*F^2-B*D*F+D^2;
 
 % Find the three roots for y values
 y = roots([c1 c2 c3 c4]);

 % Get the corresponding x values
 x = -(D+G*y+y.^2)./(F+y);

 % Results:
 x = 1.0e+02 *
   1.62669423030425
   1.40465468520115
  -0.02615082158746

 y = 1.0e+02 *
  -1.62669423100515
  -1.40455558867504
   0.02615082161392

  The x and y values are nearly the negative of one another. Note that the
second of the above x,y pairs satisfies the original equations very poorly,
which is apparently due to the fact that the divisor F+y is nearly zero for
them, and this seems to have caused excessively large computation errors.
The equation

 c1*y.^3+c2*y.^2+c3*y+c4

is reasonably close to zero for all three y values.

Roger Stafford

Subject: How to solve these equations?

From: Liqing

Date: 19 Mar, 2008 01:08:13

Message: 9 of 13

Dear Roger,

Thank you for your detailed code. I appreciate your
precious time very very much!

I wonder if I can ask a last question? When I use 'solve',
I got a y that has hundreds of thousands of characters
instead of the following 'rootof', which is so nice.

Do you know what might be the matter?


 y = RootOf((-B+C+2*G-2*F)*_Z^3 ...
           +(A-G*B-B*F+2*C*F+G^2+2*D-F^2)*_Z^2 ...
           +(2*A*F-B*G*F-B*D+C*F^2+2*G*D)*_Z ...
           +A*F^2-B*D*F+D^2)


--------------------------------------------------

"Roger Stafford"
<ellieandrogerxyzzy@mindspring.com.invalid> wrote in
message <frpil6$69i$1@fred.mathworks.com>...
> "Liqing " <liqingj@gmail.com> wrote in message <frpa78$c8s
> $1@fred.mathworks.com>...
> > Hi Roger,
> > ......
> > Sorry but could you send me the code? I still did not
figure
> > out how to get the real numbers at the end.
> > ......
> > LQ
> -------
> Here are my matlab computations. Bear in mind that my
matlab version is
> very old, so the notation used for the Symbolic Toolbox
expressions will be
> somewhat different from those in newer versions.
>
> First I write the equations to be solved in this form
for simplicity:
>
> f1 = 'A+B*x+C*y+x^2-y^2';
> f2 = 'D+F*x+G*y+x*y+y^2';
> [x,y]=solve(f1,f2);
>
> The solution to y is this:
>
> y = RootOf((-B+C+2*G-2*F)*_Z^3 ...
> +(A-G*B-B*F+2*C*F+G^2+2*D-F^2)*_Z^2 ...
> +(2*A*F-B*G*F-B*D+C*F^2+2*G*D)*_Z ...
> +A*F^2-B*D*F+D^2)
>
> which indicates the cubic coefficients involved in
finding y. From pretty(x) (or
> scratchings on the back of an envelope) can be deduced
this equality:
>
> x = -(D+G*y+y^2)/(F+y)
>
> Now set about computing these expressions numerically:
>
> % The given constants
> a =13.9497;
> b = 1.8350e+003;
> c = 140.4557;
> d =5.8077e-009;
> k1_2 =9.3378e-007;
> k2_2 =6.6098e-010;
>
> % Evaluate the six above coefficients in the two
equations
> A = b*d-k1_2*a;
> B = d+b+k1_2;
> C = b-d;
> D = c*d-k2_2*b;
> F = c-k2_2;
> G = c+d+k2_2;
>
> % The cubic coefficients from the "RootOf" expression
above
> c1 = -B+C+2*G-2*F;
> c2 = A-G*B-B*F+2*C*F+G^2+2*D-F^2;
> c3 = 2*A*F-B*G*F-B*D+C*F^2+2*G*D;
> c4 = A*F^2-B*D*F+D^2;
>
> % Find the three roots for y values
> y = roots([c1 c2 c3 c4]);
>
> % Get the corresponding x values
> x = -(D+G*y+y.^2)./(F+y);
>
> % Results:
> x = 1.0e+02 *
> 1.62669423030425
> 1.40465468520115
> -0.02615082158746
>
> y = 1.0e+02 *
> -1.62669423100515
> -1.40455558867504
> 0.02615082161392
>
> The x and y values are nearly the negative of one
another. Note that the
> second of the above x,y pairs satisfies the original
equations very poorly,
> which is apparently due to the fact that the divisor F+y
is nearly zero for
> them, and this seems to have caused excessively large
computation errors.
> The equation
>
> c1*y.^3+c2*y.^2+c3*y+c4
>
> is reasonably close to zero for all three y values.
>
> Roger Stafford
>
>

Subject: How to solve these equations?

From: Roger Stafford

Date: 19 Mar, 2008 02:54:02

Message: 10 of 13

"Liqing " <liqingj@gmail.com> wrote in message <frpp1t$o3f
$1@fred.mathworks.com>...
> Dear Roger,
>
> Thank you for your detailed code. I appreciate your
> precious time very very much!
>
> I wonder if I can ask a last question? When I use 'solve',
> I got a y that has hundreds of thousands of characters
> instead of the following 'rootof', which is so nice.
>
> Do you know what might be the matter?
----------
  If I use the original form of f1 and f2

 f1 = '(b*d-k1_2*a)+(d+b+k1_2)*x+(b-d)*y+x^2-y^2';
 f2 = '(c*d-k2_2*b)+(c-k2_2)*x+(c+d+k2_2)*y+x*y+y^2';
 [x,y]=solve(f1,f2);

I get 1785 characters for x and 351 for y from 'solve', whereas if I use

 f1 = 'A+B*x+C*y+x^2-y^2';
 f2 = 'D+F*x+G*y+x*y+y^2';
 [x,y]=solve(f1,f2);

I get only 350 for x and 113 for y. You could get shorter strings by using
only single characters in place of k1_2 and k2_2. I don't know how you
managed to get hundreds of thousands of characters, though.

  Maybe 'solve' has gotten dumber. No, I'm just kidding. My own version
refuses to handle any expression string, even intermediate ones, longer than
8192 which sometimes makes life a little difficult.

Roger Stafford

Subject: How to solve these equations?

From: Roger Stafford

Date: 19 Mar, 2008 08:57:36

Message: 11 of 13

"Liqing " <liqingj@gmail.com> wrote in message <frpp1t$o3f
$1@fred.mathworks.com>...
> .....
> I wonder if I can ask a last question? When I use 'solve',
> I got a y that has hundreds of thousands of characters
> instead of the following 'rootof', which is so nice.
>
> Do you know what might be the matter?
> .....
-------
  One other point occurs to me. You should always specify the variables you
are solving for in 'solve' to make sure no mistake is made in that respect.
Therefore in your case you should really write:

 [x,y] = solve(f1,f2,'x,y'); % Specify 'x' and 'y' as the variables to be solved

In the absence of such specification, 'solve' has its own method of deciding
which are the variables to be solved for, based on the particular letters of the
alphabet used, and that is something you don't want to depend on.

Roger Stafford

Subject: How to solve these equations?

From: Steven Lord

Date: 19 Mar, 2008 13:56:41

Message: 12 of 13


"Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid> wrote in
message news:frqki0$fcc$1@fred.mathworks.com...
> "Liqing " <liqingj@gmail.com> wrote in message <frpp1t$o3f
> $1@fred.mathworks.com>...
>> .....
>> I wonder if I can ask a last question? When I use 'solve',
>> I got a y that has hundreds of thousands of characters
>> instead of the following 'rootof', which is so nice.
>>
>> Do you know what might be the matter?
>> .....
> -------
> One other point occurs to me. You should always specify the variables
> you
> are solving for in 'solve' to make sure no mistake is made in that
> respect.
> Therefore in your case you should really write:
>
> [x,y] = solve(f1,f2,'x,y'); % Specify 'x' and 'y' as the variables to be
> solved
>
> In the absence of such specification, 'solve' has its own method of
> deciding
> which are the variables to be solved for, based on the particular letters
> of the
> alphabet used, and that is something you don't want to depend on.

In addition, for the OP I recommend using:

S = solve(f1, f2, 'x', 'y');

S will be a struct array containing two fields, x and y. This way there's
no potential confusion about which output received the solution for which of
the variables. Roger's solution prevents this potential confusion as well,
but with his method if you change the number or ordering of the variables
inside the SOLVE call, you have to remember to change the number or ordering
of the outputs as well.

--
Steve Lord
slord@mathworks.com

Subject: How to solve these equations?

From: Liqing

Date: 19 Mar, 2008 14:04:01

Message: 13 of 13

Thank you!

I still got the same answer. "Output truncated. Text
exceeds maximum line length of 25,000 characters for Command
Window display."

Below is my code:

f1 = 'A + B*x + C*y + x^2 - y^2';
f2 = 'D + F*x + G*y + x*y + y^2';

[x,y] = solve(f1,f2,'x,y');
%---------------------------------

And below is my results:


y =
                                  

Tags for this Thread

No tags are associated with this thread.

What are tags?

A tag is like a keyword or category label associated with each thread. Tags make it easier for you to find threads of interest.

Anyone can tag a thread. Tags are public and visible to everyone.

Contact us