az_dd_overlap.c

并行解法器,功能强大

    (*data_org)[AZ_matrix_type] = AZ_VBR_MATRIX;

    /* move the information to the back of the send_ptr */

    for (i = j - 1; i >= 0; i--) send_ptr[total_to_be_sent+i-j] = send_ptr[i];
    tj = total_to_be_sent - j;
  }
  else (*data_org)[AZ_matrix_type] = AZ_MSR_MATRIX;

  total_to_be_sent = j;

  /* replace global indices by local indices */

  ii = 0;
  firstone  = AZ_send_length;
  current   = firstone;

  while ((current-firstone < num_recv_neighbors) && ((*data_org)[current] == 0))
    current++;

  bins = (int *) BV_ALLOC((N_update / 4 + 10)*sizeof(int));
  if (bins == NULL) {
    (void) fprintf(stderr, "%sERROR: Not enough dynamic space.\n", yo);
    exit(-1);
  }

  AZ_init_quick_find(update, N_update, &shift, bins);
  for (i = tj; i < total_to_be_sent + tj; i++) {

    /*
     * j = AZ_find_index(send_ptr[i], update, N_update);
     */

    j = AZ_quick_find(send_ptr[i], update, N_update, shift, bins);
    if (j == -1) {
      (void) fprintf(stderr, "%sERROR: A point (%d) requested from processor "
                     "%d\nwas not found among the update ele of proc %d\n", yo,
                     send_ptr[i], (*data_org)[AZ_neighbors+current-firstone],
                     proc);
      exit(-1);
    }

    send_ptr[ii] = j;
    ii +=  1;
  }

  BV_FREE((char *) bins);
  BV_FREE((char *) recv_list);
  BV_FREE((char *) new_external);
  BV_FREE((char *) new_extern_proc);

  (*data_org)[AZ_N_neigh] = num_send_neighbors;

} /* PAZ_set_message_info */


