Skip to content

Commit

Permalink
Refactoring.
Browse files Browse the repository at this point in the history
  • Loading branch information
@CamelliaDPG
CamelliaDPG committed Apr 2, 2020
1 parent be512fd commit fb98ff9
Showing 1 changed file with 57 additions and 35 deletions.
92 changes: 57 additions & 35 deletions packages/intrepid2/src/Discretization/Basis/Intrepid2_IntegratedLegendreBasis_HGRAD_TET.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -318,23 +318,33 @@ namespace Intrepid2
auto jacobi_alpha = subview(jacobi_values1_at_point, alphaOrdinal, ALL);
Polynomials::integratedJacobiValues(jacobi_alpha, alpha, polyOrder_-3, lambda[3], jacobiScaling);
}
for (int i=2; i<polyOrder_-1; i++)
const int min_i = 2;
const int min_j = 1;
const int min_k = 1;
const int min_ij = min_i + min_j;
const int min_ijk = min_ij + min_k;
int localInteriorBasisOrdinal = 0;
for (int totalPolyOrder_ijk=min_ijk; totalPolyOrder_ijk <= polyOrder_; totalPolyOrder_ijk++)
{
for (int j=1; i+j<polyOrder_; j++)
int localFaceBasisOrdinal = 0;
for (int totalPolyOrder_ij=min_ij; totalPolyOrder_ij <= totalPolyOrder_ijk-min_j; totalPolyOrder_ij++)
{
const int faceBasisOrdinalInFace = localTetFaceBasisOrdinalESEAS(polyOrder_,i,j);
const int faceBasisOrdinal = numEdges*num1DEdgeFunctions + numVertices + faceBasisOrdinalInFace;
const auto & faceValue = output_(faceBasisOrdinal,pointOrdinal);
const int alphaOrdinal = (i+j)-3;

for (int k=1; i+j+k<=polyOrder_; k++)
for (int i=2; i <= totalPolyOrder_ij-min_j; i++)
{
const int fieldOrdinal = fieldOrdinalOffset + localTetInteriorBasisOrdinalESEAS(polyOrder_,i,j,k);
const int j = totalPolyOrder_ij - i;
const int k = totalPolyOrder_ijk - totalPolyOrder_ij;
const int faceBasisOrdinal = numEdges*num1DEdgeFunctions + numVertices + localFaceBasisOrdinal;
const auto & faceValue = output_(faceBasisOrdinal,pointOrdinal);
const int alphaOrdinal = (i+j)-3;
localFaceBasisOrdinal++;

const int fieldOrdinal = fieldOrdinalOffset + localInteriorBasisOrdinal;
const auto & jacobiValue = jacobi_values1_at_point(alphaOrdinal,k);
output_(fieldOrdinal,pointOrdinal) = faceValue * jacobiValue;
} // end k loop
} // end j loop
} // end i loop
localInteriorBasisOrdinal++;
} // end i loop
} // end totalPolyOrder_ij loop
} // end totalPolyOrder_ijk loop
fieldOrdinalOffset += numInteriorBasisFunctions;
} // end OPERATOR_VALUE
break;
Expand Down Expand Up @@ -561,39 +571,51 @@ namespace Intrepid2
Polynomials::shiftedScaledJacobiValues(P, alpha, polyOrder_-3, lambda[3], jacobiScaling);
}

for (int i=2; i<polyOrder_-1; i++)
const int min_i = 2;
const int min_j = 1;
const int min_k = 1;
const int min_ij = min_i + min_j;
const int min_ijk = min_ij + min_k;
int localInteriorBasisOrdinal = 0;
for (int totalPolyOrder_ijk=min_ijk; totalPolyOrder_ijk <= polyOrder_; totalPolyOrder_ijk++)
{
for (int j=1; i+j<polyOrder_; j++)
int localFaceBasisOrdinal = 0;
for (int totalPolyOrder_ij=min_ij; totalPolyOrder_ij <= totalPolyOrder_ijk-min_j; totalPolyOrder_ij++)
{
// interior functions use basis values belonging to the first face, 012
const int faceBasisOrdinalInFace = localTetFaceBasisOrdinalESEAS(polyOrder_,i,j);
const int faceBasisOrdinal = numEdges*num1DEdgeFunctions + numVertices + faceBasisOrdinalInFace;
const auto & faceValue_dx = output_(faceBasisOrdinal,pointOrdinal,0);
const auto & faceValue_dy = output_(faceBasisOrdinal,pointOrdinal,1);
const auto & faceValue_dz = output_(faceBasisOrdinal,pointOrdinal,2);

// determine faceValue (on face 0)
OutputScalar faceValue;
for (int i=2; i <= totalPolyOrder_ij-min_j; i++)
{
const auto & edgeValue = legendre_values1_at_point(i);
const int alphaOrdinal = i-2;
const auto & jacobiValue = jacobi_values3_at_point(alphaOrdinal,j);
faceValue = edgeValue * jacobiValue;
}
const int j = totalPolyOrder_ij - i;
const int k = totalPolyOrder_ijk - totalPolyOrder_ij;
// interior functions use basis values belonging to the first face, 012
const int faceBasisOrdinal = numEdges*num1DEdgeFunctions + numVertices + localFaceBasisOrdinal;

const auto & faceValue_dx = output_(faceBasisOrdinal,pointOrdinal,0);
const auto & faceValue_dy = output_(faceBasisOrdinal,pointOrdinal,1);
const auto & faceValue_dz = output_(faceBasisOrdinal,pointOrdinal,2);

// determine faceValue (on face 0)
OutputScalar faceValue;
{
const auto & edgeValue = legendre_values1_at_point(i);
const int alphaOrdinal = i-2;
const auto & jacobiValue = jacobi_values3_at_point(alphaOrdinal,j);
faceValue = edgeValue * jacobiValue;
}
localFaceBasisOrdinal++;

const int alphaOrdinal = (i+j)-3;
const int alphaOrdinal = (i+j)-3;

for (int k=1; i+j+k<=polyOrder_; k++)
{
const int fieldOrdinal = fieldOrdinalOffset + localTetInteriorBasisOrdinalESEAS(polyOrder_,i,j,k);
const int fieldOrdinal = fieldOrdinalOffset + localInteriorBasisOrdinal;
const auto & integratedJacobiValue = L_alpha(alphaOrdinal,k);
const auto & jacobiValue = P_alpha(alphaOrdinal,k-1);
output_(fieldOrdinal,pointOrdinal,0) = integratedJacobiValue * faceValue_dx + faceValue * jacobiValue * lambda_dx[3];
output_(fieldOrdinal,pointOrdinal,1) = integratedJacobiValue * faceValue_dy + faceValue * jacobiValue * lambda_dy[3];
output_(fieldOrdinal,pointOrdinal,2) = integratedJacobiValue * faceValue_dz + faceValue * jacobiValue * lambda_dz[3];
} // end k loop
} // end j loop
} // end i loop

localInteriorBasisOrdinal++;
} // end i loop
} // end totalPolyOrder_ij loop
} // end totalPolyOrder_ijk loop
fieldOrdinalOffset += numInteriorBasisFunctions;
}
break;
Expand Down

0 comments on commit fb98ff9

Please sign in to comment.