本文整理汇总了C++中ArrayScalar类的典型用法代码示例。如果您正苦于以下问题:C++ ArrayScalar类的具体用法?C++ ArrayScalar怎么用?C++ ArrayScalar使用的例子?那么, 这里精选的类代码示例或许可以为您提供帮助。
在下文中一共展示了ArrayScalar类的15个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的C++代码示例。
示例1: dZ
void TabulatorTri<Scalar,ArrayScalar,derivOrder>::tabulate( ArrayScalar &outputValues ,
const int deg ,
const ArrayScalar &z )
{
const int np = z.dimension(0);
const int card = outputValues.dimension(0);
FieldContainer<Sacado::Fad::DFad<Scalar> > dZ( z.dimension(0) , z.dimension(1) );
for (int i=0;i<np;i++) {
for (int j=0;j<2;j++) {
dZ(i,j) = Sacado::Fad::DFad<Scalar>( z(i,j) );
dZ(i,j).diff(j,2);
}
}
FieldContainer<Sacado::Fad::DFad<Scalar> > dResult(card,np,derivOrder+1);
TabulatorTri<Sacado::Fad::DFad<Scalar>,FieldContainer<Sacado::Fad::DFad<Scalar> >,derivOrder-1>::tabulate(dResult ,
deg ,
dZ );
for (int i=0;i<card;i++) {
for (int j=0;j<np;j++) {
outputValues(i,j,0) = dResult(i,j,0).dx(0);
for (unsigned k=0;k<derivOrder;k++) {
outputValues(i,j,k+1) = dResult(i,j,k).dx(1);
}
}
}
return;
}
示例2: cell
void Basis_HGRAD_QUAD_C2_FEM<Scalar, ArrayScalar>::getDofCoords(ArrayScalar & DofCoords) const {
#ifdef HAVE_INTREPID2_DEBUG
// Verify rank of output array.
TEUCHOS_TEST_FOR_EXCEPTION( !(DofCoords.rank() == 2), std::invalid_argument,
">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_C2_FEM::getDofCoords) rank = 2 required for DofCoords array");
// Verify 0th dimension of output array.
TEUCHOS_TEST_FOR_EXCEPTION( !( static_cast<index_type>(DofCoords.dimension(0)) == static_cast<index_type>(this -> basisCardinality_) ), std::invalid_argument,
">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_C2_FEM::getDofCoords) mismatch in number of DoF and 0th dimension of DofCoords array");
// Verify 1st dimension of output array.
TEUCHOS_TEST_FOR_EXCEPTION( !( static_cast<index_type>(DofCoords.dimension(1)) == static_cast<index_type>(this -> basisCellTopology_.getDimension()) ), std::invalid_argument,
">>> ERROR: (Intrepid2::Basis_HGRAD_QUAD_C2_FEM::getDofCoords) incorrect reference cell (1st) dimension in DofCoords array");
#endif
DofCoords(0,0) = -1.0; DofCoords(0,1) = -1.0;
DofCoords(1,0) = 1.0; DofCoords(1,1) = -1.0;
DofCoords(2,0) = 1.0; DofCoords(2,1) = 1.0;
DofCoords(3,0) = -1.0; DofCoords(3,1) = 1.0;
DofCoords(4,0) = 0.0; DofCoords(4,1) = -1.0;
DofCoords(5,0) = 1.0; DofCoords(5,1) = 0.0;
DofCoords(6,0) = 0.0; DofCoords(6,1) = 1.0;
DofCoords(7,0) = -1.0; DofCoords(7,1) = 0.0;
DofCoords(8,0) = 0.0; DofCoords(8,1) = 0.0;
}
示例3: ptsx_
Basis_HGRAD_QUAD_Cn_FEM<Scalar,ArrayScalar>::Basis_HGRAD_QUAD_Cn_FEM( const int orderx , const int ordery,
const ArrayScalar &pts_x ,
const ArrayScalar &pts_y ):
ptsx_( pts_x.dimension(0) , 1 ) ,
ptsy_( pts_y.dimension(0) , 1 )
{
Array<Array<RCP<Basis<Scalar,ArrayScalar> > > > bases(1);
bases[0].resize(2);
bases[0][0] = Teuchos::rcp( new Basis_HGRAD_LINE_Cn_FEM<Scalar,ArrayScalar>( orderx , pts_x ) );
bases[0][1] = Teuchos::rcp( new Basis_HGRAD_LINE_Cn_FEM<Scalar,ArrayScalar>( ordery , pts_y ) );
this->setBases( bases );
this->basisCardinality_ = (orderx+1)*(ordery+1);
if (orderx > ordery) {
this->basisDegree_ = orderx;
}
else {
this->basisDegree_ = ordery;
}
this -> basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Quadrilateral<4> >() );
this -> basisType_ = BASIS_FEM_FIAT;
this -> basisCoordinates_ = COORDINATES_CARTESIAN;
this -> basisTagsAreSet_ = false;
for (int i=0;i<pts_x.dimension(0);i++)
{
ptsx_(i,0) = pts_x(i,0);
}
for (int i=0;i<pts_y.dimension(0);i++)
{
ptsy_(i,0) = pts_y(i,0);
}
}
示例4: dZ
void TabulatorTet<Scalar,ArrayScalar,1>::tabulate( ArrayScalar &outputValues ,
const int deg ,
const ArrayScalar &z )
{
const int np = z.dimension(0);
const int card = outputValues.dimension(0);
FieldContainer<Sacado::Fad::DFad<Scalar> > dZ( z.dimension(0) , z.dimension(1) );
for (int i=0;i<np;i++) {
for (int j=0;j<3;j++) {
dZ(i,j) = Sacado::Fad::DFad<Scalar>( z(i,j) );
dZ(i,j).diff(j,3);
}
}
FieldContainer<Sacado::Fad::DFad<Scalar> > dResult(card,np);
TabulatorTet<Sacado::Fad::DFad<Scalar>,FieldContainer<Sacado::Fad::DFad<Scalar> >,0>::tabulate( dResult ,
deg ,
dZ );
for (int i=0;i<card;i++) {
for (int j=0;j<np;j++) {
for (int k=0;k<3;k++) {
outputValues(i,j,k) = dResult(i,j).dx(k);
}
}
}
return;
}
示例5: cell
void Basis_HVOL_TRI_C0_FEM<Scalar, ArrayScalar>::getDofCoords(ArrayScalar & DofCoords) const {
#ifdef HAVE_INTREPID_DEBUG
// Verify rank of output array.
TEUCHOS_TEST_FOR_EXCEPTION( !(DofCoords.rank() == 2), std::invalid_argument,
">>> ERROR: (Intrepid::Basis_HVOL_TRI_C0_FEM::getDofCoords) rank = 2 required for DofCoords array");
// Verify 0th dimension of output array.
TEUCHOS_TEST_FOR_EXCEPTION( !( DofCoords.dimension(0) == this -> basisCardinality_ ), std::invalid_argument,
">>> ERROR: (Intrepid::Basis_HVOL_TRI_C0_FEM::getDofCoords) mismatch in number of DoF and 0th dimension of DofCoords array");
// Verify 1st dimension of output array.
TEUCHOS_TEST_FOR_EXCEPTION( !( DofCoords.dimension(1) == (int)(this -> basisCellTopology_.getDimension()) ), std::invalid_argument,
">>> ERROR: (Intrepid::Basis_HVOL_TRI_C0_FEM::getDofCoords) incorrect reference cell (1st) dimension in DofCoords array");
DofCoords(0,0) = 0.0; DofCoords(0,1) = 0.0;
#endif
}
示例6: pts
Basis_HGRAD_LINE_Hermite_FEM<Scalar,ArrayScalar>::Basis_HGRAD_LINE_Hermite_FEM( const ArrayScalar &pts) :
latticePts_( pts.dimension(0), 1 ) {
int n = pts.dimension(0);
this->basisCardinality_ = 2*n;
this->basisDegree_ = 2*n-1;
this->basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Line<2> >() );
this->basisType_ = BASIS_FEM_DEFAULT;
this->basisCoordinates_ = COORDINATES_CARTESIAN;
this->basisTagsAreSet_ = false;
for( int i=0; i<n-1; ++i ) {
TEUCHOS_TEST_FOR_EXCEPTION( pts(i,0) >= pts(i+1,0), std::runtime_error ,
"Intrepid::Basis_HGRAD_LINE_Hermite_FEM Illegal points given to constructor" );
}
// copy points int latticePts, correcting endpoints if needed
if (std::abs(pts(0,0)+1.0) < INTREPID_TOL) {
latticePts_(0,0) = -1.0;
}
else {
latticePts_(0,0) = pts(0,0);
}
for (int i=1;i<n-1;i++) {
latticePts_(i,0) = pts(i,0);
}
if (std::abs(pts(n-1,0)-1.0) < INTREPID_TOL) {
latticePts_(n-1,0) = 1.0;
}
else {
latticePts_(n-1,0) = pts(n-1,0);
}
setupVandermonde();
} // Constructor with points given
示例7: P
void Basis_HGRAD_LINE_Hermite_FEM<Scalar,ArrayScalar>::recurrence( ArrayScalar &P,
ArrayScalar &Px,
const Scalar x ) const {
int n = P.dimension(0);
Scalar q = x*x-1.0;
P (0) = 1.0;
Px(0) = 0.0;
// Loop over basis indices
for( int j=0; j<n-1; ++j ) {
P (j+1) = x*P(j) + q*Px(j)/(j+1); // Compute \f$P_{j+1}(x_i)\f$
Px(j+1) = (j+1)*P(j) + x*Px(j); // Compute \f$P'_{j+1}(x_i)\f$
}
} // recurrence()
示例8: getBaseCellTopology
void Basis_Constant<Scalar, ArrayScalar>::getValues(ArrayScalar & outputValues,
const ArrayScalar & inputPoints,
const Intrepid2::EOperator operatorType) const {
// Verify arguments
#ifdef HAVE_INTREPID_DEBUG
Intrepid2::getValues_HGRAD_Args<Scalar, ArrayScalar>(outputValues,
inputPoints,
operatorType,
this -> getBaseCellTopology(),
this -> getCardinality() );
#endif
// Number of evaluation points = dim 0 of inputPoints
int dim0 = inputPoints.dimension(0);
switch (operatorType) {
case Intrepid2::OPERATOR_VALUE:
for (int i0 = 0; i0 < dim0; i0++) {
// outputValues is a rank-2 array with dimensions (basisCardinality_, dim0)
outputValues(0, i0) = 1.0;
}
break;
case Intrepid2::OPERATOR_GRAD:
case Intrepid2::OPERATOR_D1:
case Intrepid2::OPERATOR_CURL:
case Intrepid2::OPERATOR_DIV:
case Intrepid2::OPERATOR_D2:
case Intrepid2::OPERATOR_D3:
case Intrepid2::OPERATOR_D4:
case Intrepid2::OPERATOR_D5:
case Intrepid2::OPERATOR_D6:
case Intrepid2::OPERATOR_D7:
case Intrepid2::OPERATOR_D8:
case Intrepid2::OPERATOR_D9:
case Intrepid2::OPERATOR_D10:
default:
TEUCHOS_TEST_FOR_EXCEPTION( !( Intrepid2::isValidOperator(operatorType) ), std::invalid_argument,
">>> ERROR (Basis_Constant): Invalid operator type");
}
}
示例9: ptsx_
Basis_HGRAD_HEX_Cn_FEM<Scalar,ArrayScalar>::Basis_HGRAD_HEX_Cn_FEM( const int orderx ,
const int ordery ,
const int orderz ,
const ArrayScalar &pts_x ,
const ArrayScalar &pts_y ,
const ArrayScalar &pts_z ):
ptsx_( pts_x.dimension(0) , 1 ),
ptsy_( pts_y.dimension(0) , 1 ),
ptsz_( pts_z.dimension(0) , 1 )
{
for (int i=0;i<pts_x.dimension(0);i++)
{
ptsx_(i,0) = pts_x(i,0);
}
for (int i=0;i<pts_y.dimension(0);i++)
{
ptsy_(i,0) = pts_y(i,0);
}
for (int i=0;i<pts_z.dimension(0);i++)
{
ptsz_(i,0) = pts_z(i,0);
}
Array<Array<RCP<Basis<Scalar,ArrayScalar> > > > bases(1);
bases[0].resize(3);
bases[0][0] = Teuchos::rcp( new Basis_HGRAD_LINE_Cn_FEM< Scalar , ArrayScalar >( orderx , pts_x ) );
bases[0][1] = Teuchos::rcp( new Basis_HGRAD_LINE_Cn_FEM< Scalar , ArrayScalar >( ordery , pts_y ) );
bases[0][2] = Teuchos::rcp( new Basis_HGRAD_LINE_Cn_FEM< Scalar , ArrayScalar >( orderz , pts_z ) );
this->setBases( bases );
this->basisCardinality_ = (orderx+1)*(ordery+1)*(orderz+1);
if (orderx >= ordery && orderx >= orderz ) {
this->basisDegree_ = orderx;
}
else if (ordery >= orderx && ordery >= orderz) {
this->basisDegree_ = ordery;
}
else {
this->basisDegree_ = orderz;
}
this -> basisCellTopology_ = shards::CellTopology(shards::getCellTopologyData<shards::Hexahedron<8> >() );
this -> basisType_ = BASIS_FEM_FIAT;
this -> basisCoordinates_ = COORDINATES_CARTESIAN;
this -> basisTagsAreSet_ = false;
initializeTags();
this->basisTagsAreSet_ = true;
}
示例10: getBaseCellTopology
void Basis_HCURL_TRI_I1_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar & outputValues,
const ArrayScalar & inputPoints,
const EOperator operatorType) const {
// Verify arguments
#ifdef HAVE_INTREPID_DEBUG
Intrepid::getValues_HCURL_Args<Scalar, ArrayScalar>(outputValues,
inputPoints,
operatorType,
this -> getBaseCellTopology(),
this -> getCardinality() );
#endif
// Number of evaluation points = dim 0 of inputPoints
int dim0 = inputPoints.dimension(0);
// Temporaries: (x,y) coordinates of the evaluation point
Scalar x = 0.0;
Scalar y = 0.0;
switch (operatorType) {
case OPERATOR_VALUE:
for (int i0 = 0; i0 < dim0; i0++) {
x = inputPoints(i0, 0);
y = inputPoints(i0, 1);
// outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
outputValues(0, i0, 0) = 1.0 - y;
outputValues(0, i0, 1) = x;
outputValues(1, i0, 0) = -y;
outputValues(1, i0, 1) = x;
outputValues(2, i0, 0) = -y;
outputValues(2, i0, 1) = -1.0 + x;
}
break;
case OPERATOR_CURL:
for (int i0 = 0; i0 < dim0; i0++) {
x = inputPoints(i0, 0);
y = inputPoints(i0, 1);
// outputValues is a rank-2 array with dimensions (basisCardinality_, dim0)
outputValues(0, i0) = 2.0;
outputValues(1, i0) = 2.0;
outputValues(2, i0) = 2.0;
}
break;
case OPERATOR_DIV:
TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
">>> ERROR (Basis_HCURL_TRI_I1_FEM): DIV is invalid operator for HCURL Basis Functions");
break;
case OPERATOR_GRAD:
TEST_FOR_EXCEPTION( (operatorType == OPERATOR_GRAD), std::invalid_argument,
">>> ERROR (Basis_HCURL_TRI_I1_FEM): GRAD is invalid operator for HCURL Basis Functions");
break;
case OPERATOR_D1:
case OPERATOR_D2:
case OPERATOR_D3:
case OPERATOR_D4:
case OPERATOR_D5:
case OPERATOR_D6:
case OPERATOR_D7:
case OPERATOR_D8:
case OPERATOR_D9:
case OPERATOR_D10:
TEST_FOR_EXCEPTION( ( (operatorType == OPERATOR_D1) ||
(operatorType == OPERATOR_D2) ||
(operatorType == OPERATOR_D3) ||
(operatorType == OPERATOR_D4) ||
(operatorType == OPERATOR_D5) ||
(operatorType == OPERATOR_D6) ||
(operatorType == OPERATOR_D7) ||
(operatorType == OPERATOR_D8) ||
(operatorType == OPERATOR_D9) ||
(operatorType == OPERATOR_D10) ),
std::invalid_argument,
">>> ERROR (Basis_HCURL_TRI_I1_FEM): Invalid operator type");
break;
default:
TEST_FOR_EXCEPTION( ( (operatorType != OPERATOR_VALUE) &&
(operatorType != OPERATOR_GRAD) &&
(operatorType != OPERATOR_CURL) &&
(operatorType != OPERATOR_DIV) &&
(operatorType != OPERATOR_D1) &&
(operatorType != OPERATOR_D2) &&
(operatorType != OPERATOR_D3) &&
(operatorType != OPERATOR_D4) &&
(operatorType != OPERATOR_D5) &&
(operatorType != OPERATOR_D6) &&
(operatorType != OPERATOR_D7) &&
(operatorType != OPERATOR_D8) &&
(operatorType != OPERATOR_D9) &&
(operatorType != OPERATOR_D10) ),
std::invalid_argument,
//.........这里部分代码省略.........
示例11: getBaseCellTopology
void Basis_HGRAD_WEDGE_C1_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar & outputValues,
const ArrayScalar & inputPoints,
const EOperator operatorType) const {
// Verify arguments
#ifdef HAVE_INTREPID_DEBUG
Intrepid2::getValues_HGRAD_Args<Scalar, ArrayScalar>(outputValues,
inputPoints,
operatorType,
this -> getBaseCellTopology(),
this -> getCardinality() );
#endif
// Number of evaluation points = dim 0 of inputPoints
int dim0 = inputPoints.dimension(0);
// Temporaries: (x,y,z) coordinates of the evaluation point
Scalar x = 0.0;
Scalar y = 0.0;
Scalar z = 0.0;
switch (operatorType) {
case OPERATOR_VALUE:
for (int i0 = 0; i0 < dim0; i0++) {
x = inputPoints(i0, 0);
y = inputPoints(i0, 1);
z = inputPoints(i0, 2);
// outputValues is a rank-2 array with dimensions (basisCardinality_, dim0)
outputValues(0, i0) = (1.0 - x - y)*(1.0 - z)/2.0;
outputValues(1, i0) = x*(1.0 - z)/2.0;
outputValues(2, i0) = y*(1.0 - z)/2.0;
outputValues(3, i0) = (1.0 - x - y)*(1.0 + z)/2.0;
outputValues(4, i0) = x*(1.0 + z)/2.0;
outputValues(5, i0) = y*(1.0 + z)/2.0;
}
break;
case OPERATOR_GRAD:
case OPERATOR_D1:
for (int i0 = 0; i0 < dim0; i0++) {
x = inputPoints(i0,0);
y = inputPoints(i0,1);
z = inputPoints(i0, 2);
// outputValues is a rank-3 array with dimensions (basisCardinality_, dim0, spaceDim)
outputValues(0, i0, 0) = -(1.0 - z)/2.0;
outputValues(0, i0, 1) = -(1.0 - z)/2.0;
outputValues(0, i0, 2) = -(1.0 - x - y)/2.0;
outputValues(1, i0, 0) = (1.0 - z)/2.0;
outputValues(1, i0, 1) = 0.0;
outputValues(1, i0, 2) = -x/2.0;
outputValues(2, i0, 0) = 0.0;
outputValues(2, i0, 1) = (1.0 - z)/2.0;
outputValues(2, i0, 2) = -y/2.0;
outputValues(3, i0, 0) = -(1.0 + z)/2.0;
outputValues(3, i0, 1) = -(1.0 + z)/2.0;
outputValues(3, i0, 2) = (1.0 - x - y)/2.0;
outputValues(4, i0, 0) = (1.0 + z)/2.0;
outputValues(4, i0, 1) = 0.0;
outputValues(4, i0, 2) = x/2.0;
outputValues(5, i0, 0) = 0.0;
outputValues(5, i0, 1) = (1.0 + z)/2.0;
outputValues(5, i0, 2) = y/2.0;
}
break;
case OPERATOR_CURL:
TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_CURL), std::invalid_argument,
">>> ERROR (Basis_HGRAD_WEDGE_C1_FEM): CURL is invalid operator for rank-0 (scalar) functions in 3D");
break;
case OPERATOR_DIV:
TEUCHOS_TEST_FOR_EXCEPTION( (operatorType == OPERATOR_DIV), std::invalid_argument,
">>> ERROR (Basis_HGRAD_WEDGE_C1_FEM): DIV is invalid operator for rank-0 (scalar) functions in 3D");
break;
case OPERATOR_D2:
for (int i0 = 0; i0 < dim0; i0++) {
outputValues(0, i0, 0) = 0.0; outputValues(3, i0, 0) = 0.0;
outputValues(0, i0, 1) = 0.0; outputValues(3, i0, 1) = 0.0;
outputValues(0, i0, 2) = 0.5; outputValues(3, i0, 2) =-0.5;
outputValues(0, i0, 3) = 0.0; outputValues(3, i0, 3) = 0.0;
outputValues(0, i0, 4) = 0.5; outputValues(3, i0, 4) =-0.5;
outputValues(0, i0, 5) = 0.0; outputValues(3, i0, 5) = 0.0;
outputValues(1, i0, 0) = 0.0; outputValues(4, i0, 0) = 0.0;
outputValues(1, i0, 1) = 0.0; outputValues(4, i0, 1) = 0.0;
outputValues(1, i0, 2) =-0.5; outputValues(4, i0, 2) = 0.5;
outputValues(1, i0, 3) = 0.0; outputValues(4, i0, 3) = 0.0;
outputValues(1, i0, 4) = 0.0; outputValues(4, i0, 4) = 0.0;
outputValues(1, i0, 5) = 0.0; outputValues(4, i0, 5) = 0.0;
//.........这里部分代码省略.........
示例12: getBaseCellTopology
void Basis_HCURL_QUAD_In_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar & outputValues,
const ArrayScalar & inputPoints,
const EOperator operatorType) const {
// Verify arguments
#ifdef HAVE_INTREPID_DEBUG
Intrepid2::getValues_HCURL_Args<Scalar, ArrayScalar>(outputValues,
inputPoints,
operatorType,
this -> getBaseCellTopology(),
this -> getCardinality() );
#endif
// Number of evaluation points = dim 0 of inputPoints
int dim0 = inputPoints.dimension(0);
// separate out points
FieldContainer<Scalar> xPoints(dim0,1);
FieldContainer<Scalar> yPoints(dim0,1);
for (int i=0;i<dim0;i++) {
xPoints(i,0) = inputPoints(i,0);
yPoints(i,0) = inputPoints(i,1);
}
switch (operatorType) {
case OPERATOR_VALUE:
{
FieldContainer<Scalar> closedBasisValsXPts( closedBasis_.getCardinality() , dim0 );
FieldContainer<Scalar> closedBasisValsYPts( closedBasis_.getCardinality() , dim0 );
FieldContainer<Scalar> openBasisValsXPts( openBasis_.getCardinality() , dim0 );
FieldContainer<Scalar> openBasisValsYPts( openBasis_.getCardinality() , dim0 );
closedBasis_.getValues( closedBasisValsXPts , xPoints , OPERATOR_VALUE );
closedBasis_.getValues( closedBasisValsYPts , yPoints , OPERATOR_VALUE );
openBasis_.getValues( openBasisValsXPts , xPoints , OPERATOR_VALUE );
openBasis_.getValues( openBasisValsYPts , yPoints , OPERATOR_VALUE );
int bfcur = 0;
// x component bfs are (open(x) closed(y),0)
for (int j=0;j<closedBasis_.getCardinality();j++)
{
for (int i=0;i<openBasis_.getCardinality();i++)
{
for (int l=0;l<dim0;l++)
{
outputValues(bfcur,l,0) = closedBasisValsYPts(j,l) * openBasisValsXPts(i,l);
outputValues(bfcur,l,1) = 0.0;
}
bfcur++;
}
}
// y component bfs are (0,closed(x) open(y))
for (int j=0;j<openBasis_.getCardinality();j++)
{
for (int i=0;i<closedBasis_.getCardinality();i++)
{
for (int l=0;l<dim0;l++)
{
outputValues(bfcur,l,0) = 0.0;
outputValues(bfcur,l,1) = openBasisValsYPts(j,l) * closedBasisValsXPts(i,l);
}
bfcur++;
}
}
}
break;
case OPERATOR_CURL:
{
FieldContainer<Scalar> closedBasisDerivsXPts( closedBasis_.getCardinality() , dim0 , 1 );
FieldContainer<Scalar> closedBasisDerivsYPts( closedBasis_.getCardinality() , dim0 , 1 );
FieldContainer<Scalar> openBasisValsXPts( openBasis_.getCardinality() , dim0 );
FieldContainer<Scalar> openBasisValsYPts( openBasis_.getCardinality() , dim0 );
closedBasis_.getValues( closedBasisDerivsXPts , xPoints , OPERATOR_D1 );
closedBasis_.getValues( closedBasisDerivsYPts , yPoints , OPERATOR_D1 );
openBasis_.getValues( openBasisValsXPts , xPoints , OPERATOR_VALUE );
openBasis_.getValues( openBasisValsYPts , yPoints , OPERATOR_VALUE );
int bfcur = 0;
// x component basis functions first
for (int j=0;j<closedBasis_.getCardinality();j++)
{
for (int i=0;i<openBasis_.getCardinality();i++)
{
for (int l=0;l<dim0;l++) {
outputValues(bfcur,l) = -closedBasisDerivsYPts(j,l,0) * openBasisValsXPts(i,l);
}
bfcur++;
}
}
// now y component basis functions
for (int j=0;j<openBasis_.getCardinality();j++)
{
for (int i=0;i<closedBasis_.getCardinality();i++)
{
for (int l=0;l<dim0;l++)
//.........这里部分代码省略.........
示例13: getBaseCellTopology
void Basis_HGRAD_LINE_Cn_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar & outputValues,
const ArrayScalar & inputPoints,
const EOperator operatorType) const {
// Verify arguments
#ifdef HAVE_INTREPID2_DEBUG
Intrepid2::getValues_HGRAD_Args<Scalar, ArrayScalar>(outputValues,
inputPoints,
operatorType,
this -> getBaseCellTopology(),
this -> getCardinality() );
#endif
const int numPts = inputPoints.dimension(0);
const int numBf = this->getCardinality();
try {
switch (operatorType) {
case OPERATOR_VALUE:
{
ArrayScalar phisCur( numBf , numPts );
Phis_.getValues( phisCur , inputPoints , operatorType );
for (int i=0;i<outputValues.dimension(0);i++) {
for (int j=0;j<outputValues.dimension(1);j++) {
outputValues(i,j) = 0.0;
for (int k=0;k<this->getCardinality();k++) {
outputValues(i,j) += this->Vinv_(k,i) * phisCur(k,j);
}
}
}
}
break;
case OPERATOR_GRAD:
case OPERATOR_D1:
case OPERATOR_D2:
case OPERATOR_D3:
case OPERATOR_D4:
case OPERATOR_D5:
case OPERATOR_D6:
case OPERATOR_D7:
case OPERATOR_D8:
case OPERATOR_D9:
case OPERATOR_D10:
{
const int dkcard =
(operatorType == OPERATOR_GRAD)? getDkCardinality(OPERATOR_D1,1): getDkCardinality(operatorType,1);
ArrayScalar phisCur( numBf , numPts , dkcard );
Phis_.getValues( phisCur , inputPoints , operatorType );
for (int i=0;i<outputValues.dimension(0);i++) {
for (int j=0;j<outputValues.dimension(1);j++) {
for (int k=0;k<outputValues.dimension(2);k++) {
outputValues(i,j,k) = 0.0;
for (int l=0;l<this->getCardinality();l++) {
outputValues(i,j,k) += this->Vinv_(l,i) * phisCur(l,j,k);
}
}
}
}
}
break;
default:
TEUCHOS_TEST_FOR_EXCEPTION( true , std::invalid_argument,
">>> ERROR (Basis_HGRAD_LINE_Cn_FEM): Operator type not implemented" );
break;
}
}
catch (std::invalid_argument &exception){
TEUCHOS_TEST_FOR_EXCEPTION( true , std::invalid_argument,
">>> ERROR (Basis_HGRAD_LINE_Cn_FEM): Operator failed");
}
}
示例14: outputValues
void TabulatorTri<Scalar,ArrayScalar,0>::tabulate(ArrayScalar &outputValues ,
const int deg ,
const ArrayScalar &z )
{
const int np = z.dimension( 0 );
// each point needs to be transformed from Pavel's element
// z(i,0) --> (2.0 * z(i,0) - 1.0)
// z(i,1) --> (2.0 * z(i,1) - 1.0)
// set up constant term
int idx_cur = TabulatorTri<Scalar,ArrayScalar,0>::idx(0,0);
int idx_curp1,idx_curm1;
// set D^{0,0} = 1.0
for (int i=0;i<np;i++) {
outputValues(idx_cur,i) = Scalar( 1.0 ) + z(i,0) - z(i,0) + z(i,1) - z(i,1);
}
if (deg > 0) {
Teuchos::Array<Scalar> f1(np),f2(np),f3(np);
for (int i=0;i<np;i++) {
f1[i] = 0.5 * (1.0+2.0*(2.0*z(i,0)-1.0)+(2.0*z(i,1)-1.0));
f2[i] = 0.5 * (1.0-(2.0*z(i,1)-1.0));
f3[i] = f2[i] * f2[i];
}
// set D^{1,0} = f1
idx_cur = TabulatorTri<Scalar,ArrayScalar,0>::idx(1,0);
for (int i=0;i<np;i++) {
outputValues(idx_cur,i) = f1[i];
}
// recurrence in p
for (int p=1;p<deg;p++) {
idx_cur = TabulatorTri<Scalar,ArrayScalar,0>::idx(p,0);
idx_curp1 = TabulatorTri<Scalar,ArrayScalar,0>::idx(p+1,0);
idx_curm1 = TabulatorTri<Scalar,ArrayScalar,0>::idx(p-1,0);
Scalar a = (2.0*p+1.0)/(1.0+p);
Scalar b = p / (p+1.0);
for (int i=0;i<np;i++) {
outputValues(idx_curp1,i) = a * f1[i] * outputValues(idx_cur,i)
- b * f3[i] * outputValues(idx_curm1,i);
}
}
// D^{p,1}
for (int p=0;p<deg;p++) {
int idxp0 = TabulatorTri<Scalar,ArrayScalar,0>::idx(p,0);
int idxp1 = TabulatorTri<Scalar,ArrayScalar,0>::idx(p,1);
for (int i=0;i<np;i++) {
outputValues(idxp1,i) = outputValues(idxp0,i)
*0.5*(1.0+2.0*p+(3.0+2.0*p)*(2.0*z(i,1)-1.0));
}
}
// recurrence in q
for (int p=0;p<deg-1;p++) {
for (int q=1;q<deg-p;q++) {
int idxpqp1=TabulatorTri<Scalar,ArrayScalar,0>::idx(p,q+1);
int idxpq=TabulatorTri<Scalar,ArrayScalar,0>::idx(p,q);
int idxpqm1=TabulatorTri<Scalar,ArrayScalar,0>::idx(p,q-1);
Scalar a,b,c;
TabulatorTri<Scalar,ArrayScalar,0>::jrc((Scalar)(2*p+1),(Scalar)0,q,a,b,c);
for (int i=0;i<np;i++) {
outputValues(idxpqp1,i)
= (a*(2.0*z(i,1)-1.0)+b)*outputValues(idxpq,i)
- c*outputValues(idxpqm1,i);
}
}
}
}
// orthogonalize
for (int p=0;p<=deg;p++) {
for (int q=0;q<=deg-p;q++) {
for (int i=0;i<np;i++) {
outputValues(TabulatorTri<Scalar,ArrayScalar,0>::idx(p,q),i) *= sqrt( (p+0.5)*(p+q+1.0));
}
}
}
return;
}
示例15: getBaseCellTopology
void Basis_HGRAD_LINE_Hermite_FEM<Scalar, ArrayScalar>::getValues(ArrayScalar & outputValues,
const ArrayScalar & inputPoints,
const EOperator operatorType) const {
// Verify arguments
#ifdef HAVE_INTREPID_DEBUG
Intrepid::getValues_HGRAD_Args<Scalar, ArrayScalar>(outputValues,
inputPoints,
operatorType,
this -> getBaseCellTopology(),
this -> getCardinality() );
#endif
// Number of evaluation points = dim 0 of inputPoints
int nPts = inputPoints.dimension(0);
int nBf = this->getCardinality();
int n = nBf/2;
// Legendre polynomials and their derivatives evaluated on inputPoints
SerialDenseMatrix legendre(nBf,nPts);
// Hermite interpolants evaluated on inputPoints
SerialDenseMatrix hermite(nBf,nPts);
ArrayScalar P (nBf);
ArrayScalar Px(nBf);
int derivative_order;
int derivative_case = static_cast<int>(operatorType);
if( derivative_case == 0 ) {
derivative_order = 0;
}
else if( derivative_case > 0 && derivative_case < 5 ) {
derivative_order = 1;
}
else {
derivative_order = derivative_case - 3;
}
try {
// GRAD,CURL,DIV, and D1 are all the first derivative
switch (operatorType) {
case OPERATOR_VALUE:
{
for( int i=0; i<nPts; ++i ) {
recurrence( P, Px, inputPoints(i,0) );
for( int j=0; j<nBf; ++j ) {
legendre(j,i) = P(j);
}
}
break;
}
case OPERATOR_GRAD:
case OPERATOR_DIV:
case OPERATOR_CURL:
case OPERATOR_D1:
{
for( int i=0; i<nPts; ++i ) {
recurrence( P, Px, inputPoints(i,0) );
for( int j=0; j<nBf; ++j ) {
legendre(j,i) = Px(j);
}
}
break;
}
case OPERATOR_D2:
case OPERATOR_D3:
case OPERATOR_D4:
case OPERATOR_D5:
case OPERATOR_D6:
case OPERATOR_D7:
case OPERATOR_D8:
case OPERATOR_D9:
case OPERATOR_D10:
{
for( int i=0; i<nPts; ++i ) {
legendre_d( P, Px, derivative_order, inputPoints(i,0));
for( int j=0; j<nBf; ++j ) {
legendre(j,i) = Px(j);
}
}
break;
}
default:
TEUCHOS_TEST_FOR_EXCEPTION( !( Intrepid::isValidOperator(operatorType) ), std::invalid_argument,
">>> ERROR (Basis_HGRAD_LINE_Hermite_FEM): Invalid operator type");
} // switch(operatorType)
}
catch (std::invalid_argument &exception){
TEUCHOS_TEST_FOR_EXCEPTION( true , std::invalid_argument,
">>> ERROR (Basis_HGRAD_LINE_Hermite_FEM): Operator failed");
}
if( !isFactored_ ) {
solver_.factorWithEquilibration(true);
solver_.factor();
}
//.........这里部分代码省略.........