az_tools.c

并行解法器,功能强大

/*====================================================================
 * ------------------------
 * | CVS File Information |
 * ------------------------
 *
 * $RCSfile: az_tools.c,v $
 *
 * $Author: tuminaro $
 *
 * $Date: 2000/06/02 16:48:32 $
 *
 * $Revision: 1.54 $
 *
 * $Name:  $
 *====================================================================*/
#ifndef lint
static char rcsid[] = "$Id: az_tools.c,v 1.54 2000/06/02 16:48:32 tuminaro Exp $";
#endif


/*******************************************************************************
 * Copyright 1995, Sandia Corporation.  The United States Government retains a *
 * nonexclusive license in this software as prescribed in AL 88-1 and AL 91-7. *
 * Export of this program may require a license from the United States         *
 * Government.                                                                 *
 ******************************************************************************/

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <malloc.h>
#include <float.h>
#include "az_aztec.h"

extern int AZ_sys_msg_type;
extern int AZ_using_fortran;

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

void AZ_order_ele(int update_index[], int extern_index[], int *internal,
                  int *border, int N_update, int bpntr[], int bindx[],
                  int extern_proc[], int N_external, int option, int m_type)

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

  Find an ordering for the elements in 'update' and 'external' (if option ==
  AZ_EXTERNS, we only order 'external'). We first order 'update' by numbering
  the internal elements first followed by the border elements. On output
  update_index[] contains the numbers 0->N_update-1 ordered such that
  update_index[k] < update_index[j] if update[k] is internal and update[j] is
  border. Next, we order the external elements such that elements that are
  updated by the same processor are contiguous.

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

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

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

  update_index:    On output, update_index[i] gives the local numbering of
                   global point 'update[i]'. See User's Guide.

  extern_index:    On output, extern_index[i] gives the local numbering of
                   global point 'external[i]'. See User's Guide.

  internal:        On output, equals the number of internal elements on this
                   processor.

  border:          On output, equals the number of border elements on this
                   processor.

  N_update:        Number of elements updated on this processor.

  bpntr, bindx:    Nonzero matrix information (see User's Guide).
                   Note: for MSR arrays the user invokes this routine by calling
                   AZ_order_ele(...,bindx,bindx, ....). That is both 'bpntr' and
                   'bindx' point to 'bindx'.

  extern_proc:     extern_proc[i] is updating processor of external[i].

  N_external:      Number of external elements on this processor.

  option:          option = AZ_ALL ==> order internal, border and extern ele.
                   option = AZ_EXTERNS ==> order only external elements.
  m_type:          On input, m_type = AZ_MSR_MATRIX
                     or      m_type = AZ_VBR_MATRIX

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

{

  /* local variables */

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

  int  i, j, lilstatus, count;
  int  *t_ptr;
  char *yo = "AZ_order_ele: ";

  *internal = 0;
  *border   = 0;

  /*
   * Sort through update elements. Classify them as either internal or
   * border. Give an index for the internal elements. Temporarily give a
   * negative index for border elements.
   */

  if (option == AZ_ALL) {
    if (m_type == AZ_MSR_MATRIX)      t_ptr = bindx;
    else if (m_type == AZ_VBR_MATRIX) t_ptr = bpntr;
    else {
      (void) fprintf(stderr, "%sERROR: Unknown matrix type (%d)\n", yo, m_type);
      exit(1);
    }

    for (i = 0; i < N_update; i++) {
      lilstatus = 0;

      for (j = t_ptr[i]; j < t_ptr[i + 1]; j++) {
        if (bindx[j] >= N_update) {
          update_index[i] = -*border - 1;

          /*
           * We encode the border elements as negative numbers so that we can
           * pull them out later.
           */

          (*border)++;
          lilstatus = 1;
          break;
        }
      }

      if (lilstatus == 0) {
        update_index[i] = *internal;
        (*internal)++;
      }
    }

    /* Use negative indices to number the border elements */

    for (i = 0; i < N_update; i++) {
      if (update_index[i] < 0)
        update_index[i] = *internal - update_index[i] - 1;
    }
  }

  else {
    for (i = 0; i < N_update; i++) update_index[i] = i;
    *internal = 0;
    *border = N_update;
  }

  /*
   * Sift through the external elements. For each newly encountered external
   * point assign it the next index in the sequence. Then look for other
   * external elements who are update by the same node and assign them the next
   * set of index numbers in the sequence (ie. elements updated by the same node
   * have consecutive indices).
   */

  count = N_update;
  for (i = 0; i < N_external; i++) extern_index[i] = -1;

  for (i = 0; i < N_external; i++) {
    if (extern_index[i] == -1) {
      extern_index[i] = count++;

      for (j = i + 1; j < N_external; j++) {
        if (extern_proc[j] == extern_proc[i]) extern_index[j] = count++;
      }
    }
  }

} /* AZ_order_ele */

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

void AZ_invorder_vec(double vector[], int data_org[], int update_index[], 
	int rpntr[], double new[])
{

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

  Reorder a vector (could be right hand side or solution vector) which 
  conforms to a matrix reordered via AZ_transform or AZ_reorder_matrix
  so that it conforms to the unAZ_transformed matrix.

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

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

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

  vector:          On input, a vector of length 'length'. On output,
                   'vector' is reordered to be consistant with 
                   update_index[]. 

  data_org:        Array use to specifiy communication information. See User's
                   Guide.

  update_index:    update_index[i] gives the local numbering of global point
                   'update[i]'.

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

    int i,j, ii,length, current;

    length = data_org[AZ_N_int_blk] + data_org[AZ_N_bord_blk];

    if (data_org[AZ_matrix_type] == AZ_MSR_MATRIX) {
      for (i = 0 ; i < length ; i++ ) new[i] = vector[ update_index[i] ] ;
    }
    else if (data_org[AZ_matrix_type] == AZ_VBR_MATRIX) {
      current = 0;
      for (i = 0 ; i < length ; i++ ) {
        ii = update_index[i];
        for (j = rpntr[ii] ; j < rpntr[ii+1] ; j++ ) { 
          new[current++] = vector[j] ;
        }
      }
    }
    else {
      (void) fprintf(stderr,"Error: Unknown matrix type (%d) in AZ_reorder_vec\n",
		     data_org[AZ_matrix_type]);
      exit(-1);
    }

}
/******************************************************************************/
/******************************************************************************/
/******************************************************************************/

void AZ_reorder_vec(double vector[], int data_org[], int update_index[], 
	int rpntr[])
{

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

  Reorder the vector (could be right hand side or solution vector) so
  that it conforms to a matrix reordered via AZ_transform or AZ_reorder_matrix.

  IMPORTANT: This routine assumes that update_index[] contains two sequencies of
  numbers that are ordered but intertwined. For example,

  update_index:  4 5 0 6 1 2 3 7

  seq 1 =>    0   1 2 3

  seq 2 =>4 5   6 7

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

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

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

  vector:          On input, a vector of length 'length'. On output,
                   'vector' is reordered to be consistant with 
                   update_index[]. 

  data_org:        Array use to specifiy communication information. See User's
                   Guide.

  update_index:    update_index[i] gives the local numbering of global point
                   'update[i]'.

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

    int i,ii,length, *temp;

    length = data_org[AZ_N_int_blk] + data_org[AZ_N_bord_blk];
    if (data_org[AZ_matrix_type] == AZ_MSR_MATRIX) 
      AZ_sortqlists((char *) vector, 0, update_index, length, sizeof(double), length);
    else if (data_org[AZ_matrix_type] == AZ_VBR_MATRIX) {
       temp = (int *) AZ_allocate( (length+1)*sizeof(int));
       if (temp == NULL) {
          (void) fprintf(stderr, "Out of memory in AZ_reorder_vec()\n");
          exit(-1);
       }
       for (i = 0 ; i < length ; i++ ) {
         ii = update_index[i];
         temp[i] = rpntr[ii+1]-rpntr[ii];
       }
       AZ_sortqlists((char *)vector,temp,update_index,rpntr[length],
                     sizeof(double),length);
       AZ_free(temp);
    }
    else {
      (void) fprintf(stderr,"Error: Unknown matrix type (%d) in AZ_reorder_vec\n",
		     data_org[AZ_matrix_type]);
      exit(-1);
    }

}


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

void AZ_reorder_matrix(int N_update, int bindx[], double val[],
                       int update_index[], int extern_index[], int indx[],
                       int rnptr[], int bnptr[], int N_external,
                       int cnptr[], int option, int mat_type)

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

  Reorder the matrix so that it corresponds to the new ordering given by
  'update_index' and 'extern_index'.

  IMPORTANT: This routine assumes that update_index[] contains two sequencies of
  numbers that are ordered but intertwined. For example,

  update_index:  4 5 0 6 1 2 3 7

  seq 1 =>    0   1 2 3

  seq 2 =>4 5   6 7

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

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

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

  N_update:        Number of elements updated on this processor.

  val,
  rnptr,
  bindx,
  indx,
  bnptr,
  cnptr:           Sparse matrix arrays. See User's Guide.
                   On input, the matrix corresponds to the initial ordering
                   (e.g. row i corresponds to global row update[i]).
                   On output, the matrix rows and columns are renumbered to
                   correspond to the ordering given by 'update_index' and
                   'extern_index'. (e.g. global row update[i] now appears
                   as row update_index[i] in the matrix).

  update_index:    update_index[i] gives the local numbering of global point
                   'update[i]'.

  extern_index:    extern_index[i] gives the local numbering of global point
                   'external[i]'.

  N_external:      Number of external elements on this processor.

  option:

  mat_type:        Indicates whether this is an MSR (= AZ_MSR_MATRIX) or a
                   VBR (= AZ_VBR_MATRIX).

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

{

  /* local variables */

  int   start, end;
  int   val_length, indx_length;
  int  *temp;
  int   i, j;
  char *yo = "AZ_reorder_matrix: ";

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

  if (mat_type == AZ_MSR_MATRIX) {
    start = N_update+1;      /* first nonzero offdiag */
    end   = bindx[N_update]; /* last nonzero          */
  }
  else if (mat_type == AZ_VBR_MATRIX) {
    start = 0;               /* first nonzero */
    end   = bnptr[N_update]; /* last nonzero  */

    /* reorder cnptr[] */

    /* 1) convert the cnptr array to give the blk size */

    AZ_convert_ptrs_to_values(cnptr, N_update + N_external);

    /* 2) order the internal blocks.  NOTE: AZ_sortqlists() can only be used for
     * the internal elements as it expects the list to correspond to 2 ordered
     * sequences that are intermixed.