az_util.c

并行解法器,功能强大

                        message_recv_length, message_send_length);
    write_t += AZ_second() - start_t;

    start_t = AZ_second();

    temp1                      = data_org[AZ_N_bord_blks];
    temp2                      = data_org[AZ_N_border];
    data_org[AZ_N_bord_blks] = 0;
    data_org[AZ_N_border] = 0;

    Amat->matvec(x, r, Amat, proc_config);


    data_org[AZ_N_bord_blks] = temp1;
    data_org[AZ_N_border] = temp2;

    sparax_overlap_internal_comp_t += AZ_second() - start_t;

    start_t = AZ_second();
    AZ_read_local_info(data_org, message_recv_add, message_recv_length);
    read_t += AZ_second() - start_t;

    start_t = AZ_second();

    /* compute boundary portion of the sparse matrix - vector product */

    bpntr_index = bpntr[Num_Internal_Blks];

    temp1                      = data_org[AZ_N_int_blks];
    temp2                      = data_org[AZ_N_internal];
    data_org[AZ_N_int_blks]    = data_org[AZ_N_bord_blks];
    data_org[AZ_N_internal]    = data_org[AZ_N_border];
    data_org[AZ_N_bord_blks]   = 0;
    data_org[AZ_N_border]      = 0;

/*
    Amat2.rpntr = &(Amat->rpntr[Num_Internal_Blks]);
    Amat2.cpntr =   Amat->cpntr;
    Amat2.bpntr = &(Amat->bpntr[Num_Internal_Blks]);
    Amat2.bindx = &(Amat->bindx[bpntr_index]);
    Amat2.matvec      = Amat->matvec;
    Amat2.matrix_type = Amat->matrix_type;
    Amat2.data_org = &(Amat->data_org);
    Amat2.val   = &(Amat->val[bpntr_index]);
    Amat2.indx  = &(Amat->indx[bpntr_index]);
*/

    Amat->matvec(x, r, Amat, proc_config);

