1 const char help[] = "Tests PetscDTBaryToIndex(), PetscDTIndexToBary(), PetscDTIndexToGradedOrder() and PetscDTGradedOrderToIndex()"; 2 3 #include <petsc/private/petscimpl.h> 4 #include <petsc/private/dtimpl.h> 5 #include <petsc/private/petscfeimpl.h> 6 7 int main(int argc, char **argv) { 8 PetscInt d, n, maxdim = 4; 9 PetscInt *btupprev, *btup; 10 PetscInt *gtup; 11 12 PetscFunctionBeginUser; 13 PetscCall(PetscInitialize(&argc, &argv, NULL, help)); 14 PetscCall(PetscMalloc3(maxdim + 1, &btup, maxdim + 1, &btupprev, maxdim, >up)); 15 for (d = 0; d <= maxdim; d++) { 16 for (n = 0; n <= d + 2; n++) { 17 PetscInt j, k, Nk, kchk; 18 19 PetscCall(PetscDTBinomialInt(d + n, d, &Nk)); 20 for (k = 0; k < Nk; k++) { 21 PetscInt sum; 22 23 PetscCall(PetscDTIndexToBary(d + 1, n, k, btup)); 24 for (j = 0, sum = 0; j < d + 1; j++) { 25 PetscCheck(btup[j] >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTIndexToBary, d = %" PetscInt_FMT ", n = %" PetscInt_FMT ", k = %" PetscInt_FMT " negative entry", d, n, k); 26 sum += btup[j]; 27 } 28 PetscCheck(sum == n, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTIndexToBary, d = %" PetscInt_FMT ", n = %" PetscInt_FMT ", k = %" PetscInt_FMT " incorrect sum", d, n, k); 29 PetscCall(PetscDTBaryToIndex(d + 1, n, btup, &kchk)); 30 PetscCheck(kchk == k, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTBaryToIndex, d = %" PetscInt_FMT ", n = %" PetscInt_FMT ", k = %" PetscInt_FMT " mismatch", d, n, k); 31 if (k) { 32 j = d; 33 while (j >= 0 && btup[j] == btupprev[j]) j--; 34 PetscCheck(j >= 0, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTIndexToBary, d = %" PetscInt_FMT ", n = %" PetscInt_FMT ", k = %" PetscInt_FMT " equal to previous", d, n, k); 35 PetscCheck(btup[j] >= btupprev[j], PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTIndexToBary, d = %" PetscInt_FMT ", n = %" PetscInt_FMT ", k = %" PetscInt_FMT " less to previous", d, n, k); 36 } else PetscCall(PetscArraycpy(btupprev, btup, d + 1)); 37 PetscCall(PetscDTIndexToGradedOrder(d, Nk - 1 - k, gtup)); 38 PetscCall(PetscDTGradedOrderToIndex(d, gtup, &kchk)); 39 PetscCheck(kchk == Nk - 1 - k, PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTGradedOrderToIndex, d = %" PetscInt_FMT ", n = %" PetscInt_FMT ", k = %" PetscInt_FMT " mismatch", d, n, Nk - 1 - k); 40 for (j = 0; j < d; j++) { PetscCheck(gtup[j] == btup[d - 1 - j], PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTIndexToGradedOrder, d = %" PetscInt_FMT ", n = %" PetscInt_FMT ", k = %" PetscInt_FMT " incorrect", d, n, Nk - 1 - k); } 41 } 42 } 43 } 44 PetscCall(PetscFree3(btup, btupprev, gtup)); 45 PetscCall(PetscFinalize()); 46 return 0; 47 } 48 49 /*TEST 50 51 test: 52 53 TEST*/ 54