az_check.c

来自「并行解法器,功能强大」· C语言 代码 · 共 889 行 · 第 1/2 页

  sprintf(yo, "AZ_check_options: ");
  switch (options[AZ_solver]) {

  case AZ_cg:
    break;
  case AZ_analyze:
    break;
  case AZ_GMRESR:
    break;
  case AZ_gmres:
    break;
  case AZ_cgs:
    break;
  case AZ_tfqmr:
    break;
  case AZ_bicgstab:
    break;
  case AZ_symmlq:
    break;
  case AZ_fixed_pt:
    break;
  case AZ_lu:
    if (az_nprocs != 1) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sERROR: LU not implemented in parallel."
                       "\n       Try domain decompostion with LU "
                       "preconditioning.\n\n", yo);
      }
      return 0;
    }

/*
    if ((data_org[AZ_matrix_type] != AZ_MSR_MATRIX) &&
        (data_org[AZ_matrix_type] != AZ_VBR_MATRIX)) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "ERROR: LU factorization can only be "
                       "used with MSR or VBR matrices.\n"
                       "       data_org[AZ_matrix_type] = %d\n\n",
                       data_org[AZ_matrix_type]);
      }
      return 0;
    }
*/
    break;

  default:
    if (az_proc == 0) {
      (void) fprintf(stderr, "%sERROR: options[AZ_solver] has improper value "
                     "(%d)\n\n", yo, options[AZ_solver]);
    }
    return 0;
  }

  switch (options[AZ_precond]) {

  case AZ_none:
    break;
  case AZ_Neumann:
    break;
  case AZ_ls:
    if (options[AZ_poly_ord] > AZ_MAX_POLY_ORDER) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sERROR: Exceeds highest order least_squares "
                       "polynomial available: %d vs %d\n\n", yo,
                       options[AZ_poly_ord], AZ_MAX_POLY_ORDER);
      }
      return 0;
    }
    break;

  case AZ_Jacobi:
    break;

  case AZ_sym_GS:
    if (data_org[AZ_matrix_type] != AZ_MSR_MATRIX) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sERROR: sym GS preconditioning can only "
                       "be used with MSR matrices.\n"
                       "       data_org[AZ_matrix_type] = %d\n\n", yo,
                       data_org[AZ_matrix_type]);
      }
      return 0;
    }
/*
    if (precond->Pmat!=Amat) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sERROR: sym GS preconditioning can only %s",
                       yo,"be used with Pmat=Amat .\n");
      }
      return 0;
    }
*/



    break;

  case AZ_bilu:
/* Begin Aztec 2.1 mheroux mod */
  case AZ_bilu_ifp:
/* End Aztec 2.1 mheroux mod */
     if (options[AZ_reorder]) {
        options[AZ_reorder] = 0;
        if ((options[AZ_output] != AZ_none) && (az_proc  == 0)) {
           printf("\t\t***** Reordering not implemented for Block ILU.\n");
           printf("\t\t***** Continuing without reordering\n");
        }
     }

    if (data_org[AZ_matrix_type] != AZ_VBR_MATRIX) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sERROR: Block ILU can only be used on VBR "
                       "matrices\n       data_org[AZ_matrix_type] = %d\n", yo,
                       data_org[AZ_matrix_type]);
        (void) fprintf(stderr, "       options[AZ_precond]      = %d\n\n",
                       options[AZ_precond]);
      }
      return 0;
    }
    if ( (options[AZ_solver] == AZ_cg) && (options[AZ_overlap] > 0) &&
         (options[AZ_type_overlap] != AZ_symmetric)) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: Preconditioned matrix may not be"
                       " symmetric (due to overlap).\n\n", yo);
      }
    }
    break;
  case AZ_multilevel:
    if (az_proc == 0) 
      printf("Are you sure you want the multilevel preconditioner\n");
    break;
  case AZ_dom_decomp:
