La page de J-P. Merlet

J-P. Merlet home page

La page Présentation de HEPHAISTOS

HEPHAISTOS home page

La page "Présentation" de l'INRIA

INRIA home page

Subsections

• 6-body
• Ferraris
• Wilkinson polynomial

Examples

We present here some examples of Maple file for solving systems of equations.

6-body

The problem is here to find the solutions of a set of 6 equations (eq1 to eq6). Two simplification procedures are used together with the 3B method. The HessianSolve procedure will allow to get exact results.

 
eq1:=3*(b-d)*(B-D)+B+D-6*F+4:
eq2:=3*(b-d)*(B+D-2*F)+5*(B-D):
eq3:=3*(b-d)^2-6*(b+d)+4*f+3:
eq4:=B^2*b^3-1:
eq5:=D^2*d^3-1:
eq6:=F^2*f^3-1:
with(ALIAS):
`ALIAS/lib`:="../Test/Lib":
`ALIAS/profil`:="../Profil":
`ALIAS/single_bisection`:=2:
`ALIAS/3B`:=0:
`ALIAS/Full3B`:=1:
`ALIAS/Full3B_Change`:=0.1:
`ALIAS/Delta3B`:=0.01:
`ALIAS/Max3B`:=10^16:
`ALIAS/Use_Simp_3B`:=0:
`ALIAS/debug`:=1:
`ALIAS/storage_mode`:=10:
`ALIAS/full_simplex`:=9:
`ALIAS/min_diam_simplex`:=1000:
`ALIAS/diam_simplex`:=0.1:
u:=[0,10^5]:
SimplexConsistency([eq1,eq2,eq3],[B,D,F,b,d,f],"Simp1");
HullConsistency([eq1,eq2,eq3,eq4,eq5,eq6],[B,D,F,b,d,f],"Simp2");
CatSimp("Simp","Simp1","Simp2",0.01);


HessianSolve([eq1,eq2,eq3,eq4,eq5,eq6],[B,D,F,b,d,f],
	[u,u,u,u,u,u],"Simp");

Ferraris

Here we deal with a non polynomial system of 2 equations with a very large initial search domain. An additional problem is that the second equation involves the interval evaluation of the exponential of a variable that may leads to an overflow. Hence we use the mechanism described in section 2.1.4 to avoid the evaluation of the exponential if the upper bound of the first variable is larger than 1e4.

 
eq1:= -0.25/Pi*x2 - 0.5*x1  + 0.5*sin(x1*x2);
eq2:= exp(1)/Pi*x2 - 2*exp(1)*x1  + (1 - 0.25/Pi)*(exp(2*x1) - exp(1));

EQ:=[evalf(eq1),evalf(eq2)]:
VAR:=[x1,x2]:
INIT:=[[-10^8,10^8],[-10^8,10^8]]:

with(ALIAS):
`ALIAS/debug`:=1:
`ALIAS/single_bisection`:=2:
`ALIAS/storage_mode`:=40:
`ALIAS/lib`:="../Test/Lib":
`ALIAS/profil`:="../Profil":
`ALIAS/3B`:=1:
`ALIAS/Max3B`:=5*10^9:
`ALIAS/Delta3B`:=0.0001:
`ALIAS/Full3B`:=0:
`ALIAS/Full3B_Change`:=0.0001:
`ALIAS/maxgradient`:=20:
`ALIAS/maxkraw`:=1:
`ALIAS/maxnewton`:=1:

`ALIAS/user_FINIT`:="INTERVAL U;":
ALIAS_F:=proc(fid,i)
if i=2 then
fprintf(fid,"if(Sup(v_IS(1))>1.e4)V(2)=INTERVAL(-1.e6,1.e12);else\n"):
fi:
RETURN(0):
end:

HullConsistency(EQ,VAR,"Simp",0.1);
HessianSolve(EQ,VAR,INIT,"Simp"):

To manage the possible overflow problem instead of designing the ALIAS_F procedure by hand we may have also used the procedures described in the section 2.1.5 that create automatically the ALIAS_F and ALIAS_J procedures.

 
with(ALIAS):
with(linalg):
`ALIAS/gpp_linux`:="/usr/local/gcc-2.95/bin/g++":
`ALIAS/debug`:=1:
`ALIAS/single_bisection`:=1:
`ALIAS/storage_mode`:=40:
`ALIAS/lib`:="../Test/Lib":
`ALIAS/profil`:="../Profil":
`ALIAS/3B`:=1:
`ALIAS/Max3B`:=5*10^9:
`ALIAS/Delta3B`:=0.0001:
`ALIAS/Full3B`:=0:
`ALIAS/Full3B_Change`:=0.0001:
`ALIAS/maxgradient`:=20:
`ALIAS/maxkraw`:=1:
`ALIAS/maxnewton`:=1:
`ALIAS/optimized`:=0:

J:=jacobian(EQ,VAR):

Verify_Problem_Expression(EQ,VAR):
Verify_Problem_ExpressionJ(eval(J),VAR):

HullConsistency(EQ,VAR,"Simp1",0.1);
HessianSolve(EQ,VAR,INIT,"Simp1"):

Wilkinson polynomial

The Wilkinson polynomial of order $n$ is obtained as the product of the terms $(x-i)$ for i from 1 to $n$. Its expanded form is difficult to solve as soon as $n$ is larger than 12. If eq is the polynomial a procedure for solving it is:

 
EQ:=[eq]:
VAR:=[x]:
INIT:=[[-100,100]]:

`ALIAS/fepsilon`:=1e-9:
`ALIAS/epsilon`:=1e-14:
`ALIAS/storage_mode`:=10:
`ALIAS/3B`:=1:
`ALIAS/Max3B`:=10000:
`ALIAS/Delta3B`:=1e-9:
`ALIAS/optimized`:=0:
`ALIAS/rand`:=500:
`ALIAS/eps_inflation`:=1e-10:
`ALIAS/dist`:=1e-16:

Rouche(EQ,"p",VAR,"rouche"):
WeylFilter(EQ,VAR,VAR,["automatic"],["symbolic"],"SimpWeyl");
DeflationUP(eq,x,"F","J","simp_quo_new",debug);

CatSimp("SIMP","rouche","NO_3B","simp_quo_new","SimpWeyl",0.1);

SOL:=GradientSolve(EQ,VAR,INIT,"SIMP");

To filter the boxes we use WeylFilter, DeflationUP and to find the solutions we use Rouche.

J-P. Merlet home page

La page Présentation de HEPHAISTOS

HEPHAISTOS home page

La page "Présentation" de l'INRIA

INRIA home page