1c4762a1bSJed Brown static char help[] = "Tests 1D cell-based discretization tools.\n\n";
2c4762a1bSJed Brown
3c4762a1bSJed Brown #include <petscdt.h>
4c4762a1bSJed Brown #include <petscviewer.h>
5c4762a1bSJed Brown
main(int argc,char ** argv)6d71ae5a4SJacob Faibussowitsch int main(int argc, char **argv)
7d71ae5a4SJacob Faibussowitsch {
8c4762a1bSJed Brown PetscInt i, j, degrees[1000], ndegrees, nsrc_points, ntarget_points;
9c4762a1bSJed Brown PetscReal src_points[1000], target_points[1000], *R;
10c4762a1bSJed Brown PetscBool flg;
11c4762a1bSJed Brown
12327415f7SBarry Smith PetscFunctionBeginUser;
13*c8025a54SPierre Jolivet PetscCall(PetscInitialize(&argc, &argv, NULL, help));
14d0609cedSBarry Smith PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "Discretization tools test options", NULL);
15c4762a1bSJed Brown {
16c4762a1bSJed Brown ndegrees = 1000;
17c4762a1bSJed Brown degrees[0] = 1;
18c4762a1bSJed Brown degrees[1] = 2;
19c4762a1bSJed Brown degrees[2] = 3;
209566063dSJacob Faibussowitsch PetscCall(PetscOptionsIntArray("-degrees", "list of max degrees to evaluate", "", degrees, &ndegrees, &flg));
21c4762a1bSJed Brown if (!flg) ndegrees = 3;
22c4762a1bSJed Brown
23c4762a1bSJed Brown nsrc_points = 1000;
24c4762a1bSJed Brown src_points[0] = -1.;
25c4762a1bSJed Brown src_points[1] = 0.;
26c4762a1bSJed Brown src_points[2] = 1.;
279566063dSJacob Faibussowitsch PetscCall(PetscOptionsRealArray("-src_points", "list of points defining intervals on which to integrate", "", src_points, &nsrc_points, &flg));
28c4762a1bSJed Brown if (!flg) nsrc_points = 3;
29c4762a1bSJed Brown
30c4762a1bSJed Brown ntarget_points = 1000;
31c4762a1bSJed Brown target_points[0] = -1.;
32c4762a1bSJed Brown target_points[1] = 0.;
33c4762a1bSJed Brown target_points[2] = 1.;
349566063dSJacob Faibussowitsch PetscCall(PetscOptionsRealArray("-target_points", "list of points defining intervals on which to integrate", "", target_points, &ntarget_points, &flg));
35c4762a1bSJed Brown if (!flg) ntarget_points = 3;
36c4762a1bSJed Brown }
37d0609cedSBarry Smith PetscOptionsEnd();
38c4762a1bSJed Brown
399566063dSJacob Faibussowitsch PetscCall(PetscMalloc1((nsrc_points - 1) * (ntarget_points - 1), &R));
40c4762a1bSJed Brown for (i = 0; i < ndegrees; i++) {
419566063dSJacob Faibussowitsch PetscCall(PetscDTReconstructPoly(degrees[i], nsrc_points - 1, src_points, ntarget_points - 1, target_points, R));
42c4762a1bSJed Brown for (j = 0; j < (ntarget_points - 1) * (nsrc_points - 1); j++) { /* Truncate to zero for nicer output */
43c4762a1bSJed Brown if (PetscAbs(R[j]) < 10 * PETSC_MACHINE_EPSILON) R[j] = 0;
44c4762a1bSJed Brown }
45c4762a1bSJed Brown for (j = 0; j < ntarget_points - 1; j++) {
4663a3b9bcSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Degree %" PetscInt_FMT " target interval (%g,%g)\n", degrees[i], (double)target_points[j], (double)target_points[j + 1]));
479566063dSJacob Faibussowitsch PetscCall(PetscRealView(nsrc_points - 1, R + j * (nsrc_points - 1), PETSC_VIEWER_STDOUT_WORLD));
48c4762a1bSJed Brown }
49c4762a1bSJed Brown }
509566063dSJacob Faibussowitsch PetscCall(PetscFree(R));
519566063dSJacob Faibussowitsch PetscCall(PetscFinalize());
52b122ec5aSJacob Faibussowitsch return 0;
53c4762a1bSJed Brown }
54c4762a1bSJed Brown
55c4762a1bSJed Brown /*TEST
56c4762a1bSJed Brown test:
57c4762a1bSJed Brown suffix: 1
58c4762a1bSJed Brown args: -degrees 1,2,3 -target_points -0.3,0,.2 -src_points -1,-.3,0,.2,1
59c4762a1bSJed Brown TEST*/
60