Komplex
Development
|
Include file for Aztec Komplex Library. More...
#include <Komplex_config.h>
#include "az_aztec.h"
Classes | |
struct | AZ_KOMPLEX_STRUCT |
Macros | |
#define | AZK_True_Complex 1 |
#define | AZK_Komplex_Same_Structure 2 |
#define | AZK_Komplex_General 3 |
#define | AZK_Komplex_No_Copy 4 |
Typedefs | |
typedef struct AZ_KOMPLEX_STRUCT | AZ_KOMPLEX |
Functions | |
void | AZK_create_linsys_c2k (double *xc, double *bc, int *options, double *params, int *proc_config, AZ_MATRIX *Amat_complex, double **x, double **b, AZ_MATRIX **Amat_komplex) |
Create Komplex System from Complex System. More... | |
void | AZK_create_linsys_g2k (double *xr, double *xi, double *br, double *bi, int *options, double *params, int *proc_config, double c0r, double c0i, AZ_MATRIX *Amat_mat0, double c1r, double c1i, AZ_MATRIX *Amat_mat1, double **x, double **b, AZ_MATRIX **Amat_komplex) |
Create Komplex System from General System. More... | |
void | AZK_create_linsys_ri2k (double *xr, double *xi, double *br, double *bi, int *options, double *params, int *proc_config, AZ_MATRIX *Amat_real, double *val_imag, double **x, double **b, AZ_MATRIX **Amat_komplex) |
Create Komplex System from Real and Imaginary Parts. More... | |
void | AZK_create_linsys_no_copy (double *xr, double *xi, double *br, double *bi, int *options, double *params, int *proc_config, AZ_MATRIX *Amat_real, AZ_MATRIX *Amat_imag, double **x, double **b, AZ_MATRIX **Amat) |
void | AZK_create_matrix_g2k (int options[], double params[], int proc_config[], double c0r, double c0i, AZ_MATRIX *Amat_mat0, double c1r, double c1i, AZ_MATRIX *Amat_mat1, AZ_MATRIX **Amat_komplex) |
Create Komplex Matrix from General Matrix. More... | |
void | AZK_create_matrix_g2k_fill_entry (int nrow, int ncol, double c11, double c12, double *mat0v, double c21, double c22, double *mat1v, double *komplex) |
void | AZK_create_matrix_c2k (int options[], double params[], int proc_config[], AZ_MATRIX *Amat_complex, AZ_MATRIX **Amat_komplex) |
Create Komplex matrix from Complex matrix. More... | |
void | AZK_create_matrix_c2k_fill_entry (int nrow, int ncol, double *cur_complex, double *cur_komplex) |
void | AZK_create_matrix_ri2k (int options[], double params[], int proc_config[], AZ_MATRIX *Amat_real, double *val_imag, AZ_MATRIX **Amat_komplex) |
Create Komplex Matrix from Real and Imaginary Parts. More... | |
void | AZK_create_matrix_ri2k_fill_entry (int nrow, int ncol, double *realv, double *imagv, double *komplex) |
void | AZK_create_precon (int *options, double *params, int *proc_config, double *x, double *b, AZ_MATRIX *Amat, AZ_PRECOND **Prec) |
Create a Preconditioner for a Komplex matrix. More... | |
void | AZK_create_vector_c2k (int *options, double *params, int *proc_config, AZ_MATRIX *Amat_komplex, double *vc, double **vk) |
Create Komplex vector from Complex vector. More... | |
void | AZK_create_vector_g2k (int *options, double *params, int *proc_config, AZ_MATRIX *Amat_komplex, double *vr, double *vi, double **vk) |
Create Komplex vector from Real and Imaginary Parts. More... | |
void | AZK_create_vector_ri2k (int *options, double *params, int *proc_config, AZ_MATRIX *Amat_komplex, double *vr, double *vi, double **vk) |
Create Komplex vector from Real and Imaginary Parts. More... | |
void | AZK_destroy_linsys (int *options, double *params, int *proc_config, double **x, double **b, AZ_MATRIX **Amat_komplex) |
Destroy a Komplex System. More... | |
void | AZK_destroy_matrix (int options[], double params[], int proc_config[], AZ_MATRIX **Amat_komplex) |
Destroy a Komplex Matrix. More... | |
void | AZK_destroy_precon (int *options, double *params, int *proc_config, AZ_MATRIX *Amat, AZ_PRECOND **Prec) |
Destroy a Komplex preconditioner. More... | |
void | AZK_destroy_vector (int *options, double *params, int *proc_config, AZ_MATRIX *Amat_komplex, double **vk) |
Destroy a Komplex vector. More... | |
void | AZK_extract_solution_k2c (int *options, double *params, int *proc_config, AZ_MATRIX *Amat_komplex, AZ_PRECOND *Prec, double *vk, double *vc) |
Extract a Complex vector from a Komplex vector. More... | |
void | AZK_extract_solution_k2g (int *options, double *params, int *proc_config, AZ_MATRIX *Amat_komplex, AZ_PRECOND *Prec, double *vk, double *vr, double *vi) |
void | AZK_extract_solution_k2ri (int *options, double *params, int *proc_config, AZ_MATRIX *Amat_komplex, AZ_PRECOND *Prec, double *vk, double *vr, double *vi) |
Extract real/imaginary parts of a complex vector from a Komplex vector. More... | |
void | AZK_iterate (double *xr, double *xi, double *br, double *bi, int *options, double *params, double *status, int *proc_config, AZ_MATRIX *Amat_real, AZ_MATRIX *Amat_imag) |
void | AZK_matvec_no_copy (double *x, double *y, AZ_MATRIX *Amat, int proc_config[]) |
void | AZK_permute_ri (int *options, double *params, int *proc_config, double *b, AZ_MATRIX *Amat_komplex) |
Permute a Komplex system for better numerical stability. More... | |
void | AZK_precon (double x[], int options[], int proc_config[], double params[], AZ_MATRIX *Amat, AZ_PRECOND *Prec) |
double | AZK_residual_norm_no_copy (double *xr, double *xi, double *br, double *bi, int *options, double *params, int *proc_config, AZ_MATRIX *Amat_real, AZ_MATRIX *Amat_imag) |
double | AZK_residual_norm (double *xk, double *bk, int *options, double *params, int *proc_config, AZ_MATRIX *Amat_komplex) |
Include file for Aztec Komplex Library.
The file, along with az_aztec.h must be included in every file that uses the Komplex library.
KOMPLEX is an add-on module to AZTEC that allows users to solve complex-valued linear systems.
KOMPLEX solves a complex-valued linear system Ax = b by solving an equivalent real-valued system of twice the dimension. Specifically, writing in terms of real and imaginary parts, we have
or by separating into real and imaginary equations we have
which is a real-valued system of twice the size. If we find xr and xi, we can form the solution to the original system as x = xr +i*xi.
KOMPLEX accept user linear systems in three forms with either global or local index values.
1) The first form is true complex. The user passes in an MSR or VBR format matrix where the values are stored like Fortran complex numbers. Thus, the values array is of type double that is twice as long as the number of complex values. Each complex entry is stored with real part followed by imaginary part (as in Fortran).
2) The second form stores real and imaginary parts separately, but the pattern for each is identical. Thus only the values of the imaginary part are passed to the creation routines.
3) The third form accepts two real-valued matrices with no assumption about the structure of the matrices. Each matrix is multiplied by a user-supplied complex constant. This is the most general form.
Each of the above forms supports a global or local index set. By this we mean that the index values (stored in bindx) refer to the global problem indices, or the local indices (for example after calling AZ_transform).
void AZK_create_linsys_c2k | ( | double * | xc, |
double * | bc, | ||
int * | options, | ||
double * | params, | ||
int * | proc_config, | ||
AZ_MATRIX * | Amat_complex, | ||
double ** | x, | ||
double ** | b, | ||
AZ_MATRIX ** | Amat_komplex | ||
) |
Create Komplex System from Complex System.
Transforms a complex-valued system
Amat_complex * xc = bc
where double precision arrays hold the complex values of Amat_complex, xc and bc in Fortran complex format, i.e., if dimension of complex system is N then xc is of length 2*N and the first complex value is stored with the real part in xc[0] and the imaginary part in xc[1] and so on.
xc | (In) Contains the complex initial guess/solution vector with the real/imag parts interleaved as in Fortran complex format. |
bc | (In) RHS in Fortran complex format. |
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_complex | (In) An AZ_MATRIX structure where Amat_complex->val contain the values of the complex matrix in Fortran complex format. |
x | (Out) Komplex version of initial guess and solution. |
b | (Out) Komplex version of RHS. |
Amat_komplex | (Out) Komplex version of matrix stored as an AZ_MATRIX structure. |
References AZK_create_matrix_c2k(), AZK_create_vector_c2k(), and AZK_permute_ri().
void AZK_create_linsys_g2k | ( | double * | xr, |
double * | xi, | ||
double * | br, | ||
double * | bi, | ||
int * | options, | ||
double * | params, | ||
int * | proc_config, | ||
double | c0r, | ||
double | c0i, | ||
AZ_MATRIX * | Amat_mat0, | ||
double | c1r, | ||
double | c1i, | ||
AZ_MATRIX * | Amat_mat1, | ||
double ** | x, | ||
double ** | b, | ||
AZ_MATRIX ** | Amat_komplex | ||
) |
Create Komplex System from General System.
Transforms a complex-valued system
(c0r+i*c0i)*A0 +(c1r+i*c1i)*A1) * (xr+i*xi) = (br+i*bi)
to a Komplex system.
xr | (In) Real part of initial guess. |
xi | (In) Imaginary part of initial guess. |
br | (In) Real part of right hand side of linear system. |
bi | (In) Imaginary part of right hand side of linear system. |
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
c0r | (In) Real part of constant to be multiplied with first matrix. |
c0i | (In) Imaginary part of constant to be multiplied with first matrix. |
c1r | (In) Real part of constant to be multiplied with second matrix. |
c1i | (In) Imaginary part of constant to be multiplied with second matrix. |
Amat_mat0 | (In) AZ_MATRIX object containing first real-valued matrix. |
Amat_mat1 | (In) AZ_MATRIX object containing second real-valued matrix. |
x | (Out) Komplex version of initial guess and solution. |
b | (Out) Komplex version of RHS. |
Amat_komplex | (Out) Komplex version of matrix stored as an AZ_MATRIX structure. |
References AZK_create_matrix_g2k(), AZK_create_vector_g2k(), and AZK_permute_ri().
void AZK_create_linsys_ri2k | ( | double * | xr, |
double * | xi, | ||
double * | br, | ||
double * | bi, | ||
int * | options, | ||
double * | params, | ||
int * | proc_config, | ||
AZ_MATRIX * | Amat_real, | ||
double * | val_imag, | ||
double ** | x, | ||
double ** | b, | ||
AZ_MATRIX ** | Amat_komplex | ||
) |
Create Komplex System from Real and Imaginary Parts.
Transforms a complex-valued system
(Ar +i*Ai) * (xr + i*xi) = (br + i*bi)
where double precision arrays hold the real and imaginary parts separately. The pattern of the imaginary part matches the real part. Thus no structure for the imaginary part is passed in.
xr | (In) Real part of initial guess. |
xi | (In) Imaginary part of initial guess. |
br | (In) Real part of right hand side of linear system. |
bi | (In) Imaginary part of right hand side of linear system. |
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_real | (In) AZ_MATRIX object containing real matrix. |
val_imag | (In) Double arrya containing the values ONLY for imaginary matrix. |
x | (Out) Komplex version of initial guess and solution. |
b | (Out) Komplex version of RHS. |
Amat_komplex | (Out) Komplex version of matrix stored as an AZ_MATRIX structure. |
References AZK_create_matrix_ri2k(), AZK_create_vector_ri2k(), and AZK_permute_ri().
void AZK_create_matrix_c2k | ( | int | options[], |
double | params[], | ||
int | proc_config[], | ||
AZ_MATRIX * | Amat_complex, | ||
AZ_MATRIX ** | Amat_komplex | ||
) |
Create Komplex matrix from Complex matrix.
Transforms a complex-valued matrix where double precision array hold the complex values of Amat_complex in Fortran complex format.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_complex | (In) An AZ_MATRIX structure where Amat_complex->val contain the values of the complex matrix in Fortran complex format. |
Amat_komplex | (Out) Komplex version of matrix stored as an AZ_MATRIX structure. |
Referenced by AZK_create_linsys_c2k().
void AZK_create_matrix_g2k | ( | int | options[], |
double | params[], | ||
int | proc_config[], | ||
double | c0r, | ||
double | c0i, | ||
AZ_MATRIX * | Amat_mat0, | ||
double | c1r, | ||
double | c1i, | ||
AZ_MATRIX * | Amat_mat1, | ||
AZ_MATRIX ** | Amat_komplex | ||
) |
Create Komplex Matrix from General Matrix.
Transforms a complex-valued Matrix
(c0r+i*c0i)*A0 +(c1r+i*c1i)*A1)
to a Komplex matrix.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
c0r | (In) Real part of constant to be multiplied with first matrix. |
c0i | (In) Imaginary part of constant to be multiplied with first matrix. |
Amat_mat0 | (In) AZ_MATRIX object containing first real-valued matrix. |
c1r | (In) Real part of constant to be multiplied with second matrix. |
c1i | (In) Imaginary part of constant to be multiplied with second matrix. |
Amat_mat1 | (In) AZ_MATRIX object containing second real-valued matrix. |
Amat_komplex | (Out) Komplex version of matrix stored as an AZ_MATRIX structure. |
Referenced by AZK_create_linsys_g2k().
void AZK_create_matrix_ri2k | ( | int | options[], |
double | params[], | ||
int | proc_config[], | ||
AZ_MATRIX * | Amat_real, | ||
double * | val_imag, | ||
AZ_MATRIX ** | Amat_komplex | ||
) |
Create Komplex Matrix from Real and Imaginary Parts.
Transforms a complex-valued matrix
(Ar +i*Ai)
where double precision arrays hold the real and imaginary parts separately. The pattern of the imaginary part matches the real part. Thus no structure for the imaginary part is passed in.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_real | (In) AZ_MATRIX object containing real matrix. |
val_imag | (In) Double arrya containing the values ONLY for imaginary matrix. |
Amat_komplex | (Out) Komplex version of matrix stored as an AZ_MATRIX structure. |
Referenced by AZK_create_linsys_ri2k().
void AZK_create_precon | ( | int * | options, |
double * | params, | ||
int * | proc_config, | ||
double * | x, | ||
double * | b, | ||
AZ_MATRIX * | Amat, | ||
AZ_PRECOND ** | Prec | ||
) |
Create a Preconditioner for a Komplex matrix.
Constructs a preconditioner for a Komplex matrix Amat. All preconditioning options available in Aztec are supported.
options | (In) Determines specific preconditioner method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
x | (In/Out) Komplex version of initial guess and solution. May be modified depending on preconditioner options. |
b | (In/Out) Komplex version of RHS. May be modified depending on preconditioner options. |
Amat | (In) Komplex version of matrix stored as an AZ_MATRIX structure. |
Prec | (Out) Preconditioner for Amat stored as an AZ_PRECOND structure. |
void AZK_create_vector_c2k | ( | int * | options, |
double * | params, | ||
int * | proc_config, | ||
AZ_MATRIX * | Amat_komplex, | ||
double * | vc, | ||
double ** | vk | ||
) |
Create Komplex vector from Complex vector.
Transforms a complex-valued vector vc to a real vector where vc in Fortran complex format, i.e., if dimension of complex system is N then vc is of length 2*N and the first complex value is stored with the real part in vc[0] and the imaginary part in vc[1] and so on.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_komplex | (In) Komplex version of matrix stored as an AZ_MATRIX structure. |
vc | (In) Contains a complex vector with the real/imag parts interleaved as in Fortran complex format. |
vk | (Out) Komplex version of vc. |
Referenced by AZK_create_linsys_c2k().
void AZK_create_vector_g2k | ( | int * | options, |
double * | params, | ||
int * | proc_config, | ||
AZ_MATRIX * | Amat_komplex, | ||
double * | vr, | ||
double * | vi, | ||
double ** | vk | ||
) |
Create Komplex vector from Real and Imaginary Parts.
Transforms a complex vector where double precision arrays hold the real and imaginary parts separately.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_komplex | (Out) Komplex version of matrix stored as an AZ_MATRIX structure. |
vr | (In) Real part of input vector. |
vi | (In) Imaginary part of input vector. |
vk | (Out) Komplex version of input vector. |
References AZK_create_vector_ri2k().
Referenced by AZK_create_linsys_g2k().
void AZK_create_vector_ri2k | ( | int * | options, |
double * | params, | ||
int * | proc_config, | ||
AZ_MATRIX * | Amat_komplex, | ||
double * | vr, | ||
double * | vi, | ||
double ** | vk | ||
) |
Create Komplex vector from Real and Imaginary Parts.
Transforms a complex vector where double precision arrays hold the real and imaginary parts separately.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_komplex | (Out) Komplex version of matrix stored as an AZ_MATRIX structure. |
vr | (In) Real part of input vector. |
vi | (In) Imaginary part of input vector. |
vk | (Out) Komplex version of input vector. |
Referenced by AZK_create_linsys_ri2k(), and AZK_create_vector_g2k().
void AZK_destroy_linsys | ( | int * | options, |
double * | params, | ||
int * | proc_config, | ||
double ** | x, | ||
double ** | b, | ||
AZ_MATRIX ** | Amat_komplex | ||
) |
Destroy a Komplex System.
Destroys a komplex system created by any of the AZK_create_linsys functions. Deletes any memory allocated by creation routine.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
x | (Out) Deleted komplex version of solution. Remember to call AZK_extract_solution_[k2c,g2k,ri2k] before calling this routine. |
b | (Out) Deleted komplex version of RHS. |
Amat_komplex | (Out) Deleted komplex version of matrix stored as an AZ_MATRIX structure. |
References AZK_destroy_matrix(), and AZK_destroy_vector().
void AZK_destroy_matrix | ( | int | options[], |
double | params[], | ||
int | proc_config[], | ||
AZ_MATRIX ** | Amat_komplex | ||
) |
Destroy a Komplex Matrix.
Destroys a komplex matrix created by any of the AZK_create_matrix functions. Deletes any memory allocated by creation routine.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_komplex | (Out) Deleted komplex version of matrix stored as an AZ_MATRIX structure. |
Referenced by AZK_destroy_linsys().
void AZK_destroy_precon | ( | int * | options, |
double * | params, | ||
int * | proc_config, | ||
AZ_MATRIX * | Amat, | ||
AZ_PRECOND ** | Prec | ||
) |
Destroy a Komplex preconditioner.
Destroys a komplex preconditioner created by the AZK_create_preconditioner function. Deletes any memory allocated by creation routine.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat | (In) Komplex version of matrix stored as an AZ_MATRIX structure. |
Prec | (Out) Deleted komplex version of preconditioner stored as an AZ_PRECOND structure. |
void AZK_destroy_vector | ( | int * | options, |
double * | params, | ||
int * | proc_config, | ||
AZ_MATRIX * | Amat_komplex, | ||
double ** | vk | ||
) |
Destroy a Komplex vector.
Destroys a komplex vector created by any of the AZK_create_vector functions. Deletes any memory allocated by creation routine.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_komplex | (In) Komplex version of matrix stored as an AZ_MATRIX structure. |
vk | (Out) Deleted komplex version of a vector. Remember to call AZK_extract_solution_[k2c,g2k,ri2k] before calling this routine. |
Referenced by AZK_destroy_linsys().
void AZK_extract_solution_k2c | ( | int * | options, |
double * | params, | ||
int * | proc_config, | ||
AZ_MATRIX * | Amat_komplex, | ||
AZ_PRECOND * | Prec, | ||
double * | vk, | ||
double * | vc | ||
) |
Extract a Complex vector from a Komplex vector.
Transforms a komplex vector to a complex vector.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_komplex | (In) Komplex version of matrix stored as an AZ_MATRIX structure. |
Prec | (In) Preconditioner for Amat stored as an AZ_PRECOND structure. |
vk | (In) Komplex version of vector. |
vc | (Out) Contains a complex vector with the real/imag parts interleaved as in Fortran complex format. Note that the user must allocate sufficient storage for results. |
void AZK_extract_solution_k2ri | ( | int * | options, |
double * | params, | ||
int * | proc_config, | ||
AZ_MATRIX * | Amat_komplex, | ||
AZ_PRECOND * | Prec, | ||
double * | vk, | ||
double * | vr, | ||
double * | vi | ||
) |
Extract real/imaginary parts of a complex vector from a Komplex vector.
Transforms a komplex vector to real and imaginary parts.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
Amat_komplex | (In) Komplex version of matrix stored as an AZ_MATRIX structure. |
Prec | (In) Preconditioner for Amat stored as an AZ_PRECOND structure. |
vk | (In) Komplex version of vector. |
vc | (Out) Contains a complex vector with the real/imag parts interleaved as in Fortran complex format. Note that the user must allocate sufficient storage for results. |
void AZK_permute_ri | ( | int * | options, |
double * | params, | ||
int * | proc_config, | ||
double * | b, | ||
AZ_MATRIX * | Amat_komplex | ||
) |
Permute a Komplex system for better numerical stability.
An alternative to the standard Komplex formulation is to permute the block rows so that the imaginary part is on the main diagonal. For example:
This action may be desirable, or necessary in situations where the real part has small or zero diagonal entries. This routine looks at each real/imaginary pair and, based on a heuristic may swap the real and imaginary parts. This action does not affect the sparsity pattern, but only the mapping from the complex (or real/imaginary) mapping to the komplex mapping, and back.
options | (In) Determines specific solution method and other parameters. |
params | (In) Drop tolerance and convergence tolerance info. |
proc_config | (In) Machine configuration. proc_config[AZ_node] is the node number. proc_config[AZ_N_procs] is the number of processors. |
b | (Out) Komplex version of RHS, possibly permuted. |
Amat_komplex | (Out) Komplex version of matrix stored as an AZ_MATRIX structure, possibly permuted. |
Referenced by AZK_create_linsys_c2k(), AZK_create_linsys_g2k(), and AZK_create_linsys_ri2k().