solve non linear equations symbolically

조회 수: 19 (최근 30일)
Dhiraj R
Dhiraj R 2016년 12월 1일
편집: Walter Roberson 2017년 8월 25일
I have 4 variables & 3 non linear equations. Objective is to get 4th non linear equations expressed in 2 variables.
Variables are
a b M phi
Equations are
x/sqrt(x^2+y^2) = -1 + 1/M * cos(phi)
y/sqrt(x^2+y^2) = 1/M * sin(phi)
arg(diff((x- jy)/sqrt(x^2+y^2)),a) = pi/2 - phi
where,
x = cos(a) - a*b*sin(a) - 1
y = a*b*cos(a) + sin(a)
Objective: Express variable 'b' in terms of 'M'
  댓글 수: 5
Walter Roberson
Walter Roberson 2016년 12월 1일
Can we constrain a, b, phi, M to real-valued? If though you use imaginary components in the expression?
Dhiraj R
Dhiraj R 2016년 12월 2일
Yes a, b, phi, M are real valued.

댓글을 달려면 로그인하십시오.

채택된 답변

Walter Roberson
Walter Roberson 2016년 12월 2일
My investigations indicate that if there are symbolic solutions at all, then they are impractical to reach.
It is practical to solve the first two equations to get b and phi, leaving the third equation in terms of a and M . However, solving the resulting third equation (after substitution) for a is not practical (if it is possible at all.)
I have had a symbolic session running for the last hour just trying to evaluate that third equation at a = -Pi/5
My investigations show that there are an infinite number of matching real-valued a values for any M value for which there is a real-valued solution at all -- the arg() makes the solution either periodic or quasi-periodic (that is, same shape but possibly the distance changes.)
There do not happen to be any real-valued solutions for positive a . If you rewrite the equation from the form arg(X)=Y to be the expression arg(X)-Y and look for the zeros, then for a number of positive a values, there are both positive and negative values (depending on M), but there are no 0s because the positive values are discontinuous from the negative values.
There are solutions for some negative a values. There is an area near 0 in negative a where the above-mentioned expression has both positive and negative values but is discontinuous. By a = -Pi/4 there are real-valued solutions in M. My attempt to evaluate at a = -Pi/5 was probing to see where the discontinuity boundary is, more precisely. Guess I should have used floating point...
If you consider M values instead of a values, then there is a range of M values on both sides of 0 for which there are no real-valued solutions. I have not yet characterized how far that range extends.
  댓글 수: 4
Walter Roberson
Walter Roberson 2016년 12월 2일
Notice:
"Non-linear Eq. (40) can be solved easily by Newton iterative method when Ms is given between 1.2 and 2.0 with step value of 0.05 respectively."
Newton iterative method is a numeric method, not a symbolic method.
Dhiraj R
Dhiraj R 2016년 12월 4일
편집: Dhiraj R 2016년 12월 4일
I tried to solve the problem using Newton iterative method, but inverse of Jacobian matrix fails. I wrote the following code
clc
clear all
syms a b phi
syms x(a) y(a)
x(a) = cos(a) - a*b*sin(a) - 1;
y(a) = a*b*cos(a) + sin(a);
temp = diff((x(a)- j*y(a))/sqrt(x(a)^2+y(a)^2), a);
f1 = x(a)/sqrt(x(a)^2+y(a)^2) + 1 - 1/1.2 * cos(phi);
f2 = y(a)/sqrt(x(a)^2+y(a)^2) - 1/1.2 * sin(phi);
f3 = angle(temp) - pi/2 + phi;
F = [f1; f2; f3]
J = jacobian(F)
% assume initial value a = 1, b = 3, phi = pi/4
Ja = subs(J,a,1);
Jb = subs(Ja,b,3);
Jc = subs(Jb,phi,pi/4);
invJ = inv(Jc)
Inverse of Jacobian Matrix failed here.

댓글을 달려면 로그인하십시오.

추가 답변 (1개)

