1c4762a1bSJed Brown static char help[] = "Tests PetscOptionsGetScalar(), PetscOptionsScalarArray() for complex numbers\n";
2c4762a1bSJed Brown
3c4762a1bSJed Brown #include <petscsys.h>
4c4762a1bSJed Brown
main(int argc,char ** argv)5d71ae5a4SJacob Faibussowitsch int main(int argc, char **argv)
6d71ae5a4SJacob Faibussowitsch {
7d0609cedSBarry Smith PetscInt n, i;
8c4762a1bSJed Brown PetscScalar a, array[10];
9c4762a1bSJed Brown PetscReal rarray[10];
10c4762a1bSJed Brown
11327415f7SBarry Smith PetscFunctionBeginUser;
12*c8025a54SPierre Jolivet PetscCall(PetscInitialize(&argc, &argv, NULL, help));
139566063dSJacob Faibussowitsch PetscCall(PetscOptionsGetScalar(NULL, NULL, "-a", &a, NULL));
149566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_SELF, "Scalar a = %g + %gi\n", (double)PetscRealPart(a), (double)PetscImaginaryPart(a)));
15c4762a1bSJed Brown
16d0609cedSBarry Smith PetscOptionsBegin(PETSC_COMM_WORLD, NULL, "test options", NULL);
17c4762a1bSJed Brown n = 10; /* max num of input values */
189566063dSJacob Faibussowitsch PetscCall(PetscOptionsRealArray("-rarray", "Input a real array", "ex14.c", rarray, &n, NULL));
19c4762a1bSJed Brown if (n) {
209566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_SELF, "Real rarray of length %" PetscInt_FMT "\n", n));
2148a46eb9SPierre Jolivet for (i = 0; i < n; i++) PetscCall(PetscPrintf(PETSC_COMM_SELF, " %g,\n", (double)rarray[i]));
22c4762a1bSJed Brown }
23c4762a1bSJed Brown
24c4762a1bSJed Brown n = 10; /* max num of input values */
259566063dSJacob Faibussowitsch PetscCall(PetscOptionsScalarArray("-array", "Input a scalar array", "ex14.c", array, &n, NULL));
26c4762a1bSJed Brown if (n) {
279566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_SELF, "Scalar rarray of length %" PetscInt_FMT "\n", n));
28c4762a1bSJed Brown for (i = 0; i < n; i++) {
29c4762a1bSJed Brown if (PetscImaginaryPart(array[i]) < 0.0) {
309566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_SELF, " %g - %gi\n", (double)PetscRealPart(array[i]), (double)PetscAbsReal(PetscImaginaryPart(array[i]))));
31c4762a1bSJed Brown } else {
329566063dSJacob Faibussowitsch PetscCall(PetscPrintf(PETSC_COMM_SELF, " %g + %gi\n", (double)PetscRealPart(array[i]), (double)PetscImaginaryPart(array[i])));
33c4762a1bSJed Brown }
34c4762a1bSJed Brown }
35c4762a1bSJed Brown }
36d0609cedSBarry Smith PetscOptionsEnd();
379566063dSJacob Faibussowitsch PetscCall(PetscFinalize());
38b122ec5aSJacob Faibussowitsch return 0;
39c4762a1bSJed Brown }
40c4762a1bSJed Brown
41c4762a1bSJed Brown /*TEST
42c4762a1bSJed Brown
43c4762a1bSJed Brown test:
44c4762a1bSJed Brown requires: complex
45c4762a1bSJed Brown args: -array 1.0,-2-3i,4.5+6.2i,4.5,6.8+4i,i,-i,-1.2i -rarray 1,2,3 -a 1.5+2.1i
46c4762a1bSJed Brown
47c4762a1bSJed Brown TEST*/
48