/* Begin Aztec 2.1 mheroux mod */
    if ((options[AZ_subdomain_solve]==AZ_bilu ||
        options[AZ_subdomain_solve]==AZ_bilu_ifp)&&(options[AZ_reorder])){
/* End Aztec 2.1 mheroux mod */
        options[AZ_reorder] = 0;
        if ((options[AZ_output] != AZ_none) && (az_proc  == 0)) {
           printf("\t\t***** Reordering not implemented for Block ILU.\n");
           printf("\t\t***** Continuing without reordering\n");
        }
    }
   if ( (options[AZ_solver] == AZ_cg) && 
          ((options[AZ_subdomain_solve] == AZ_ilu) ||
           (options[AZ_subdomain_solve] == AZ_rilu) ||
           (options[AZ_subdomain_solve] == AZ_bilu_ifp) ||
           (options[AZ_subdomain_solve] == AZ_ilut))) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: Preconditioned matrix may not be"
                       " symmetric.\n\n", yo);
      }
    }

    if ( (options[AZ_subdomain_solve] != AZ_lu  ) &&
         (options[AZ_subdomain_solve] != AZ_ilu ) &&
         (options[AZ_subdomain_solve] != AZ_icc ) &&
         (options[AZ_subdomain_solve] != AZ_rilu) &&
         (options[AZ_subdomain_solve] != AZ_ilut) &&
/* Begin Aztec 2.1 mheroux mod */
         (options[AZ_subdomain_solve] != AZ_bilu_ifp) &&
/* End Aztec 2.1 mheroux mod */
         (options[AZ_subdomain_solve] != AZ_bilu) ) {
       if (options[AZ_subdomain_solve] >= AZ_SOLVER_PARAMS) {
          if (az_proc == 0) {
             (void) fprintf(stderr, "%sERROR: options[AZ_subdomain_solve]"
                    " has improper value = %d\n\n", yo, 
                    options[AZ_subdomain_solve]);
           }
           return 0; 
        }
        else {
           AZ_recover_sol_params(options[AZ_subdomain_solve], 
			         &sub_options,&sub_params,&sub_status,
                                 &sub_matrix,&sub_precond,&sub_scaling);
           if (!AZ_check_options(sub_options, az_proc, data_org, az_nprocs,
				    sub_params, sub_matrix, sub_precond)) 
	      return 0;
        }
    }
#ifdef eigen
  case AZ_slu:
#endif
  case AZ_lu:
  case AZ_ilu:
  case AZ_icc:
  case AZ_rilu:
  case AZ_ilut:
    if (options[AZ_overlap] < AZ_diag) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sERROR: Negative overlapping not allowed\n",
                       yo);
      }
       return 0;
    }
    if ( (options[AZ_solver] == AZ_cg) && (options[AZ_overlap] > 0) &&
         (options[AZ_type_overlap] != AZ_symmetric)) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: Preconditioned matrix may not be"
                       " symmetric (due to overlap).\n\n", yo);
      }
    }

    break;

  case AZ_smoother:
    break;
  case AZ_user_precond:
    break;

  default:
    if (options[AZ_precond] >= AZ_SOLVER_PARAMS) {
    if (az_proc == 0) {
      (void) fprintf(stderr, "%sERROR: options[AZ_precond] has improper "
                     "value = %d\n\n", yo, options[AZ_precond]);
    }
    return 0; }
    else {
      AZ_recover_sol_params(options[AZ_precond], &sub_options, &sub_params,
                            &sub_status, &sub_matrix, &sub_precond, &sub_scaling);
      if (!AZ_check_options(sub_options, az_proc, data_org, az_nprocs,
                     sub_params, sub_matrix, sub_precond)) return 0;
    }
  }

  switch (options[AZ_scaling]) {

  case AZ_none:
    break;
  case AZ_Jacobi:
    if ((data_org[AZ_matrix_type] != AZ_MSR_MATRIX) &&
        (data_org[AZ_matrix_type] != AZ_VBR_MATRIX)) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: Jacobi scaling can only be "
                       "with MSR or VBR\n                 matrices. "
                       "Turning off scaling.\n\n", yo);
      }
    }

    if (data_org[AZ_matrix_type] == AZ_VBR_MATRIX) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: Jacobi scaling for VBR matrices "
                       "is not implemented. Substituting\n"
                       "         block Jacobi instead.\n\n", yo);
      }
    }

    if (options[AZ_solver] == AZ_cg) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: Jacobi scaling may make matrix "
                       "unsymmetric for CG.\n\n", yo);
      }
    }
    break;

  case AZ_BJacobi:
    if ((data_org[AZ_matrix_type] != AZ_MSR_MATRIX) &&
        (data_org[AZ_matrix_type] != AZ_VBR_MATRIX)) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: block Jacobi scaling can only be "
                       "used with MSR or VBR\n                 matrices."
                       "Turning off scaling.\n\n", yo);
      }
    }

    if (data_org[AZ_matrix_type] == AZ_MSR_MATRIX) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: Block Jacobi for MSR matrices is "
                       "not implemented. Substituting \n"
                       "         Jacobi instead.\n\n", yo);
      }
    }

    if (options[AZ_solver] == AZ_cg) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: Jacobi scaling may make matrix "
                       "unsymmetric for CG.\n\n", yo);
      }
    }
    break;

  case AZ_row_sum:
    if ((data_org[AZ_matrix_type] != AZ_MSR_MATRIX) &&
        (data_org[AZ_matrix_type] != AZ_VBR_MATRIX)) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: left row-sum scaling can only be "
                       "with MSR or VBR\n                  matrices. "
                       "Turning off.\n\n", yo);
      }
    }

    if (options[AZ_solver] == AZ_cg) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: row sum scaling may make matrix "
                       "unsymmetric for CG.\n\n", yo);
      }
    }
    break;

  case AZ_sym_diag:
    if ((data_org[AZ_matrix_type] != AZ_MSR_MATRIX) &&
        (data_org[AZ_matrix_type] != AZ_VBR_MATRIX)) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: sym diag scaling can only be "
                       "used with MSR or VBR\n                 matrices. "
                       "Turning off.\n\n", yo);
      }
    }
    break;


  case AZ_sym_BJacobi:
    if (data_org[AZ_matrix_type] != AZ_VBR_MATRIX) {
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: sym block diag. scaling can only be "
                       "used with VBR\n                 matrices."
                       " Turning off.\n\n", yo);
      }
    }
    break;

  case AZ_sym_row_sum:
    if (data_org[AZ_matrix_type] != AZ_MSR_MATRIX){
      if (az_proc == 0) {
        (void) fprintf(stderr, "%sWARNING: sym row scaling can only be "
                       "used with MSR matrices.\n"
                       "                     Turning off.\n\n", yo);
      }
    }
    break;

  default:
    if (az_proc == 0) {
      (void) fprintf(stderr, "%sERROR: scaling flag %d not implemented.\n\n",
                     yo, options[AZ_scaling]);
    }
    return 0;
  }

  /* check the the norm used */
  if (((options[AZ_conv] == AZ_sol) || (options[AZ_conv] == AZ_Anorm)) &&
      ((data_org[AZ_matrix_type] != AZ_MSR_MATRIX) &&
       (data_org[AZ_matrix_type] != AZ_VBR_MATRIX))) {
        if (az_proc == 0) {
        (void) fprintf(stderr, "%sERROR: The norm of the matrix can only be "
                       "computed  with MSR or VBR matrices.\n"
                       "       data_org[AZ_matrix_type] = %d\n\n", yo,
                       data_org[AZ_matrix_type]);
      }
      return 0;
  }



  if (((options[AZ_conv] == 4) || (options[AZ_conv] == 3)) &&
      ((options[AZ_solver] == AZ_gmres) || (options[AZ_solver] == AZ_tfqmr))){
    if (az_proc == 0) {
      (void) fprintf(stderr, "%sWARNING: This convergence option requires "
                     "slightly more work for this solver\n"
                     "         in order to compute the residual.\n\n", yo);
    }
  }

  /* check the weights */

  if (options[AZ_conv] == 4) {
    for (i = 0; i < data_org[AZ_N_internal] + data_org[AZ_N_border]; i++) {
      if (params[AZ_weights+i] <= 0.0) {
        (void) fprintf(stderr, "%sWARNING: A weight vector component is <= 0, "
                       "check params[AZ_WEIGHTS]\n\n", yo);
        return 0;
      }
    }
  }

  return 1;

} /* AZ_check_options */

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

