Solving Equation Depending on Parameter

2 vues (au cours des 30 derniers jours)
onsagerian
onsagerian le 24 Juil 2018
Commenté : Star Strider le 24 Juil 2018
Hello,
An equation in a reference paper has been simply translated into a Matlab code (see below) for solving it, and it must have particular values of x depending on the parameter "gamma", but I got an error. Would you check the equation and answer me what I need to make a correction?
format long e
T=0.1;
k=1.0; %k(1) default value is 1.0
k_1=0.2; %k(-1) default value is 2.0 (modified: 1.0)
k_2=0.00001; %k(-1*) default value is 0.1 (Note: k_2=0.00001 higher error rate for n=1 compared to n>1)
%k_2=[0.1:0.1:1.0];
a_p=0.0001; %k(p) default value is 0.0001
k_p=a_p/10; %k(-p) Assuming this ratio is the constant (a_p/10)
a_q=0.00015; %k(q) default value is 0.00015
k_q=0.00002; %k(-q) default value is 0.00002
theta=0.01; %default value is 0.01
m1=0.0000005; %default value is 0.0000005
m2=0.0000002; %default value is 0.0000002
W=0.00002; %default value is 0.002
n=6;
syms x
gamma=(a_p*k*k_2)/(k_p*m1*k_1);
equation=gamma-[1+(theta-x)/(1-theta)*[(1-theta)/(x^(1/(n+1))-theta)]^(n+1)]==0;
sol=solve(equation,x);
  1 commentaire
Aquatris
Aquatris le 24 Juil 2018
Modifié(e) : Aquatris le 24 Juil 2018
You should put your code in a nice format. I cannot say which parts are comment which parts are not. For instance, you cannot do
syms x gamma=(a_p*k*k_2)/(k_p*m1*k_1);
in a single line. It should be
syms x
gamma = (a_p*k*k_2)/(k_p*m1*k_1);
Is it due to your copy paste or do you actually use it this way cannot be determined from your question.

Connectez-vous pour commenter.

Réponse acceptée

Star Strider
Star Strider le 24 Juil 2018
There appears to be one valid value for ‘x’:
x =
0.2014042781988740154717685102689
The reason is that although there are 7 roots to the equation for ‘x’, it is the only one that meets the criteria set out in the ‘conditions’ result, that being:
-0.4487989505128276054946633404685 < angle(z2) & angle(z2) <= 0.4487989505128276054946633404685
This results from instructing solve to return the conditions as well:
sol = solve(equation,x, 'ReturnConditions',1);
that returns a polynomial in ‘z2’ with roots:
z2 =
0.2014042781988740154717685102689
- 0.16246739411681022588655095122335 + 0.083052533086836812066634941174803i
- 0.16246739411681022588655095122335 - 0.083052533086836812066634941174804i
- 0.032602888099564513824045864450354 + 0.1866174981726704349497552909463i
- 0.032602888099564513824045864450354 - 0.1866174981726704349497552909463i
0.12933518938326112554046904728674 - 0.14965579028283080032788468715005i
0.12933518938326112554046904728674 + 0.14965579028283080032788468715005i
the angles of these roots being:
angle_z2 =
0
2.6690291408928938742331335149207
-2.6690291408928938742331335149207
1.7437551108769301908411899136693
-1.7437551108769301908411899136693
-0.85810582106481765317999373280298
0.85810582106481765317999373280298
with the only value meeting the ‘angle’ requirement being the real root.
Also, please highlight your code, then press the [{}Code] button to format it properly. This makes it easier for everyone to read it and work with it. (I did that for you this time.)
  2 commentaires
onsagerian
onsagerian le 24 Juil 2018
I appreciate your help, and also thank you for letting me know how to make the format in a proper manner using the button "{} code". But,I still have a question. I used the format sol = solve(equation,x, 'ReturnConditions',1); to get proper solutions for the equation, but there was no return value. No error message, but not response. It is the only thing I made a correction. Do you know why this happened?
Star Strider
Star Strider le 24 Juil 2018
My pleasure.
I have no idea. (I am using R2018a.)
The ‘return value’ when I ran it was only ‘z2^7’. I had to use the ‘conditions’ field to get anything, and I posted everything I got.

Connectez-vous pour commenter.

Plus de réponses (0)

Catégories

En savoir plus sur Particle & Nuclear Physics dans Help Center et File Exchange

Community Treasure Hunt

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

Start Hunting!

Translated by