Solve multiple non-linear equations with vector variables

Hi all,

How can I solve multiple equations with vector variables? Say I have two vectors X and Y: X=[x1,x2], Y =[y1,y2], and two equations: X.^2+Y.^2=A, X.^2-Y.^2=B, where A=[20,5], and B =[12,3]. How can I solve this problem using "fsolve"?

In the real case, my equations are more complicated and I have 50,000 rows for vectors X and Y. Instead of looping each row and solve X(n)^2+Y(n)^2=A(n),X(n)^2-Y(n)^2=B(n), I wonder if there is a more effecient way. Thanks!

[James Tursa]

You have all the A and B values and you are trying to find the X and Y values? If so, why not just use a little algebra to solve for them directly?

[Yang Li]

Hi James, sorry for not making it clear. Here I would like to simplify my question, so I gave 2 by 1 vectors as an example, and the equations are also simplified. The real-world problem is that I got 50,000 sample points, so the X,Y, A,and B are are vectors with 50,000 rows. Instead of computing the multivariate nonlinear equations 50,000 times, I would like to know how to do it more efficiently in Matlab. Any idea will be appreciated.

[Torsten]

If MATLAB's "solve" does not return a symbolic solution for the master-system (in this case x^2+y^2 = a, x^2-y^2=b), you will have to solve your 50000 systems one after the other.

[Yang Li]

Hi Torsten, Thanks for your comment. Yes I didn't find other more efficient way, so I tried parellel computing in the end and it saved me a bunch of time.

[James Tursa]

I still don't get it. If you have these equations:

x.^2 + y.^2 = A

x.^2 - y.^2 = B

Just add them to get

2*x.^2 = A+B

and solve for x.^2 and y.^2

x.^2 = (A+B)/2

y.^2 = x.^2 - B = (A-B)/2

From that you can sqrt( ) to get x and y, and maybe pick signs as appropriate to your problem.

What am I missing here? If this isn't your actual problem with actual equations, then please post your actual problem with the actual equations.

[Yang Li]

My real equations are two physical equations:

eqn1 = C*X+C*D*(((1-X-Y)^2)/(1-D*X))==A;

eqn2 = C*(1-X-Y)/(1-D*X)==B.

Here A,B,C,D are know vectors with 50,000 rows, representing the states of four variables of 50,000 samples. For each A(i), B(i), C(i), and D(i), the corresponding X(i) and Y(i) can be solved as you demonstrated . But now I have 50,000 independent samples and looping from the first sample to the last is really time-consuming. That's why I post this question and would like to know if there is anyway more efficient.

[Torsten]

But it's not difficult to solve this 2-equation system for X and Y. I thought your equations were much harder.

And once you have solved for X and Y, you don't need to loop, but you can instantly insert the complete 50000 element vectors A,B,C and D to get back the 50000 element vectors X and Y.

[Yang Li]

Hi Torsten, Matt gives detailed explanation of how to reduce the complicated equations to a linear system. I accepted his answer. Thank you very much for your help.

Answer 1

[Matt J]

But it's not difficult to solve this 2-equation system for X and Y. I thought your equations were much harder.

And in fact, the equations can be reduced to a 2x2 linear system. When D=0, the equations reduce to linear equations in X and Y, with a simple solution.

X=A./C;
Y=1-((A+B)./C);

When D is not 0, you can make the change of variables P=(1-X-Y) and Q=(1-D*X). The equations then become,

eqn1 = C/D*(1-Q) +  C*D*(P^2/Q) ==A; 
eqn2 = C*P/Q==B

The second equation can be used to simplify the second term in the first equation,

 eqn1 = C/D*(1-Q) +  D*B*P == A; 

which is a linear eqaution and the second equation can also be rearranged as linear,

 eqn2 = C*P-B*Q==0

Simplifying everything leads to the linear matrix equations

 [ D*B  -C/D;
    C   -B   ]*[P;Q] =  [ A-C/D; 0]

whos analytical (and vectorized) solution is,

 d=C.^2./D-D*B.^2; %determinant
 
 P = -B.*(A-C./D)./d;
 Q = -C.*(A-C./D)./d;
 
 X=(1-Q)./D;
 Y=1-X-P;

[Yang Li]

This helps a lot. Thanks Matt!

Answer 2

[Matt J]

 A=[20,5]; B =[12,3]; 

 XY0=ones(2); %initial guess
 
[XY,fval]=fsolve(@(XY) Equations(XY, A,B)   , XY0);

Equation solved.

fsolve completed because the vector of function values is near zero
as measured by the value of the function tolerance, and
the problem appears regular as measured by the gradient.

X=XY(:,1), Y=XY(:,2),fval

X = 2×1

4.0000 2.0000

Y = 2×1

2.0000 1.0000

fval = 4×1

1.0e+-12 * 0.3730 0.3713 -0.3730 0.3713

 
function F=Equations(XY, A,B)

  X=XY(:,1); Y=XY(:,2);
  
  F=[X.^2+Y.^2-A(:); X.^2-Y.^2-B(:)];

end