#include <stdlib.h>
#include "anamod.h"
#include "anamodsalsamodules.h"
#include "petscerror.h"
#include "petscmat.h"
#include "src/mat/impls/aij/seq/aij.h"
Go to the source code of this file.
Defines | |
#define | RC_NORM(i, j, v, w) |
#define | ASSERT3(v, s, i, j) if (!(v)) SETERRQ3(1,"Violation of <%s> @ (%d,%d): %d",s,i,j); |
#define | ASSERT(v, s, i, j, t) if (!(v)) SETERRQ4(1,"Violation of <%s> @ (%d,%d): %d",s,i,j,t); |
Functions | |
static PetscErrorCode | computetrace (Mat A) |
static PetscErrorCode | Trace (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | TraceAbs (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | norm1 (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | normF (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | normInf (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | compute_dd (Mat A, PetscTruth *flg) |
static PetscErrorCode | DiagonalDominance (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | MatSymmPartNormInf_seq (Mat A, PetscReal *sn_r, PetscReal *an_r, PetscReal *srn_r, PetscReal *arn_r, int *nun_r, PetscTruth *done) |
static PetscErrorCode | MatSymmPartNormInf (Mat A, PetscTruth *done) |
static PetscErrorCode | SymmetrySNorm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | SymmetryANorm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | SymmetryFSNorm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | SymmetryFANorm (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
static PetscErrorCode | NUnstruct (AnaModNumericalProblem prob, AnalysisItem *rv, int *lv, PetscTruth *flg) |
PetscErrorCode | RegisterSimpleModules () |
PetscErrorCode | DeRegisterSimpleModules (void) |
PetscErrorCode RegisterSimpleModules();
Compute these elements with
ComputeQuantity("simple",<element>,A,(AnalysisItem*)&res,&flg);
Available elements are:
Definition in file simple.c.
#define ASSERT | ( | v, | |||
s, | |||||
i, | |||||
j, | |||||
t | ) | if (!(v)) SETERRQ4(1,"Violation of <%s> @ (%d,%d): %d",s,i,j,t); |
Referenced by MatSymmPartNormInf_seq().
#define ASSERT3 | ( | v, | |||
s, | |||||
i, | |||||
j | ) | if (!(v)) SETERRQ3(1,"Violation of <%s> @ (%d,%d): %d",s,i,j); |
Referenced by MatSymmPartNormInf_seq().
#define RC_NORM | ( | i, | |||
j, | |||||
v, | |||||
w | ) |
Value:
{ \ PetscReal sum = PetscAbsScalar(v+w)/2, dif = PetscAbsScalar(v-w)/2; \ ASSERT3(!(i==j && v!=w),"same elements on diagonal",i,j); \ sronorm[i] += sum; if (i!=j) sronorm[j] += sum; \ aronorm[i] += dif; if (i!=j) aronorm[j] += dif; \ if (i==j) snorm += sum*sum; \ else {snorm += 2*sum*sum; anorm += 2*dif*dif;} \ }
Referenced by MatSymmPartNormInf_seq().
static PetscErrorCode compute_dd | ( | Mat | A, | |
PetscTruth * | flg | |||
) | [static] |
Compute the diagonal dominance of a matrix: minimum value of diagonal element (not absolute!) minus off-diagonal elements absolute.
Definition at line 246 of file simple.c.
References GetDataID(), HASTOEXIST, and id.
Referenced by DiagonalDominance().
00247 { 00248 PetscReal r,gr; 00249 MPI_Comm comm; 00250 int first,last,m,row; PetscErrorCode ierr; 00251 00252 PetscFunctionBegin; 00253 ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr); 00254 ierr = MatGetOwnershipRange(A,&first,&last); CHKERRQ(ierr); 00255 m = last-first; 00256 r = 0; 00257 for (row=first; row<last; row++) { 00258 int col,ncols; const int *cols; 00259 const PetscScalar *vals; PetscReal rr=0,d=0; 00260 ierr = MatGetRow(A,row,&ncols,&cols,&vals); CHKERRQ(ierr); 00261 for (col=0; col<ncols; col++) { 00262 PetscReal v = PetscRealPart(vals[col]); 00263 if (cols[col]==row) 00264 {d = v; rr += v;} else rr -= v; 00265 } 00266 r = PetscMin(r,rr); 00267 ierr = MatRestoreRow(A,row,&ncols,&cols,&vals); CHKERRQ(ierr); 00268 } 00269 MPI_Allreduce((void*)&r,(void*)&gr,1,MPIU_REAL,MPI_MIN,comm); 00270 { 00271 int id; PetscTruth has; 00272 ierr = GetDataID("simple","diagonal-dominance",&id,&has); CHKERRQ(ierr); 00273 HASTOEXIST(has); 00274 ierr = PetscObjectComposedDataSetReal 00275 ((PetscObject)A,id,gr); CHKERRQ(ierr); 00276 } 00277 *flg = PETSC_TRUE; 00278 PetscFunctionReturn(0); 00279 }
static PetscErrorCode computetrace | ( | Mat | A | ) | [static] |
This is the computational routine for trace calculation. It computes both the trace and the absolute trace, that is, using absolute values of matrix elements.
A Petsc vector is created and destroyed.
Definition at line 63 of file simple.c.
References GetDataID(), HASTOEXIST, and id.
Referenced by Trace(), and TraceAbs().
00064 { 00065 MPI_Comm comm; Vec D; PetscScalar *d; 00066 PetscTruth has; PetscScalar vsum,gvsum; PetscReal asum,gasum; 00067 int local_size,i,id; 00068 PetscErrorCode ierr; 00069 00070 PetscFunctionBegin; 00071 00072 ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr); 00073 ierr = MatGetLocalSize(A,&local_size,PETSC_NULL); CHKERRQ(ierr); 00074 ierr = VecCreateMPI(comm,local_size,PETSC_DECIDE,&D); CHKERRQ(ierr); 00075 ierr = MatGetDiagonal(A,D); CHKERRQ(ierr); 00076 ierr = VecGetArray(D,&d); CHKERRQ(ierr); 00077 vsum = 0; asum = 0; 00078 for (i=0; i<local_size; i++) { 00079 vsum += d[i]; asum += PetscAbsScalar(d[i]); 00080 } 00081 ierr = VecRestoreArray(D,&d); CHKERRQ(ierr); 00082 ierr = VecDestroy(D); CHKERRQ(ierr); 00083 ierr = MPI_Allreduce 00084 ((void*)&vsum,(void*)&gvsum,1,MPIU_SCALAR,MPI_SUM,comm); CHKERRQ(ierr); 00085 ierr = MPI_Allreduce 00086 ((void*)&asum,(void*)&gasum,1,MPI_DOUBLE,MPI_SUM,comm); CHKERRQ(ierr); 00087 00088 ierr = GetDataID("simple","trace",&id,&has); CHKERRQ(ierr); 00089 HASTOEXIST(has); 00090 ierr = PetscObjectComposedDataSetScalar 00091 ((PetscObject)A,id,gvsum); CHKERRQ(ierr); 00092 ierr = GetDataID("simple","trace-abs",&id,&has); CHKERRQ(ierr); 00093 HASTOEXIST(has); 00094 ierr = PetscObjectComposedDataSetReal 00095 ((PetscObject)A,id,gasum); CHKERRQ(ierr); 00096 00097 PetscFunctionReturn(0); 00098 }
PetscErrorCode DeRegisterSimpleModules | ( | void | ) |
static PetscErrorCode DiagonalDominance | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Definition at line 282 of file simple.c.
References compute_dd(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSimpleModules().
00283 { 00284 Mat A = (Mat)prob; 00285 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00286 PetscFunctionBegin; 00287 00288 ierr = GetDataID("simple","diagonal-dominance",&id,&has); CHKERRQ(ierr); 00289 HASTOEXIST(has); 00290 ierr = PetscObjectComposedDataGetReal 00291 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00292 if (*flg) rv->r = v; 00293 else { 00294 ierr = compute_dd(A,flg); CHKERRQ(ierr); 00295 if (*flg) { 00296 ierr = PetscObjectComposedDataGetReal 00297 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00298 if (*flg) rv->r = v; 00299 } 00300 } 00301 00302 PetscFunctionReturn(0); 00303 }
static PetscErrorCode MatSymmPartNormInf | ( | Mat | A, | |
PetscTruth * | done | |||
) | [static] |
This is the computational routine for all quantities relating to symmetric and anti-symmetric parts of the matrix.
This computation can only be done in single-processor mode. If the commandline option (see Commandline options) "-anamod-force sequential"
is specified, the matrix is gathered on processor zero to compute this quantity.
See MatSymmPartNormInf_seq() for the actual computation.
Definition at line 432 of file simple.c.
References AnaModGetSequentialMatrix(), GetDataID(), HASTOEXIST, id, and MatSymmPartNormInf_seq().
Referenced by NUnstruct(), SymmetryANorm(), SymmetryFANorm(), SymmetryFSNorm(), and SymmetrySNorm().
00433 { 00434 MPI_Comm comm; Mat Awork; int id,nun; 00435 PetscReal snorm,anorm,sronorm,aronorm; 00436 PetscTruth has,alloc_work,do_work,any_work; PetscErrorCode ierr; 00437 00438 PetscFunctionBegin; 00439 00440 ierr = PetscObjectGetComm((PetscObject)A,&comm); CHKERRQ(ierr); 00441 ierr = AnaModGetSequentialMatrix 00442 (A,&Awork,&alloc_work,&do_work,&any_work); CHKERRQ(ierr); 00443 if (!any_work) { 00444 *done = PETSC_FALSE; PetscFunctionReturn(0); 00445 } 00446 if (do_work) { 00447 ierr = MatSymmPartNormInf_seq 00448 (Awork,&snorm,&anorm,&sronorm,&aronorm,&nun,done); CHKERRQ(ierr); 00449 } 00450 if (alloc_work) { 00451 ierr = MatDestroy(Awork); CHKERRQ(ierr); 00452 } 00453 00454 MPI_Bcast((void*)done,1,MPI_INT,0,comm); 00455 MPI_Bcast((void*)&snorm,1,MPIU_REAL,0,comm); 00456 MPI_Bcast((void*)&anorm,1,MPIU_REAL,0,comm); 00457 MPI_Bcast((void*)&sronorm,1,MPIU_REAL,0,comm); 00458 MPI_Bcast((void*)&aronorm,1,MPIU_REAL,0,comm); 00459 MPI_Bcast((void*)&nun,1,MPI_INT,0,comm); 00460 00461 /* infinity norms of symmetric / anti-symmetric part */ 00462 ierr = GetDataID("simple","symmetry-snorm",&id,&has); CHKERRQ(ierr); 00463 HASTOEXIST(has); 00464 ierr = PetscObjectComposedDataSetReal 00465 ((PetscObject)A,id,sronorm); CHKERRQ(ierr); 00466 00467 ierr = GetDataID("simple","symmetry-anorm",&id,&has); CHKERRQ(ierr); 00468 HASTOEXIST(has); 00469 ierr = PetscObjectComposedDataSetReal 00470 ((PetscObject)A,id,aronorm); CHKERRQ(ierr); 00471 00472 /* frobenius norms of symmetric / anti-symmetric part */ 00473 ierr = GetDataID("simple","symmetry-fsnorm",&id,&has); CHKERRQ(ierr); 00474 HASTOEXIST(has); 00475 ierr = PetscObjectComposedDataSetReal 00476 ((PetscObject)A,id,snorm); CHKERRQ(ierr); 00477 00478 ierr = GetDataID("simple","symmetry-fanorm",&id,&has); CHKERRQ(ierr); 00479 HASTOEXIST(has); 00480 ierr = PetscObjectComposedDataSetReal 00481 ((PetscObject)A,id,anorm); CHKERRQ(ierr); 00482 00483 ierr = GetDataID("structure","n-struct-unsymm",&id,&has); CHKERRQ(ierr); 00484 HASTOEXIST(has); 00485 ierr = PetscObjectComposedDataSetInt 00486 ((PetscObject)A,id,nun); CHKERRQ(ierr); 00487 00488 PetscFunctionReturn(0); 00489 }
static PetscErrorCode MatSymmPartNormInf_seq | ( | Mat | A, | |
PetscReal * | sn_r, | |||
PetscReal * | an_r, | |||
PetscReal * | srn_r, | |||
PetscReal * | arn_r, | |||
int * | nun_r, | |||
PetscTruth * | done | |||
) | [static] |
This is where the real work of MatSymmPartNormInf() is done
A few temporary arrays of the size of a vector are allocated and freed. Otherwise this is just a whole bunch of pointer chasing.
Definition at line 318 of file simple.c.
References ASSERT, ASSERT3, and RC_NORM.
Referenced by MatSymmPartNormInf().
00320 { 00321 Mat_SeqAIJ *aij; 00322 int *adx,*bdx,*aii,*bii,*aptr,*bptr,ma,na,i,nun; PetscScalar *va,*vb; 00323 PetscReal snorm,anorm,*sronorm,*aronorm; 00324 PetscErrorCode ierr; 00325 PetscFunctionBegin; 00326 00327 aij = (Mat_SeqAIJ *) A->data; 00328 ierr = MatGetLocalSize(A,&ma,&na); CHKERRQ(ierr); 00329 00330 /* two arrays of pointers into the matrix */ 00331 aii = bii = aij->i; 00332 adx = bdx = aij->j; 00333 va = vb = aij->a; 00334 ierr = PetscMalloc(ma*sizeof(int),&aptr); CHKERRQ(ierr); 00335 ierr = PetscMalloc(ma*sizeof(int),&bptr); CHKERRQ(ierr); 00336 for (i=0; i<ma; i++) aptr[i] = bptr[i] = aii[i]; 00337 00338 /* temp storage for the inf norm */ 00339 ierr = PetscMalloc(ma*sizeof(PetscReal),&sronorm); CHKERRQ(ierr); 00340 ierr = PetscMalloc(ma*sizeof(PetscReal),&aronorm); CHKERRQ(ierr); 00341 ierr = PetscMemzero(sronorm,ma*sizeof(PetscReal)); CHKERRQ(ierr); 00342 ierr = PetscMemzero(aronorm,ma*sizeof(PetscReal)); CHKERRQ(ierr); 00343 00344 /* update row & frobenius norms with a(i,j)=v and a(j,i)=w */ 00345 #define RC_NORM(i,j,v,w) \ 00346 { \ 00347 PetscReal sum = PetscAbsScalar(v+w)/2, dif = PetscAbsScalar(v-w)/2; \ 00348 ASSERT3(!(i==j && v!=w),"same elements on diagonal",i,j); \ 00349 sronorm[i] += sum; if (i!=j) sronorm[j] += sum; \ 00350 aronorm[i] += dif; if (i!=j) aronorm[j] += dif; \ 00351 if (i==j) snorm += sum*sum; \ 00352 else {snorm += 2*sum*sum; anorm += 2*dif*dif;} \ 00353 } 00354 #define ASSERT3(v,s,i,j) \ 00355 if (!(v)) SETERRQ3(1,"Violation of <%s> @ (%d,%d): %d",s,i,j); 00356 #define ASSERT(v,s,i,j,t) \ 00357 if (!(v)) SETERRQ4(1,"Violation of <%s> @ (%d,%d): %d",s,i,j,t); 00358 00359 snorm = anorm = 0; nun = 0; 00360 for (i=0; i<ma; i++) { 00361 int idc,idr; PetscScalar vc,vr; 00362 /* get the a-pointer into the upper triangle */ 00363 while (aptr[i]<aii[i+1] && adx[aptr[i]]<i) aptr[i]++; 00364 ASSERT3(aptr[i]==aii[i+1]||adx[aptr[i]]>=i,"A in upper0",i,adx[aptr[i]]); 00365 /* incorporate L elements that are not structurally symmetric */ 00366 while (bptr[i]<bii[i+1] && bdx[bptr[i]]<i) { 00367 vr = PetscAbsScalar(vb[bptr[i]]); nun++; 00368 RC_NORM(i,bdx[bptr[i]],vr,0); 00369 bptr[i]++; 00370 } 00371 ASSERT3(bptr[i]==bii[i+1] || bdx[bptr[i]]>=i,"B part done",i,bdx[bptr[i]]); 00372 while (aptr[i]<aii[i+1]) { 00373 /* in row i we are at element j = idc @ adx( aptr[i] ); 00374 vc = A( i, j ) */ 00375 idc = adx[aptr[i]]; vc = va[aptr[i]]; 00376 ASSERT(idc>=i,"A in upper",i,idc,i-idc); 00377 /* in row j=idc look at elements i' = idr @ bdx(bptr[idc]); vr = B(j,i') 00378 first catch up with structurally unsymmetric elements */ 00379 idr = bdx[bptr[idc]]; vr = vb[bptr[idc]]; 00380 while (idr<i) { 00381 RC_NORM(idc,idr,vr,0); 00382 bptr[idc]++; idr = bdx[bptr[idc]]; vr = vb[bptr[idc]]; nun++; 00383 } 00384 /* now see if (i,j) == (j,i') */ 00385 if (idr==i) { 00386 /* structurally symmetric; process these elements of L and U */ 00387 RC_NORM(i,idc,vc,vr); 00388 bptr[idc]++; 00389 } else { 00390 /* element (j,i') in L corresponds to a later row in U, 00391 so we have a structurally unsymmetric element in U */ 00392 nun++; RC_NORM(i,idc,vc,0); 00393 } 00394 aptr[i]++; 00395 } 00396 ASSERT(aptr[i]==aii[i+1],"Aptr has traversed row i",i+1,adx[aptr[i]],aii[i+1]); 00397 } 00398 00399 ierr = PetscFree(aptr); CHKERRQ(ierr); 00400 ierr = PetscFree(bptr); CHKERRQ(ierr); 00401 00402 snorm = sqrt(snorm); anorm = sqrt(anorm); 00403 00404 for (i=1; i<ma; i++) { 00405 if (sronorm[i]>sronorm[0]) sronorm[0] = sronorm[i]; 00406 if (aronorm[i]>aronorm[0]) aronorm[0] = aronorm[i]; 00407 } 00408 00409 *sn_r = snorm; *an_r = anorm; *srn_r = sronorm[0]; *arn_r = aronorm[0]; 00410 *nun_r = nun; 00411 *done = PETSC_TRUE; 00412 00413 ierr = PetscFree(sronorm); CHKERRQ(ierr); 00414 ierr = PetscFree(aronorm); CHKERRQ(ierr); 00415 00416 PetscFunctionReturn(0); 00417 }
static PetscErrorCode norm1 | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute the 1-norm of a matrix.
Definition at line 159 of file simple.c.
References GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSimpleModules().
00160 { 00161 Mat A = (Mat)prob; 00162 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00163 PetscFunctionBegin; 00164 00165 ierr = GetDataID("simple","norm1",&id,&has); CHKERRQ(ierr); 00166 HASTOEXIST(has); 00167 ierr = PetscObjectComposedDataGetReal 00168 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00169 if (*flg) rv->r = v; 00170 else { 00171 ierr = MatNorm(A,NORM_1,&v); CHKERRQ(ierr); 00172 ierr = PetscObjectComposedDataSetReal 00173 ((PetscObject)A,id,v); CHKERRQ(ierr); 00174 *flg = PETSC_TRUE; 00175 rv->r = v; 00176 } 00177 00178 *flg = PETSC_TRUE; 00179 PetscFunctionReturn(0); 00180 }
static PetscErrorCode normF | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute the Frobenius norm of a matrix.
Definition at line 186 of file simple.c.
References GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSimpleModules().
00187 { 00188 Mat A = (Mat)prob; 00189 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00190 PetscFunctionBegin; 00191 00192 ierr = GetDataID("simple","normF",&id,&has); CHKERRQ(ierr); 00193 HASTOEXIST(has); 00194 ierr = PetscObjectComposedDataGetReal 00195 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00196 if (*flg) rv->r = v; 00197 else { 00198 ierr = MatNorm(A,NORM_FROBENIUS,&v); CHKERRQ(ierr); 00199 ierr = PetscObjectComposedDataSetReal 00200 ((PetscObject)A,id,v); CHKERRQ(ierr); 00201 *flg = PETSC_TRUE; 00202 rv->r = v; 00203 } 00204 00205 *flg = PETSC_TRUE; 00206 PetscFunctionReturn(0); 00207 }
static PetscErrorCode normInf | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute the infinity-norm of a matrix.
Definition at line 213 of file simple.c.
References GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSimpleModules().
00214 { 00215 Mat A = (Mat)prob; 00216 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00217 PetscFunctionBegin; 00218 00219 ierr = GetDataID("simple","normInf",&id,&has); CHKERRQ(ierr); 00220 HASTOEXIST(has); 00221 ierr = PetscObjectComposedDataGetReal 00222 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00223 if (*flg) rv->r = v; 00224 else { 00225 ierr = MatNorm(A,NORM_INFINITY,&v); CHKERRQ(ierr); 00226 ierr = PetscObjectComposedDataSetReal 00227 ((PetscObject)A,id,v); CHKERRQ(ierr); 00228 *flg = PETSC_TRUE; 00229 rv->r = v; 00230 } 00231 00232 *flg = PETSC_TRUE; 00233 PetscFunctionReturn(0); 00234 }
static PetscErrorCode NUnstruct | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute and store the number of structurally unsymmetric elements of a matrix. The actual computation is done in MatSymmPartNormInf(); see that routine for remarks on parallelism.
The "n-struct-unsymm" feature is declared to be in category "structure", rather than simple.
Definition at line 617 of file simple.c.
References GetDataID(), HASTOEXIST, AnalysisItem::i, id, and MatSymmPartNormInf().
Referenced by RegisterSimpleModules().
00618 { 00619 Mat A = (Mat)prob; 00620 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00621 PetscFunctionBegin; 00622 ierr = GetDataID("structure","n-struct-unsymm",&id,&has); CHKERRQ(ierr); 00623 HASTOEXIST(has); 00624 ierr = PetscObjectComposedDataGetInt 00625 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00626 if (*flg) rv->i = v; 00627 else { 00628 ierr = MatSymmPartNormInf(A,flg); CHKERRQ(ierr); 00629 if (*flg) { 00630 ierr = PetscObjectComposedDataGetInt 00631 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00632 if (*flg) rv->i = v; 00633 } 00634 } 00635 00636 PetscFunctionReturn(0); 00637 }
PetscErrorCode RegisterSimpleModules | ( | void | ) |
Declare norm-like modules that can be performed with simple calculations.
Definition at line 644 of file simple.c.
References ANALYSISDOUBLE, ANALYSISINTEGER, DiagonalDominance(), norm1(), normF(), normInf(), NUnstruct(), RegisterModule(), SymmetryANorm(), SymmetryFANorm(), SymmetryFSNorm(), SymmetrySNorm(), Trace(), and TraceAbs().
Referenced by AnaModRegisterSalsaModules(), and AnaModRegisterStandardModules().
00645 { 00646 PetscErrorCode ierr; 00647 PetscFunctionBegin; 00648 00649 ierr = RegisterModule 00650 ("simple","trace",ANALYSISDOUBLE,&Trace); CHKERRQ(ierr); 00651 ierr = RegisterModule 00652 ("simple","trace-abs",ANALYSISDOUBLE,&TraceAbs); CHKERRQ(ierr); 00653 00654 ierr = RegisterModule 00655 ("simple","norm1",ANALYSISDOUBLE,&norm1); CHKERRQ(ierr); 00656 ierr = RegisterModule 00657 ("simple","normInf",ANALYSISDOUBLE,&normInf); CHKERRQ(ierr); 00658 ierr = RegisterModule 00659 ("simple","normF",ANALYSISDOUBLE,&normF); CHKERRQ(ierr); 00660 00661 ierr = RegisterModule 00662 ("simple","diagonal-dominance",ANALYSISDOUBLE,&DiagonalDominance); CHKERRQ(ierr); 00663 00664 ierr = RegisterModule 00665 ("simple","symmetry-snorm",ANALYSISDOUBLE,&SymmetrySNorm); CHKERRQ(ierr); 00666 ierr = RegisterModule 00667 ("simple","symmetry-anorm",ANALYSISDOUBLE,&SymmetryANorm); CHKERRQ(ierr); 00668 ierr = RegisterModule 00669 ("simple","symmetry-fsnorm",ANALYSISDOUBLE,&SymmetryFSNorm); CHKERRQ(ierr); 00670 ierr = RegisterModule 00671 ("simple","symmetry-fanorm",ANALYSISDOUBLE,&SymmetryFANorm); CHKERRQ(ierr); 00672 00673 ierr = RegisterModule 00674 ("structure","n-struct-unsymm",ANALYSISINTEGER,&NUnstruct); CHKERRQ(ierr); 00675 00676 PetscFunctionReturn(0); 00677 }
static PetscErrorCode SymmetryANorm | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute and store the infinity norm of the antisymmetric part of the matrix. The actual computation is done in MatSymmPartNormInf(); see that routine for remarks on parallelism.
Definition at line 527 of file simple.c.
References GetDataID(), HASTOEXIST, id, MatSymmPartNormInf(), and AnalysisItem::r.
Referenced by RegisterSimpleModules().
00528 { 00529 Mat A = (Mat)prob; 00530 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00531 PetscFunctionBegin; 00532 ierr = GetDataID("simple","symmetry-anorm",&id,&has); CHKERRQ(ierr); 00533 HASTOEXIST(has); 00534 ierr = PetscObjectComposedDataGetReal 00535 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00536 if (*flg) rv->r = v; 00537 else { 00538 ierr = MatSymmPartNormInf(A,flg); CHKERRQ(ierr); 00539 if (*flg) { 00540 ierr = PetscObjectComposedDataGetReal 00541 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00542 if (*flg) rv->r = v; 00543 } 00544 } 00545 00546 PetscFunctionReturn(0); 00547 }
static PetscErrorCode SymmetryFANorm | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute and store the Frobenius norm of the antisymmetric part of the matrix. The actual computation is done in MatSymmPartNormInf(); see that routine for remarks on parallelism.
Definition at line 585 of file simple.c.
References GetDataID(), HASTOEXIST, id, MatSymmPartNormInf(), and AnalysisItem::r.
Referenced by RegisterSimpleModules().
00586 { 00587 Mat A = (Mat)prob; 00588 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00589 PetscFunctionBegin; 00590 ierr = GetDataID("simple","symmetry-fanorm",&id,&has); CHKERRQ(ierr); 00591 HASTOEXIST(has); 00592 ierr = PetscObjectComposedDataGetReal 00593 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00594 if (*flg) rv->r = v; 00595 else { 00596 ierr = MatSymmPartNormInf(A,flg); CHKERRQ(ierr); 00597 if (*flg) { 00598 ierr = PetscObjectComposedDataGetReal 00599 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00600 if (*flg) rv->r = v; 00601 } 00602 } 00603 00604 PetscFunctionReturn(0); 00605 }
static PetscErrorCode SymmetryFSNorm | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute and store the Frobenius norm of the symmetric part of the matrix. The actual computation is done in MatSymmPartNormInf(); see that routine for remarks on parallelism.
Definition at line 556 of file simple.c.
References GetDataID(), HASTOEXIST, id, MatSymmPartNormInf(), and AnalysisItem::r.
Referenced by RegisterSimpleModules().
00557 { 00558 Mat A = (Mat)prob; 00559 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00560 PetscFunctionBegin; 00561 ierr = GetDataID("simple","symmetry-fsnorm",&id,&has); CHKERRQ(ierr); 00562 HASTOEXIST(has); 00563 ierr = PetscObjectComposedDataGetReal 00564 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00565 if (*flg) rv->r = v; 00566 else { 00567 ierr = MatSymmPartNormInf(A,flg); CHKERRQ(ierr); 00568 if (*flg) { 00569 ierr = PetscObjectComposedDataGetReal 00570 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00571 if (*flg) rv->r = v; 00572 } 00573 } 00574 00575 PetscFunctionReturn(0); 00576 }
static PetscErrorCode SymmetrySNorm | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute and store the infinity norm of the symmetric part of the matrix. The actual computation is done in MatSymmPartNormInf(); see that routine for remarks on parallelism.
Definition at line 498 of file simple.c.
References GetDataID(), HASTOEXIST, id, MatSymmPartNormInf(), and AnalysisItem::r.
Referenced by RegisterSimpleModules().
00499 { 00500 Mat A = (Mat)prob; 00501 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00502 PetscFunctionBegin; 00503 ierr = GetDataID("simple","symmetry-snorm",&id,&has); CHKERRQ(ierr); 00504 HASTOEXIST(has); 00505 ierr = PetscObjectComposedDataGetReal 00506 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00507 if (*flg) rv->r = v; 00508 else { 00509 ierr = MatSymmPartNormInf(A,flg); CHKERRQ(ierr); 00510 if (*flg) { 00511 ierr = PetscObjectComposedDataGetReal 00512 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00513 if (*flg) rv->r = v; 00514 } 00515 } 00516 00517 PetscFunctionReturn(0); 00518 }
static PetscErrorCode Trace | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute and store the trace of a matrix; the computation is done in computetrace().
Definition at line 106 of file simple.c.
References computetrace(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSimpleModules().
00107 { 00108 Mat A = (Mat)prob; 00109 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00110 PetscFunctionBegin; 00111 ierr = GetDataID("simple","trace",&id,&has); CHKERRQ(ierr); 00112 HASTOEXIST(has); 00113 ierr = PetscObjectComposedDataGetScalar 00114 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00115 if (*flg) { 00116 rv->r = v; 00117 } else { 00118 ierr = computetrace(A); CHKERRQ(ierr); 00119 ierr = PetscObjectComposedDataGetScalar 00120 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00121 if (*flg) { 00122 rv->r = v; 00123 } else SETERRQ(1,"Should have computed trace by now"); 00124 } 00125 00126 PetscFunctionReturn(0); 00127 }
static PetscErrorCode TraceAbs | ( | AnaModNumericalProblem | prob, | |
AnalysisItem * | rv, | |||
int * | lv, | |||
PetscTruth * | flg | |||
) | [static] |
Compute and store the absolute trace of a matrix; the computation is done in computetrace().
Definition at line 135 of file simple.c.
References computetrace(), GetDataID(), HASTOEXIST, id, and AnalysisItem::r.
Referenced by RegisterSimpleModules().
00136 { 00137 Mat A = (Mat)prob; 00138 PetscReal v=0; PetscTruth has; int id; PetscErrorCode ierr; 00139 PetscFunctionBegin; 00140 ierr = GetDataID("simple","trace-abs",&id,&has); CHKERRQ(ierr); 00141 HASTOEXIST(has); 00142 ierr = PetscObjectComposedDataGetReal 00143 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00144 if (*flg) rv->r = v; 00145 else { 00146 ierr = computetrace(A); CHKERRQ(ierr); 00147 ierr = PetscObjectComposedDataGetReal 00148 ((PetscObject)A,id,v,*flg); CHKERRQ(ierr); 00149 if (*flg) rv->r = v; 00150 } 00151 00152 PetscFunctionReturn(0); 00153 }