1c4762a1bSJed Brown 2c4762a1bSJed Brown static char help[] = "Tests PetscOptionsGetScalar(), PetscOptionsScalarArray() for complex numbers\n"; 3c4762a1bSJed Brown 4c4762a1bSJed Brown #include <petscsys.h> 5c4762a1bSJed Brown 6c4762a1bSJed Brown int main(int argc,char **argv) 7c4762a1bSJed Brown { 8c4762a1bSJed Brown PetscInt ierr,n,i; 9c4762a1bSJed Brown PetscScalar a,array[10]; 10c4762a1bSJed Brown PetscReal rarray[10]; 11c4762a1bSJed Brown 12*b122ec5aSJacob Faibussowitsch CHKERRQ(PetscInitialize(&argc,&argv,(char*)0,help)); 135f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsGetScalar(NULL,NULL,"-a",&a,NULL)); 145f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"Scalar a = %g + %gi\n",(double)PetscRealPart(a),(double)PetscImaginaryPart(a))); 15c4762a1bSJed Brown 16c4762a1bSJed Brown ierr = PetscOptionsBegin(PETSC_COMM_WORLD,NULL,"test options",NULL);CHKERRQ(ierr); 17c4762a1bSJed Brown n = 10; /* max num of input values */ 185f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsRealArray("-rarray", "Input a real array", "ex14.c", rarray, &n, NULL)); 19c4762a1bSJed Brown if (n) { 205f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"Real rarray of length %" PetscInt_FMT "\n",n)); 21c4762a1bSJed Brown for (i=0; i<n; i++) { 225f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_SELF," %g,\n",(double)rarray[i])); 23c4762a1bSJed Brown } 24c4762a1bSJed Brown } 25c4762a1bSJed Brown 26c4762a1bSJed Brown n = 10; /* max num of input values */ 275f80ce2aSJacob Faibussowitsch CHKERRQ(PetscOptionsScalarArray("-array", "Input a scalar array", "ex14.c", array, &n, NULL)); 28c4762a1bSJed Brown if (n) { 295f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_SELF,"Scalar rarray of length %" PetscInt_FMT "\n",n)); 30c4762a1bSJed Brown for (i=0; i<n; i++) { 31c4762a1bSJed Brown if (PetscImaginaryPart(array[i]) < 0.0) { 325f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_SELF," %g - %gi\n",(double)PetscRealPart(array[i]),(double)PetscAbsReal(PetscImaginaryPart(array[i])))); 33c4762a1bSJed Brown } else { 345f80ce2aSJacob Faibussowitsch CHKERRQ(PetscPrintf(PETSC_COMM_SELF," %g + %gi\n",(double)PetscRealPart(array[i]),(double)PetscImaginaryPart(array[i]))); 35c4762a1bSJed Brown } 36c4762a1bSJed Brown } 37c4762a1bSJed Brown } 38c4762a1bSJed Brown ierr = PetscOptionsEnd();CHKERRQ(ierr); 39*b122ec5aSJacob Faibussowitsch CHKERRQ(PetscFinalize()); 40*b122ec5aSJacob Faibussowitsch return 0; 41c4762a1bSJed Brown } 42c4762a1bSJed Brown 43c4762a1bSJed Brown /*TEST 44c4762a1bSJed Brown 45c4762a1bSJed Brown test: 46c4762a1bSJed Brown requires: complex 47c4762a1bSJed 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 48c4762a1bSJed Brown 49c4762a1bSJed Brown TEST*/ 50