symprintperformance.c

数学算法的实现库。可以实现常见的线性计算。

#ifdef USE_MPILUC
    printf("M-ILUPACK, multilevel structure\n");
#else
    printf("ILUPACK,   multilevel structure\n");
#endif

    fprintf(fo,"%3d|",nlev);

    next=⪯
    nnzU=0;
    tmp=0;
    ierr=0;
    flag=0;

    for (i=1; i<=nlev; i++) {
        // fill-in LU
        printf("level %3d, block size %7d\n",i,next->LU.nr); fflush(stdout);
	l=nnzU;
	if (i<nlev || next->LU.ja!=NULL) {
	   nnzU+=next->LU.ja[next->LU.nr-1]-next->LU.ja[0]+2*next->nB;
	   k=0;
	   for (m=0; m<next->LU.nr; m++) {
	     if (next->LU.ja[next->LU.nr+1+m]>0) 
	        k+=2;
	   } // end for m
	   ierr+=k;
	   flag+=next->LU.nr;
	   printf("           2x2 pivots %5.1f%%\n",(100.0*k)/next->LU.nr); fflush(stdout);
	}
	if (i==nlev) {
	   if (next->LU.ja==NULL) {
	      printf("switched to full matrix processing\n");fflush(STDOUT);
	      tmp=-1;
	      j=next->LU.nr;
	      nnzU+=(j*(j-1))/2;
	   }
	}
	printf("  local fill-in %7d(%5.1lfav)\n",
	       nnzU-l+next->LU.nr,(1.0*(nnzU-l+next->LU.nr))/next->LU.nr);

	if (i<nlev) {
	   // fill-in F
	   nnzU+=next->F.ia[next->F.nr]-1;
	   printf("level %3d->%3d, block size (%7d,%7d)\n",i,i+1,next->LU.nr,next->F.nc);
	   printf("  local fill-in F %7d(%5.1lfav pr)\n",
		  next->F.ia[next->F.nr]-1,(1.0*(next->F.ia[next->F.nr]-1))/next->LU.nr);
	}
	next=next->next;
    }
    printf("\ntotal fill-in sum%8d(%5.1lfav)\n",
	   nnzU+n,(1.0*(nnzU))/n);
    printf("fill-in factor:      %5.1lf\n",(1.0*nnzU)/nz);
    printf("total number of sparse  2x2 pivots %5.1f%%\n",(100.0*ierr)/flag); fflush(stdout);
    ierr=flag=0;

    if (tmp) {
       // nnzU-j*(j+1)/2+n-j memory for sparse data structures
       //                    indices (weight 1/3) and values (weight 2/3)
       // j*(j+1)/2          memory for dense data, no indices (weight 2/3)
       printf("memory usage factor: %5.1lf\n",(1.0*(nnzU-(j*(j+1))/2-j))/nz
	                                     +(j*(j+1))/(3.0*nz));
       fprintf(fo,"%5.1f|",(1.0*(nnzU-(j*(j+1))/2-j))/nz
                          +(j*(j+1))/(3.0*nz));
    }
    else {
       printf("memory usage factor: %5.1lf\n",(1.0*(nnzU))/nz);
       fprintf(fo,"%5.1f|",(1.0*(nnzU))/nz);
    }
    printf("total time: %8.1le [sec]\n",  (double)secnds); 
    printf("            %8.1le [sec]\n\n",(double)ILUPACK_secnds[7]); 
    fprintf(fo,"%7.1le|",(double)secnds);

#ifdef USE_MPILUC
    printf("refined timings for M-ILUPACK multilevel factorization\n"); 
#else
    printf("refined timings for   ILUPACK multilevel factorization\n"); 
#endif
    printf("initial preprocessing:         %8.1le [sec]\n",ILUPACK_secnds[0]); 
    printf("reorderings remaining levels:  %8.1le [sec]\n",ILUPACK_secnds[1]); 
#ifdef USE_MPILUC
    printf("SYMMPILUC(sum over all levels):%8.1le [sec]\n",ILUPACK_secnds[2]); 
#else
    printf("SYMPILUC (sum over all levels):%8.1le [sec]\n",ILUPACK_secnds[2]); 
