0027275: Unused formal parameter in BSplCLib::EvalBsplineBasis
[occt.git] / src / FairCurve / FairCurve_DistributionOfTension.cxx
CommitLineData
973c2be1 1// Copyright (c) 1999-2014 OPEN CASCADE SAS
b311480e 2//
973c2be1 3// This file is part of Open CASCADE Technology software library.
b311480e 4//
d5f74e42 5// This library is free software; you can redistribute it and/or modify it under
6// the terms of the GNU Lesser General Public License version 2.1 as published
973c2be1 7// by the Free Software Foundation, with special exception defined in the file
8// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
9// distribution for complete text of the license and disclaimer of any warranty.
b311480e 10//
973c2be1 11// Alternatively, this file may be used under the terms of Open CASCADE
12// commercial license or contractual agreement.
b311480e 13
7fd59977 14// 30-01-1996 : PMN Version originale
15
0797d9d3 16#ifndef OCCT_DEBUG
7fd59977 17#define No_Standard_RangeError
18#define No_Standard_OutOfRange
19#endif
20
7fd59977 21
42cf5bc1 22#include <BSplCLib.hxx>
23#include <FairCurve_BattenLaw.hxx>
24#include <FairCurve_DistributionOfTension.hxx>
7fd59977 25#include <gp_Pnt2d.hxx>
42cf5bc1 26#include <gp_XY.hxx>
7fd59977 27#include <math_Matrix.hxx>
28#include <math_Vector.hxx>
7fd59977 29
30 FairCurve_DistributionOfTension::FairCurve_DistributionOfTension(const Standard_Integer BSplOrder,
31 const Handle(TColStd_HArray1OfReal)& FlatKnots,
32 const Handle(TColgp_HArray1OfPnt2d)& Poles,
33 const Standard_Integer DerivativeOrder,
34 const Standard_Real LengthSliding,
35 const FairCurve_BattenLaw& Law,
36 const Standard_Integer NbValAux,
37 const Standard_Boolean Uniform ) :
38 FairCurve_DistributionOfEnergy(BSplOrder,
39 FlatKnots,
40 Poles,
41 DerivativeOrder,
42 NbValAux) ,
43 MyLengthSliding (LengthSliding),
44 MyLaw (Law)
45{
46 if (Uniform) {MyLaw.Value(0.5, MyHeight);} // it used in MVC to avoid Parametrization Problemes
47 else {MyHeight = 0;}
48}
49
50
51
52Standard_Boolean FairCurve_DistributionOfTension::Value(const math_Vector& TParam, math_Vector& FTension)
53{
54 Standard_Boolean Ok = Standard_True;
55 Standard_Integer ier, ii, jj, kk;
56 gp_XY CPrim (0., 0.);
57 Standard_Integer LastGradientIndex, FirstNonZero, LastZero;
58
59 // (0.0) initialisations generales
60 FTension.Init(0.0);
61 math_Matrix Base(1, 3, 1, MyBSplOrder ); // On shouhaite utiliser la derive premieres
62 // Dans EvalBsplineBasis C' <=> DerivOrder = 2
63 // et il faut ajouter 1 rang dans la matrice Base => 3 rang
64
6143f12f 65 ier = BSplCLib::EvalBsplineBasis(1, MyBSplOrder,
7fd59977 66 MyFlatKnots->Array1(), TParam(TParam.Lower()),
67 FirstNonZero, Base );
68 if (ier != 0) return Standard_False;
69 LastZero = FirstNonZero - 1;
70 FirstNonZero = 2*LastZero+1;
71
72 // (0.1) evaluation de CPrim
73 for (ii= 1; ii<= MyBSplOrder; ii++) {
74 CPrim += Base(2, ii) * MyPoles->Value(ii+LastZero).Coord();
75 }
76
77 // (1) Evaluation de la tension locale --------------------------------
78 Standard_Real NormeCPrim = CPrim.Modulus();
79 Standard_Real Hauteur, Difference;
80
81 if (MyHeight > 0) {Hauteur = MyHeight;} // it used in MVC to avoid Parametrization Problemes
82 else {
83 Ok = MyLaw.Value (TParam(TParam.Lower()), Hauteur);
84 if (!Ok) return Ok;
85 }
86 Difference = NormeCPrim - MyLengthSliding;
87
88 FTension(FTension.Lower()) = Hauteur * pow(Difference, 2) / MyLengthSliding ;
89
90 if (MyDerivativeOrder >= 1) {
91 // (2) Evaluation du gradient de la tension locale ----------------------
92 math_Vector GradDifference (1, 2*MyBSplOrder+MyNbValAux);
93 Standard_Real Xaux, Yaux, Facteur;
94
95 Xaux = CPrim.X() / NormeCPrim;
96 Yaux = CPrim.Y() / NormeCPrim;
97 Facteur = 2 * Hauteur * Difference / MyLengthSliding;
98
99 kk = FTension.Lower() + FirstNonZero;
100 jj = 1;
101 for (ii=1; ii<= MyBSplOrder; ii++) {
102 GradDifference(jj) = Base(2, ii) * Xaux;
103 FTension(kk) = Facteur * GradDifference(jj);
104 jj +=1;
105 GradDifference(jj) = Base(2, ii) * Yaux;
106 FTension(kk+1) = Facteur * GradDifference(jj);
107 jj += 1;
108 kk += 2;
109 }
110 if (MyNbValAux == 1) {
111 LastGradientIndex = FTension.Lower() + 2*MyPoles->Length() + 1;
112 GradDifference( GradDifference.Upper()) = (1 - pow( NormeCPrim/MyLengthSliding, 2));
113 FTension(LastGradientIndex) = Hauteur * GradDifference(GradDifference.Upper());
114 }
115
116 else { LastGradientIndex = FTension.Lower() + 2*MyPoles->Length(); }
117
118
119 if (MyDerivativeOrder >= 2) {
120
121// (3) Evaluation du Hessien de la tension locale ----------------------
122
123 Standard_Real FacteurX = Difference * (1-pow(Xaux,2)) / NormeCPrim;
124 Standard_Real FacteurY = Difference * (1-pow(Yaux,2)) / NormeCPrim;
125 Standard_Real FacteurXY = - Difference * Xaux*Yaux / NormeCPrim;
126 Standard_Real Produit;
127 Standard_Integer k1, k2;
128
129 Facteur = 2 * Hauteur / MyLengthSliding;
130
131 kk = FirstNonZero;
132 k2 = LastGradientIndex + (kk-1)*kk/2;
133
134 for (ii=2; ii<= 2*MyBSplOrder; ii+=2) {
135 k1 = k2+FirstNonZero;
136 k2 = k1+kk;
137 kk += 2;
138 for (jj=2; jj< ii; jj+=2) {
139 Produit = Base(2, ii/2) * Base(2, jj/2);
140
141 FTension(k1) = Facteur * ( GradDifference(ii-1)*GradDifference(jj-1)
142 + FacteurX * Produit) ; // derivation en XiXj
143 k1++;
144 FTension(k1) = Facteur * ( GradDifference(ii)*GradDifference(jj-1)
145 + FacteurXY * Produit); // derivation en YiXj
146 k1++;
147 FTension(k2) = Facteur * ( GradDifference(ii-1)*GradDifference(jj)
148 + FacteurXY * Produit); // derivation en XiYj
149 k2++;
150 FTension(k2) = Facteur * ( GradDifference(ii)*GradDifference(jj)
151 + FacteurY * Produit); // derivation en YiYj
152 k2++;
153 }
154 // cas ou jj = ii : remplisage en triangle
155 Produit = pow (Base(2, ii/2), 2);
156
157 FTension(k1) = Facteur * ( GradDifference(ii-1)*GradDifference(ii-1)
158 + FacteurX * Produit) ; // derivation en XiXi
159 FTension(k2) = Facteur * ( GradDifference(ii)*GradDifference(ii-1)
160 + FacteurXY * Produit); // derivation en XiYi
161 k2++;
162 FTension(k2) = Facteur * ( GradDifference(ii)*GradDifference(ii)
163 + FacteurY * Produit); // derivation en YiYi
164 }
165 if (MyNbValAux == 1) {
166 FacteurX = -2*CPrim.X()*Hauteur / pow (MyLengthSliding, 2);
167 FacteurY = -2*CPrim.Y()*Hauteur / pow (MyLengthSliding, 2);
168
169 ii = LastGradientIndex-FTension.Lower();
170 kk = LastGradientIndex + (ii-1)*ii/2 + FirstNonZero;
171 for (ii=1; ii<= MyBSplOrder; ii++) {
172 FTension(kk) = FacteurX * Base(2, ii);
173 kk ++;
174 FTension(kk) = FacteurY * Base(2, ii);
175 kk ++;
176 }
177 FTension(FTension.Upper()) = 2 * Hauteur * pow (NormeCPrim/MyLengthSliding, 2)
178 / MyLengthSliding;
179 }
180 }
181 }
182
183// sortie standard
184 return Ok;
185}
186