Finding Solution to Inequality in Matlab

3 visualizzazioni (ultimi 30 giorni)
Laura Freitas
Laura Freitas il 13 Ott 2023
Commentato: Walter Roberson il 14 Ott 2023
I have a system of five equations. As below
syms x y z v w q beta rho theta lambda b tau
eqn1 = x == q*(2*rho/(1-lambda) + beta*y) + (1-q)*((1-theta)*(1-tau)/(1-lambda)+beta*z);
eqn2 = y == q*beta*v + (1-q)*(b*(1-theta)/(1-lambda)+beta*z);
eqn3 = v == q*(2*rho/(1-lambda)+ beta*y) + (1-q)*(b*(1-theta)/(1-lambda)+beta*z);
eqn4 = z == (x + w)/2;
eqn5 = w == q*beta*y + (1-q)*((1-tau)*(1-theta)\(1-lambda)+ beta*z);
I want to find the solution for x, in terms of q, beta, rho, theta, lambda, b and tau.
sol = solve([eqn1, eqn2, eqn3, eqn4, eqn5], [x, y, z, w, v]);
sol.x
After that, I want to find the solution to the following inequality
eqn6 = (1-tau)*(1-theta)/((1-beta)*(1-lambda))> sol.x;
This is where the problem begins. When I run the following command:
solution = solve (eqn6, theta);
I get the following error message:
Warning: Unable to find explicit solution. For options, see help. > In sym/solve (line 317)
In dynamic_model_8 (line 13)
I get that these equations are complex, but is there a way to get matlaf to find the result?
Thanks in advance.

Risposte (3)

Torsten
Torsten il 13 Ott 2023
Modificato: Torsten il 13 Ott 2023
I don't understand the result, but there is some output from the symbolic toolbox.
As you can see, the numerator of sol.x is a linear function of theta and the denominator does not contain theta.
Thus the whole thing boils down to determine the zero of a linear function - I'd do it using pencil and paper considering several cases (tau < 1, tau > 1, beta < 1, beta > 1, lambda < 1, lambda > 1) and combinations of these.
I changed
eqn5 = w == q*beta*y + (1-q)*((1-tau)*(1-theta)\(1-lambda)+ beta*z);
to
eqn5 = w == q*beta*y + (1-q)*((1-tau)*(1-theta)/(1-lambda)+ beta*z);
I think it was a typo.
syms x y z v w q beta rho theta lambda b tau
eqn1 = x == q*(2*rho/(1-lambda) + beta*y) + (1-q)*((1-theta)*(1-tau)/(1-lambda)+beta*z);
eqn2 = y == q*beta*v + (1-q)*(b*(1-theta)/(1-lambda)+beta*z);
eqn3 = v == q*(2*rho/(1-lambda)+ beta*y) + (1-q)*(b*(1-theta)/(1-lambda)+beta*z);
eqn4 = z == (x + w)/2;
eqn5 = w == q*beta*y + (1-q)*((1-tau)*(1-theta)/(1-lambda)+ beta*z);
sol = solve([eqn1, eqn2, eqn3, eqn4, eqn5], [x, y, z, w, v]);
sol.x
ans = 
[N,D] = numden(sol.x);
collect(N,theta)
ans = 
eqn6 = (1-tau)*(1-theta)/((1-beta)*(1-lambda)) - sol.x > 0;
solution = solve(eqn6, theta)
solution = 

John D'Errico
John D'Errico il 13 Ott 2023
Simple enough. I'll just work with eqn6 here.
syms tau beta lambda theta Solx
eqn6 = (1-tau)*(1-theta)/((1-beta)*(1-lambda))> Solx;
solve(eqn6,theta,'returnconditions',true)
ans = struct with fields:
theta: ((x + (tau - 1)/((beta - 1)*(lambda - 1)))*(beta - 1)*(lambda - 1))/(tau - 1) parameters: x conditions: Solx < x & tau ~= 1
  2 Commenti
Dyuman Joshi
Dyuman Joshi il 13 Ott 2023
Modificato: Dyuman Joshi il 13 Ott 2023
Also, as every other division operator used is a backslash or rdivide, ./, there seems to be a typo here -
% v
qn5 = w == q*beta*y + (1-q)*((1-tau)*(1-theta)\(1-lambda)+ beta*z);
Torsten
Torsten il 13 Ott 2023
sol.x is a function of theta.

Accedi per commentare.


Walter Roberson
Walter Roberson il 13 Ott 2023
Turn inequalities into equalities by adding a positive (or non-negative) variable representing how much more one side is compared to the other.
syms x y z v w q beta rho theta lambda b tau
eqn1 = x == q*(2*rho/(1-lambda) + beta*y) + (1-q)*((1-theta)*(1-tau)/(1-lambda)+beta*z);
eqn2 = y == q*beta*v + (1-q)*(b*(1-theta)/(1-lambda)+beta*z);
eqn3 = v == q*(2*rho/(1-lambda)+ beta*y) + (1-q)*(b*(1-theta)/(1-lambda)+beta*z);
eqn4 = z == (x + w)/2;
eqn5 = w == q*beta*y + (1-q)*((1-tau)*(1-theta)\(1-lambda)+ beta*z);
sol = solve([eqn1, eqn2, eqn3, eqn4, eqn5], [x, y, z, w, v]);
sol.x
ans = 
syms Excess positive
eqn6 = (1-tau)*(1-theta)/((1-beta)*(1-lambda)) - Excess == sol.x
eqn6 = 
solution = solve (eqn6, theta)
solution = 
sympref('AbbreviateOutput', 0)
ans = logical
1
collect(simplify(solution, 'steps', 50), Excess)
ans = 
  6 Commenti
Sam Chak
Sam Chak il 14 Ott 2023
Thanks @Walter Roberson. I guess that your tests further reinforce that the 'collect()' function only works on strictly polynomials, , but not on rational functions, even though they are quotients of two polynomials, .
Walter Roberson
Walter Roberson il 14 Ott 2023
If you look at coeffs4 the denominator is collected with respect to x: it just isn't as factored as for the coeffs3 case. It did not fail to collect: it just didn't fully factor after collecting, even though fully factoring is not a documented outcome in the first place.

Accedi per commentare.

Prodotti


Release

R2023b

Community Treasure Hunt

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

Start Hunting!

Translated by