/**************************************************************/
/**************************************************************/
/**************************************************************/
char *AZ_allocate_or_free(void *ptr, unsigned int input_size, int action)
{
/*
 * Simple memory allocated that tries to allocate space off the
 * back end of bindx and val. The general scheme is as follows:
 *   1) we give this routine the length of bindx and val.
 *   2) we give this routine an estimated value of how much
 *      space we will use in bindx and val
 *   3) we first try to allocate space from bindx or val if there
 *      is additional memory beyond the estimated requirements.
 *   4) if 3) is unsuccessful, we try and malloc space.
 *   5) if 4) is unsuccessful, we try to allocate space
 *      from bindx and val that might violate the estimates.
*/

    long int           size;
    static long int    b_estimated_need, v_estimated_need;
    static long int    v_freelist, b_freelist;
    static long int    b_lastused, v_lastused;
    static double      *v,*b;
    static long int    v_end, b_end;
    static long int    v_smallest_free, b_smallest_free;

    double *dptr;
    char *t_ptr;
    long int msize, where; 
    long int current, prev;

    size = (long int) input_size;

    if (action == ALLOCATE) {
       msize = size;
       size += sizeof(double) - size%sizeof(double);
       size = size/sizeof(double);
       size = size + 1;

       /* check the v free list */

       current = v_freelist;
       prev    = -1;
       while (current != -1) {
          if ( (long int) v[current-1] >= size) break;
          prev = current;
          current = (long int) v[current];
       }
       if ( current != -1) {
          if (prev == -1) v_freelist = (long int) v[current];
          else v[prev] = v[current];
          return( (char *) &(v[current])   );
       }

       /* check the b free list */

       current = b_freelist;
       prev    = -1;
       while (current != -1) {
          if ( (long int) b[current-1] >= size) break;
          prev = current;
          current = (long int) b[current];
       }
       if ( current != -1) {
          if (prev == -1) b_freelist = (long int) b[current];
          else b[prev] = b[current];
          return( (char *) &(b[current])   );
       }

       /* check (conservative) if there is  */
       /* additional space in either v or b */
     
       if (v_smallest_free - v_estimated_need >= b_smallest_free -
                                                   b_estimated_need) {
          if (v_smallest_free - size > v_estimated_need) {
             v_smallest_free    -= size;
             v[v_smallest_free]  = (double) size;
             return( (char *) &(v[v_smallest_free+1]) );
          }
       }
       else {
          if (b_smallest_free - size > b_estimated_need ) {
             b_smallest_free    -= size;
             b[b_smallest_free]  = (double) size;
             return( (char *) &(b[b_smallest_free+1]) );
          }
       }

       /* check if there is malloc space */
    
       t_ptr = (char *) AZ_allocate((unsigned int) msize);
       if (t_ptr != NULL) return(t_ptr);


       /* check if there is any space  */
       /* in either v or b             */
     
       if (v_smallest_free - v_lastused >= b_smallest_free - b_lastused){
          if (v_smallest_free - size > v_lastused ) {
             v_smallest_free -= size;
             v[v_smallest_free] = size;
             return( (char *) &(v[v_smallest_free+1]) );
          }
       }
       else {
          if (b_smallest_free - size > b_lastused) {
             b_smallest_free -= size;
             b[b_smallest_free] = size;
             return( (char *) &(b[b_smallest_free+1]) );
          }
       }
       return(NULL);

    }
    else if (action == FREE_IT) {

       dptr = (double *) ptr;

       where = OUT_OF_BOUNDS;
       if ( (dptr >= v) && (dptr < &(v[v_end]) ) ) where = IN_V;

       if ( (dptr >= b) && (dptr < &(b[b_end]) ) ) where = IN_B;

       if (where == OUT_OF_BOUNDS) {
          AZ_free(ptr);
          return(NULL);
       }

       dptr--;
       size = (int) dptr[0];

       if (where == IN_V) {
          dptr[1] = (int) v_freelist;
          dptr++;
          v_freelist = ((long int) dptr - (long int) v)/sizeof(double);
       }
       else {
          dptr[1] = (long int) b_freelist;
          dptr++;
          b_freelist = ((long int) dptr - (long int) b)/sizeof(double);
       }
    }
    else if (action == V_SET) {
        v          = (double *) ptr;
        v_freelist = -1;
        v_end      = size;
        v_smallest_free = v_end;
    }
    else if (action == B_SET) {
        b          = (double *) ptr;
        b_freelist = -1;
        size *= sizeof(int);
        size -= size%sizeof(double);
        size /= sizeof(double);
        b_end      = size;
        b_smallest_free = b_end;
    }
    else if (action == LASTUSED_SET) {
        size *= sizeof(int);
        if (size%sizeof(double) != 0) 
           size += (sizeof(double) - size%sizeof(double));
        
        b_lastused = size/sizeof(double);
        v_lastused = size/sizeof(int);
     
        if (b_lastused > b_smallest_free) {
          printf("Error: Out of space due to poor estimate of memory needed\n");
          printf("       for overlapping.\n");
          exit(1);
        }
        if (v_lastused > v_smallest_free) {
          printf("Error: Out of space due to poor estimate of memory needed\n");
          printf("       for overlapping.\n");
        }
    }
    else if (action == ESTIMATED_SET) {
        size *= sizeof(int);
        if (size%sizeof(double) != 0) 
           size += (sizeof(double) - size%sizeof(double));

        b_estimated_need = size/sizeof(double);
        v_estimated_need = size/sizeof(int);
    }
    else if (action == -43) {
       printf("v_list: ");
       current = v_freelist;
       prev    = -1;
       while (current != -1) {
          printf("(%d, %d) ",(int) current,(int) v[current-1]);
          prev = current;
          current = (long int) v[current];
       }
       printf("\n");
       printf("b_list: ");
       current = b_freelist;
       prev    = -1;
       while (current != -1) {
          printf("(%d, %d) ",(int) current,(int) b[current-1]);
          prev = current;
          current = (long int) b[current];
       }
       printf("\n\n");

    }
    else {
        printf("Error: unknown option (%d) for allocate_or_free\n",action);
        exit(1);
    }
    return((char *) NULL);

}

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

/* This subroutine uses the external node list input (ext_nodelist) to
 * construct the send-to-processor information (send_proc, send_rowcnt,
 * and send_rownum.  A sorted local node list and the auxilliary list
 * (the one obtained in AZ_sort) are also required
 *
 * Input :
 *
 *    N_ext             : the length of the external node list
 *    externs           : the actual external node list. This list must
 *                        be sorted.
 *    N_rows            : the number of local rows
 *    Rownum            : the indices of the local rows (MUST BE SORTED)
 *    proc_config       : parallel machine information
 *
 * Output :
 *
 *    send_proc         : if send_proc[i] = 1, then there is something to
 *                        send to processor i. (Note : this array should be
 *                        allocated as nprocs integers before entering)
 *    send_rowcnt       : send_rowcnt[i] holds the number of rows to be
 *                        sent to processor i. (Again, this array should be
 *                        allocated as nprocs integers before entering)
 *    send_rownum       : send_rownum[i] holds the actual row numbers (local)
 *                        of the rows to be sent to processor i. (Again, this
 *                        array should be allocated as nprocs int pointers
 *                        before entering)
 */