/*
                   &bindx[bpntr_index], &rpntr[Num_Internal_Blks], cpntr,
                   &bpntr[Num_Internal_Blks], x, &y[rpntr[Num_Internal_Blks]],
                   0, data_org);
*/

    data_org[AZ_N_bord_blks] = data_org[AZ_N_int_blks];
    data_org[AZ_N_border]    = data_org[AZ_N_internal];
    data_org[AZ_N_int_blks]  = temp1;
    data_org[AZ_N_internal]  = temp2;

    sparax_overlap_border_comp_t += AZ_second() - start_t;
  }

  overall_t = AZ_second() - start_overall_t;

  for (i = 1; i <= ntimes; i++) {

    /* time the reads in the nonoverlapped case */

    AZ_gather_mesg_info(x, data_org, message_recv_add, message_send_add,
                        message_recv_length, message_send_length);

    AZ_write_local_info(data_org, message_recv_add, message_send_add,
                        message_recv_length, message_send_length);

    start_t = AZ_second();
    AZ_read_local_info(data_org, message_recv_add, message_recv_length);
    read_nonoverlap_t += AZ_second() - start_t;
  }

  read_t_max            = AZ_gmax_double(read_t, proc_config);
  read_nonoverlap_t_max = AZ_gmax_double(read_nonoverlap_t, proc_config);
  time_overlap_max      = AZ_gmax_double(time, proc_config);

  /* dot */

  if (iout > 0 && Proc == 0) {
    (void) printf("time the individual routines for ddot\n");
  }

  start_t = AZ_second();
  for (i = 1; i <= ntimes; i++) Djunk = AZ_gdot(N, x, y, proc_config);
  total_t = AZ_second() - start_t;

  start_t = AZ_second();
  for (i = 1; i <= ntimes; i++) AZ_gsum_double(Djunk, proc_config);

  dot_comm_t = AZ_gmax_double(AZ_second() - start_t, proc_config);
  dot_comp_t = AZ_gmax_double(total_t - dot_comm_t, proc_config);

  /* daxpy */

  if (iout > 0 && Proc == 0) {
    (void) printf("time the individual routines for daxpy\n");
  }
  start_t = AZ_second();
  for (i = 1; i <= ntimes; i++) {
    Djunk = (double) i;
    daxpy_(&N, &Djunk, x, &one, y, &one);
  }
  daxpy_comp_t = AZ_gmax_double(AZ_second() - start_t, proc_config);

  if (Proc == 0) {              /* calculate and print results */
    (void) printf("\n********** sparax ***********\n\n");
    (void) printf("nonoverlapped\n");
    (void) printf("\t\tmax\t\tavg\t\tmin\n");
  }

  /* dzpax */

  if (iout > 0 && Proc == 0) {
    (void) printf("time the individual routines for dzpax\n");
  }

  z1 = (double *) AZ_allocate(N * sizeof(double));
  z2 = (double *) AZ_allocate(N * sizeof(double));
  z3 = (double *) AZ_allocate(N * sizeof(double));
  for (i = 0; i < N; i++) {
    z2[i] = (double) i;
    z3[i] = (double) 2*i;
  }

  start_t = AZ_second();
  for (i = 1; i <= ntimes; i++) {
    Djunk = 3.14 * (double) i;
    for (j = 0; j < N; j++)
      z1[j] = z2[j] + Djunk*z3[j];
    /*    dzpax_(&N, &Djunk, z3, &one, z2, &one, z1, &one);*/
  }
  dzpax_comp_t = AZ_gmax_double(AZ_second() - start_t, proc_config);

  AZ_free((void *) z1);
  AZ_free((void *) z2);
  AZ_free((void *) z3);

  if (Proc == 0) {              /* calculate and print results */
    (void) printf("\n********** sparax ***********\n\n");
    (void) printf("nonoverlapped\n");
    (void) printf("\t\tmax\t\tavg\t\tmin\n");
  }

  max = AZ_gmax_double(sprx_t, proc_config);
  avg = AZ_gavg_double(sprx_t, proc_config);
  min = AZ_gmin_double(sprx_t, proc_config);
  if (Proc == 0) (void) printf("total_time\t%e\t%e\t%e\n", max, avg, min);

  max = AZ_gmax_double(sprx_comp_t, proc_config);
  avg = AZ_gavg_double(sprx_comp_t, proc_config);
  min = AZ_gmin_double(sprx_comp_t, proc_config);
  if (Proc == 0) (void) printf("comp_time\t%e\t%e\t%e\n", max, avg, min);

  max = AZ_gmax_double(sprx_comm_t, proc_config);
  avg = AZ_gavg_double(sprx_comm_t, proc_config);
  min = AZ_gmin_double(sprx_comm_t, proc_config);
  if (Proc == 0) (void) printf("comm_time\t%e\t%e\t%e\n", max, avg, min);

  sprx_comp_t = AZ_gmax_double(sprx_comp_t, proc_config);
  sprx_t      = AZ_gmax_double(sprx_t, proc_config);

  if (Proc == 0) {
    Mflop             = (double) ntimes * (gnzeros + gnzeros) * 1.0e-6;
    Mflops_comp       = Mflop/sprx_comp_t;
    Mflops_node_comp  = Mflops_comp/(double) Num_Proc;
    Mflops_comm       = Mflop/sprx_t;
    Mflops_node_comm  = Mflops_comm/(double) Num_Proc;

    (void) printf("computation Mflops: %e \n", Mflops_comp);
    (void) printf("computation Mflops per node: %e \n", Mflops_node_comp);
    (void) printf("comp & comm Mflops: %e \n", Mflops_comm);
    (void) printf("comp & comm Mflops per node: %e \n\n", Mflops_node_comm);
  }

  /* statistics */

  if (Proc == 0) {
    (void) printf("overlapped\n\t\tmax\t\tavg\t\tmin\n");
  }

  max = AZ_gmax_double(time, proc_config);
  min = AZ_gmin_double(time, proc_config);
  avg = AZ_gavg_double(time, proc_config);
  if (Proc == 0) (void) printf("total_time\t%e\t%e\t%e\n", max, avg, min);

  max = AZ_gmax_double(gather_t, proc_config);
  min = AZ_gmin_double(gather_t, proc_config);
  avg = AZ_gavg_double(gather_t, proc_config);
  if (Proc == 0) (void) printf("gather_t\t%e\t%e\t%e\n", max, avg, min);

  max = AZ_gmax_double(write_t, proc_config);
  min = AZ_gmin_double(write_t, proc_config);
  avg = AZ_gavg_double(write_t, proc_config);
  if (Proc == 0) (void) printf("write_t \t%e\t%e\t%e\n", max, avg, min);

  max = AZ_gmax_double(sparax_overlap_internal_comp_t, proc_config);
  min = Az_gmin_double(sparax_overlap_internal_comp_t, proc_config);
  avg = AZ_gavg_double(sparax_overlap_internal_comp_t, proc_config);
  if (Proc == 0) (void) printf("internal_t\t%e\t%e\t%e\n", max, avg, min);

  max = AZ_gmax_double(read_t, proc_config);
  min = Az_gmin_double(read_t, proc_config);
  avg = AZ_gavg_double(read_t, proc_config);
  if (Proc == 0) (void) printf("read_t   \t%e\t%e\t%e\n", max, avg, min);

  max = AZ_gmax_double(sparax_overlap_border_comp_t, proc_config);
  min = AZ_gmin_double(sparax_overlap_border_comp_t, proc_config);
  avg = AZ_gavg_double(sparax_overlap_border_comp_t, proc_config);
  if (Proc == 0) (void) printf("border_t\t%e\t%e\t%e\n", max, avg, min);

  if (Proc == 0) {
    Mflops_comm      = Mflop/time_overlap_max;
    Mflops_node_comm = Mflops_comm/(double) Num_Proc;
    (void) printf("comp & comm Mflops: %e \n", Mflops_comm);
    (void) printf("comp & comm Mflops per node: %e \n\n", Mflops_node_comm);

    (void) printf("Ratio of overlapped/nonoverlapped sparax times: %e\n\n",
                  time_overlap_max/sprx_t);
    (void) printf("Ratio of overlapped/nonoverlapped read times: %e\n\n",
                  read_t_max/read_nonoverlap_t_max);
  }

  max = AZ_gmax_double(time, proc_config);
  if (time  ==  max) {
    (void) printf("max time proc\n");
    (void) printf("\nProc:%d\tNum_Neighbors: %d\n", Proc, data_org[AZ_N_neigh]);
    (void) printf("total_time: %e\n", time);
    (void) printf("gather_t  : %e\n", gather_t);
    (void) printf("write_t   : %e\n", write_t);
    (void) printf("internal_t: %e\n", sparax_overlap_internal_comp_t);
    (void) printf("read_t    : %e\n", read_t);
    (void) printf("border_t  : %e\n\n", sparax_overlap_border_comp_t);
  }

  max = AZ_gmax_double(overall_t, proc_config);
  if (overall_t  ==  max) {
    (void) printf("overall max time proc\n");
    (void) printf("\nProc:%d\tNum_Neighbors: %d\n", Proc, data_org[AZ_N_neigh]);
    (void) printf("overall total_time: %e\n", overall_t);
    (void) printf("total_time: %e\n", time);
    (void) printf("gather_t  : %e\n", gather_t);
    (void) printf("write_t   : %e\n", write_t);
    (void) printf("internal_t: %e\n", sparax_overlap_internal_comp_t);
    (void) printf("read_t    : %e\n", read_t);
    (void) printf("border_t  : %e\n\n", sparax_overlap_border_comp_t);
  }
  /*
    if (Proc == 0) {
    (void) printf("cop overlapped\n");
    (void) printf("\t\tmax\t\tavg\t\tmin\n");
    }

    max = AZ_gmax_double(sprx_overlap_cop_t, proc_config);
    min = AZ_gmin_double(sprx_overlap_cop_t, proc_config);
    avg = AZ_gavg_double(sprx_overlap_cop_t, proc_config);
    if (Proc == 0) (void) printf("total_time\t%e\t%e\t%e\n", max, avg, min);

    if (Proc == 0) {
    Mflops_comm = Mflop/max;
    Mflops_node_comm  = Mflops_comm/(double)Num_Proc;
    (void) printf("comp & comm Mflops: %e \n", Mflops_comm);
    (void) printf("comp & comm Mflops per node: %e \n\n", Mflops_node_comm);
    }
    */
  if (Proc == 0) {
    (void) printf("\n********* exchange **********\n\n");
    (void) printf("comm_time = %7.4e sec.\n", sprx_comm_t);

    (void) printf("\n************ dot ************\n\n");
    Mflop = (double) ntimes * 2.0 * ((double) gN * 1.0e-06);

    (void) printf("comp_time = %7.4e sec.\n", dot_comp_t);
    (void) printf("comm_time = %7.4e sec.\n", dot_comm_t);

    Mflops_comp      = Mflop/dot_comp_t;
    Mflops_node_comp = Mflops_comp/(double) Num_Proc;
    Mflops_comm      = Mflop/(dot_comm_t + dot_comp_t);
    Mflops_node_comm = Mflops_comm/(double) Num_Proc;

    (void) printf("computation Mflops: %e \n", Mflops_comp);
    (void) printf("computation Mflops per node: %e \n", Mflops_node_comp);
    (void) printf("comp & comm Mflops: %e \n", Mflops_comm);
    (void) printf("comp & comm Mflops per node: %e \n\n", Mflops_node_comm);

    (void) printf("\n*********** daxpy ***********\n\n");
    (void) printf("comp_time = %7.4e sec.\n", daxpy_comp_t);
    Mflops_comp      = Mflop/daxpy_comp_t;
    Mflops_node_comp = Mflops_comp/(double) Num_Proc;

    (void) printf("computation Mflops: %e \n", Mflops_comp);
    (void) printf("computation Mflops per node: %e \n", Mflops_node_comp);

    (void) printf("\n*********** dzpax ***********\n\n");
    (void) printf("comp_time = %7.4e sec.\n", dzpax_comp_t);
    Mflops_comp  = Mflop/dzpax_comp_t;
    Mflops_node_comp  = Mflops_comp/(double) Num_Proc;

    (void) printf("computation Mflops: %e \n", Mflops_comp);
    (void) printf("computation Mflops per node: %e \n", Mflops_node_comp);
  }

} /* AZ_time_kernals */

