1 // Copyright (c) 1997-1999 Matra Datavision
2 // Copyright (c) 1999-2014 OPEN CASCADE SAS
4 // This file is part of Open CASCADE Technology software library.
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License version 2.1 as published
8 // by the Free Software Foundation, with special exception defined in the file
9 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
10 // distribution for complete text of the license and disclaimer of any warranty.
12 // Alternatively, this file may be used under the terms of Open CASCADE
13 // commercial license or contractual agreement.
16 #define No_Standard_RangeError
17 #define No_Standard_OutOfRange
18 #define No_Standard_DimensionError
22 #include <math_FunctionAllRoots.hxx>
23 #include <math_FunctionRoots.hxx>
24 #include <math_FunctionSample.hxx>
25 #include <math_FunctionWithDerivative.hxx>
26 #include <Standard_NumericError.hxx>
27 #include <Standard_OutOfRange.hxx>
28 #include <StdFail_NotDone.hxx>
30 math_FunctionAllRoots::math_FunctionAllRoots (
31 math_FunctionWithDerivative& F,
32 const math_FunctionSample& S,
33 const Standard_Real EpsX, const Standard_Real EpsF,
34 const Standard_Real EpsNul) {
40 Standard_Boolean Nul, PNul, InterNul, Nuld, Nulf;
41 Standard_Real DebNul = 0., FinNul = 0.;
42 Standard_Integer Indd = 0, Indf = 0;
43 Standard_Real cst,val,valsav=0,valbid;
44 Standard_Boolean fini;
45 Standard_Integer Nbp,i;
48 F.Value(S.GetParameter(1),val);
49 PNul=Abs(val)<=EpsNul;
50 if (!PNul) {valsav=val;}
51 InterNul=Standard_False;
60 F.Value(S.GetParameter(i),val);
65 if (InterNul && (!Nul)) {
66 InterNul=Standard_False;
75 math_FunctionRoots Res1(F,S.GetParameter(i-1),S.GetParameter(i),10,
77 Standard_NumericError_Raise_if((!Res1.IsDone()) ||
79 (Res1.NbSolutions()==0), " ");
82 Indf=Res1.StateNumber(1);
85 math_FunctionRoots Res2(F,S.GetParameter(i-1),S.GetParameter(i),10,
87 Standard_NumericError_Raise_if((!Res2.IsDone()) ||
88 (Res2.IsAllNull())," ");
90 //-- || (Res2.NbSolutions()!=0), " "); lbr le 13 mai 87 (!=0 -> ==0)
91 if (Res2.NbSolutions()!=0) {
92 if (Res2.Value(1) < FinNul) {
93 FinNul = Res2.Value(1);
94 Indf = Res2.StateNumber(1);
100 else if ((!InterNul) && PNul && Nul) {
101 InterNul=Standard_True;
103 DebNul=S.GetParameter(1);
104 F.Value(DebNul,valbid);
105 Indd = F.GetStateNumber();
115 math_FunctionRoots Res1(F,S.GetParameter(i-2),S.GetParameter(i-1),10,
117 Standard_NumericError_Raise_if((!Res1.IsDone()) ||
118 (Res1.IsAllNull()) ||
119 (Res1.NbSolutions()==0), " ");
120 DebNul=Res1.Value(Res1.NbSolutions());
121 Indd = Res1.StateNumber(Res1.NbSolutions());
124 math_FunctionRoots Res3(F,S.GetParameter(i-2),S.GetParameter(i-1),10,
126 Standard_NumericError_Raise_if((!Res3.IsDone()) ||
127 (Res3.IsAllNull()), " ");
129 if (Res3.NbSolutions()!=0) {
130 if (Res3.Value(Res3.NbSolutions()) > DebNul) {
131 DebNul = Res3.Value(Res3.NbSolutions());
132 Indd = Res3.StateNumber(Res3.NbSolutions());
142 if (InterNul) { // rajouter l intervalle finissant au dernier pt
145 FinNul = S.GetParameter(Nbp);
146 F.Value(FinNul,valbid);
147 Indf = F.GetStateNumber();
153 if (pdeb.Length()==0) { // Pas d intervalle nul
155 math_FunctionRoots Res(F,S.GetParameter(1),S.GetParameter(Nbp),Nbp,
157 Standard_NumericError_Raise_if((!Res.IsDone()) ||
158 (Res.IsAllNull()), " ");
160 for (Standard_Integer j=1; j<=Res.NbSolutions(); j++) {
161 piso.Append(Res.Value(j));
162 iiso.Append(Res.StateNumber(j));
166 Standard_Integer NbpMin = 3;
167 Standard_Integer Nbrpt;
168 if (!Nuld) { // Recherche des solutions entre S.GetParameter(1)
169 // et le debut du 1er intervalle nul
171 Nbrpt=(Standard_Integer ) IntegerPart(
172 Abs((pdeb.Value(1)-S.GetParameter(1))/
173 (S.GetParameter(Nbp)-S.GetParameter(1)))*Nbp);
174 math_FunctionRoots Res(F,S.GetParameter(1),pdeb.Value(1),
175 Max(Nbrpt, NbpMin), EpsX,EpsF,0.0);
176 Standard_NumericError_Raise_if((!Res.IsDone()) ||
177 (Res.IsAllNull()), " ");
179 for (Standard_Integer j=1; j<=Res.NbSolutions(); j++) {
180 piso.Append(Res.Value(j));
181 iiso.Append(Res.StateNumber(j));
184 for (Standard_Integer k=2; k<=pdeb.Length(); k++) {
185 Nbrpt=(Standard_Integer )
186 IntegerPart(Abs((pdeb.Value(k)-pfin.Value(k-1))/
187 (S.GetParameter(Nbp)-S.GetParameter(1)))*Nbp);
188 math_FunctionRoots Res(F,pfin.Value(k-1),pdeb.Value(k),
189 Max(Nbrpt, NbpMin), EpsX,EpsF,0.0);
190 Standard_NumericError_Raise_if((!Res.IsDone()) ||
191 (Res.IsAllNull()), " ");
193 for (Standard_Integer j=1; j<=Res.NbSolutions(); j++) {
194 piso.Append(Res.Value(j));
195 iiso.Append(Res.StateNumber(j));
198 if (!Nulf) { // Recherche des solutions entre la fin du
199 // dernier intervalle nul et Value(Nbp).
201 Nbrpt=(Standard_Integer )
202 IntegerPart(Abs((S.GetParameter(Nbp)-
203 pfin.Value(pdeb.Length()))/
204 (S.GetParameter(Nbp)-S.GetParameter(1)))*Nbp);
205 math_FunctionRoots Res(F,pfin.Value(pdeb.Length()),
206 S.GetParameter(Nbp), Max(Nbrpt, NbpMin),
208 Standard_NumericError_Raise_if((!Res.IsDone()) ||
209 (Res.IsAllNull()), " ");
211 for (Standard_Integer j=1; j<=Res.NbSolutions(); j++) {
212 piso.Append(Res.Value(j));
213 iiso.Append(Res.StateNumber(j));
221 void math_FunctionAllRoots::Dump(Standard_OStream& o) const {
223 o<< "math_FunctionAllRoots ";
225 o<< " Status = Done \n";
226 o << " Number of null intervals = " << pdeb.Length() << std::endl;
227 o << " Number of points where the function is null: " << piso.Length() << std::endl;
230 o<< " Status = not Done \n";