void AZ_defaults(int options[], double params[])

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

  Routine to set default parameters for iterative solver options.

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

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

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

  options:         Determines specific solution method and other parameters.

  params:          Drop tolerance and convergence tolerance info.

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

{
  int i;

  for (i = 0 ; i < AZ_OPTIONS_SIZE; i++ ) options[i] = 0;
  for (i = 0 ; i < AZ_PARAMS_SIZE; i++ )  params[i] = 0.0;

  /* setup default options */

  options[AZ_solver]   = AZ_gmres;
  options[AZ_scaling]  = AZ_none;
  options[AZ_precond]  = AZ_none;
  options[AZ_conv]     = AZ_r0;
  options[AZ_output]   = 1;
  options[AZ_pre_calc] = AZ_calc;
  options[AZ_max_iter] = 500;
  options[AZ_poly_ord] = 3;
  options[AZ_overlap]  = AZ_none;
  options[AZ_type_overlap]  = AZ_standard;
  options[AZ_kspace]   = 30;
  options[AZ_orthog]   = AZ_classic;
  options[AZ_aux_vec]  = AZ_resid;
  options[AZ_reorder]  = 1;
  options[AZ_keep_info]= 0;
  options[AZ_subdomain_solve] = AZ_ilut;
  options[AZ_graph_fill] = 0;
  options[AZ_init_guess] = AZ_NOT_ZERO;
  options[AZ_keep_kvecs] = 0;
  options[AZ_apply_kvecs]= AZ_FALSE;
  options[AZ_orth_kvecs] = AZ_FALSE;
  options[AZ_ignore_scaling] = AZ_FALSE;
  options[AZ_check_update_size] = AZ_FALSE;
  options[AZ_extreme] = AZ_high;

  params[AZ_tol]  = 1.0e-06;
  params[AZ_drop] = 0.0;
  params[AZ_ilut_fill] = 1.;
  params[AZ_omega]= 1.;
  params[AZ_update_reduction] = 10e10;
}