| Line | Branch | Exec | Source |
|---|---|---|---|
| 1 | /* | ||
| 2 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 3 | SLEPc - Scalable Library for Eigenvalue Problem Computations | ||
| 4 | Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain | ||
| 5 | |||
| 6 | This file is part of SLEPc. | ||
| 7 | SLEPc is distributed under a 2-clause BSD license (see LICENSE). | ||
| 8 | - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - | ||
| 9 | */ | ||
| 10 | |||
| 11 | static char help[] = "Test DSHEP with block size larger than one.\n\n"; | ||
| 12 | |||
| 13 | #include <slepcds.h> | ||
| 14 | |||
| 15 | 10 | int main(int argc,char **argv) | |
| 16 | { | ||
| 17 | 10 | DS ds; | |
| 18 | 10 | SlepcSC sc; | |
| 19 | 10 | PetscScalar *A,*eig; | |
| 20 | 10 | PetscInt i,j,n,ld,bs,maxbw=3,nblks=8; | |
| 21 | 10 | PetscViewer viewer; | |
| 22 | 10 | PetscBool verbose; | |
| 23 | |||
| 24 |
1/2✓ Branch 0 taken 1 times.
✗ Branch 1 not taken.
|
10 | PetscFunctionBeginUser; |
| 25 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(SlepcInitialize(&argc,&argv,NULL,help)); |
| 26 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-maxbw",&maxbw,NULL)); |
| 27 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(PetscOptionsGetInt(NULL,NULL,"-nblks",&nblks,NULL)); |
| 28 | 10 | n = maxbw*nblks; | |
| 29 | 10 | bs = maxbw; | |
| 30 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Solve a block HEP Dense System - dimension %" PetscInt_FMT " (bandwidth=%" PetscInt_FMT ", blocks=%" PetscInt_FMT ").\n",n,maxbw,nblks)); |
| 31 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(PetscOptionsHasName(NULL,NULL,"-verbose",&verbose)); |
| 32 | |||
| 33 | /* Create DS object */ | ||
| 34 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSCreate(PETSC_COMM_WORLD,&ds)); |
| 35 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSSetType(ds,DSHEP)); |
| 36 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSSetMethod(ds,3)); /* Select block divide-and-conquer */ |
| 37 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSSetBlockSize(ds,bs)); |
| 38 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSSetFromOptions(ds)); |
| 39 | 10 | ld = n; | |
| 40 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSAllocate(ds,ld)); |
| 41 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSSetDimensions(ds,n,0,0)); |
| 42 | |||
| 43 | /* Set up viewer */ | ||
| 44 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(PetscViewerASCIIGetStdout(PETSC_COMM_WORLD,&viewer)); |
| 45 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_INFO_DETAIL)); |
| 46 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSView(ds,viewer)); |
| 47 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(PetscViewerPopFormat(viewer)); |
| 48 |
1/8✗ Branch 0 not taken.
✓ Branch 1 taken 5 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
10 | if (verbose) PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_ASCII_MATLAB)); |
| 49 | |||
| 50 | /* Fill with a symmetric band Toeplitz matrix */ | ||
| 51 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSGetArray(ds,DS_MAT_A,&A)); |
| 52 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
145 | for (i=0;i<n;i++) A[i+i*ld]=2.0; |
| 53 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
40 | for (j=1;j<=bs;j++) { |
| 54 |
2/2✓ Branch 0 taken 5 times.
✓ Branch 1 taken 5 times.
|
375 | for (i=0;i<n-j;i++) { A[i+(i+j)*ld]=1.0; A[(i+j)+i*ld]=1.0; } |
| 55 | } | ||
| 56 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSRestoreArray(ds,DS_MAT_A,&A)); |
| 57 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSSetState(ds,DS_STATE_RAW)); |
| 58 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 5 times.
|
10 | if (verbose) { |
| 59 | ✗ | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Initial - - - - - - - - -\n")); | |
| 60 | ✗ | PetscCall(DSView(ds,viewer)); | |
| 61 | } | ||
| 62 | |||
| 63 | /* Solve */ | ||
| 64 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(PetscMalloc1(n,&eig)); |
| 65 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSGetSlepcSC(ds,&sc)); |
| 66 | 10 | sc->comparison = SlepcCompareSmallestReal; | |
| 67 | 10 | sc->comparisonctx = NULL; | |
| 68 | 10 | sc->map = NULL; | |
| 69 | 10 | sc->mapobj = NULL; | |
| 70 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSSolve(ds,eig,NULL)); |
| 71 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSSort(ds,eig,NULL,NULL,NULL,NULL)); |
| 72 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 5 times.
|
10 | if (verbose) { |
| 73 | ✗ | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"After solve - - - - - - - - -\n")); | |
| 74 | ✗ | PetscCall(DSView(ds,viewer)); | |
| 75 | } | ||
| 76 | |||
| 77 | /* Print eigenvalues */ | ||
| 78 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Computed eigenvalues =\n")); |
| 79 |
7/8✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 5 times.
✓ Branch 3 taken 4 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
✓ Branch 6 taken 1 times.
✓ Branch 7 taken 1 times.
|
145 | for (i=0;i<n;i++) PetscCall(PetscViewerASCIIPrintf(viewer," %.5f\n",(double)PetscRealPart(eig[i]))); |
| 80 | |||
| 81 |
5/8✓ Branch 0 taken 5 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
|
10 | PetscCall(PetscFree(eig)); |
| 82 |
4/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 4 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
|
10 | PetscCall(DSDestroy(&ds)); |
| 83 |
3/6✓ Branch 0 taken 1 times.
✓ Branch 1 taken 1 times.
✓ Branch 2 taken 1 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
10 | PetscCall(SlepcFinalize()); |
| 84 | return 0; | ||
| 85 | } | ||
| 86 | |||
| 87 | /*TEST | ||
| 88 | |||
| 89 | test: | ||
| 90 | suffix: 1 | ||
| 91 | requires: !complex !single | ||
| 92 | |||
| 93 | test: | ||
| 94 | suffix: 2 | ||
| 95 | args: -nblks 1 | ||
| 96 | requires: !complex | ||
| 97 | |||
| 98 | TEST*/ | ||
| 99 |