Alternative to add assumptions to vpasolve

4 visualizaciones (últimos 30 días)
vaggelis papasot
vaggelis papasot el 19 de Feb. de 2022
Comentada: vaggelis papasot el 22 de Feb. de 2022
I need to numerically solve the following system of equations, in order to calculate the parameters w1,w2,w3,w4,a1,a2,a3,a4,mu1,mu2,mu3, and mu4. Also, I need the constraint that w1+w2+w3+w4=1. I know that vpasolve does not use assumptions and that one can only define ranges for the parameters that need to be calculated. I have the following question. Is there a way to define the constaint w1+w2+w3+w4=1 to get a numerical solution with vpasolve?
C1=0.9945;
C2=1;
C3=1.0162;
C4=1.0431;
C5=1.0809;
C6=1.13;
C7=1.1911;
C8=1.2650;
C9=1.3530;
C10=1.4564;
C11=1.5769;
C12=1.7165;
syms w1 a1 mu1 w2 a2 mu2 w3 a3 mu3 w4 a4 mu4
eq1=w1*( gamma( mu1+(1/a1) )/( ( mu1^(1/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(1/a2))/( ( mu2^(1/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(1/a3) )/( ( mu3^(1/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(1/a4) )/( ( mu4^(1/a4) )*gamma(mu4) ) )==C1;
eq2=w1*( gamma( mu1+(2/a1) )/( ( mu1^(2/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(2/a2))/( ( mu2^(2/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(2/a3) )/( ( mu3^(2/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(2/a4) )/( ( mu4^(2/a4) )*gamma(mu4) ) )==C2;
eq3=w1*( gamma( mu1+(3/a1) )/( ( mu1^(3/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(3/a2))/( ( mu2^(3/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(3/a3) )/( ( mu3^(3/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(3/a4) )/( ( mu4^(3/a4) )*gamma(mu4) ) )==C3;
eq4=w1*( gamma( mu1+(4/a1) )/( ( mu1^(4/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(4/a2))/( ( mu2^(4/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(4/a3) )/( ( mu3^(4/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(4/a4) )/( ( mu4^(4/a4) )*gamma(mu4) ) )==C4;
eq5=w1*( gamma( mu1+(5/a1) )/( ( mu1^(5/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(5/a2))/( ( mu2^(5/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(5/a3) )/( ( mu3^(5/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(5/a4) )/( ( mu4^(5/a4) )*gamma(mu4) ) )==C5;
eq6=w1*( gamma( mu1+(6/a1) )/( ( mu1^(6/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(6/a2))/( ( mu2^(6/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(6/a3) )/( ( mu3^(6/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(6/a4) )/( ( mu4^(6/a4) )*gamma(mu4) ) )==C6;
eq7=w1*( gamma( mu1+(7/a1) )/( ( mu1^(7/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(7/a2))/( ( mu2^(7/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(7/a3) )/( ( mu3^(7/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(7/a4) )/( ( mu4^(7/a4) )*gamma(mu4) ) )==C7;
eq8=w1*( gamma( mu1+(8/a1) )/( ( mu1^(8/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(8/a2))/( ( mu2^(8/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(8/a3) )/( ( mu3^(8/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(8/a4) )/( ( mu4^(8/a4) )*gamma(mu4) ) )==C8;
eq9=w1*( gamma( mu1+(9/a1) )/( ( mu1^(9/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(9/a2))/( ( mu2^(9/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(9/a3) )/( ( mu3^(9/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(9/a4) )/( ( mu4^(9/a4) )*gamma(mu4) ) )==C9;
eq10=w1*( gamma( mu1+(10/a1) )/( ( mu1^(10/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(10/a2))/( ( mu2^(10/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(10/a3) )/( ( mu3^(10/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(10/a4) )/( ( mu4^(10/a4) )*gamma(mu4) ) )==C10;
eq11=w1*( gamma( mu1+(11/a1) )/( ( mu1^(11/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(11/a2))/( ( mu2^(11/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(11/a3) )/( ( mu3^(11/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(11/a4) )/( ( mu4^(11/a4) )*gamma(mu4) ) )==C11;
eq12=w1*( gamma( mu1+(12/a1) )/( ( mu1^(12/a1) )*gamma(mu1) ) )...
+w2*( gamma(mu2+(12/a2))/( ( mu2^(12/a2) )*gamma(mu2) ) )...
+w3*( gamma( mu3+(12/a3) )/( ( mu3^(12/a3) )*gamma(mu3) ) )...
+w4*( gamma( mu1+(12/a4) )/( ( mu4^(12/a4) )*gamma(mu4) ) )==C12;
eqn=[eq1 eq2 eq3 eq4 eq5 eq6 eq7 eq8 eq9 eq10 eq11 eq12];
vars=[w1,a1,mu1,w2,a2,mu2,w3,a3,mu3,w4,a4,mu4];
range=[0 1; 0 Inf; 0 Inf; 0 1; 0 Inf; 0 Inf; 0 1; 0 Inf;0 Inf; 0 1; 0 Inf; 0 Inf];
S=vpasolve(eqn,vars,range,'Random',true);

Respuesta aceptada

Torsten
Torsten el 19 de Feb. de 2022
I don't know if your system of equations has a solution, but you can enforce sum(w(i))=1 by replacing
wi by wi/(w1+w2+w3+w4) (i=1,...,4)
  1 comentario
vaggelis papasot
vaggelis papasot el 22 de Feb. de 2022
Thank you @Torsten. This is a good workaround for this problem. However, the system now yields no solution.

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Gamma Functions en Help Center y File Exchange.

Productos


Versión

R2020a

Community Treasure Hunt

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

Start Hunting!

Translated by