#endif
    printf("SYMILUC (if used):             %8.1le [sec]\n",ILUPACK_secnds[3]); 
    printf("SPTRF, LAPACK (if used):       %8.1le [sec]\n",ILUPACK_secnds[4]); 
    printf("remaining parts:               %8.1le [sec]\n\n",MAX(0.0,(double)secnds
	                                                   -ILUPACK_secnds[0]
							   -ILUPACK_secnds[1]
							   -ILUPACK_secnds[2]
							   -ILUPACK_secnds[3]
							   -ILUPACK_secnds[4]));

    fflush(STDOUT);

    /*
    next=&PRE;
    for (i=1; i<=nlev; i++) {
      printf("%d,%d\n",next->n,next->nB);
      for (j=0; j<next->n; j++)
	printf("%8d",next->p[j]);
      printf("\n");
      for (j=0; j<next->n; j++)
	printf("%8d",next->invq[j]);
      printf("\n");
        next=next->next;
    }
    fflush(STDOUT);

    printf("multilevel structure\n");
    fflush(STDOUT);
    printf("number of levels: %d\n",nlev);
    fflush(STDOUT);
    next=&PRE;
    for (i=1; i<=nlev; i++) {
      printf("total size %d\n",next->n);
      fflush(STDOUT);
      printf("leading block %d\n",next->nB);
      fflush(STDOUT);

      printf("row permutation\n");
      for (j=0; j<next->n; j++)
	printf("%4d",next->p[j]);
      printf("\n");
      fflush(STDOUT);
      printf("inverse column permutation\n");
      for (j=0; j<next->n; j++)
	printf("%4d",next->invq[j]);
      printf("\n");
      fflush(STDOUT);

      if (nlev==1) {
	printf("row scaling\n");
	for (j=0; j<next->n; j++)
	  printf("%12.4e",next->rowscal[j]);
	printf("\n");
	fflush(STDOUT);
	printf("column scaling\n");
	for (j=0; j<next->n; j++)
	  printf("%12.4e",next->colscal[j]);
	printf("\n");
	fflush(STDOUT);
      }

      printf("(2,1) block E (%d,%d)\n",next->E.nr,next->E.nc);
      fflush(STDOUT);
      for (k=0; k<next->E.nr; k++) {
	printf("%3d: ",k+1);
	for (j=next->E.ia[k]-1; j<next->E.ia[k+1]-1;j++)
	  fprintf(STDOUT,"%12d",next->E.ja[j]);
	printf("\n");
	fflush(STDOUT);
	printf("     ");
	for (j=next->E.ia[k]-1; j<next->E.ia[k+1]-1;j++)
	  fprintf(STDOUT,"%12.4e",next->E.a[j]);
	printf("\n");
	fflush(STDOUT);
      }
      printf("(1,2) block F (%d,%d)\n",next->F.nr,next->F.nc);
      fflush(STDOUT);
      for (k=0; k<next->F.nr; k++) {
	printf("%3d: ",k+1);
	for (j=next->F.ia[k]-1; j<next->F.ia[k+1]-1;j++)
	  fprintf(STDOUT,"%12d",next->F.ja[j]);
	printf("\n");
	fflush(STDOUT);
	printf("     ");
	for (j=next->F.ia[k]-1; j<next->F.ia[k+1]-1;j++)
	  fprintf(STDOUT,"%12.4e",next->F.a[j]);
	printf("\n");
	fflush(STDOUT);
      }

      printf("ILU...\n");
      printf("Diagonal part\n");
      for (k=0; k<next->LU.nr; k++) {
	fprintf(STDOUT,"%12.4e",next->LU.a[k]);
      }
      printf("\n");
      printf("L part\n");
      for (k=0; k<next->LU.nr; k++) {
	printf("col %3d: ",k+1);
	for (j=next->LU.ja[k]-1; j<next->LU.ia[k]-1;j++)
	  fprintf(STDOUT,"%12d",next->LU.ja[j]);
	printf("\n");
	fflush(STDOUT);
	printf("         ");
	for (j=next->LU.ja[k]-1; j<next->LU.ia[k]-1;j++)
	  fprintf(STDOUT,"%12.4e",next->LU.a[j]);
	printf("\n");
	fflush(STDOUT);
      }
      printf("U part\n");
      for (k=0; k<next->LU.nr; k++) {
	printf("row %3d: ",k+1);
	for (j=next->LU.ia[k]-1; j<next->LU.ja[k+1]-1;j++)
	  fprintf(STDOUT,"%12d",next->LU.ja[j]);
	printf("\n");
	fflush(STDOUT);
	printf("         ");
	for (j=next->LU.ia[k]-1; j<next->LU.ja[k+1]-1;j++)
	  fprintf(STDOUT,"%12.4e",next->LU.a[j]);
	printf("\n");
	fflush(STDOUT);
      }

      next=next->next;
    }
    */