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 { 9 PetscInt d, n, maxdim = 4; 10 PetscInt *btupprev, *btup; 11 PetscInt *gtup; 12 PetscErrorCode ierr; 13 14 ierr = PetscInitialize(&argc, &argv, NULL, help);if (ierr) return ierr; 15 CHKERRQ(PetscMalloc3(maxdim + 1, &btup, maxdim + 1, &btupprev, maxdim, >up)); 16 for (d = 0; d <= maxdim; d++) { 17 for (n = 0; n <= d + 2; n++) { 18 PetscInt j, k, Nk, kchk; 19 20 CHKERRQ(PetscDTBinomialInt(d + n, d, &Nk)); 21 for (k = 0; k < Nk; k++) { 22 PetscInt sum; 23 24 CHKERRQ(PetscDTIndexToBary(d + 1, n, k, btup)); 25 for (j = 0, sum = 0; j < d + 1; j++) { 26 PetscCheckFalse(btup[j] < 0,PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTIndexToBary, d = %D, n = %D, k = %D negative entry", d, n, k); 27 sum += btup[j]; 28 } 29 PetscCheckFalse(sum != n,PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTIndexToBary, d = %D, n = %D, k = %D incorrect sum", d, n, k); 30 CHKERRQ(PetscDTBaryToIndex(d + 1, n, btup, &kchk)); 31 PetscCheckFalse(kchk != k,PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTBaryToIndex, d = %D, n = %D, k = %D mismatch", d, n, k); 32 if (k) { 33 j = d; 34 while (j >= 0 && btup[j] == btupprev[j]) j--; 35 PetscCheckFalse(j < 0,PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTIndexToBary, d = %D, n = %D, k = %D equal to previous", d, n, k); 36 PetscCheckFalse(btup[j] < btupprev[j],PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTIndexToBary, d = %D, n = %D, k = %D less to previous", d, n, k); 37 } else { 38 CHKERRQ(PetscArraycpy(btupprev, btup, d + 1)); 39 } 40 CHKERRQ(PetscDTIndexToGradedOrder(d, Nk - 1 - k, gtup)); 41 CHKERRQ(PetscDTGradedOrderToIndex(d, gtup, &kchk)); 42 PetscCheckFalse(kchk != Nk - 1 - k,PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTGradedOrderToIndex, d = %D, n = %D, k = %D mismatch", d, n, Nk - 1 - k); 43 for (j = 0; j < d; j++) { 44 PetscCheckFalse(gtup[j] != btup[d - 1 - j],PETSC_COMM_SELF, PETSC_ERR_PLIB, "PetscDTIndexToGradedOrder, d = %D, n = %D, k = %D incorrect", d, n, Nk - 1 - k); 45 } 46 } 47 } 48 } 49 CHKERRQ(PetscFree3(btup, btupprev, gtup)); 50 ierr = PetscFinalize(); 51 return ierr; 52 } 53 54 /*TEST 55 56 test: 57 58 TEST*/ 59