| 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 | #include <slepc/private/slepcimpl.h> /*I "slepcsys.h" I*/ | ||
| 12 | |||
| 13 | 598 | static PetscErrorCode MatCreateTile_Seq(PetscScalar a,Mat A,PetscScalar b,Mat B,PetscScalar c,Mat C,PetscScalar d,Mat D,Mat G) | |
| 14 | { | ||
| 15 | 598 | PetscInt i,j,M1,M2,N1,N2,ncols,*scols; | |
| 16 | 598 | PetscScalar *svals,*buf; | |
| 17 | 598 | const PetscInt *cols; | |
| 18 | 598 | const PetscScalar *vals; | |
| 19 | |||
| 20 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
598 | PetscFunctionBegin; |
| 21 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
598 | PetscCall(MatGetSize(A,&M1,&N1)); |
| 22 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
598 | PetscCall(MatGetSize(D,&M2,&N2)); |
| 23 | |||
| 24 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
598 | PetscCall(PetscMalloc2(PetscMax(N1,N2),&buf,PetscMax(N1,N2),&scols)); |
| 25 | /* Transfer A */ | ||
| 26 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
598 | if (a!=0.0) { |
| 27 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
34724 | for (i=0;i<M1;i++) { |
| 28 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
34220 | PetscCall(MatGetRow(A,i,&ncols,&cols,&vals)); |
| 29 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
34220 | if (a!=1.0) { |
| 30 | 2840 | svals=buf; | |
| 31 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
13880 | for (j=0;j<ncols;j++) svals[j] = vals[j]*a; |
| 32 | 31380 | } else svals=(PetscScalar*)vals; | |
| 33 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
34220 | PetscCall(MatSetValues(G,1,&i,ncols,cols,svals,INSERT_VALUES)); |
| 34 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
34220 | PetscCall(MatRestoreRow(A,i,&ncols,&cols,&vals)); |
| 35 | } | ||
| 36 | } | ||
| 37 | /* Transfer B */ | ||
| 38 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
598 | if (b!=0.0) { |
| 39 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
34104 | for (i=0;i<M1;i++) { |
| 40 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
33620 | PetscCall(MatGetRow(B,i,&ncols,&cols,&vals)); |
| 41 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
33620 | if (b!=1.0) { |
| 42 | 2840 | svals=buf; | |
| 43 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
11540 | for (j=0;j<ncols;j++) svals[j] = vals[j]*b; |
| 44 | 30780 | } else svals=(PetscScalar*)vals; | |
| 45 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
300882 | for (j=0;j<ncols;j++) scols[j] = cols[j]+N1; |
| 46 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
33620 | PetscCall(MatSetValues(G,1,&i,ncols,scols,svals,INSERT_VALUES)); |
| 47 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
33620 | PetscCall(MatRestoreRow(B,i,&ncols,&cols,&vals)); |
| 48 | } | ||
| 49 | } | ||
| 50 | /* Transfer C */ | ||
| 51 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
598 | if (c!=0.0) { |
| 52 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
24545 | for (i=0;i<M2;i++) { |
| 53 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
24071 | PetscCall(MatGetRow(C,i,&ncols,&cols,&vals)); |
| 54 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
24071 | if (c!=1.0) { |
| 55 | 5360 | svals=buf; | |
| 56 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
25972 | for (j=0;j<ncols;j++) svals[j] = vals[j]*c; |
| 57 | 18711 | } else svals=(PetscScalar*)vals; | |
| 58 | 24071 | j = i+M1; | |
| 59 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
24071 | PetscCall(MatSetValues(G,1,&j,ncols,cols,svals,INSERT_VALUES)); |
| 60 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
24071 | PetscCall(MatRestoreRow(C,i,&ncols,&cols,&vals)); |
| 61 | } | ||
| 62 | } | ||
| 63 | /* Transfer D */ | ||
| 64 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
598 | if (d!=0.0) { |
| 65 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
30419 | for (i=0;i<M2;i++) { |
| 66 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
29871 | PetscCall(MatGetRow(D,i,&ncols,&cols,&vals)); |
| 67 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
29871 | if (d!=1.0) { |
| 68 | 8240 | svals=buf; | |
| 69 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
24392 | for (j=0;j<ncols;j++) svals[j] = vals[j]*d; |
| 70 | 21631 | } else svals=(PetscScalar*)vals; | |
| 71 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
119842 | for (j=0;j<ncols;j++) scols[j] = cols[j]+N1; |
| 72 | 29871 | j = i+M1; | |
| 73 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
29871 | PetscCall(MatSetValues(G,1,&j,ncols,scols,svals,INSERT_VALUES)); |
| 74 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
29871 | PetscCall(MatRestoreRow(D,i,&ncols,&cols,&vals)); |
| 75 | } | ||
| 76 | } | ||
| 77 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
598 | PetscCall(PetscFree2(buf,scols)); |
| 78 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
116 | PetscFunctionReturn(PETSC_SUCCESS); |
| 79 | } | ||
| 80 | |||
| 81 | 60 | static PetscErrorCode MatCreateTile_MPI(PetscScalar a,Mat A,PetscScalar b,Mat B,PetscScalar c,Mat C,PetscScalar d,Mat D,Mat G) | |
| 82 | { | ||
| 83 | 60 | PetscMPIInt np; | |
| 84 | 60 | PetscInt p,i,j,N1,N2,m1,m2,*map1,*map2; | |
| 85 | 60 | PetscInt ncols,*scols,start,gstart; | |
| 86 | 60 | PetscScalar *svals,*buf; | |
| 87 | 60 | const PetscInt *cols,*mapptr1,*mapptr2; | |
| 88 | 60 | const PetscScalar *vals; | |
| 89 | |||
| 90 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
60 | PetscFunctionBegin; |
| 91 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatGetSize(A,NULL,&N1)); |
| 92 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatGetLocalSize(A,&m1,NULL)); |
| 93 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatGetSize(D,NULL,&N2)); |
| 94 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatGetLocalSize(D,&m2,NULL)); |
| 95 | |||
| 96 | /* Create mappings */ | ||
| 97 |
14/28✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
|
60 | PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)G),&np)); |
| 98 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatGetOwnershipRangesColumn(A,&mapptr1)); |
| 99 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatGetOwnershipRangesColumn(B,&mapptr2)); |
| 100 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(PetscMalloc4(PetscMax(N1,N2),&buf,PetscMax(N1,N2),&scols,N1,&map1,N2,&map2)); |
| 101 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
180 | for (p=0;p<np;p++) { |
| 102 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1520 | for (i=mapptr1[p];i<mapptr1[p+1];i++) map1[i] = i+mapptr2[p]; |
| 103 | } | ||
| 104 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
180 | for (p=0;p<np;p++) { |
| 105 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1520 | for (i=mapptr2[p];i<mapptr2[p+1];i++) map2[i] = i+mapptr1[p+1]; |
| 106 | } | ||
| 107 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatGetOwnershipRange(G,&gstart,NULL)); |
| 108 | |||
| 109 | /* Transfer A */ | ||
| 110 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
60 | if (a!=0.0) { |
| 111 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(MatGetOwnershipRange(A,&start,NULL)); |
| 112 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
440 | for (i=0;i<m1;i++) { |
| 113 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
400 | PetscCall(MatGetRow(A,i+start,&ncols,&cols,&vals)); |
| 114 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
400 | if (a!=1.0) { |
| 115 | 100 | svals=buf; | |
| 116 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
380 | for (j=0;j<ncols;j++) svals[j] = vals[j]*a; |
| 117 | 300 | } else svals=(PetscScalar*)vals; | |
| 118 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
980 | for (j=0;j<ncols;j++) scols[j] = map1[cols[j]]; |
| 119 | 400 | j = gstart+i; | |
| 120 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
400 | PetscCall(MatSetValues(G,1,&j,ncols,scols,svals,INSERT_VALUES)); |
| 121 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
400 | PetscCall(MatRestoreRow(A,i+start,&ncols,&cols,&vals)); |
| 122 | } | ||
| 123 | } | ||
| 124 | /* Transfer B */ | ||
| 125 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
60 | if (b!=0.0) { |
| 126 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
40 | PetscCall(MatGetOwnershipRange(B,&start,NULL)); |
| 127 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
440 | for (i=0;i<m1;i++) { |
| 128 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
400 | PetscCall(MatGetRow(B,i+start,&ncols,&cols,&vals)); |
| 129 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
400 | if (b!=1.0) { |
| 130 | 100 | svals=buf; | |
| 131 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
200 | for (j=0;j<ncols;j++) svals[j] = vals[j]*b; |
| 132 | 300 | } else svals=(PetscScalar*)vals; | |
| 133 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
800 | for (j=0;j<ncols;j++) scols[j] = map2[cols[j]]; |
| 134 | 400 | j = gstart+i; | |
| 135 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
400 | PetscCall(MatSetValues(G,1,&j,ncols,scols,svals,INSERT_VALUES)); |
| 136 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
400 | PetscCall(MatRestoreRow(B,i+start,&ncols,&cols,&vals)); |
| 137 | } | ||
| 138 | } | ||
| 139 | /* Transfer C */ | ||
| 140 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
60 | if (c!=0.0) { |
| 141 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
20 | PetscCall(MatGetOwnershipRange(C,&start,NULL)); |
| 142 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
320 | for (i=0;i<m2;i++) { |
| 143 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
300 | PetscCall(MatGetRow(C,i+start,&ncols,&cols,&vals)); |
| 144 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
300 | if (c!=1.0) { |
| 145 | 300 | svals=buf; | |
| 146 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1180 | for (j=0;j<ncols;j++) svals[j] = vals[j]*c; |
| 147 | ✗ | } else svals=(PetscScalar*)vals; | |
| 148 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1180 | for (j=0;j<ncols;j++) scols[j] = map1[cols[j]]; |
| 149 | 300 | j = gstart+m1+i; | |
| 150 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
300 | PetscCall(MatSetValues(G,1,&j,ncols,scols,svals,INSERT_VALUES)); |
| 151 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
300 | PetscCall(MatRestoreRow(C,i+start,&ncols,&cols,&vals)); |
| 152 | } | ||
| 153 | } | ||
| 154 | /* Transfer D */ | ||
| 155 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
60 | if (d!=0.0) { |
| 156 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(MatGetOwnershipRange(D,&start,NULL)); |
| 157 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
760 | for (i=0;i<m2;i++) { |
| 158 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
700 | PetscCall(MatGetRow(D,i+start,&ncols,&cols,&vals)); |
| 159 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
700 | if (d!=1.0) { |
| 160 | 400 | svals=buf; | |
| 161 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
1560 | for (j=0;j<ncols;j++) svals[j] = vals[j]*d; |
| 162 | 300 | } else svals=(PetscScalar*)vals; | |
| 163 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
2160 | for (j=0;j<ncols;j++) scols[j] = map2[cols[j]]; |
| 164 | 700 | j = gstart+m1+i; | |
| 165 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
700 | PetscCall(MatSetValues(G,1,&j,ncols,scols,svals,INSERT_VALUES)); |
| 166 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
700 | PetscCall(MatRestoreRow(D,i+start,&ncols,&cols,&vals)); |
| 167 | } | ||
| 168 | } | ||
| 169 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
60 | PetscCall(PetscFree4(buf,scols,map1,map2)); |
| 170 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
12 | PetscFunctionReturn(PETSC_SUCCESS); |
| 171 | } | ||
| 172 | |||
| 173 | /*@ | ||
| 174 | MatCreateTile - Explicitly build a matrix from four blocks. | ||
| 175 | |||
| 176 | Collective | ||
| 177 | |||
| 178 | Input Parameters: | ||
| 179 | + A - matrix for top-left block | ||
| 180 | . a - scaling factor for block A | ||
| 181 | . B - matrix for top-right block | ||
| 182 | . b - scaling factor for block B | ||
| 183 | . C - matrix for bottom-left block | ||
| 184 | . c - scaling factor for block C | ||
| 185 | . D - matrix for bottom-right block | ||
| 186 | - d - scaling factor for block D | ||
| 187 | |||
| 188 | Output Parameter: | ||
| 189 | . G - the resulting matrix | ||
| 190 | |||
| 191 | Notes: | ||
| 192 | The resulting matrix is built as | ||
| 193 | $$G = \begin{bmatrix} aA & bB\\ cC & dD \end{bmatrix}.$$ | ||
| 194 | |||
| 195 | In the case of a parallel matrix, a permuted version of `G` is returned. The permutation | ||
| 196 | is a perfect shuffle such that the local parts of `A`, `B`, `C`, `D` remain in the local part of | ||
| 197 | `G` for the same process. | ||
| 198 | |||
| 199 | Matrix `G` must be destroyed by the user. | ||
| 200 | |||
| 201 | The blocks can be of different type. They can be either `MATCONSTANTDIAGONAL`, or a standard | ||
| 202 | type such as `MATAIJ`, or any other type provided that it supports the `MatGetRow()` operation. | ||
| 203 | The type of the output matrix will be the same as the first block that is not | ||
| 204 | `MATCONSTANTDIAGONAL` (checked in the `A`,`B`,`C`,`D` order). | ||
| 205 | |||
| 206 | Level: developer | ||
| 207 | |||
| 208 | .seealso: `MatCreateNest()` | ||
| 209 | @*/ | ||
| 210 | 658 | PetscErrorCode MatCreateTile(PetscScalar a,Mat A,PetscScalar b,Mat B,PetscScalar c,Mat C,PetscScalar d,Mat D,Mat *G) | |
| 211 | { | ||
| 212 | 658 | PetscInt i,k,M1,M2,N1,N2,M,N,m1,m2,n1,n2,m,n,bs; | |
| 213 | 658 | PetscBool diag[4]; | |
| 214 | 658 | Mat block[4] = {A,B,C,D}; | |
| 215 | 658 | MatType type[4]; | |
| 216 | 658 | PetscMPIInt size; | |
| 217 | |||
| 218 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
658 | PetscFunctionBegin; |
| 219 |
3/16✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
658 | PetscValidHeaderSpecific(A,MAT_CLASSID,2); |
| 220 |
3/16✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
658 | PetscValidHeaderSpecific(B,MAT_CLASSID,4); |
| 221 |
3/16✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
658 | PetscValidHeaderSpecific(C,MAT_CLASSID,6); |
| 222 |
3/16✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
658 | PetscValidHeaderSpecific(D,MAT_CLASSID,8); |
| 223 |
13/32✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
|
658 | PetscCheckSameTypeAndComm(A,2,B,4); |
| 224 |
13/32✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
|
658 | PetscCheckSameTypeAndComm(A,2,C,6); |
| 225 |
13/32✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
|
658 | PetscCheckSameTypeAndComm(A,2,D,8); |
| 226 |
30/68✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✓ Branch 22 taken 2 times.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✓ Branch 26 taken 2 times.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✓ Branch 30 taken 2 times.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✓ Branch 33 taken 2 times.
✓ Branch 34 taken 2 times.
✗ Branch 35 not taken.
✓ Branch 36 taken 2 times.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✓ Branch 39 taken 2 times.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✓ Branch 43 taken 2 times.
✗ Branch 44 not taken.
✗ Branch 45 not taken.
✗ Branch 46 not taken.
✓ Branch 47 taken 2 times.
✓ Branch 48 taken 2 times.
✗ Branch 49 not taken.
✗ Branch 50 not taken.
✓ Branch 51 taken 2 times.
✓ Branch 52 taken 2 times.
✗ Branch 53 not taken.
✓ Branch 54 taken 2 times.
✗ Branch 55 not taken.
✗ Branch 56 not taken.
✓ Branch 57 taken 2 times.
✗ Branch 58 not taken.
✗ Branch 59 not taken.
✓ Branch 60 taken 2 times.
✗ Branch 61 not taken.
✓ Branch 62 taken 2 times.
✗ Branch 63 not taken.
✗ Branch 64 not taken.
✓ Branch 65 taken 2 times.
✗ Branch 66 not taken.
✗ Branch 67 not taken.
|
658 | PetscValidLogicalCollectiveScalar(A,a,1); |
| 227 |
30/68✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✓ Branch 22 taken 2 times.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✓ Branch 26 taken 2 times.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✓ Branch 30 taken 2 times.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✓ Branch 33 taken 2 times.
✓ Branch 34 taken 2 times.
✗ Branch 35 not taken.
✓ Branch 36 taken 2 times.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✓ Branch 39 taken 2 times.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✓ Branch 43 taken 2 times.
✗ Branch 44 not taken.
✗ Branch 45 not taken.
✗ Branch 46 not taken.
✓ Branch 47 taken 2 times.
✓ Branch 48 taken 2 times.
✗ Branch 49 not taken.
✗ Branch 50 not taken.
✓ Branch 51 taken 2 times.
✓ Branch 52 taken 2 times.
✗ Branch 53 not taken.
✓ Branch 54 taken 2 times.
✗ Branch 55 not taken.
✗ Branch 56 not taken.
✓ Branch 57 taken 2 times.
✗ Branch 58 not taken.
✗ Branch 59 not taken.
✓ Branch 60 taken 2 times.
✗ Branch 61 not taken.
✓ Branch 62 taken 2 times.
✗ Branch 63 not taken.
✗ Branch 64 not taken.
✓ Branch 65 taken 2 times.
✗ Branch 66 not taken.
✗ Branch 67 not taken.
|
658 | PetscValidLogicalCollectiveScalar(A,b,3); |
| 228 |
30/68✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✓ Branch 22 taken 2 times.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✓ Branch 26 taken 2 times.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✓ Branch 30 taken 2 times.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✓ Branch 33 taken 2 times.
✓ Branch 34 taken 2 times.
✗ Branch 35 not taken.
✓ Branch 36 taken 2 times.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✓ Branch 39 taken 2 times.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✓ Branch 43 taken 2 times.
✗ Branch 44 not taken.
✗ Branch 45 not taken.
✗ Branch 46 not taken.
✓ Branch 47 taken 2 times.
✓ Branch 48 taken 2 times.
✗ Branch 49 not taken.
✗ Branch 50 not taken.
✓ Branch 51 taken 2 times.
✓ Branch 52 taken 2 times.
✗ Branch 53 not taken.
✓ Branch 54 taken 2 times.
✗ Branch 55 not taken.
✗ Branch 56 not taken.
✓ Branch 57 taken 2 times.
✗ Branch 58 not taken.
✗ Branch 59 not taken.
✓ Branch 60 taken 2 times.
✗ Branch 61 not taken.
✓ Branch 62 taken 2 times.
✗ Branch 63 not taken.
✗ Branch 64 not taken.
✓ Branch 65 taken 2 times.
✗ Branch 66 not taken.
✗ Branch 67 not taken.
|
658 | PetscValidLogicalCollectiveScalar(A,c,5); |
| 229 |
30/68✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 2 times.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✓ Branch 22 taken 2 times.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✓ Branch 26 taken 2 times.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✓ Branch 29 taken 2 times.
✓ Branch 30 taken 2 times.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✓ Branch 33 taken 2 times.
✓ Branch 34 taken 2 times.
✗ Branch 35 not taken.
✓ Branch 36 taken 2 times.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✓ Branch 39 taken 2 times.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✓ Branch 43 taken 2 times.
✗ Branch 44 not taken.
✗ Branch 45 not taken.
✗ Branch 46 not taken.
✓ Branch 47 taken 2 times.
✓ Branch 48 taken 2 times.
✗ Branch 49 not taken.
✗ Branch 50 not taken.
✓ Branch 51 taken 2 times.
✓ Branch 52 taken 2 times.
✗ Branch 53 not taken.
✓ Branch 54 taken 2 times.
✗ Branch 55 not taken.
✗ Branch 56 not taken.
✓ Branch 57 taken 2 times.
✗ Branch 58 not taken.
✗ Branch 59 not taken.
✓ Branch 60 taken 2 times.
✗ Branch 61 not taken.
✓ Branch 62 taken 2 times.
✗ Branch 63 not taken.
✗ Branch 64 not taken.
✓ Branch 65 taken 2 times.
✗ Branch 66 not taken.
✗ Branch 67 not taken.
|
658 | PetscValidLogicalCollectiveScalar(A,d,7); |
| 230 |
2/8✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
658 | PetscAssertPointer(G,9); |
| 231 | |||
| 232 | /* check row 1 */ | ||
| 233 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetSize(A,&M1,NULL)); |
| 234 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetLocalSize(A,&m1,NULL)); |
| 235 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetSize(B,&M,NULL)); |
| 236 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetLocalSize(B,&m,NULL)); |
| 237 |
2/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
658 | PetscCheck(M==M1 && m==m1,PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_INCOMP,"Incompatible dimensions"); |
| 238 | /* check row 2 */ | ||
| 239 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetSize(C,&M2,NULL)); |
| 240 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetLocalSize(C,&m2,NULL)); |
| 241 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetSize(D,&M,NULL)); |
| 242 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetLocalSize(D,&m,NULL)); |
| 243 |
2/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
658 | PetscCheck(M==M2 && m==m2,PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_INCOMP,"Incompatible dimensions"); |
| 244 | /* check column 1 */ | ||
| 245 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetSize(A,NULL,&N1)); |
| 246 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetLocalSize(A,NULL,&n1)); |
| 247 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetSize(C,NULL,&N)); |
| 248 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetLocalSize(C,NULL,&n)); |
| 249 |
2/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
658 | PetscCheck(N==N1 && n==n1,PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_INCOMP,"Incompatible dimensions"); |
| 250 | /* check column 2 */ | ||
| 251 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetSize(B,NULL,&N2)); |
| 252 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetLocalSize(B,NULL,&n2)); |
| 253 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetSize(D,NULL,&N)); |
| 254 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetLocalSize(D,NULL,&n)); |
| 255 |
2/6✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 10 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
|
658 | PetscCheck(N==N2 && n==n2,PetscObjectComm((PetscObject)A),PETSC_ERR_ARG_INCOMP,"Incompatible dimensions"); |
| 256 | |||
| 257 | /* check matrix types */ | ||
| 258 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
3290 | for (i=0;i<4;i++) { |
| 259 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2632 | PetscCall(MatGetType(block[i],&type[i])); |
| 260 |
3/6✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2632 | PetscCall(PetscStrcmp(type[i],MATCONSTANTDIAGONAL,&diag[i])); |
| 261 | } | ||
| 262 |
3/4✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✗ Branch 3 not taken.
|
793 | for (k=0;k<4;k++) if (!diag[k]) break; |
| 263 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
658 | PetscCheck(k<4,PetscObjectComm((PetscObject)A),PETSC_ERR_SUP,"Not implemented for 4 diagonal blocks"); |
| 264 | |||
| 265 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatGetBlockSize(block[k],&bs)); |
| 266 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatCreate(PetscObjectComm((PetscObject)block[k]),G)); |
| 267 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatSetSizes(*G,m1+m2,n1+n2,M1+M2,N1+N2)); |
| 268 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatSetType(*G,type[k])); |
| 269 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatSetBlockSize(*G,bs)); |
| 270 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatSetUp(*G)); |
| 271 | |||
| 272 |
14/28✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
|
658 | PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)*G),&size)); |
| 273 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
658 | if (size>1) PetscCall(MatCreateTile_MPI(a,A,b,B,c,C,d,D,*G)); |
| 274 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
598 | else PetscCall(MatCreateTile_Seq(a,A,b,B,c,C,d,D,*G)); |
| 275 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatAssemblyBegin(*G,MAT_FINAL_ASSEMBLY)); |
| 276 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
658 | PetscCall(MatAssemblyEnd(*G,MAT_FINAL_ASSEMBLY)); |
| 277 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
128 | PetscFunctionReturn(PETSC_SUCCESS); |
| 278 | } | ||
| 279 | |||
| 280 | /*@ | ||
| 281 | MatCreateVecsEmpty - Get vector(s) compatible with the matrix, i.e., with the same | ||
| 282 | parallel layout, but without internal array. | ||
| 283 | |||
| 284 | Collective | ||
| 285 | |||
| 286 | Input Parameter: | ||
| 287 | . mat - the matrix | ||
| 288 | |||
| 289 | Output Parameters: | ||
| 290 | + right - (optional) vector that the matrix can be multiplied against | ||
| 291 | - left - (optional) vector that the matrix vector product can be stored in | ||
| 292 | |||
| 293 | Note: | ||
| 294 | This is similar to `MatCreateVecs()`, but the new vectors do not have an internal | ||
| 295 | array, so the intended usage is with `VecPlaceArray()`. | ||
| 296 | |||
| 297 | Level: developer | ||
| 298 | |||
| 299 | .seealso: `MatCreateVecs()`, `VecDuplicateEmpty()` | ||
| 300 | @*/ | ||
| 301 | 35750 | PetscErrorCode MatCreateVecsEmpty(Mat mat,Vec *right,Vec *left) | |
| 302 | { | ||
| 303 | 35750 | PetscBool standard,cuda=PETSC_FALSE,hip=PETSC_FALSE,skip=PETSC_FALSE; | |
| 304 | 35750 | PetscInt M,N,mloc,nloc,rbs,cbs; | |
| 305 | 35750 | PetscMPIInt size; | |
| 306 | 35750 | Vec v; | |
| 307 | |||
| 308 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
35750 | PetscFunctionBegin; |
| 309 |
3/16✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
35750 | PetscValidHeaderSpecific(mat,MAT_CLASSID,1); |
| 310 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
35750 | PetscValidType(mat,1); |
| 311 | |||
| 312 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
35750 | PetscCall(PetscObjectTypeCompareAny((PetscObject)mat,&standard,MATSEQAIJ,MATMPIAIJ,MATSEQBAIJ,MATMPIBAIJ,MATSEQSBAIJ,MATMPISBAIJ,MATSEQDENSE,MATMPIDENSE,"")); |
| 313 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
35750 | PetscCall(PetscObjectTypeCompareAny((PetscObject)mat,&cuda,MATSEQAIJCUSPARSE,MATMPIAIJCUSPARSE,"")); |
| 314 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
35750 | PetscCall(PetscObjectTypeCompareAny((PetscObject)mat,&hip,MATSEQAIJHIPSPARSE,MATMPIAIJHIPSPARSE,"")); |
| 315 |
6/6✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 10 times.
✓ Branch 5 taken 2 times.
|
35750 | if (!standard && !cuda && !hip) { |
| 316 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2822 | PetscCall(MatCreateVecs(mat,right,left)); |
| 317 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
2822 | v = right? *right: *left; |
| 318 |
1/2✓ Branch 0 taken 10 times.
✗ Branch 1 not taken.
|
2822 | if (v) { |
| 319 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2822 | PetscCall(PetscObjectTypeCompareAny((PetscObject)v,&standard,VECSEQ,VECMPI,"")); |
| 320 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2822 | PetscCall(PetscObjectTypeCompareAny((PetscObject)v,&cuda,VECSEQCUDA,VECMPICUDA,"")); |
| 321 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
2822 | PetscCall(PetscObjectTypeCompareAny((PetscObject)v,&hip,VECSEQHIP,VECMPIHIP,"")); |
| 322 | } | ||
| 323 |
5/6✓ Branch 0 taken 4 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
|
2822 | if (!standard && !cuda && !hip) skip = PETSC_TRUE; |
| 324 | else { | ||
| 325 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
2822 | if (right) PetscCall(VecDestroy(right)); |
| 326 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
2822 | if (left) PetscCall(VecDestroy(left)); |
| 327 | } | ||
| 328 | } | ||
| 329 | if (!skip) { | ||
| 330 |
14/28✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
✓ Branch 12 taken 2 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 2 times.
✓ Branch 16 taken 2 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 2 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 2 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 2 times.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
|
35750 | PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)mat),&size)); |
| 331 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
35750 | PetscCall(MatGetLocalSize(mat,&mloc,&nloc)); |
| 332 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
35750 | PetscCall(MatGetSize(mat,&M,&N)); |
| 333 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
35750 | PetscCall(MatGetBlockSizes(mat,&rbs,&cbs)); |
| 334 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
35750 | if (right) { |
| 335 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 2 times.
|
19657 | if (cuda) { |
| 336 | #if defined(PETSC_HAVE_CUDA) | ||
| 337 |
3/4✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
|
273 | if (size>1) PetscCall(VecCreateMPICUDAWithArray(PetscObjectComm((PetscObject)mat),cbs,nloc,N,NULL,right)); |
| 338 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
|
261 | else PetscCall(VecCreateSeqCUDAWithArray(PetscObjectComm((PetscObject)mat),cbs,N,NULL,right)); |
| 339 | #endif | ||
| 340 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 2 times.
|
19384 | } else if (hip) { |
| 341 | #if defined(PETSC_HAVE_HIP) | ||
| 342 |
3/4✓ Branch 0 taken 2 times.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
|
273 | if (size>1) PetscCall(VecCreateMPIHIPWithArray(PetscObjectComm((PetscObject)mat),cbs,nloc,N,NULL,right)); |
| 343 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
|
261 | else PetscCall(VecCreateSeqHIPWithArray(PetscObjectComm((PetscObject)mat),cbs,N,NULL,right)); |
| 344 | #endif | ||
| 345 | } else { | ||
| 346 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
19111 | if (size>1) PetscCall(VecCreateMPIWithArray(PetscObjectComm((PetscObject)mat),cbs,nloc,N,NULL,right)); |
| 347 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
16234 | else PetscCall(VecCreateSeqWithArray(PetscObjectComm((PetscObject)mat),cbs,N,NULL,right)); |
| 348 | } | ||
| 349 | } | ||
| 350 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
35750 | if (left) { |
| 351 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 2 times.
|
17814 | if (cuda) { |
| 352 | #if defined(PETSC_HAVE_CUDA) | ||
| 353 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
210 | if (size>1) PetscCall(VecCreateMPICUDAWithArray(PetscObjectComm((PetscObject)mat),rbs,mloc,M,NULL,left)); |
| 354 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
|
210 | else PetscCall(VecCreateSeqCUDAWithArray(PetscObjectComm((PetscObject)mat),rbs,M,NULL,left)); |
| 355 | #endif | ||
| 356 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 2 times.
|
17604 | } else if (hip) { |
| 357 | #if defined(PETSC_HAVE_HIP) | ||
| 358 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
210 | if (size>1) PetscCall(VecCreateMPIHIPWithArray(PetscObjectComm((PetscObject)mat),rbs,mloc,M,NULL,left)); |
| 359 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
|
210 | else PetscCall(VecCreateSeqHIPWithArray(PetscObjectComm((PetscObject)mat),rbs,M,NULL,left)); |
| 360 | #endif | ||
| 361 | } else { | ||
| 362 |
6/8✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
17394 | if (size>1) PetscCall(VecCreateMPIWithArray(PetscObjectComm((PetscObject)mat),rbs,mloc,M,NULL,left)); |
| 363 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
15778 | else PetscCall(VecCreateSeqWithArray(PetscObjectComm((PetscObject)mat),rbs,M,NULL,left)); |
| 364 | } | ||
| 365 | } | ||
| 366 | } | ||
| 367 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
7068 | PetscFunctionReturn(PETSC_SUCCESS); |
| 368 | } | ||
| 369 | |||
| 370 | /*@ | ||
| 371 | MatNormEstimate - Estimate the 2-norm of a matrix. | ||
| 372 | |||
| 373 | Collective | ||
| 374 | |||
| 375 | Input Parameters: | ||
| 376 | + A - the matrix | ||
| 377 | . vrn - random vector with normally distributed entries (can be `NULL`) | ||
| 378 | - w - workspace vector (can be `NULL`) | ||
| 379 | |||
| 380 | Output Parameter: | ||
| 381 | . nrm - the norm estimate | ||
| 382 | |||
| 383 | Notes: | ||
| 384 | Does not need access to the matrix entries, just performs a matrix-vector product. | ||
| 385 | Based on work by I. Ipsen and coworkers <https://ipsen.math.ncsu.edu/ps/slides_ima.pdf>. | ||
| 386 | |||
| 387 | The input vector `vrn` must have unit 2-norm. | ||
| 388 | If `vrn` is `NULL`, then it is created internally and filled with `VecSetRandomNormal()`. | ||
| 389 | |||
| 390 | Level: developer | ||
| 391 | |||
| 392 | .seealso: `VecSetRandomNormal()` | ||
| 393 | @*/ | ||
| 394 | 693 | PetscErrorCode MatNormEstimate(Mat A,Vec vrn,Vec w,PetscReal *nrm) | |
| 395 | { | ||
| 396 | 693 | PetscInt n; | |
| 397 | 693 | Vec vv=NULL,ww=NULL; | |
| 398 | |||
| 399 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
693 | PetscFunctionBegin; |
| 400 |
3/16✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
|
693 | PetscValidHeaderSpecific(A,MAT_CLASSID,1); |
| 401 |
1/4✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
|
693 | PetscValidType(A,1); |
| 402 |
3/14✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
693 | if (vrn) PetscValidHeaderSpecific(vrn,VEC_CLASSID,2); |
| 403 |
3/14✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
|
693 | if (w) PetscValidHeaderSpecific(w,VEC_CLASSID,3); |
| 404 |
2/8✗ Branch 0 not taken.
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
|
693 | PetscAssertPointer(nrm,4); |
| 405 | |||
| 406 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
693 | if (!vrn) { |
| 407 | ✗ | PetscCall(MatCreateVecs(A,&vv,NULL)); | |
| 408 | ✗ | vrn = vv; | |
| 409 | ✗ | PetscCall(VecSetRandomNormal(vv,NULL,NULL,NULL)); | |
| 410 | ✗ | PetscCall(VecNormalize(vv,NULL)); | |
| 411 | } | ||
| 412 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 10 times.
|
693 | if (!w) { |
| 413 | ✗ | PetscCall(MatCreateVecs(A,&ww,NULL)); | |
| 414 | ✗ | w = ww; | |
| 415 | } | ||
| 416 | |||
| 417 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
693 | PetscCall(MatGetSize(A,&n,NULL)); |
| 418 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
693 | PetscCall(MatMult(A,vrn,w)); |
| 419 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
693 | PetscCall(VecNorm(w,NORM_2,nrm)); |
| 420 | 693 | *nrm *= PetscSqrtReal((PetscReal)n); | |
| 421 | |||
| 422 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
693 | PetscCall(VecDestroy(&vv)); |
| 423 |
4/6✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
|
693 | PetscCall(VecDestroy(&ww)); |
| 424 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
135 | PetscFunctionReturn(PETSC_SUCCESS); |
| 425 | } | ||
| 426 |