1
0
mirror of https://git.dev.opencascade.org/repos/occt.git synced 2025-04-07 18:30:55 +03:00
occt/src/PLib/PLib_HermitJacobi.cxx
abv 92efcf78a6 0026936: Drawbacks of inlining in new type system in OCCT 7.0 -- automatic
Automatic restore of IMPLEMENT_STANDARD_RTTIEXT macro (upgrade -rtti)
2015-12-04 14:15:06 +03:00

310 lines
10 KiB
C++

// Created on: 1997-10-22
// Created by: Sergey SOKOLOV
// Copyright (c) 1997-1999 Matra Datavision
// Copyright (c) 1999-2014 OPEN CASCADE SAS
//
// This file is part of Open CASCADE Technology software library.
//
// This library is free software; you can redistribute it and/or modify it under
// the terms of the GNU Lesser General Public License version 2.1 as published
// by the Free Software Foundation, with special exception defined in the file
// OCCT_LGPL_EXCEPTION.txt. Consult the file LICENSE_LGPL_21.txt included in OCCT
// distribution for complete text of the license and disclaimer of any warranty.
//
// Alternatively, this file may be used under the terms of Open CASCADE
// commercial license or contractual agreement.
#include <NCollection_LocalArray.hxx>
#include <PLib.hxx>
#include <PLib_HermitJacobi.hxx>
#include <PLib_JacobiPolynomial.hxx>
#include <Standard_ConstructionError.hxx>
#include <Standard_Type.hxx>
#include <TColStd_HArray1OfReal.hxx>
IMPLEMENT_STANDARD_RTTIEXT(PLib_HermitJacobi,PLib_Base)
//=======================================================================
//function : PLib_HermitJacobi
//purpose :
//=======================================================================
PLib_HermitJacobi::PLib_HermitJacobi(const Standard_Integer WorkDegree,
const GeomAbs_Shape ConstraintOrder) :
myH(1,2*(PLib::NivConstr(ConstraintOrder)+1),
1,2*(PLib::NivConstr(ConstraintOrder)+1)),
myWCoeff(1,2*(PLib::NivConstr(ConstraintOrder)+1)+1)
{
Standard_Integer NivConstr = PLib::NivConstr(ConstraintOrder);
PLib::HermiteCoefficients(-1.,1.,NivConstr,NivConstr,myH);
myJacobi = new PLib_JacobiPolynomial (WorkDegree,ConstraintOrder);
myWCoeff.Init(0.);
myWCoeff(1) = 1.;
switch(NivConstr) {
case 0: myWCoeff(3) = -1.; break;
case 1: myWCoeff(3) = -2.; myWCoeff(5) = 1.; break;
case 2: myWCoeff(3) = -3.; myWCoeff(5) = 3.; myWCoeff(7) = -1.; break;
}
}
//=======================================================================
//function : MaxError
//purpose :
//=======================================================================
Standard_Real PLib_HermitJacobi::MaxError(const Standard_Integer Dimension,
Standard_Real& HermJacCoeff,
const Standard_Integer NewDegree) const
{
return myJacobi->MaxError(Dimension,HermJacCoeff,NewDegree);
}
//=======================================================================
//function : ReduceDegree
//purpose :
//=======================================================================
void PLib_HermitJacobi::ReduceDegree(const Standard_Integer Dimension,
const Standard_Integer MaxDegree,
const Standard_Real Tol,
Standard_Real& HermJacCoeff,
Standard_Integer& NewDegree,
Standard_Real& MaxError) const
{
myJacobi->ReduceDegree(Dimension,MaxDegree,Tol,
HermJacCoeff,NewDegree,MaxError);
}
//=======================================================================
//function : AverageError
//purpose :
//=======================================================================
Standard_Real PLib_HermitJacobi::AverageError(const Standard_Integer Dimension,
Standard_Real& HermJacCoeff,
const Standard_Integer NewDegree) const
{
return myJacobi->AverageError(Dimension,HermJacCoeff,NewDegree);
}
//=======================================================================
//function : ToCoefficients
//purpose :
//=======================================================================
void PLib_HermitJacobi::ToCoefficients(const Standard_Integer Dimension,
const Standard_Integer Degree,
const TColStd_Array1OfReal& HermJacCoeff,
TColStd_Array1OfReal& Coefficients) const
{
Standard_Integer i,k,idim,i1,i2;
Standard_Real h1, h2;
Standard_Integer NivConstr = this->NivConstr(),
DegreeH = 2*NivConstr+1;
Standard_Integer ibegHJC = HermJacCoeff.Lower(), kdim;
TColStd_Array1OfReal AuxCoeff(0,(Degree+1)*Dimension-1);
AuxCoeff.Init(0.);
for (k=0; k<=DegreeH; k++) {
kdim = k*Dimension;
for (i=0; i<=NivConstr; i++) {
h1 = myH(i+1, k+1);
h2 = myH(i+NivConstr+2,k+1);
i1 = ibegHJC + i*Dimension;
i2 = ibegHJC + (i+NivConstr+1)*Dimension;
for (idim=0; idim<Dimension; idim++) {
AuxCoeff(idim + kdim) += HermJacCoeff(i1 + idim) * h1 +
HermJacCoeff(i2 + idim) * h2;
}
}
}
kdim = (Degree+1)*Dimension;
for (k=(DegreeH+1)*Dimension; k<kdim; k++) {
AuxCoeff(k) = HermJacCoeff(ibegHJC + k);
}
if(Degree > DegreeH)
myJacobi->ToCoefficients(Dimension,Degree,AuxCoeff,Coefficients);
else {
Standard_Integer ibegC = Coefficients.Lower();
kdim = (Degree+1)*Dimension;
for(k=0; k < kdim; k++)
Coefficients(ibegC+k) = AuxCoeff(k);
}
}
//=======================================================================
//function : D0123
//purpose : common part of D0,D1,D2,D3 (FORTRAN subroutine MPOBAS)
//=======================================================================
void PLib_HermitJacobi::D0123(const Standard_Integer NDeriv,
const Standard_Real U,
TColStd_Array1OfReal& BasisValue,
TColStd_Array1OfReal& BasisD1,
TColStd_Array1OfReal& BasisD2,
TColStd_Array1OfReal& BasisD3)
{
NCollection_LocalArray<Standard_Real> jac0 (4 * 20);
NCollection_LocalArray<Standard_Real> jac1 (4 * 20);
NCollection_LocalArray<Standard_Real> jac2 (4 * 20);
NCollection_LocalArray<Standard_Real> jac3 (4 * 20);
NCollection_LocalArray<Standard_Real> wvalues (4);
Standard_Integer i, j;
Standard_Integer NivConstr = this->NivConstr(),
WorkDegree = this->WorkDegree(),
DegreeH = 2*NivConstr+1;
Standard_Integer ibeg0 = BasisValue.Lower(),
ibeg1 = BasisD1.Lower(),
ibeg2 = BasisD2.Lower(),
ibeg3 = BasisD3.Lower();
Standard_Integer JacDegree = WorkDegree-DegreeH-1;
Standard_Real W0;
TColStd_Array1OfReal JacValue0(jac0[0], 0, Max(0,JacDegree));
TColStd_Array1OfReal WValues(wvalues[0],0,NDeriv);
WValues.Init(0.);
// Evaluation des polynomes d'hermite
math_Matrix HermitValues(0,DegreeH, 0, NDeriv, 0.);
if(NDeriv == 0)
for (i=0; i<=DegreeH; i++) {
PLib::NoDerivativeEvalPolynomial(U,DegreeH,1, DegreeH,
myH(i+1,1), HermitValues(i,0));
}
else
for (i=0; i<=DegreeH; i++) {
PLib::EvalPolynomial(U,NDeriv,DegreeH,1,
myH(i+1,1), HermitValues(i,0));
}
// Evaluation des polynomes de Jaccobi
if(JacDegree >= 0) {
switch (NDeriv) {
case 0 :
myJacobi->D0(U,JacValue0);
break;
case 1 :
{
TColStd_Array1OfReal JacValue1(jac1[0], 0, JacDegree);
myJacobi->D1(U,JacValue0,JacValue1);
break;
}
case 2 :
{
TColStd_Array1OfReal JacValue1(jac1[0], 0, JacDegree);
TColStd_Array1OfReal JacValue2(jac2[0], 0, JacDegree);
myJacobi->D2(U,JacValue0,JacValue1,JacValue2);
break;
}
case 3 :
{
TColStd_Array1OfReal JacValue1(jac1[0], 0, JacDegree);
TColStd_Array1OfReal JacValue2(jac2[0], 0, JacDegree);
TColStd_Array1OfReal JacValue3(jac3[0], 0, JacDegree);
myJacobi->D3(U,JacValue0,JacValue1,JacValue2,JacValue3);
}
}
// Evaluation de W(t)
if(NDeriv == 0)
PLib::NoDerivativeEvalPolynomial(U,DegreeH+1,1,DegreeH+1,myWCoeff(1), WValues(0));
else
PLib::EvalPolynomial(U,NDeriv,DegreeH+1,1,myWCoeff(1), WValues(0));
}
// Evaluation a l'ordre 0
for (i=0; i<=DegreeH; i++) {
BasisValue(ibeg0+i) = HermitValues(i,0);
}
W0 = WValues(0);
for (i=DegreeH+1, j=0; i<=WorkDegree; i++, j++) {
BasisValue(ibeg0+i) = W0 * jac0[j];
}
// Evaluation a l'ordre 1
if (NDeriv >= 1) {
Standard_Real W1=WValues(1);
for (i=0; i<=DegreeH; i++) {
BasisD1(ibeg1+i) = HermitValues(i,1);
}
for (i=DegreeH+1, j=0; i<=WorkDegree; i++, j++) {
BasisD1(ibeg1+i) = W0 * jac1[j] +
W1 * jac0[j];
}
// Evaluation a l'ordre 2
if (NDeriv >= 2) {
Standard_Real W2=WValues(2);
for (i=0; i<=DegreeH; i++) {
BasisD2(ibeg2+i) = HermitValues(i,2);
}
for (i=DegreeH+1, j=0; i<=WorkDegree; i++, j++) {
BasisD2(ibeg2+i) =
W0 * jac2[j] + 2 * W1 * jac1[j] + W2 * jac0[j];
}
// Evaluation a l'ordre 3
if (NDeriv == 3) {
Standard_Real W3 = WValues(3);
for (i=0; i<=DegreeH; i++) {
BasisD3(ibeg3+i) = HermitValues(i,3);
}
for (i=DegreeH+1, j=0; i<=WorkDegree; i++,j++) {
BasisD3(ibeg3+i) = W0*jac3[j] + W3*jac0[j]
+ 3*(W1*jac2[j] + W2*jac1[j]);
}
}
}
}
}
//=======================================================================
//function : D0
//purpose :
//=======================================================================
void PLib_HermitJacobi::D0(const Standard_Real U, TColStd_Array1OfReal& BasisValue)
{
D0123(0,U,BasisValue,BasisValue,BasisValue,BasisValue);
}
//=======================================================================
//function : D1
//purpose :
//=======================================================================
void PLib_HermitJacobi::D1(const Standard_Real U,
TColStd_Array1OfReal& BasisValue, TColStd_Array1OfReal& BasisD1)
{
D0123(1,U,BasisValue,BasisD1,BasisD1,BasisD1);
}
//=======================================================================
//function : D2
//purpose :
//=======================================================================
void PLib_HermitJacobi::D2(const Standard_Real U, TColStd_Array1OfReal& BasisValue,
TColStd_Array1OfReal& BasisD1,TColStd_Array1OfReal& BasisD2)
{
D0123(2,U,BasisValue,BasisD1,BasisD2,BasisD2);
}
//=======================================================================
//function : D3
//purpose :
//=======================================================================
void PLib_HermitJacobi::D3(const Standard_Real U,
TColStd_Array1OfReal& BasisValue, TColStd_Array1OfReal& BasisD1,
TColStd_Array1OfReal& BasisD2, TColStd_Array1OfReal& BasisD3)
{
D0123(3,U,BasisValue,BasisD1,BasisD2,BasisD3);
}