#endif

/******************************************************************************/
/******************************************************************************/
/******************************************************************************/

int AZ_find_index(int key, int list[], int length)

/*******************************************************************************

  Find 'key' in 'list' and return the index number.

  Author:          Ray Tuminaro, SNL, 1422
  =======

  Return code:     int, -1 = key not found, i = list[i] = key
  ============

  Parameter list:
  ===============

  key:             Element to be search for in list.

  list:            List to be searched.

  length:          Length of list.

*******************************************************************************/

{

  /* local variables */

  int start, end;
  int mid;

  /**************************** execution begins ******************************/

  if (length == 0) return -1;

  start = 0;
  end   = length - 1;

  while (end - start > 1) {
    mid = (start + end) / 2;
    if (list[mid] < key) start = mid;
    else end = mid;
  }

  if (list[start] == key) return start;
  if (list[end] == key)   return end;
  return -1;

} /* AZ_find_index */

int AZ_exit(int input)
{
#ifdef AZ_MPI
   MPI_Finalize();
#endif
   exit(input);
   return(1);
}

/******************************************************************************/
/******************************************************************************/
/******************************************************************************/

void AZ_free_memory(int label)

/*******************************************************************************

  Clear all internal memory that has been allocated by previous calls to
  AZ_manage_memory(*,*,label,*,*). This memory is typically includes
  preconditioner and scaling information.

  Author:          Ray S. Tuminaro, SNL, 1422
  =======

  Return code:     void
  ============

  Parameter list:
  ===============

  label:           An integer associated with memory that was allocated with
                   AZ_manage_memory(). On ouput, all memory allocated with this
                   integer will be freed.

*******************************************************************************/

{
  (void) AZ_manage_memory((int) NULL, AZ_CLEAR, label, (char *) NULL,