#define proc(i)  (externs[i]/max_per_proc)
void AZ_setup_sendlist(int N_ext, int externs[], int send_proc[],
		       int send_rowcnt[], int *send_rownum[], 
                       int proc_config[], int max_per_proc, int N_rows,
                       int Rownum[])
{
   int node, nprocs;
   int N_send, N_recv;                     /* Number of messages processor */
                                           /* must send and recv           */
   MPI_AZRequest request[AZ_MAX_NEIGHBORS];  /* Message handle               */
   int type,st;
   int length;                             /* message send length          */
   int *temp1,*temp2;
   int i, j, start, processor, *ttt;


   nprocs    = proc_config[AZ_N_procs];
   node      = proc_config[AZ_node];
   temp1     = send_proc;
   temp2     = send_rowcnt;

   /* determine how many messages I need to send  */
   /* to update external rows of other processors */

   for (i = 0 ; i < nprocs; i++) temp1[i] = 0;
   for (i = 0 ; i < N_ext ; i++) temp1[proc(i)] = 1;
   AZ_gsum_vec_int(temp1,temp2, nprocs,proc_config);
   N_send = temp1[node];
   
   /* exchange information so that each processor knows: */
   /*    1) processor to which it must send rows         */
   /*    2) number of rows to be sent to each processor  */

   type            =AZ_sys_msg_type;
   AZ_sys_msg_type =(AZ_sys_msg_type+1-AZ_MSG_TYPE) % AZ_NUM_MSGS + AZ_MSG_TYPE;

   for (i = 0 ; i < N_send; i++ ) {
      send_proc[i] = -1;
      (void) mdwrap_iread((void *) &(send_rowcnt[i]), sizeof(int), 
                           &(send_proc[i]), &type, &(request[i]));
   }

   N_recv = 0;
   length = 1;
   for (i = 1 ; i < N_ext ; i++ ) {
      if (proc(i) != proc(i-1)) {
         (void) mdwrap_write((void *) &length, sizeof(int),proc(i-1),type,&st);
         length = 0;
         N_recv++;
      }
      length++;
   }
   if (N_ext != 0) {
      (void) mdwrap_write((void *)&length, sizeof(int),proc(N_ext-1),type,&st);
      N_recv++;
   }

   for (i = 0 ; i < N_send; i++ ) {
       (void) mdwrap_wait((void *) &(send_rowcnt[i]), sizeof(int), 
                           &(send_proc[i]), &type, &st, &(request[i]));
   }

   AZ_sort(send_proc,N_send,send_rowcnt,NULL);   /* sort the processor info */

   /* now receive the actual row numbers */

   type            =AZ_sys_msg_type;
   AZ_sys_msg_type =(AZ_sys_msg_type+1-AZ_MSG_TYPE) % AZ_NUM_MSGS + AZ_MSG_TYPE;


   for (i = 0 ; i < N_send; i++ ) {
      send_rownum[i] = (int *) BV_ALLOC((send_rowcnt[i]+1) * sizeof(int));
      (void) mdwrap_iread((void *) send_rownum[i], send_rowcnt[i]*sizeof(int), 
                           &(send_proc[i]), &type, &(request[i]));
   }
   start = 0;
   length = 1;
   for (i = 1 ; i < N_ext ; i++ ) {
      if (proc(i) != proc(i-1)) {
         (void) mdwrap_write((void *) &(externs[start]), 
                               length*sizeof(int),proc(i-1),type,&st);
         start += length;
         length = 0;
      }
      length++;
   }
   if (N_ext != 0) {
      (void) mdwrap_write((void *) &(externs[start]), length*sizeof(int),
                           proc(N_ext-1),type,&st);
   }

   for (i = 0 ; i < N_send; i++ ) {
       (void) mdwrap_wait((void *) send_rownum[i], send_rowcnt[i]*sizeof(int), 
                           &(send_proc[i]), &type, &st, &(request[i]));
   }


   for (i = N_send ; i < nprocs ; i++ ) {
      send_proc[i] = 0;
      send_rowcnt[i]  = 0;
   }

   for (i = N_send-1 ; i >= 0 ; i-- ) {
      processor               = send_proc[i];
      j                       = send_rowcnt[i];
      ttt                     = send_rownum[i];
      if (i < processor) {send_proc[i] = 0;       /* The 'if' statement is  */
                          send_rowcnt[i] = 0;     /* needed.                */
                          send_rownum[i] = NULL;}
      send_rowcnt[processor]  = send_rowcnt[i];
      send_rowcnt[processor]  = j;
      send_rownum[processor]  = ttt;

      for (j = 0 ; j < send_rowcnt[processor] ; j++ ) 
         send_rownum[processor][j] =PAZ_sorted_search(send_rownum[processor][j],
                                                      N_rows,Rownum);

      send_proc[processor]    = 1;
   }
}