1 static char help[] = "Tests VecMDot(),VecDot(),VecMTDot(), and VecTDot()\n"; 2 3 #include <petscvec.h> 4 5 int main(int argc, char **argv) 6 { 7 PetscErrorCode ierr; 8 Vec *V,t; 9 PetscInt i,j,reps,n=15,k=6; 10 PetscRandom rctx; 11 PetscScalar *val_dot,*val_mdot,*tval_dot,*tval_mdot; 12 13 ierr = PetscInitialize(&argc,&argv,(char*)0,help);if (ierr) return ierr; 14 ierr = PetscOptionsGetInt(NULL,NULL,"-n",&n,NULL);CHKERRQ(ierr); 15 ierr = PetscOptionsGetInt(NULL,NULL,"-k",&k,NULL);CHKERRQ(ierr); 16 ierr = PetscPrintf(PETSC_COMM_WORLD,"Test with %D random vectors of length %D",k,n);CHKERRQ(ierr); 17 ierr = PetscPrintf(PETSC_COMM_WORLD,"\n",k,n);CHKERRQ(ierr); 18 ierr = PetscRandomCreate(PETSC_COMM_WORLD,&rctx);CHKERRQ(ierr); 19 ierr = PetscRandomSetFromOptions(rctx);CHKERRQ(ierr); 20 #if defined(PETSC_USE_COMPLEX) 21 ierr = PetscRandomSetInterval(rctx,-1.+4.*PETSC_i,1.+5.*PETSC_i);CHKERRQ(ierr); 22 #else 23 ierr = PetscRandomSetInterval(rctx,-1.,1.);CHKERRQ(ierr); 24 #endif 25 ierr = VecCreate(PETSC_COMM_WORLD,&t);CHKERRQ(ierr); 26 ierr = VecSetSizes(t,n,PETSC_DECIDE);CHKERRQ(ierr); 27 ierr = VecSetFromOptions(t);CHKERRQ(ierr); 28 ierr = VecDuplicateVecs(t,k,&V);CHKERRQ(ierr); 29 ierr = VecSetRandom(t,rctx);CHKERRQ(ierr); 30 ierr = VecViewFromOptions(t,NULL,"-t_view");CHKERRQ(ierr); 31 ierr = PetscMalloc1(k,&val_dot);CHKERRQ(ierr); 32 ierr = PetscMalloc1(k,&val_mdot);CHKERRQ(ierr); 33 ierr = PetscMalloc1(k,&tval_dot);CHKERRQ(ierr); 34 ierr = PetscMalloc1(k,&tval_mdot);CHKERRQ(ierr); 35 for (i=0; i<k; i++) { ierr = VecSetRandom(V[i],rctx);CHKERRQ(ierr); } 36 for (reps=0; reps<20; reps++) { 37 for (i=1; i<k; i++) { 38 ierr = VecMDot(t,i,V,val_mdot);CHKERRQ(ierr); 39 ierr = VecMTDot(t,i,V,tval_mdot);CHKERRQ(ierr); 40 for (j=0;j<i;j++) { 41 ierr = VecDot(t,V[j],&val_dot[j]);CHKERRQ(ierr); 42 ierr = VecTDot(t,V[j],&tval_dot[j]);CHKERRQ(ierr); 43 } 44 /* Check result */ 45 for (j=0;j<i;j++) { 46 if (PetscAbsScalar(val_mdot[j] - val_dot[j])/PetscAbsScalar(val_dot[j]) > 1e-5) { 47 ierr = PetscPrintf(PETSC_COMM_WORLD, "[TEST FAILED] i=%D, j=%D, val_mdot[j]=%g, val_dot[j]=%g\n",i,j,(double)PetscAbsScalar(val_mdot[j]), (double)PetscAbsScalar(val_dot[j]));CHKERRQ(ierr); 48 break; 49 } 50 if (PetscAbsScalar(tval_mdot[j] - tval_dot[j])/PetscAbsScalar(tval_dot[j]) > 1e-5) { 51 ierr = PetscPrintf(PETSC_COMM_WORLD, "[TEST FAILED] i=%D, j=%D, tval_mdot[j]=%g, tval_dot[j]=%g\n",i,j,(double)PetscAbsScalar(tval_mdot[j]), (double)PetscAbsScalar(tval_dot[j]));CHKERRQ(ierr); 52 break; 53 } 54 } 55 } 56 } 57 ierr = PetscPrintf(PETSC_COMM_WORLD,"Test completed successfully!\n",k,n);CHKERRQ(ierr); 58 ierr = PetscFree(val_dot);CHKERRQ(ierr); 59 ierr = PetscFree(val_mdot);CHKERRQ(ierr); 60 ierr = PetscFree(tval_dot);CHKERRQ(ierr); 61 ierr = PetscFree(tval_mdot);CHKERRQ(ierr); 62 ierr = VecDestroyVecs(k,&V);CHKERRQ(ierr); 63 ierr = VecDestroy(&t);CHKERRQ(ierr); 64 ierr = PetscRandomDestroy(&rctx);CHKERRQ(ierr); 65 ierr = PetscFinalize(); 66 return ierr; 67 } 68 69 /*TEST 70 71 test: 72 73 testset: 74 output_file: output/ex43_1.out 75 76 test: 77 suffix: cuda 78 args: -vec_type cuda -random_type curand 79 requires: cuda 80 81 test: 82 suffix: kokkos 83 args: -vec_type kokkos 84 requires: kokkos_kernels 85 86 TEST*/ 87