1*c4762a1bSJed Brown static char help[] = "Tests PetscIsCloseAtTol() routine.\n"; 2*c4762a1bSJed Brown 3*c4762a1bSJed Brown #include <petscsys.h> 4*c4762a1bSJed Brown 5*c4762a1bSJed Brown PETSC_INTERN PetscReal zero; 6*c4762a1bSJed Brown PetscReal zero = 0; 7*c4762a1bSJed Brown 8*c4762a1bSJed Brown #define CALL(call) do { \ 9*c4762a1bSJed Brown PetscErrorCode _ierr; \ 10*c4762a1bSJed Brown _ierr = PetscPrintf(PETSC_COMM_WORLD,"%s -> %s\n",#call,(call)?"True":"False");CHKERRQ(_ierr); \ 11*c4762a1bSJed Brown } while(0); 12*c4762a1bSJed Brown 13*c4762a1bSJed Brown int main(int argc, char **argv) { 14*c4762a1bSJed Brown 15*c4762a1bSJed Brown PetscReal eps = PETSC_MACHINE_EPSILON; 16*c4762a1bSJed Brown PetscReal neg_zero = PetscRealConstant(-0.0); 17*c4762a1bSJed Brown PetscReal pos_zero = PetscRealConstant(+0.0); 18*c4762a1bSJed Brown PetscReal neg_one = PetscRealConstant(-1.0); 19*c4762a1bSJed Brown PetscReal pos_one = PetscRealConstant(+1.0); 20*c4762a1bSJed Brown PetscReal neg_inf = neg_one/zero; /* -inf */ 21*c4762a1bSJed Brown PetscReal pos_inf = pos_one/zero; /* +inf */ 22*c4762a1bSJed Brown PetscReal x_nan = zero/zero; /* NaN */ 23*c4762a1bSJed Brown 24*c4762a1bSJed Brown PetscErrorCode ierr; 25*c4762a1bSJed Brown ierr = PetscInitialize(&argc,&argv,NULL,help);if (ierr) return ierr; 26*c4762a1bSJed Brown 27*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_zero,neg_zero,0,0)); 28*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one,pos_one,0,0)); 29*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one,neg_one,0,0)); 30*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one,neg_one,0,2)); 31*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one,neg_one,2,0)); 32*c4762a1bSJed Brown 33*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one+eps,pos_one,0,0)); 34*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one-eps,pos_one,0,0)); 35*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one+eps,pos_one,0,0)); 36*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one-eps,pos_one,0,0)); 37*c4762a1bSJed Brown 38*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one+eps,pos_one,0,eps)); 39*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one-eps,pos_one,0,eps)); 40*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one+eps,pos_one,eps,0)); 41*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one-eps,pos_one,eps,0)); 42*c4762a1bSJed Brown 43*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one+2*eps,pos_one,eps,0)); 44*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one-2*eps,pos_one,eps,0)); 45*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one+2*eps,pos_one,0,eps)); 46*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_one-2*eps,pos_one,0,eps)); 47*c4762a1bSJed Brown 48*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(neg_inf,neg_zero,2,2)); 49*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(neg_inf,pos_zero,2,2)); 50*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(neg_inf,neg_one,2,2)); 51*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(neg_inf,pos_one,2,2)); 52*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(neg_inf,neg_inf,2,2)); 53*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(neg_inf,pos_inf,2,2)); 54*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(neg_inf,x_nan,2,2)); 55*c4762a1bSJed Brown 56*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_inf,neg_zero,2,2)); 57*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_inf,pos_zero,2,2)); 58*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_inf,neg_one,2,2)); 59*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_inf,pos_one,2,2)); 60*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_inf,neg_inf,2,2)); 61*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_inf,pos_inf,2,2)); 62*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(pos_inf,x_nan,2,2)); 63*c4762a1bSJed Brown 64*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(x_nan,neg_zero,2,2)); 65*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(x_nan,pos_zero,2,2)); 66*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(x_nan,neg_one,2,2)); 67*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(x_nan,pos_one,2,2)); 68*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(x_nan,neg_inf,2,2)); 69*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(x_nan,pos_inf,2,2)); 70*c4762a1bSJed Brown CALL(PetscIsCloseAtTol(x_nan,x_nan,2,2)); 71*c4762a1bSJed Brown 72*c4762a1bSJed Brown ierr = PetscFinalize(); 73*c4762a1bSJed Brown return ierr; 74*c4762a1bSJed Brown } 75*c4762a1bSJed Brown 76*c4762a1bSJed Brown 77*c4762a1bSJed Brown /*TEST 78*c4762a1bSJed Brown 79*c4762a1bSJed Brown test: 80*c4762a1bSJed Brown output_file: output/ex39.out 81*c4762a1bSJed Brown 82*c4762a1bSJed Brown TEST*/ 83