Juan Ponce
Juan Ponce 2017년 8월 25일
편집: Walter Roberson 2017년 8월 25일
Hi Dhiraj R,
I obtain the following solutions, only with b and phi as variables.
% Variables Modelo
syms a b phi real
% a:Variable asociada a w*l ()
% b:Variable asociada a n2/l ()
% phi:Posicion angular (rad)
% n2:Constante de tiempo del filtro (s)
% Parametros Simbolicos
syms M real
% M:Sensitividad maxima Ms ()
% Ecuaciones auxiliares
x=cos(a)-a*b*sin(a)-1;
y=a*b*cos(a)+sin(a);
L=(x-j*y)/(x^2+y^2)^(1/2);
% Calculo Derivada Parcial
dL=diff(L,a);
% Calculo Argumento o Angulo
arg=atan(dL);
% Construccion Ecuaciones
% Nodos o Mallas, Fuerzas, Torques, Conservacion Energias
eq1=x/(x^2+y^2)^(1/2)-(-1+1/M*cos(phi)); % x/(x^2+y^2)^(1/2)=-1+1/M*cos(phi)
eq2=y/(x^2+y^2)^(1/2)-(1/M*sin(phi)); %y/(x^2+y^2)^(1/2)=1/M*sin(phi)
eq3=arg-(pi/2-phi); %arg=pi/2-phi
%%%%%%%%%%%%%%%%
F =
1 - (a*b*sin(a) - cos(a) + 1)/((sin(a) + a*b*cos(a))^2 + (a*b*sin(a) - cos(a) + 1)^2)^(1/2) - cos(phi)/M
(sin(a) + a*b*cos(a))/((sin(a) + a*b*cos(a))^2 + (a*b*sin(a) - cos(a) + 1)^2)^(1/2) - sin(phi)/M
phi - pi/2 - atan((sin(a) + i*cos(a) + b*sin(a) + a*b*cos(a) + b*i*cos(a) - a*b*i*sin(a))/((sin(a) + a*b*cos(a))^2 + (a*b*sin(a) - cos(a) + 1)^2)^(1/2) - ((2*(sin(a) + a*b*cos(a))*(cos(a) + b*cos(a) - a*b*sin(a)) + 2*(a*b*sin(a) - cos(a) + 1)*(sin(a) + b*sin(a) + a*b*cos(a)))*(i*sin(a) - cos(a) + a*b*sin(a) + a*b*i*cos(a) + 1))/(2*((sin(a) + a*b*cos(a))^2 + (a*b*sin(a) - cos(a) + 1)^2)^(3/2)))
X =
b
phi
mxd_nuevo =
[ -(2*sin(a - 4*pi*k) + 2*sin(a + 4*pi*k) + 4*sin(a) + 8*cos(2*pi*k)*(2*M*cos(2*pi*k) - cos(2*k*pi)^2)^(1/2) - 4*cos(2*pi*k - a)*(2*M*cos(2*pi*k) - cos(2*k*pi)^2)^(1/2) + 2*cos(2*pi*k - 2*a)*(2*M*cos(2*pi*k) - cos(2*k*pi)^2)^(1/2) - 2*cos(2*a - 2*pi*k)*(2*M*cos(2*pi*k) - cos(2*k*pi)^2)^(1/2) + cos(a)*(8*M*(2*M*cos(2*pi*k) - cos(2*k*pi)^2)^(1/2) + 8*M^2*sin(a)) - 4*cos(a + 2*pi*k)*(2*M*cos(2*pi*k) - cos(2*k*pi)^2)^(1/2) - 8*M*(2*M*cos(2*pi*k) - cos(2*k*pi)^2)^(1/2) - 8*M*sin(a - 2*pi*k) - 8*M*sin(a + 2*pi*k))/(a*cos(4*pi*k - 2*a) - a*cos(2*a - 4*pi*k) + 8*a*cos(2*pi*k)^2 - 16*M*a*cos(2*pi*k) + 8*M^2*a*cos(a)^2 - 4*M*a*cos(2*pi*k - 2*a) + 4*M*a*cos(2*a - 2*pi*k)), pi - acos(M*(sin(a) - 1))]
[ (sin(a - 4*pi*k) + sin(a + 4*pi*k) - sin(2*pi*k)*(4*(2*sin(a/2)^2 - 1)*(M^2 - sin(2*k*pi)^2)^(1/2) + 4*(M^2 - sin(2*k*pi)^2)^(1/2)) - sin(a)*(4*M^2*(2*sin(a/2)^2 - 1) + 2))/(4*a*(sin(a)^2 - 1)*M^2 + 4*a + 4*a*(sin(2*pi*k)^2 - 1)), pi - acos(M*(sin(a) - 1))]
Note: Each row of mxd_nuevo representa a solution bi,phi (i=1,2)
I dont know if this is right.
Good luck!!
Juan Ponce

카테고리

Help CenterFile Exchange에서 Numeric Solvers에 대해 자세히 알아보기

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by