0025748: Parallel version of progress indicator
[occt.git] / src / Plate / Plate_Plate.cxx
1 // Created on: 1995-10-19
2 // Created by: Andre LIEUTIER
3 // Copyright (c) 1995-1999 Matra Datavision
4 // Copyright (c) 1999-2014 OPEN CASCADE SAS
5 //
6 // This file is part of Open CASCADE Technology software library.
7 //
8 // This library is free software; you can redistribute it and/or modify it under
9 // the terms of the GNU Lesser General Public License version 2.1 as published
10 // by the Free Software Foundation, with special exception defined in the file
11 // OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
12 // distribution for complete text of the license and disclaimer of any warranty.
13 //
14 // Alternatively, this file may be used under the terms of Open CASCADE
15 // commercial license or contractual agreement.
16
17 #include <gp_XY.hxx>
18 #include <gp_XYZ.hxx>
19 #include <math_Gauss.hxx>
20 #include <math_Matrix.hxx>
21 #include <math_Vector.hxx>
22 #include <Plate_FreeGtoCConstraint.hxx>
23 #include <Plate_GlobalTranslationConstraint.hxx>
24 #include <Plate_GtoCConstraint.hxx>
25 #include <Plate_LinearScalarConstraint.hxx>
26 #include <Plate_LinearXYZConstraint.hxx>
27 #include <Plate_LineConstraint.hxx>
28 #include <Plate_PinpointConstraint.hxx>
29 #include <Plate_PlaneConstraint.hxx>
30 #include <Plate_Plate.hxx>
31 #include <Plate_SampledCurveConstraint.hxx>
32 #include <Standard_ErrorHandler.hxx>
33
34 //=======================================================================
35 //function : Plate_Plate
36 //purpose  : 
37 //=======================================================================
38 Plate_Plate::Plate_Plate()
39 : order(0), n_el(0), n_dim(0),
40   solution(0),points(0),deru(0),derv(0),
41   OK(Standard_False),maxConstraintOrder(0),
42   Uold (1.e20),
43   Vold (1.e20),
44   U2 (0.0),
45   R (0.0),
46   L (0.0)
47 {
48   PolynomialPartOnly = Standard_False;
49   memset (ddu, 0, sizeof (ddu));
50   memset (ddv, 0, sizeof (ddv));
51 }
52
53 //=======================================================================
54 //function : Plate_Plate
55 //purpose  : 
56 //=======================================================================
57
58 Plate_Plate::Plate_Plate(const Plate_Plate& Ref)
59 : order(Ref.order),n_el(Ref.n_el),n_dim(Ref.n_dim),
60   solution(0),points(0),deru(0),derv(0),
61   OK (Ref.OK),
62   Uold (1.e20),
63   Vold (1.e20),
64   U2 (0.0),
65   R (0.0),
66   L (0.0)
67 {
68   Standard_Integer i;
69   if (Ref.OK) {
70     if (n_dim >0 && Ref.solution != 0) {
71       solution = new gp_XYZ[n_dim];
72       for(i=0; i<n_dim ;i++) {
73         Solution(i) = Ref.Solution(i);
74       }
75     }
76     
77     if (n_el >0) {
78       if (Ref.points != 0) {
79         points = new gp_XY[n_el];
80         for(i=0; i<n_el;i++) {
81           Points(i) = Ref.Points(i);
82         }
83       }
84       
85       if (Ref.deru != 0) {
86         deru = new Standard_Integer[n_el] ;
87         for (i = 0 ; i < n_el ; i++) {
88           Deru(i) = Ref.Deru(i);
89         }
90       }
91       
92       if (Ref.derv != 0) {
93         derv = new Standard_Integer[n_el] ;
94         for (i = 0 ; i < n_el ; i++) {
95           Derv(i) = Ref.Derv(i);
96         }
97       }
98     }
99   }
100
101   myConstraints = Ref.myConstraints;
102   myLXYZConstraints = Ref.myLXYZConstraints;
103   myLScalarConstraints = Ref.myLScalarConstraints;
104   maxConstraintOrder = Ref.maxConstraintOrder;
105   PolynomialPartOnly = Ref.PolynomialPartOnly;
106   for (i=0; i<10;i++) {
107     ddu[i]=Ref.ddu[i];
108     ddv[i]=Ref.ddv[i];
109   }
110 }
111 //=======================================================================
112 //function : Copy
113 //purpose  : 
114 //=======================================================================
115  Plate_Plate& Plate_Plate::Copy(const Plate_Plate& Ref) 
116 {
117   Init();
118    order = Ref.order;
119    n_el = Ref.n_el;
120    n_dim = Ref.n_dim;
121    OK = Ref.OK;
122   Standard_Integer i;
123   if (Ref.OK) {
124     if (n_dim >0 && Ref.solution != 0) {
125       solution = new gp_XYZ[n_dim];
126       for(i=0; i<n_dim ;i++) {
127         Solution(i) = Ref.Solution(i);
128       }
129     }
130     
131     if (n_el >0) {
132       if (Ref.points != 0) {
133         points = new gp_XY[n_el];
134         for(i=0; i<n_el;i++) {
135           Points(i) = Ref.Points(i);
136         }
137       }
138       
139       if (Ref.deru != 0) {
140         deru = new Standard_Integer[n_el] ;
141         for (i = 0 ; i < n_el ; i++) {
142           Deru(i) = Ref.Deru(i);
143         }
144       }
145       
146       if (Ref.derv != 0) {
147         derv = new Standard_Integer[n_el] ;
148         for (i = 0 ; i < n_el ; i++) {
149           Derv(i) = Ref.Derv(i);
150         }
151       }
152     }
153   }
154
155   myConstraints = Ref.myConstraints;
156   myLXYZConstraints = Ref.myLXYZConstraints;
157   myLScalarConstraints = Ref.myLScalarConstraints;
158   maxConstraintOrder = Ref.maxConstraintOrder;
159   PolynomialPartOnly = Ref.PolynomialPartOnly;
160
161    for (i=0; i<10;i++) {
162     ddu[i]=Ref.ddu[i];
163     ddv[i]=Ref.ddv[i];
164   }
165   return *this;
166 }
167 //=======================================================================
168 //function : Load
169 //purpose  : 
170 //=======================================================================
171
172 void Plate_Plate::Load(const Plate_PinpointConstraint& PConst)
173 {
174   OK = Standard_False;
175   n_el++;
176   myConstraints.Append(PConst);
177   Standard_Integer OrdreConst = PConst.Idu() + PConst.Idv();
178   if(maxConstraintOrder<OrdreConst) maxConstraintOrder = OrdreConst;
179 }
180
181  void Plate_Plate::Load(const Plate_LinearXYZConstraint& LXYZConst) 
182 {
183   OK = Standard_False;
184   n_el += LXYZConst.Coeff().RowLength();
185
186   myLXYZConstraints.Append(LXYZConst);
187   for(Standard_Integer j=1;j <= LXYZConst.GetPPC().Length() ; j++)
188     {
189       Standard_Integer OrdreConst = LXYZConst.GetPPC()(j).Idu() + LXYZConst.GetPPC()(j).Idv();
190       if(maxConstraintOrder<OrdreConst) maxConstraintOrder = OrdreConst;
191   }
192 }
193
194  void Plate_Plate::Load(const Plate_LinearScalarConstraint& LScalarConst) 
195 {
196   OK = Standard_False;
197   n_el += LScalarConst.Coeff().RowLength();
198   myLScalarConstraints.Append(LScalarConst);
199   for(Standard_Integer j=1;j <= LScalarConst.GetPPC().Length() ; j++)
200     {
201       Standard_Integer OrdreConst = LScalarConst.GetPPC()(j).Idu() + LScalarConst.GetPPC()(j).Idv();
202       if(maxConstraintOrder<OrdreConst) maxConstraintOrder = OrdreConst;
203   }
204 }
205
206 void Plate_Plate::Load(const Plate_LineConstraint& LConst)
207 {
208   Load(LConst.LSC());
209 }
210
211 void Plate_Plate::Load(const Plate_PlaneConstraint& PConst)
212 {
213   Load(PConst.LSC());
214 }
215
216 void Plate_Plate::Load(const Plate_SampledCurveConstraint& SCConst)
217 {
218   Load(SCConst.LXYZC());
219 }
220
221 void Plate_Plate::Load(const Plate_GtoCConstraint& GtoCConst)
222 {
223   for(Standard_Integer i=0;i< GtoCConst.nb_PPC();i++) 
224     Load(GtoCConst.GetPPC(i));
225 }
226
227 void Plate_Plate::Load(const Plate_FreeGtoCConstraint& FGtoCConst) 
228 {
229   Standard_Integer i ;
230   for( i=0;i< FGtoCConst.nb_PPC();i++) 
231     Load(FGtoCConst.GetPPC(i));
232   for(i=0;i< FGtoCConst.nb_LSC();i++) 
233     Load(FGtoCConst.LSC(i));
234 }
235
236 void Plate_Plate::Load(const Plate_GlobalTranslationConstraint& GTConst) 
237 {
238   Load(GTConst.LXYZC());
239 }
240
241 //=======================================================================
242 //function : SolveTI
243 //purpose  : to solve the set of constraints
244 //=======================================================================
245
246 void Plate_Plate::SolveTI(const Standard_Integer ord, 
247                           const Standard_Real anisotropie, 
248                           const Message_ProgressRange& theProgress)
249 {
250   Standard_Integer IterationNumber=0;
251   OK = Standard_False;
252   order = ord;
253   if(ord <=1) return;
254   if(ord > 9) return;
255   if(n_el<1) return;
256   if(anisotropie < 1.e-6) return;
257   if(anisotropie > 1.e+6) return;
258
259 // computation of the bounding box of the 2d PPconstraints
260   Standard_Real xmin,xmax,ymin,ymax;
261   UVBox(xmin,xmax,ymin,ymax);
262
263   Standard_Real du = 0.5*(xmax - xmin);
264   if(anisotropie >1.) du *= anisotropie;
265   if(du < 1.e-10) return;
266   ddu[0] = 1;
267   Standard_Integer i ;
268   for( i=1;i<=9;i++) ddu[i] = ddu[i-1] / du;
269
270   Standard_Real dv = 0.5*(ymax - ymin);
271   if(anisotropie <1.) dv /= anisotropie;
272   if(dv < 1.e-10) return;
273   ddv[0] = 1;
274   for(i=1;i<=9;i++) ddv[i] = ddv[i-1] / dv;
275
276   if(myLScalarConstraints.IsEmpty())
277     {
278       if(myLXYZConstraints.IsEmpty())
279         SolveTI1(IterationNumber, theProgress);
280       else
281         SolveTI2(IterationNumber, theProgress);
282     }
283   else
284     SolveTI3(IterationNumber, theProgress);
285
286 }
287
288 //=======================================================================
289 //function : SolveTI1
290 //purpose  : to solve the set of constraints in the easiest case,
291 //           only PinPointConstraints are loaded
292 //=======================================================================
293
294 void Plate_Plate::SolveTI1(const Standard_Integer IterationNumber,
295                            const Message_ProgressRange& theProgress)
296 {
297 // computation of square matrix members
298
299
300   n_dim = n_el + order*(order+1)/2;
301   math_Matrix mat(0, n_dim-1, 0, n_dim-1, 0.);
302   
303   delete [] (gp_XY*)points;
304   points = new gp_XY[n_el];
305   Standard_Integer i ;
306   for( i=0; i<n_el;i++)  Points(i) = myConstraints(i+1).Pnt2d();
307
308   delete [] (Standard_Integer*)deru;
309   deru = new Standard_Integer[n_el];
310   for(i=0; i<n_el;i++)  Deru(i) = myConstraints(i+1).Idu();
311
312   delete [] (Standard_Integer*)derv;
313   derv = new Standard_Integer[n_el];
314   for(i=0; i<n_el;i++)  Derv(i) = myConstraints(i+1).Idv();
315
316   for(i=0; i<n_el;i++) {
317     for(Standard_Integer j=0;j<i;j++) {
318       Standard_Real signe = 1;
319       if ( ((Deru(j)+Derv(j))%2) == 1) signe = -1;
320       Standard_Integer iu =  Deru(i) + Deru(j);
321       Standard_Integer iv =  Derv(i) + Derv(j);
322       mat(i,j) = signe * SolEm(Points(i) - Points(j),iu,iv);
323     }
324   }
325   
326   i = n_el;
327   for(Standard_Integer iu = 0; iu< order; iu++) {
328     for(Standard_Integer iv =0; iu+iv < order; iv++) {
329       for(Standard_Integer j=0;j<n_el;j++) {
330         Standard_Integer idu = Deru(j);
331         Standard_Integer idv = Derv(j);
332         mat(i,j) = Polm (Points(j), iu, iv, idu, idv);
333       }
334       i++;
335     }
336   }
337
338   for(i=0;i<n_dim;i++) {
339     for(Standard_Integer j = i+1; j<n_dim;j++) {
340       mat(i,j) = mat(j,i);
341     }
342   }
343
344 // initialisation of the Gauss algorithm
345   Standard_Real pivot_max = 1.e-12;
346   OK = Standard_True;     
347
348   Message_ProgressScope aScope (theProgress, NULL, 10);
349   math_Gauss algo_gauss(mat,pivot_max, aScope.Next (7));
350
351   if (aScope.UserBreak())
352   {
353     OK = Standard_False;
354     return;
355   }
356
357   if(!algo_gauss.IsDone()) {
358     Standard_Integer nbm = order*(order+1)/2;
359     for(i=n_el;i<n_el+nbm;i++) {
360       mat(i,i) = 1.e-8;
361     }
362     pivot_max = 1.e-18;
363
364     math_Gauss thealgo(mat,pivot_max, aScope.Next (3));
365
366     if (aScope.UserBreak())
367     {
368       OK = Standard_False;
369       return;
370     }
371     algo_gauss = thealgo;
372     OK = algo_gauss.IsDone();
373   }
374
375   if (OK) {
376 //   computation of the linear system solution for the X, Y and Z coordinates
377     math_Vector sec_member( 0, n_dim-1, 0.);
378     math_Vector sol(0,n_dim-1);
379
380     delete [] (gp_XYZ*) solution;
381     solution = new gp_XYZ[n_dim];
382     
383     for(Standard_Integer icoor=1; icoor<=3;icoor++) {
384       for(i=0;i<n_el;i++) {
385         sec_member(i) = myConstraints(i+1).Value().Coord(icoor);
386       }
387       algo_gauss.Solve(sec_member, sol);
388       //alr iteration pour affiner la solution
389       {
390         math_Vector sol1(0,n_dim-1);
391         math_Vector sec_member1(0,n_dim-1);
392         for(i=1;i<=IterationNumber;i++)
393           {
394             sec_member1 = sec_member - mat*sol;
395             algo_gauss.Solve(sec_member1, sol1);
396             sol += sol1;
397           }
398       }
399       //finalr
400  
401       for(i=0;i<n_dim;i++) {
402         Solution(i).SetCoord (icoor, sol(i));
403       }
404     }
405   }
406 }
407
408 //=======================================================================
409 //function : SolveTI2
410 //purpose  : to solve the set of constraints in the medium case,
411 //           LinearXYZ constraints are provided but no LinearScalar one
412 //=======================================================================
413
414 void Plate_Plate::SolveTI2(const Standard_Integer IterationNumber,
415                            const Message_ProgressRange& theProgress)
416 {
417 // computation of square matrix members
418
419   Standard_Integer nCC1 = myConstraints.Length();
420   Standard_Integer nCC2 = 0;
421   Standard_Integer i ;
422   for( i = 1; i<= myLXYZConstraints.Length(); i++)
423     nCC2 += myLXYZConstraints(i).Coeff().ColLength();
424
425   Standard_Integer n_dimat = nCC1 + nCC2 + order*(order+1)/2;
426
427   
428   delete [] (gp_XY*)points;
429   points = new gp_XY[n_el];
430   delete [] (Standard_Integer*)deru;
431   deru = new Standard_Integer[n_el];
432   delete [] (Standard_Integer*)derv;
433   derv = new Standard_Integer[n_el];
434
435
436   for(i=0; i< nCC1;i++)
437     {
438       Points(i) = myConstraints(i+1).Pnt2d();
439       Deru(i) = myConstraints(i+1).Idu();
440       Derv(i) = myConstraints(i+1).Idv();
441     }
442
443   Standard_Integer k = nCC1;
444   for( i = 1; i<= myLXYZConstraints.Length(); i++)
445     for(Standard_Integer j=1;j <= myLXYZConstraints(i).GetPPC().Length() ; j++)
446       {
447         Points(k) = myLXYZConstraints(i).GetPPC()(j).Pnt2d();
448         Deru(k) =  myLXYZConstraints(i).GetPPC()(j).Idu();
449         Derv(k) =  myLXYZConstraints(i).GetPPC()(j).Idv();
450         k++;
451       }
452
453   math_Matrix mat(0, n_dimat-1, 0, n_dimat-1, 0.);
454
455   fillXYZmatrix(mat,0,0,nCC1,nCC2);
456
457
458 // initialisation of the Gauss algorithm
459   Standard_Real pivot_max = 1.e-12;
460   OK = Standard_True;      // ************ JHH
461
462   Message_ProgressScope aScope (theProgress, NULL, 10);
463   math_Gauss algo_gauss(mat,pivot_max, aScope.Next (7));
464   
465   if (aScope.UserBreak())
466   {
467     OK = Standard_False;
468     return;
469   }
470
471   if(!algo_gauss.IsDone()) {
472     for(i=nCC1+nCC2;i<n_dimat;i++) {
473       mat(i,i) = 1.e-8;
474     }
475     pivot_max = 1.e-18;
476
477     math_Gauss thealgo1(mat,pivot_max, aScope.Next (3));
478
479     if (aScope.UserBreak())
480     {
481       OK = Standard_False;
482       return;
483     }
484     algo_gauss = thealgo1;
485     OK = algo_gauss.IsDone();
486   }
487
488   if (OK) {
489 //   computation of the linear system solution for the X, Y and Z coordinates
490     math_Vector sec_member( 0, n_dimat-1, 0.);
491     math_Vector sol(0,n_dimat-1);
492
493     delete [] (gp_XYZ*) solution;
494     n_dim = n_el+order*(order+1)/2;
495     solution = new gp_XYZ[n_dim];
496     
497     for(Standard_Integer icoor=1; icoor<=3;icoor++) {
498       for(i=0;i<nCC1;i++) {
499         sec_member(i) = myConstraints(i+1).Value().Coord(icoor);
500       }
501
502       k = nCC1;
503       for(i = 1; i<= myLXYZConstraints.Length(); i++) {
504         for(Standard_Integer irow =1; irow <= myLXYZConstraints(i).Coeff().ColLength(); irow++) {
505           for(Standard_Integer icol=1; icol<=myLXYZConstraints(i).Coeff().RowLength();icol++)
506             sec_member(k) += myLXYZConstraints(i).Coeff()(irow,icol) 
507               *  myLXYZConstraints(i).GetPPC()(icol).Value().Coord(icoor);
508           k++;
509         }
510       }
511
512       algo_gauss.Solve(sec_member, sol);
513       //alr iteration pour affiner la solution
514       {
515         math_Vector sol1(0,n_dimat-1);
516         math_Vector sec_member1(0,n_dimat-1);
517         for(i=1;i<=IterationNumber;i++)
518           {
519             sec_member1 = sec_member - mat*sol;
520             algo_gauss.Solve(sec_member1, sol1);
521             sol += sol1;
522         }
523       }
524       //finalr
525
526       for(i=0;i<nCC1;i++) Solution(i).SetCoord (icoor, sol(i));
527
528       Standard_Integer kSolution = nCC1;
529       Standard_Integer ksol = nCC1;
530
531       for(i = 1; i<= myLXYZConstraints.Length(); i++) {
532         for(Standard_Integer icol=1; icol<=myLXYZConstraints(i).Coeff().RowLength();icol++){
533           Standard_Real vsol = 0;
534           for(Standard_Integer irow =1; irow <= myLXYZConstraints(i).Coeff().ColLength(); irow++)
535             vsol += myLXYZConstraints(i).Coeff()(irow,icol)*sol(ksol+irow-1);
536           Solution(kSolution).SetCoord (icoor, vsol);
537           kSolution++;
538         }
539         ksol += myLXYZConstraints(i).Coeff().ColLength();
540       }
541
542       for(i=0;i<order*(order+1)/2; i++) {
543         Solution(n_el+i).SetCoord (icoor, sol(ksol+i));
544       }
545     }
546   }
547 }
548
549 //=======================================================================
550 //function : SolveTI3
551 //purpose  : to solve the set of constraints in the most general situation
552 //=======================================================================
553
554 void Plate_Plate::SolveTI3(const Standard_Integer IterationNumber,
555                            const Message_ProgressRange& theProgress)
556 {
557 // computation of square matrix members
558
559   Standard_Integer nCC1 = myConstraints.Length();
560
561   Standard_Integer nCC2 = 0;
562   Standard_Integer i ;
563   for( i = 1; i<= myLXYZConstraints.Length(); i++)
564     nCC2 += myLXYZConstraints(i).Coeff().ColLength();
565
566   Standard_Integer nCC3 = 0;
567   for(i = 1; i<= myLScalarConstraints.Length(); i++)
568     nCC3 += myLScalarConstraints(i).Coeff().ColLength();
569
570   Standard_Integer nbm = order*(order+1)/2;
571   Standard_Integer n_dimsousmat = nCC1 + nCC2 + nbm ;
572   Standard_Integer n_dimat =3*n_dimsousmat + nCC3;
573
574   
575   delete [] (gp_XY*)points;
576   points = new gp_XY[n_el];
577   delete [] (Standard_Integer*)deru;
578   deru = new Standard_Integer[n_el];
579   delete [] (Standard_Integer*)derv;
580   derv = new Standard_Integer[n_el];
581
582
583   for(i=0; i< nCC1;i++)
584     {
585       Points(i) = myConstraints(i+1).Pnt2d();
586       Deru(i) = myConstraints(i+1).Idu();
587       Derv(i) = myConstraints(i+1).Idv();
588     }
589
590   Standard_Integer k = nCC1;
591   for(i = 1; i<= myLXYZConstraints.Length(); i++)
592     for(Standard_Integer j=1;j <= myLXYZConstraints(i).GetPPC().Length() ; j++)
593       {
594         Points(k) = myLXYZConstraints(i).GetPPC()(j).Pnt2d();
595         Deru(k) =  myLXYZConstraints(i).GetPPC()(j).Idu();
596         Derv(k) =  myLXYZConstraints(i).GetPPC()(j).Idv();
597         k++;
598       }
599   Standard_Integer nPPC2 = k;
600   for(i = 1; i<= myLScalarConstraints.Length(); i++)
601     for(Standard_Integer j=1;j <= myLScalarConstraints(i).GetPPC().Length() ; j++)
602       {
603         Points(k) = myLScalarConstraints(i).GetPPC()(j).Pnt2d();
604         Deru(k) =  myLScalarConstraints(i).GetPPC()(j).Idu();
605         Derv(k) =  myLScalarConstraints(i).GetPPC()(j).Idv();
606         k++;
607       }
608
609   math_Matrix mat(0, n_dimat-1, 0, n_dimat-1, 0.);
610
611   fillXYZmatrix(mat,0,0,nCC1,nCC2);
612   fillXYZmatrix(mat,n_dimsousmat,n_dimsousmat,nCC1,nCC2);
613   fillXYZmatrix(mat,2*n_dimsousmat,2*n_dimsousmat,nCC1,nCC2);
614
615   k = 3*n_dimsousmat;
616   Standard_Integer kppc = nPPC2;
617   Standard_Integer j ;
618   for(i = 1; i<= myLScalarConstraints.Length(); i++) {
619     for( j=0;j<nCC1;j++){
620
621       math_Vector vmat(1,myLScalarConstraints(i).GetPPC().Length());
622
623       for(Standard_Integer ippc=1;ippc <= myLScalarConstraints(i).GetPPC().Length() ; ippc++) {
624         Standard_Real signe = 1;
625         if ( ((Deru(j)+Derv(j))%2) == 1) signe = -1;
626         Standard_Integer iu =  Deru(kppc+ippc-1) + Deru(j);
627         Standard_Integer iv =  Derv(kppc+ippc-1) + Derv(j);
628         vmat(ippc) =  signe * SolEm(Points(kppc+ippc-1) - Points(j),iu,iv);
629       }
630
631       for(Standard_Integer irow=1;irow <= myLScalarConstraints(i).Coeff().ColLength() ; irow++)
632         for(Standard_Integer icol=1;icol <= myLScalarConstraints(i).Coeff().RowLength() ; icol++){
633           mat(k+irow-1,j) += myLScalarConstraints(i).Coeff()(irow,icol).X()*vmat(icol);
634           mat(k+irow-1,n_dimsousmat+j) += myLScalarConstraints(i).Coeff()(irow,icol).Y()*vmat(icol);
635           mat(k+irow-1,2*n_dimsousmat+j) += myLScalarConstraints(i).Coeff()(irow,icol).Z()*vmat(icol);
636         }
637     }
638
639     Standard_Integer k2 = nCC1;
640     Standard_Integer kppc2 = nCC1;
641     Standard_Integer i2 ;
642     for( i2 = 1; i2<=myLXYZConstraints.Length() ; i2++){
643
644       math_Matrix tmpmat(1,myLScalarConstraints(i).GetPPC().Length(),1,myLXYZConstraints(i2).GetPPC().Length() );
645
646       for(Standard_Integer ippc=1;ippc <= myLScalarConstraints(i).GetPPC().Length() ; ippc++)
647         for(Standard_Integer ippc2=1;ippc2 <= myLXYZConstraints(i2).GetPPC().Length() ; ippc2++){
648           Standard_Real signe = 1;
649           if ( ((Deru(kppc2+ippc2-1)+Derv(kppc2+ippc2-1))%2) == 1) signe = -1;
650           Standard_Integer iu =  Deru(kppc+ippc-1) + Deru(kppc2+ippc2-1);
651           Standard_Integer iv =  Derv(kppc+ippc-1) + Derv(kppc2+ippc2-1);
652           tmpmat(ippc,ippc2) = signe * SolEm(Points(kppc+ippc-1) - Points(kppc2+ippc2-1),iu,iv);
653         }
654
655       for(Standard_Integer irow=1;irow <= myLScalarConstraints(i).Coeff().ColLength() ; irow++)
656         for(Standard_Integer irow2=1;irow2 <= myLXYZConstraints(i2).Coeff().ColLength() ; irow2++)
657           for(Standard_Integer icol=1;icol <= myLScalarConstraints(i).Coeff().RowLength() ; icol++)
658             for(Standard_Integer icol2=1;icol2 <= myLXYZConstraints(i2).Coeff().RowLength() ; icol2++){
659               mat(k+irow-1,k2+irow2-1) += 
660                 myLScalarConstraints(i).Coeff()(irow,icol).X()*myLXYZConstraints(i2).Coeff()(irow2,icol2)*tmpmat(icol,icol2);
661               mat(k+irow-1,n_dimsousmat+k2+irow2-1) += 
662                 myLScalarConstraints(i).Coeff()(irow,icol).Y()*myLXYZConstraints(i2).Coeff()(irow2,icol2)*tmpmat(icol,icol2);
663               mat(k+irow-1,2*n_dimsousmat+k2+irow2-1) += 
664                 myLScalarConstraints(i).Coeff()(irow,icol).Z()*myLXYZConstraints(i2).Coeff()(irow2,icol2)*tmpmat(icol,icol2);
665             }
666
667       k2 += myLXYZConstraints(i2).Coeff().ColLength();
668       kppc2 += myLXYZConstraints(i2).Coeff().RowLength();
669     }
670
671
672
673     j = nCC1+nCC2;
674     for(Standard_Integer iu = 0; iu< order; iu++)
675       for(Standard_Integer iv =0; iu+iv < order; iv++) {
676
677         math_Vector vmat(1,myLScalarConstraints(i).GetPPC().Length());
678         for(Standard_Integer ippc=1;ippc <= myLScalarConstraints(i).GetPPC().Length() ; ippc++){
679           Standard_Integer idu = Deru(kppc+ippc-1);
680           Standard_Integer idv = Derv(kppc+ippc-1);
681           vmat(ippc) = Polm (Points(kppc+ippc-1),iu,iv,idu,idv);
682         }
683
684         for(Standard_Integer irow=1;irow <= myLScalarConstraints(i).Coeff().ColLength() ; irow++)
685           for(Standard_Integer icol=1;icol <= myLScalarConstraints(i).Coeff().RowLength() ; icol++){
686             mat(k+irow-1,j) += myLScalarConstraints(i).Coeff()(irow,icol).X()*vmat(icol);
687             mat(k+irow-1,n_dimsousmat+j) += myLScalarConstraints(i).Coeff()(irow,icol).Y()*vmat(icol);
688             mat(k+irow-1,2*n_dimsousmat+j) += myLScalarConstraints(i).Coeff()(irow,icol).Z()*vmat(icol);
689           }
690
691       j++;
692     }
693
694
695     k2 = 3*n_dimsousmat;
696     kppc2 = nPPC2;
697     for(i2 = 1; i2<=i ; i2++){
698
699       math_Matrix tmpmat(1,myLScalarConstraints(i).GetPPC().Length(),1,myLScalarConstraints(i2).GetPPC().Length() );
700
701       for(Standard_Integer ippc=1;ippc <= myLScalarConstraints(i).GetPPC().Length() ; ippc++)
702         for(Standard_Integer ippc2=1;ippc2 <= myLScalarConstraints(i2).GetPPC().Length() ; ippc2++){
703           Standard_Real signe = 1;
704           if ( ((Deru(kppc2+ippc2-1)+Derv(kppc2+ippc2-1))%2) == 1) signe = -1;
705           Standard_Integer a_iu =  Deru(kppc+ippc-1) + Deru(kppc2+ippc2-1);
706           Standard_Integer iv =  Derv(kppc+ippc-1) + Derv(kppc2+ippc2-1);
707           tmpmat(ippc,ippc2) = signe * SolEm(Points(kppc+ippc-1) - Points(kppc2+ippc2-1),a_iu,iv);
708         }
709
710       for(Standard_Integer irow=1;irow <= myLScalarConstraints(i).Coeff().ColLength() ; irow++)
711         for(Standard_Integer irow2=1;irow2 <= myLScalarConstraints(i2).Coeff().ColLength() ; irow2++)
712           for(Standard_Integer icol=1;icol <= myLScalarConstraints(i).Coeff().RowLength() ; icol++)
713             for(Standard_Integer icol2=1;icol2 <= myLScalarConstraints(i2).Coeff().RowLength() ; icol2++){
714               mat(k+irow-1,k2+irow2-1) += 
715                 myLScalarConstraints(i).Coeff()(irow,icol)*myLScalarConstraints(i2).Coeff()(irow2,icol2)*tmpmat(icol,icol2);
716             }
717
718       k2 += myLScalarConstraints(i2).Coeff().ColLength();
719       kppc2 += myLScalarConstraints(i2).Coeff().RowLength();
720     }
721
722     k += myLScalarConstraints(i).Coeff().ColLength();
723     kppc += myLScalarConstraints(i).Coeff().RowLength();
724   }
725
726   for( j=3*n_dimsousmat;j<n_dimat;j++)
727     for(i=0;i<j;i++)
728       mat(i,j)= mat(j,i);
729
730
731
732 // initialisation of the Gauss algorithm
733   Standard_Real pivot_max = 1.e-12;
734   OK = Standard_True;      // ************ JHH
735
736   Message_ProgressScope aScope (theProgress, NULL, 10);
737   math_Gauss algo_gauss(mat,pivot_max, aScope.Next (7));
738   
739   if (aScope.UserBreak())
740   {
741     OK = Standard_False;
742     return;
743   }
744
745   if(!algo_gauss.IsDone()) {
746     for(i=nCC1+nCC2;i<nCC1+nCC2+nbm;i++) {
747       mat(i,i) = 1.e-8;
748       mat(n_dimsousmat+i,n_dimsousmat+i) = 1.e-8;
749       mat(2*n_dimsousmat+i,2*n_dimsousmat+i) = 1.e-8;
750     }
751     pivot_max = 1.e-18;
752
753     math_Gauss thealgo2(mat,pivot_max, aScope.Next (3));
754
755     if (aScope.UserBreak())
756     {
757       OK = Standard_False;
758       return;
759     }
760     algo_gauss = thealgo2;
761     OK = algo_gauss.IsDone();
762   }
763
764   if (OK) {
765 //   computation of the linear system solution for the X, Y and Z coordinates
766     math_Vector sec_member( 0, n_dimat-1, 0.);
767     math_Vector sol(0,n_dimat-1);
768
769     delete [] (gp_XYZ*) solution;
770     n_dim = n_el+order*(order+1)/2;
771     solution = new gp_XYZ[n_dim];
772     
773     Standard_Integer icoor ;
774     for( icoor=1; icoor<=3;icoor++){
775       for(i=0;i<nCC1;i++)
776         sec_member((icoor-1)*n_dimsousmat+i) = myConstraints(i+1).Value().Coord(icoor);
777    
778
779       k = nCC1;
780       for(i = 1; i<= myLXYZConstraints.Length(); i++)
781         for(Standard_Integer irow =1; irow <= myLXYZConstraints(i).Coeff().ColLength(); irow++) {
782           for(Standard_Integer icol=1; icol<=myLXYZConstraints(i).Coeff().RowLength();icol++)
783             sec_member((icoor-1)*n_dimsousmat+k) += myLXYZConstraints(i).Coeff()(irow,icol) 
784               *  myLXYZConstraints(i).GetPPC()(icol).Value().Coord(icoor);
785           k++;
786         }
787     }
788     k = 3*n_dimsousmat;
789     for(i = 1; i<= myLScalarConstraints.Length(); i++)
790       for(Standard_Integer irow =1; irow <= myLScalarConstraints(i).Coeff().ColLength(); irow++) {
791         for(Standard_Integer icol=1; icol<=myLScalarConstraints(i).Coeff().RowLength();icol++)
792           sec_member(k) += myLScalarConstraints(i).Coeff()(irow,icol) 
793             *  myLScalarConstraints(i).GetPPC()(icol).Value();
794         k++;
795       }
796     
797     algo_gauss.Solve(sec_member, sol);
798     // iteration to refine the solution
799     {
800       math_Vector sol1(0,n_dimat-1);
801       math_Vector sec_member1(0,n_dimat-1);
802       for(i=1;i<=IterationNumber;i++)
803         {
804           sec_member1 = sec_member - mat*sol;
805           algo_gauss.Solve(sec_member1, sol1);
806           sol += sol1;
807         }
808     }
809     
810     for(icoor=1; icoor<=3;icoor++){
811       for(i=0;i<nCC1;i++) Solution(i).SetCoord (icoor, sol((icoor-1)*n_dimsousmat+i));
812
813       Standard_Integer kSolution = nCC1;
814       Standard_Integer ksol = nCC1;
815
816       for(i = 1; i<= myLXYZConstraints.Length(); i++) {
817         for(Standard_Integer icol=1; icol<=myLXYZConstraints(i).Coeff().RowLength();icol++){
818           Standard_Real vsol = 0;
819           for(Standard_Integer irow =1; irow <= myLXYZConstraints(i).Coeff().ColLength(); irow++)
820             vsol += myLXYZConstraints(i).Coeff()(irow,icol)*sol((icoor-1)*n_dimsousmat+ksol+irow-1);
821           Solution(kSolution).SetCoord (icoor, vsol);
822           kSolution++;
823         }
824         ksol += myLXYZConstraints(i).Coeff().ColLength();
825       }
826
827       ksol = nCC1+nCC2;
828       for(i=0;i<order*(order+1)/2; i++) {
829         Solution(n_el+i).SetCoord (icoor, sol((icoor-1)*n_dimsousmat+ksol+i));
830       }
831     }
832
833     Standard_Integer ksol = 3*n_dimsousmat;
834     Standard_Integer kSolution = nPPC2;
835     for(i = 1; i<= myLScalarConstraints.Length(); i++) {
836       for(Standard_Integer icol=1; icol<=myLScalarConstraints(i).Coeff().RowLength();icol++){
837         gp_XYZ Vsol(0.,0.,0.);
838         for(Standard_Integer irow =1; irow <= myLScalarConstraints(i).Coeff().ColLength(); irow++)
839           Vsol += myLScalarConstraints(i).Coeff()(irow,icol)*sol(ksol+irow-1);
840         Solution(kSolution) = Vsol;
841         kSolution++;
842       }
843       ksol += myLScalarConstraints(i).Coeff().ColLength();
844     }
845   }
846 }
847
848 //=======================================================================
849 //function : fillXYZmatrix
850 //purpose  : 
851 //=======================================================================
852 void Plate_Plate::fillXYZmatrix(math_Matrix &mat, 
853                                 const Standard_Integer i0,
854                                 const Standard_Integer j0,
855                                 const Standard_Integer ncc1,
856                                 const Standard_Integer ncc2) const
857 {
858   Standard_Integer i,j ;
859   for( i=0; i<ncc1;i++) {
860     for( j=0;j<i;j++) {
861       Standard_Real signe = 1;
862       if ( ((Deru(j)+Derv(j))%2) == 1) signe = -1;
863       Standard_Integer iu =  Deru(i) + Deru(j);
864       Standard_Integer iv =  Derv(i) + Derv(j);
865       mat(i0+i,j0+j) = signe * SolEm(Points(i) - Points(j),iu,iv);
866     }
867   }
868   
869   Standard_Integer k = ncc1;
870   Standard_Integer kppc = ncc1;
871   for( i = 1; i<= myLXYZConstraints.Length(); i++){
872
873     for(Standard_Integer a_j=0; a_j < ncc1; a_j++){
874
875       math_Vector vmat(1,myLXYZConstraints(i).GetPPC().Length());
876
877       for(Standard_Integer ippc=1;ippc <= myLXYZConstraints(i).GetPPC().Length() ; ippc++) {
878         Standard_Real signe = 1;
879         if ( ((Deru(a_j)+Derv(a_j))%2) == 1) signe = -1;
880         Standard_Integer iu =  Deru(kppc+ippc-1) + Deru(a_j);
881         Standard_Integer iv =  Derv(kppc+ippc-1) + Derv(a_j);
882         vmat(ippc) =  signe * SolEm(Points(kppc+ippc-1) - Points(a_j),iu,iv);
883       }
884
885       for(Standard_Integer irow=1;irow <= myLXYZConstraints(i).Coeff().ColLength() ; irow++)
886         for(Standard_Integer icol=1;icol <= myLXYZConstraints(i).Coeff().RowLength() ; icol++)
887           mat(i0+k+irow-1,j0+a_j) += myLXYZConstraints(i).Coeff()(irow,icol)*vmat(icol);
888     }
889
890     Standard_Integer k2 = ncc1;
891     Standard_Integer kppc2 = ncc1;
892     for(Standard_Integer i2 = 1; i2<= i; i2++){
893
894       math_Matrix tmpmat(1,myLXYZConstraints(i).GetPPC().Length(),1,myLXYZConstraints(i2).GetPPC().Length() );
895
896       for(Standard_Integer ippc=1;ippc <= myLXYZConstraints(i).GetPPC().Length() ; ippc++)
897         for(Standard_Integer ippc2=1;ippc2 <= myLXYZConstraints(i2).GetPPC().Length() ; ippc2++){
898           Standard_Real signe = 1;
899           if ( ((Deru(kppc2+ippc2-1)+Derv(kppc2+ippc2-1))%2) == 1) signe = -1;
900           Standard_Integer iu =  Deru(kppc+ippc-1) + Deru(kppc2+ippc2-1);
901           Standard_Integer iv =  Derv(kppc+ippc-1) + Derv(kppc2+ippc2-1);
902           tmpmat(ippc,ippc2) = signe * SolEm(Points(kppc+ippc-1) - Points(kppc2+ippc2-1),iu,iv);
903         }
904
905       for(Standard_Integer irow=1;irow <= myLXYZConstraints(i).Coeff().ColLength() ; irow++)
906         for(Standard_Integer irow2=1;irow2 <= myLXYZConstraints(i2).Coeff().ColLength() ; irow2++)
907           for(Standard_Integer icol=1;icol <= myLXYZConstraints(i).Coeff().RowLength() ; icol++)
908             for(Standard_Integer icol2=1;icol2 <= myLXYZConstraints(i2).Coeff().RowLength() ; icol2++)
909               mat(i0+k+irow-1,j0+k2+irow2-1) += 
910                 myLXYZConstraints(i).Coeff()(irow,icol)*myLXYZConstraints(i2).Coeff()(irow2,icol2)*tmpmat(icol,icol2);
911       
912
913       k2 += myLXYZConstraints(i2).Coeff().ColLength();
914       kppc2 += myLXYZConstraints(i2).Coeff().RowLength();
915     }
916
917     k += myLXYZConstraints(i).Coeff().ColLength();
918     kppc += myLXYZConstraints(i).Coeff().RowLength();
919   }
920
921
922
923
924   i = ncc1+ncc2;
925   for(Standard_Integer iu = 0; iu< order; iu++)
926     for(Standard_Integer iv =0; iu+iv < order; iv++) {
927       for(Standard_Integer a_j=0; a_j < ncc1; a_j++) {
928         Standard_Integer idu = Deru(a_j);
929         Standard_Integer idv = Derv(a_j);
930         mat(i0+i,j0+a_j) = Polm (Points(a_j), iu, iv, idu, idv);
931       }
932
933       Standard_Integer k2 = ncc1;
934       Standard_Integer kppc2 = ncc1;
935       for(Standard_Integer i2 = 1; i2<= myLXYZConstraints.Length(); i2++){
936         math_Vector vmat(1,myLXYZConstraints(i2).GetPPC().Length());
937         for(Standard_Integer ippc2=1;ippc2 <= myLXYZConstraints(i2).GetPPC().Length() ; ippc2++){
938           Standard_Integer idu = Deru(kppc2+ippc2-1);
939           Standard_Integer idv = Derv(kppc2+ippc2-1);
940           vmat(ippc2) = Polm (Points(kppc2+ippc2-1),iu,iv,idu,idv);
941         }
942
943         for(Standard_Integer irow2=1;irow2 <= myLXYZConstraints(i2).Coeff().ColLength() ; irow2++)
944           for(Standard_Integer icol2=1;icol2 <= myLXYZConstraints(i2).Coeff().RowLength() ; icol2++)
945             mat(i0+i,j0+k2+irow2-1) += myLXYZConstraints(i2).Coeff()(irow2,icol2)*vmat(icol2);
946
947         k2 += myLXYZConstraints(i2).Coeff().ColLength();
948         kppc2 += myLXYZConstraints(i2).Coeff().RowLength();
949         }
950
951       i++;
952     }
953
954   Standard_Integer n_dimat = ncc1 + ncc2 + order*(order+1)/2;
955
956   for(i=0;i<n_dimat;i++) {
957     for(Standard_Integer a_j = i+1; a_j < n_dimat; a_j++) {
958       mat(i0+i,j0+a_j) = mat(i0+a_j,j0+i);
959     }
960   }
961
962 }
963
964 //=======================================================================
965 //function : IsDone
966 //purpose  : 
967 //=======================================================================
968
969 Standard_Boolean Plate_Plate::IsDone() const 
970 {
971   return OK;
972 }
973
974
975 //=======================================================================
976 //function : destroy
977 //purpose  : 
978 //=======================================================================
979
980 void Plate_Plate::destroy()
981 {
982   Init();
983 }
984
985 //=======================================================================
986 //function : Init
987 //purpose  : 
988 //=======================================================================
989
990 void Plate_Plate::Init()
991
992   myConstraints.Clear();
993   myLXYZConstraints.Clear();
994   myLScalarConstraints.Clear();
995
996
997   delete [] (gp_XYZ*)solution;
998   solution = 0;
999
1000   delete [] (gp_XY*)points;
1001   points = 0;
1002
1003   delete [] (Standard_Integer*)deru;
1004   deru = 0;
1005
1006   delete [] (Standard_Integer*)derv;
1007   derv = 0;
1008
1009   order = 0;
1010   n_el = 0;
1011   n_dim = 0;
1012   OK = Standard_True;
1013   maxConstraintOrder=0;
1014 }
1015
1016 //=======================================================================
1017 //function : Evaluate
1018 //purpose  : 
1019 //=======================================================================
1020
1021 gp_XYZ Plate_Plate::Evaluate(const gp_XY& point2d) const 
1022 {
1023   if(solution == 0) return gp_XYZ(0,0,0);
1024   if(!OK) return gp_XYZ(0,0,0);
1025
1026   gp_XYZ valeur(0,0,0);
1027
1028   if(!PolynomialPartOnly)
1029     {
1030       for(Standard_Integer i=0; i<n_el;i++)
1031         {
1032           Standard_Real signe = 1;
1033           if ( ((Deru(i)+Derv(i))%2) == 1) signe = -1;
1034           valeur += Solution(i) *  (signe*SolEm(point2d - Points(i), Deru(i), Derv(i))) ;
1035         }
1036     }
1037   Standard_Integer i = n_el;
1038   for(Standard_Integer idu = 0; idu< order; idu++)
1039   for(Standard_Integer idv =0; idu+idv < order; idv++)
1040   {
1041     valeur += Solution(i) * Polm( point2d, idu,idv,0,0);
1042     i++;
1043   }
1044   return valeur;
1045 }
1046
1047 //=======================================================================
1048 //function : EvaluateDerivative
1049 //purpose  : 
1050 //=======================================================================
1051
1052 gp_XYZ Plate_Plate::EvaluateDerivative(const gp_XY& point2d, const Standard_Integer iu, const Standard_Integer iv) const 
1053 {
1054   if(solution == 0) return gp_XYZ(0,0,0);
1055   if(!OK) return gp_XYZ(0,0,0);
1056
1057   gp_XYZ valeur(0,0,0);
1058   if(!PolynomialPartOnly)
1059     {
1060       for(Standard_Integer i=0; i<n_el;i++)
1061         {
1062           Standard_Real signe = 1;
1063           if ( ((Deru(i)+Derv(i))%2) == 1) signe = -1;
1064           valeur += Solution(i) *  (signe*SolEm(point2d - Points(i), Deru(i)+iu, Derv(i)+iv)) ;
1065         }
1066     }
1067   Standard_Integer i = n_el;
1068   for(Standard_Integer idu = 0; idu< order; idu++)
1069   for(Standard_Integer idv =0; idu+idv < order; idv++)
1070   {
1071     valeur += Solution(i) * Polm( point2d, idu,idv, iu,iv);
1072     i++;
1073   }
1074   return valeur;
1075 }
1076 //=======================================================================
1077 //function : Plate_Plate::CoefPol
1078 //purpose  : give back the array of power basis coefficient of
1079 // the polynomial part of the Plate function
1080 //=======================================================================
1081
1082  void Plate_Plate::CoefPol(Handle(TColgp_HArray2OfXYZ)& Coefs) const
1083 {
1084   Coefs = new TColgp_HArray2OfXYZ(0,order-1,0,order-1,gp_XYZ(0.,0.,0.));
1085   Standard_Integer i = n_el;
1086   for(Standard_Integer iu = 0; iu< order; iu++)
1087   for(Standard_Integer iv =0; iu+iv < order; iv++)
1088   {
1089     Coefs->ChangeValue(iu,iv) = Solution(i)*ddu[iu]*ddv[iv];
1090     //Coefs->ChangeValue(idu,idv) = Solution(i);
1091     // it is necessary to reset this line if one remove factors in method Polm.
1092     i++;
1093   }
1094   
1095 }
1096 //=======================================================================
1097 //function : Plate_Plate::Continuity
1098 //purpose  : give back the continuity order of the Plate function
1099 //=======================================================================
1100
1101  Standard_Integer Plate_Plate::Continuity() const
1102 {
1103   return 2*order - 3 - maxConstraintOrder;
1104 }
1105
1106 //=======================================================================
1107 //function : Plate_Plate::SolEm
1108 //purpose  : compute the (iu,iv)th derivative of the fundamental solution
1109 // of Laplcian at the power order
1110 //=======================================================================
1111
1112
1113 Standard_Real Plate_Plate::SolEm(const gp_XY& point2d, const Standard_Integer iu, const Standard_Integer iv) const 
1114 {
1115   Plate_Plate* aThis = const_cast<Plate_Plate*>(this);
1116   Standard_Real U,V;
1117   Standard_Integer IU,IV;
1118
1119   if(iv>iu)
1120   {
1121     // SolEm is symmetric in (u<->v) : we swap u and v if iv>iu
1122     // to avoid some code 
1123     IU = iv;
1124     IV = iu;
1125     U = point2d.Y() *ddv[1];
1126     V = point2d.X() *ddu[1];
1127   }
1128   else
1129   {
1130     IU = iu;
1131     IV = iv;
1132     U = point2d.X() *ddu[1];
1133     V = point2d.Y() *ddv[1];
1134   }
1135   
1136    if((U==Uold)&&(V==Vold) )
1137      {
1138        if (R<1.e-20) return 0;
1139      }
1140   else
1141     {
1142         aThis->Uold = U;
1143         aThis->Vold = V;
1144         aThis->U2 = U*U;
1145         aThis->R = U2+V*V;
1146         if (R<1.e-20) return 0;
1147         aThis->L = log(R);
1148       }
1149   Standard_Real DUV = 0;
1150
1151   Standard_Integer m = order;
1152   Standard_Integer mm1 = m-1;
1153   Standard_Real &r = aThis->R;
1154
1155
1156   //Standard_Real pr = pow(R, mm1 - IU - IV);
1157   // this expression takes a lot of time 
1158   //(does not take into account a small integer value of the exponent)
1159   //
1160
1161   Standard_Integer expo =  mm1 - IU - IV;
1162   Standard_Real pr;
1163   if(expo<0)
1164         {
1165         pr = R;
1166         for(Standard_Integer i=1;i<-expo;i++) pr *= R;
1167         pr = 1./pr;
1168         }
1169   else if(expo>0)
1170         {
1171         pr = R;
1172         for(Standard_Integer i=1;i<expo;i++) pr *= R;
1173         }
1174   else pr = 1.;
1175         
1176
1177   switch (IU)
1178   {
1179   case 0:
1180     switch (IV)
1181     {
1182     case 0:
1183       {
1184       DUV = pr*L;
1185       }
1186       break;
1187
1188     default:
1189       break;
1190     }
1191   break;
1192
1193   case 1:
1194     switch (IV)
1195     {
1196     case 0:
1197       {
1198       DUV = 2*pr*U*(1+L*mm1);
1199       }
1200       break;
1201
1202     case 1:
1203       {
1204       Standard_Real m2 = m*m;
1205       //DUV = 4*pr*U*V*(-3+2*L+2*m-3*L*m+L*m2);
1206       DUV = 4*pr*U*V*((2*m-3)+(m2-3*m+2)*L);
1207       }
1208       break;
1209
1210     default:
1211       break;
1212     }
1213     break;
1214
1215   case 2:
1216     switch (IV)
1217     {
1218     case 0:
1219       {
1220       Standard_Real m2 = m*m;
1221       DUV = 2*pr*(R-L*R+L*m*R-6*U2+4*L*U2+4*m*U2-6*L*m*U2+2*L*m2*U2);
1222       }
1223       break;
1224
1225     case 1:
1226       {
1227       Standard_Real m2 = m*m;
1228       Standard_Real m3 = m2*m;
1229       DUV = -3*R+2*L*R+2*m*R-3*L*m*R+L*m2*R+22*U2-12*L*U2-24*m*U2+22*L*m*U2+6*m2*U2-12*L*m2*U2+2*L*m3*U2;
1230       DUV = DUV * 4* pr*V;
1231       }
1232       break;
1233
1234     case 2:
1235       {
1236       Standard_Real m2 = m*m;
1237       Standard_Real m3 = m2*m;
1238       Standard_Real m4 = m2*m2;
1239       Standard_Real V2 = V*V;
1240       Standard_Real R2 = R*R;
1241       DUV = -3*R2+2*L*R2+2*m*R2-3*L*m*R2+L*m2*R2+22*R*U2-12*L*R*U2-24*m*R*U2+22*L*m*R*U2+6*m2*R*U2-12*L*m2*R*U2;
1242       DUV += 2*L*m3*R*U2+22*R*V2-12*L*R*V2-24*m*R*V2+22*L*m*R*V2+6*m2*R*V2-12*L*m2*R*V2+2*L*m3*R*V2-200*U2*V2+96*L*U2*V2;
1243       DUV += 280*m*U2*V2-200*L*m*U2*V2-120*m2*U2*V2+140*L*m2*U2*V2+16*m3*U2*V2-40*L*m3*U2*V2+4*L*m4*U2*V2;
1244       DUV = 4*pr*DUV;
1245       }
1246       break;
1247
1248     default:
1249       break;
1250     }
1251     break;
1252
1253   case 3:
1254     switch (IV)
1255     {
1256     case 0:
1257       {
1258       Standard_Real m2 = m*m;
1259       Standard_Real m3 = m2*m;
1260       DUV = -9*R+6*L*R+6*m*R-9*L*m*R+3*L*m2*R+22*U2-12*L*U2-24*m*U2+22*L*m*U2+6*m2*U2-12*L*m2*U2+2*L*m3*U2;
1261       DUV = DUV * 4* pr*U;
1262       }
1263       break;
1264
1265     case 1:      
1266       {
1267       Standard_Real m2 = m*m;
1268       Standard_Real m3 = m2*m;
1269       Standard_Real m4 = m2*m2;
1270       DUV = 33*R-18*L*R-36*m*R+33*L*m*R+9*m2*R-18*L*m2*R+3*L*m3*R-100*U2+48*L*U2+140*m*U2-100*L*m*U2-60*m2*U2+70*L*m2*U2;
1271       DUV += 8*m3*U2-20*L*m3*U2+2*L*m4*U2;
1272       DUV = 8*pr*U*V*DUV;
1273       }
1274       break;
1275
1276     case 2:      
1277       {
1278       Standard_Real m2 = m*m;
1279       Standard_Real m3 = m2*m;
1280       Standard_Real m4 = m2*m2;
1281       Standard_Real m5 = m4*m;
1282       Standard_Real ru2 = R*U2;
1283       Standard_Real v2 = V*V;
1284       Standard_Real rv2 = R*v2;
1285       Standard_Real u2v2 = v2*U2;
1286       Standard_Real r2 = r*r;
1287
1288       // copy-paste the mathematics 
1289           DUV = 
1290      -100*ru2 + 48*L*ru2 + 140*m*ru2 - 100*L*m*ru2 - 60*m2*ru2 + 70*L*m2*ru2 + 8*m3*ru2 - 
1291      20*L*m3*ru2 + 2*L*m4*ru2 - 300*rv2 + 144*L*rv2 + 420*m*rv2 - 300*L*m*rv2 - 180*m2*rv2 + 210*L*m2*rv2 + 
1292      24*m3*rv2 - 60*L*m3*rv2 + 6*L*m4*rv2 + 33*r2 - 18*L*r2 - 36*m*r2 + 33*L*m*r2 + 9*m2*r2 - 18*L*m2*r2 + 
1293      3*L*m3*r2 + 1096*u2v2 - 480*L*u2v2 - 1800*m*u2v2 + 1096*L*m*u2v2 + 1020*m2*u2v2 - 900*L*m2*u2v2 - 
1294      240*m3*u2v2 + 340*L*m3*u2v2 + 20*m4*u2v2 - 60*L*m4*u2v2 + 4*L*m5*u2v2;
1295
1296       DUV = 8*pr*U*DUV;
1297       }
1298       break;
1299
1300     case 3:      
1301       {
1302       Standard_Real m2 = m*m;
1303       Standard_Real m3 = m2*m;
1304       Standard_Real m4 = m2*m2;
1305       Standard_Real m5 = m3*m2;
1306       Standard_Real m6 = m3*m3;
1307       Standard_Real ru2 = r*U2;
1308       Standard_Real v2 = V*V;
1309       Standard_Real rv2 = R*v2;
1310       Standard_Real u2v2 = v2*U2;
1311       Standard_Real r2 = r*r;
1312
1313      // copy-paste the mathematics  
1314       DUV = 
1315                 1644*ru2 - 720*L*ru2 - 2700*m*ru2 + 1644*L*m*ru2 + 1530*m2*ru2 - 1350*L*m2*ru2 - 
1316      360*m3*ru2 + 510*L*m3*ru2 + 30*m4*ru2 - 90*L*m4*ru2 + 6*L*m5*ru2 + 1644*rv2 - 720*L*rv2 - 2700*m*rv2 + 
1317      1644*L*m*rv2 + 1530*m2*rv2 - 1350*L*m2*rv2 - 360*m3*rv2 + 510*L*m3*rv2 + 30*m4*rv2 - 90*L*m4*rv2 + 
1318      6*L*m5*rv2 - 450*r2 + 216*L*r2 + 630*m*r2 - 450*L*m*r2 - 270*m2*r2 + 315*L*m2*r2 + 36*m3*r2 - 90*L*m3*r2 + 
1319      9*L*m4*r2 - 7056*u2v2 + 2880*L*u2v2 + 12992*m*u2v2 - 7056*L*m*u2v2 - 8820*m2*u2v2 + 6496*L*m2*u2v2 + 
1320      2800*m3*u2v2 - 2940*L*m3*u2v2 - 420*m4*u2v2 + 700*L*m4*u2v2 + 24*m5*u2v2 - 84*L*m5*u2v2 + 4*L*m6*u2v2;
1321          
1322       DUV = 16*pr*U*V*DUV;
1323       }
1324       break;
1325
1326     default:
1327       break;
1328     }
1329   break;
1330
1331   case 4:
1332     switch (IV)
1333     {
1334     case 0:     
1335       { 
1336       Standard_Real m2 = m*m;
1337       Standard_Real m3 = m2*m;
1338       Standard_Real m4 = m2*m2;
1339       Standard_Real U4 = U2*U2;
1340       Standard_Real R2 = R*R;
1341       DUV = -9*R2+6*L*R2+6*m*R2-9*L*m*R2+3*L*m2*R2+132*R*U2-72*L*R*U2-144*m*R*U2+132*L*m*R*U2+36*m2*R*U2-72*L*m2*R*U2;
1342       DUV += 12*L*m3*R*U2-200*U4+96*L*U4+280*m*U4-200*L*m*U4-120*m2*U4+140*L*m2*U4+16*m3*U4-40*L*m3*U4+4*L*m4*U4;
1343       DUV = 4*pr*DUV;
1344       }
1345       break;
1346
1347    case 1:     
1348       { 
1349       Standard_Real m2 = m*m;
1350       Standard_Real m3 = m2*m;
1351       Standard_Real m4 = m2*m2;
1352       Standard_Real m5 = m2*m3;
1353       Standard_Real u4 = U2*U2;
1354       Standard_Real ru2 = R*U2;
1355       Standard_Real r2 = R*R;
1356  
1357           // copy-paste the mathematics  
1358           DUV = 
1359                 -600*ru2 + 288*L*ru2 + 840*m*ru2 - 600*L*m*ru2 - 360*m2*ru2 + 420*L*m2*ru2 + 48*m3*ru2 - 
1360      120*L*m3*ru2 + 12*L*m4*ru2 + 33*r2 - 18*L*r2 - 36*m*r2 + 33*L*m*r2 + 9*m2*r2 - 18*L*m2*r2 + 3*L*m3*r2 + 
1361      1096*u4 - 480*L*u4 - 1800*m*u4 + 1096*L*m*u4 + 1020*m2*u4 - 900*L*m2*u4 - 240*m3*u4 + 340*L*m3*u4 + 20*m4*u4 - 
1362      60*L*m4*u4 + 4*L*m5*u4;  
1363
1364       DUV = 8*pr*V*DUV;
1365       }
1366       break;
1367
1368   case 2:     
1369       { 
1370       Standard_Real m2 = m*m;
1371       Standard_Real m3 = m2*m;
1372       Standard_Real m4 = m2*m2;
1373       Standard_Real m5 = m2*m3;
1374       Standard_Real m6 = m3*m3;
1375       Standard_Real u4 = U2*U2;
1376       Standard_Real r2 = r*r;
1377       Standard_Real r3 = r2*r;
1378       Standard_Real v2 = V*V;
1379       Standard_Real u2v2 = v2*U2;
1380       Standard_Real ru2v2 = R*u2v2;
1381       Standard_Real u4v2 = u4*v2;
1382       Standard_Real r2u2 = r2*U2;
1383       Standard_Real ru4 = r*u4;
1384       Standard_Real r2v2 = r2*v2;
1385
1386           // copy-paste the mathematics  
1387           DUV = 
1388                   6576*ru2v2 - 2880*L*ru2v2 - 10800*m*ru2v2 + 6576*L*m*ru2v2 + 6120*m2*ru2v2 - 5400*L*m2*ru2v2 - 
1389      1440*m3*ru2v2 + 2040*L*m3*ru2v2 + 120*m4*ru2v2 - 360*L*m4*ru2v2 + 24*L*m5*ru2v2 + 1096*ru4 - 480*L*ru4 - 
1390      1800*m*ru4 + 1096*L*m*ru4 + 1020*m2*ru4 - 900*L*m2*ru4 - 240*m3*ru4 + 340*L*m3*ru4 + 20*m4*ru4 - 60*L*m4*ru4 + 
1391      4*L*m5*ru4 - 600*r2u2 + 288*L*r2u2 + 840*m*r2u2 - 600*L*m*r2u2 - 360*m2*r2u2 + 420*L*m2*r2u2 + 48*m3*r2u2 - 
1392      120*L*m3*r2u2 + 12*L*m4*r2u2 - 300*r2v2 + 144*L*r2v2 + 420*m*r2v2 - 300*L*m*r2v2 - 180*m2*r2v2 + 210*L*m2*r2v2 + 
1393      24*m3*r2v2 - 60*L*m3*r2v2 + 6*L*m4*r2v2 + 33*r3 - 18*L*r3 - 36*m*r3 + 33*L*m*r3 + 9*m2*r3 - 18*L*m2*r3 + 
1394      3*L*m3*r3 - 14112*u4v2 + 5760*L*u4v2 + 25984*m*u4v2 - 14112*L*m*u4v2 - 17640*m2*u4v2 + 12992*L*m2*u4v2 + 
1395      5600*m3*u4v2 - 5880*L*m3*u4v2 - 840*m4*u4v2 + 1400*L*m4*u4v2 + 48*m5*u4v2 - 168*L*m5*u4v2 + 8*L*m6*u4v2;
1396
1397       DUV = 8*pr*DUV;
1398       }
1399       break;
1400
1401  case 3:     
1402       { 
1403       Standard_Real m2 = m*m;
1404       Standard_Real m3 = m2*m;
1405       Standard_Real m4 = m2*m2;
1406       Standard_Real m5 = m2*m3;
1407       Standard_Real m6 = m3*m3;
1408       Standard_Real m7 = m3*m4;
1409       Standard_Real u4 = U2*U2;
1410       Standard_Real r2 = r*r;
1411       Standard_Real r3 = r2*r;
1412       Standard_Real v2 = V*V;
1413       Standard_Real u2v2 = v2*U2;
1414       Standard_Real ru2v2 = R*u2v2;
1415       Standard_Real u4v2 = u4*v2;
1416       Standard_Real r2u2 = r2*U2;
1417       Standard_Real r2v2 = r2*v2;
1418       Standard_Real ru4 = r*u4;
1419
1420           // copy-paste the mathematics  
1421       DUV = 
1422                   -42336*ru2v2 + 17280*L*ru2v2 + 77952*m*ru2v2 - 42336*L*m*ru2v2 - 52920*m2*ru2v2 + 
1423      38976*L*m2*ru2v2 + 16800*m3*ru2v2 - 17640*L*m3*ru2v2 - 2520*m4*ru2v2 + 4200*L*m4*ru2v2 + 144*m5*ru2v2 - 
1424      504*L*m5*ru2v2 + 24*L*m6*ru2v2 - 21168*ru4 + 8640*L*ru4 + 38976*m*ru4 - 21168*L*m*ru4 - 26460*m2*ru4 + 
1425      19488*L*m2*ru4 + 8400*m3*ru4 - 8820*L*m3*ru4 - 1260*m4*ru4 + 2100*L*m4*ru4 + 72*m5*ru4 - 252*L*m5*ru4 + 
1426      12*L*m6*ru4 + 9864*r2u2 - 4320*L*r2u2 - 16200*m*r2u2 + 9864*L*m*r2u2 + 9180*m2*r2u2 - 8100*L*m2*r2u2 - 
1427      2160*m3*r2u2 + 3060*L*m3*r2u2 + 180*m4*r2u2 - 540*L*m4*r2u2 + 36*L*m5*r2u2 + 1644*r2v2 - 720*L*r2v2 - 
1428      2700*m*r2v2 + 1644*L*m*r2v2 + 1530*m2*r2v2 - 1350*L*m2*r2v2 - 360*m3*r2v2 + 510*L*m3*r2v2 + 30*m4*r2v2 - 
1429      90*L*m4*r2v2 + 6*L*m5*r2v2 - 450*r3 + 216*L*r3 + 630*m*r3 - 450*L*m*r3 - 270*m2*r3 + 315*L*m2*r3 + 36*m3*r3 - 
1430      90*L*m3*r3 + 9*L*m4*r3 + 104544*u4v2 - 40320*L*u4v2 - 210112*m*u4v2 + 104544*L*m*u4v2 + 162456*m2*u4v2 - 
1431      105056*L*m2*u4v2 - 62720*m3*u4v2 + 54152*L*m3*u4v2 + 12880*m4*u4v2 - 15680*L*m4*u4v2 - 1344*m5*u4v2 + 
1432      2576*L*m5*u4v2 + 56*m6*u4v2 - 224*L*m6*u4v2 + 8*L*m7*u4v2;
1433
1434       DUV = 16*pr*V*DUV;
1435       }
1436       break;
1437
1438    default:
1439       break;
1440     }
1441   break;
1442
1443   case 5:
1444     switch (IV)
1445     {
1446     case 0:     
1447       { 
1448       Standard_Real m2 = m*m;
1449       Standard_Real m3 = m2*m;
1450       Standard_Real m4 = m2*m2;
1451       Standard_Real m5 = m2*m3;
1452       Standard_Real u4 = U2*U2;
1453       Standard_Real r2 = R*R;
1454       Standard_Real ru2 = R*U2;
1455
1456           // copy-paste the mathematics  
1457       DUV = 
1458      -1000*ru2 + 480*L*ru2 + 1400*m*ru2 - 1000*L*m*ru2 - 600*m2*ru2 + 700*L*m2*ru2 + 80*m3*ru2 - 
1459      200*L*m3*ru2 + 20*L*m4*ru2 + 165*r2 - 90*L*r2 - 180*m*r2 + 165*L*m*r2 + 45*m2*r2 - 90*L*m2*r2 + 15*L*m3*r2 + 
1460      1096*u4 - 480*L*u4 - 1800*m*u4 + 1096*L*m*u4 + 1020*m2*u4 - 900*L*m2*u4 - 240*m3*u4 + 340*L*m3*u4 + 20*m4*u4 - 
1461      60*L*m4*u4 + 4*L*m5*u4;
1462           
1463           DUV = 8*pr*U*DUV;
1464       }
1465       break;
1466
1467    case 1:     
1468       { 
1469       Standard_Real m2 = m*m;
1470       Standard_Real m3 = m2*m;
1471       Standard_Real m4 = m2*m2;
1472       Standard_Real m5 = m2*m3;
1473       Standard_Real m6 = m3*m3;
1474       Standard_Real u4 = U2*U2;
1475       Standard_Real ru2 = r*U2;
1476       Standard_Real r2 = r*r;
1477
1478
1479           // copy-paste the mathematics 
1480      DUV = 
1481      5480*ru2 - 2400*L*ru2 - 9000*m*ru2 + 5480*L*m*ru2 + 5100*m2*ru2 - 4500*L*m2*ru2 - 1200*m3*ru2 + 
1482      1700*L*m3*ru2 + 100*m4*ru2 - 300*L*m4*ru2 + 20*L*m5*ru2 - 750*r2 + 360*L*r2 + 1050*m*r2 - 750*L*m*r2 - 
1483      450*m2*r2 + 525*L*m2*r2 + 60*m3*r2 - 150*L*m3*r2 + 15*L*m4*r2 - 7056*u4 + 2880*L*u4 + 12992*m*u4 - 7056*L*m*u4 - 
1484      8820*m2*u4 + 6496*L*m2*u4 + 2800*m3*u4 - 2940*L*m3*u4 - 420*m4*u4 + 700*L*m4*u4 + 24*m5*u4 - 84*L*m5*u4 + 
1485      4*L*m6*u4;          
1486
1487       DUV = 16*pr*U*V*DUV;
1488       }
1489       break;
1490
1491   case 2:     
1492       { 
1493       Standard_Real m2 = m*m;
1494       Standard_Real m3 = m2*m;
1495       Standard_Real m4 = m2*m2;
1496       Standard_Real m5 = m2*m3;
1497       Standard_Real m6 = m3*m3;
1498       Standard_Real m7 = m3*m4;
1499       Standard_Real u4 = U2*U2;
1500       Standard_Real r2 = r*r;
1501       Standard_Real r3 = r2*r;
1502       Standard_Real v2 = V*V;
1503       Standard_Real u2v2 = v2*U2;
1504       Standard_Real ru2v2 = R*u2v2;
1505       Standard_Real u4v2 = u4*v2;
1506       Standard_Real r2u2 = r2*U2;
1507       Standard_Real r2v2 = r2*v2;
1508       Standard_Real ru4 = r*u4;
1509
1510           // copy-paste the mathematics  
1511           DUV = 
1512                   
1513      -70560*ru2v2 + 28800*L*ru2v2 + 129920*m*ru2v2 - 70560*L*m*ru2v2 - 88200*m2*ru2v2 + 
1514      64960*L*m2*ru2v2 + 28000*m3*ru2v2 - 29400*L*m3*ru2v2 - 4200*m4*ru2v2 + 7000*L*m4*ru2v2 + 240*m5*ru2v2 - 
1515      840*L*m5*ru2v2 + 40*L*m6*ru2v2 - 7056*ru4 + 2880*L*ru4 + 12992*m*ru4 - 7056*L*m*ru4 - 8820*m2*ru4 + 
1516      6496*L*m2*ru4 + 2800*m3*ru4 - 2940*L*m3*ru4 - 420*m4*ru4 + 700*L*m4*ru4 + 24*m5*ru4 - 84*L*m5*ru4 + 4*L*m6*ru4 + 
1517      5480*r2u2 - 2400*L*r2u2 - 9000*m*r2u2 + 5480*L*m*r2u2 + 5100*m2*r2u2 - 4500*L*m2*r2u2 - 1200*m3*r2u2 + 
1518      1700*L*m3*r2u2 + 100*m4*r2u2 - 300*L*m4*r2u2 + 20*L*m5*r2u2 + 8220*r2v2 - 3600*L*r2v2 - 13500*m*r2v2 + 
1519      8220*L*m*r2v2 + 7650*m2*r2v2 - 6750*L*m2*r2v2 - 1800*m3*r2v2 + 2550*L*m3*r2v2 + 150*m4*r2v2 - 450*L*m4*r2v2 + 
1520      30*L*m5*r2v2 - 750*r3 + 360*L*r3 + 1050*m*r3 - 750*L*m*r3 - 450*m2*r3 + 525*L*m2*r3 + 60*m3*r3 - 150*L*m3*r3 + 
1521      15*L*m4*r3 + 104544*u4v2 - 40320*L*u4v2 - 210112*m*u4v2 + 104544*L*m*u4v2 + 162456*m2*u4v2 - 105056*L*m2*u4v2 - 
1522      62720*m3*u4v2 + 54152*L*m3*u4v2 + 12880*m4*u4v2 - 15680*L*m4*u4v2 - 1344*m5*u4v2 + 2576*L*m5*u4v2 + 56*m6*u4v2 - 
1523      224*L*m6*u4v2 + 8*L*m7*u4v2;
1524           
1525      DUV = 16*pr*U*DUV;
1526       }
1527       break;
1528
1529    default:
1530       break;
1531     }
1532   break;
1533
1534   case 6:
1535     switch (IV)
1536     {
1537     case 0:     
1538       { 
1539       Standard_Real m2 = m*m;
1540       Standard_Real m3 = m2*m;
1541       Standard_Real m4 = m2*m2;
1542       Standard_Real m5 = m2*m3;
1543       Standard_Real m6 = m3*m3;
1544       Standard_Real u4 = U2*U2;
1545       Standard_Real u6 = U2*u4;
1546       Standard_Real r2 = r*r;
1547       Standard_Real r3 = r2*r;
1548       Standard_Real r2u2 = r2*U2;
1549       Standard_Real ru4 = r*u4;
1550  
1551           // copy-paste the mathematics 
1552       DUV = 
1553                 16440*ru4 - 7200*L*ru4 - 27000*m*ru4 + 16440*L*m*ru4 + 15300*m2*ru4 - 13500*L*m2*ru4 - 
1554      3600*m3*ru4 + 5100*L*m3*ru4 + 300*m4*ru4 - 900*L*m4*ru4 + 60*L*m5*ru4 - 4500*r2u2 + 2160*L*r2u2 + 6300*m*r2u2 - 
1555      4500*L*m*r2u2 - 2700*m2*r2u2 + 3150*L*m2*r2u2 + 360*m3*r2u2 - 900*L*m3*r2u2 + 90*L*m4*r2u2 + 165*r3 - 90*L*r3 - 
1556      180*m*r3 + 165*L*m*r3 + 45*m2*r3 - 90*L*m2*r3 + 15*L*m3*r3 - 14112*u6 + 5760*L*u6 + 25984*m*u6 - 14112*L*m*u6 - 
1557      17640*m2*u6 + 12992*L*m2*u6 + 5600*m3*u6 - 5880*L*m3*u6 - 840*m4*u6 + 1400*L*m4*u6 + 48*m5*u6 - 168*L*m5*u6 + 
1558      8*L*m6*u6;
1559
1560       DUV = 8*pr*DUV;
1561       }
1562       break;
1563
1564    default:
1565       break;
1566     }
1567   break;
1568
1569   default:
1570   break;
1571   }
1572   
1573   return DUV * ddu[iu]*ddv[iv];
1574
1575 }
1576
1577
1578 //=======================================================================
1579 //function : UVBox
1580 //purpose  : 
1581 //=======================================================================
1582
1583 void Plate_Plate::UVBox(Standard_Real& UMin, Standard_Real& UMax,
1584                         Standard_Real& VMin, Standard_Real& VMax) const 
1585 {
1586   Standard_Integer i ;
1587   const Standard_Real Bmin = 1.e-3;
1588   UMin =  myConstraints(1).Pnt2d().X();
1589   UMax =  UMin;
1590   VMin =  myConstraints(1).Pnt2d().Y();
1591   VMax =  VMin;
1592
1593   for( i=2; i<=myConstraints.Length();i++)
1594   {
1595     Standard_Real x = myConstraints(i).Pnt2d().X();
1596     if(x<UMin) UMin = x;
1597     if(x>UMax) UMax = x;
1598     Standard_Real y = myConstraints(i).Pnt2d().Y();
1599     if(y<VMin) VMin = y;
1600     if(y>VMax) VMax = y; 
1601   }
1602
1603   for(i=1; i<=myLXYZConstraints.Length();i++)
1604     for(Standard_Integer j=1;j<= myLXYZConstraints(i).GetPPC().Length(); j++)
1605   {
1606     Standard_Real x = myLXYZConstraints(i).GetPPC()(j).Pnt2d().X();
1607     if(x<UMin) UMin = x;
1608     if(x>UMax) UMax = x;
1609     Standard_Real y = myLXYZConstraints(i).GetPPC()(j).Pnt2d().Y();
1610     if(y<VMin) VMin = y;
1611     if(y>VMax) VMax = y; 
1612   }
1613
1614   for(i=1; i<=myLScalarConstraints.Length();i++)
1615     for(Standard_Integer j=1;j<= myLScalarConstraints(i).GetPPC().Length(); j++)
1616   {
1617     Standard_Real x = myLScalarConstraints(i).GetPPC()(j).Pnt2d().X();
1618     if(x<UMin) UMin = x;
1619     if(x>UMax) UMax = x;
1620     Standard_Real y = myLScalarConstraints(i).GetPPC()(j).Pnt2d().Y();
1621     if(y<VMin) VMin = y;
1622     if(y>VMax) VMax = y; 
1623   }
1624  
1625
1626   if(UMax-UMin < Bmin)
1627     {
1628       Standard_Real UM = 0.5*(UMin+UMax);
1629       UMin = UM - 0.5*Bmin;
1630       UMax = UM + 0.5*Bmin;
1631     }
1632   if(VMax-VMin < Bmin)
1633     {
1634       Standard_Real VM = 0.5*(VMin+VMax);
1635       VMin = VM - 0.5*Bmin;
1636       VMax = VM + 0.5*Bmin;
1637     }
1638 }
1639
1640 //=======================================================================
1641 //function : UVConstraints
1642 //purpose  : 
1643 //=======================================================================
1644
1645 void Plate_Plate::UVConstraints(TColgp_SequenceOfXY& Seq) const 
1646 {
1647   for (Standard_Integer i=1;i<=myConstraints.Length();i++) {
1648     if ((myConstraints.Value(i).Idu()==0) && (myConstraints.Value(i).Idv()==0))
1649       Seq.Append((myConstraints.Value(i)).Pnt2d());
1650   }
1651 }
1652 //=======================================================================
1653
1654 void Plate_Plate::SetPolynomialPartOnly(const Standard_Boolean PPOnly)
1655 {
1656   PolynomialPartOnly = PPOnly;
1657 }