SuperLU Distributed 8.2.1
Distributed memory sparse direct solver
pzgstrs.c File Reference

Solves a system of distributed linear equations A*X = B with a general N-by-N matrix A using the LU factors computed previously. More...

#include <math.h>
#include "superlu_zdefs.h"
Include dependency graph for pzgstrs.c:

Macros

#define CACHELINE   64 /* bytes, Xeon Phi KNL, Cori haswell, Edision */
 

Functions

int_t pzReDistribute_B_to_X (doublecomplex *B, int_t m_loc, int nrhs, int_t ldb, int_t fst_row, int_t *ilsum, doublecomplex *x, zScalePermstruct_t *ScalePermstruct, Glu_persist_t *Glu_persist, gridinfo_t *grid, zSOLVEstruct_t *SOLVEstruct)
 
int_t pzReDistribute_X_to_B (int_t n, doublecomplex *B, int_t m_loc, int_t ldb, int_t fst_row, int_t nrhs, doublecomplex *x, int_t *ilsum, zScalePermstruct_t *ScalePermstruct, Glu_persist_t *Glu_persist, gridinfo_t *grid, zSOLVEstruct_t *SOLVEstruct)
 
void pzCompute_Diag_Inv (int_t n, zLUstruct_t *LUstruct, gridinfo_t *grid, SuperLUStat_t *stat, int *info)
 
void pzgstrs (superlu_dist_options_t *options, int_t n, zLUstruct_t *LUstruct, zScalePermstruct_t *ScalePermstruct, gridinfo_t *grid, doublecomplex *B, int_t m_loc, int_t fst_row, int_t ldb, int nrhs, zSOLVEstruct_t *SOLVEstruct, SuperLUStat_t *stat, int *info)
 

Detailed Description

Solves a system of distributed linear equations A*X = B with a general N-by-N matrix A using the LU factors computed previously.

Copyright (c) 2003, The Regents of the University of California, through Lawrence Berkeley National Laboratory (subject to receipt of any required approvals from U.S. Dept. of Energy)

All rights reserved.

The source code is distributed under BSD license, see the file License.txt at the top-level directory.

-- Distributed SuperLU routine (version 6.1) --
Lawrence Berkeley National Lab, Univ. of California Berkeley.
October 15, 2008
September 18, 2018  version 6.0
February 8, 2019  version 6.1.1

Macro Definition Documentation

◆ CACHELINE

#define CACHELINE   64 /* bytes, Xeon Phi KNL, Cori haswell, Edision */

Function Documentation

◆ pzCompute_Diag_Inv()

void pzCompute_Diag_Inv ( int_t  n,
zLUstruct_t LUstruct,
gridinfo_t grid,
SuperLUStat_t stat,
int *  info 
)
Purpose
=======
  Compute the inverse of the diagonal blocks of the L and U
  triangular matrices.

◆ pzgstrs()

void pzgstrs ( superlu_dist_options_t options,
int_t  n,
zLUstruct_t LUstruct,
zScalePermstruct_t ScalePermstruct,
gridinfo_t grid,
doublecomplex B,
int_t  m_loc,
int_t  fst_row,
int_t  ldb,
int  nrhs,
zSOLVEstruct_t SOLVEstruct,
SuperLUStat_t stat,
int *  info 
)
Purpose
=======

PZGSTRS solves a system of distributed linear equations
A*X = B with a general N-by-N matrix A using the LU factorization
computed by PZGSTRF.
If the equilibration, and row and column permutations were performed,
the LU factorization was performed for A1 where
    A1 = Pc*Pr*diag(R)*A*diag(C)*Pc^T = L*U
and the linear system solved is
    A1 * Y = Pc*Pr*B1, where B was overwritten by B1 = diag(R)*B, and
the permutation to B1 by Pc*Pr is applied internally in this routine.

Arguments
=========

options (input) superlu_dist_options_t*
        The structure defines the input parameters to control
        how the LU decomposition and triangular solve are performed.

n      (input) int (global)
       The order of the system of linear equations.

LUstruct (input) zLUstruct_t*
       The distributed data structures storing L and U factors.
       The L and U factors are obtained from PZGSTRF for
       the possibly scaled and permuted matrix A.
       See superlu_zdefs.h for the definition of 'zLUstruct_t'.
       A may be scaled and permuted into A1, so that
       A1 = Pc*Pr*diag(R)*A*diag(C)*Pc^T = L*U

grid   (input) gridinfo_t*
       The 2D process mesh. It contains the MPI communicator, the number
       of process rows (NPROW), the number of process columns (NPCOL),
       and my process rank. It is an input argument to all the
       parallel routines.
       Grid can be initialized by subroutine SUPERLU_GRIDINIT.
       See superlu_defs.h for the definition of 'gridinfo_t'.

B      (input/output) doublecomplex*
       On entry, the distributed right-hand side matrix of the possibly
       equilibrated system. That is, B may be overwritten by diag(R)*B.
       On exit, the distributed solution matrix Y of the possibly
       equilibrated system if info = 0, where Y = Pc*diag(C)^(-1)*X,
       and X is the solution of the original system.

m_loc  (input) int (local)
       The local row dimension of matrix B.

fst_row (input) int (global)
       The row number of B's first row in the global matrix.

ldb    (input) int (local)
       The leading dimension of matrix B.

nrhs   (input) int (global)
       Number of right-hand sides.

SOLVEstruct (input) zSOLVEstruct_t* (global)
       Contains the information for the communication during the
       solution phase.

stat   (output) SuperLUStat_t*
       Record the statistics about the triangular solves.
       See util.h for the definition of 'SuperLUStat_t'.

info   (output) int*
       = 0: successful exit
    < 0: if info = -i, the i-th argument had an illegal value
Here is the call graph for this function:

◆ pzReDistribute_B_to_X()

int_t pzReDistribute_B_to_X ( doublecomplex B,
int_t  m_loc,
int  nrhs,
int_t  ldb,
int_t  fst_row,
int_t ilsum,
doublecomplex x,
zScalePermstruct_t ScalePermstruct,
Glu_persist_t Glu_persist,
gridinfo_t grid,
zSOLVEstruct_t SOLVEstruct 
)
Purpose
=======
  Re-distribute B on the diagonal processes of the 2D process mesh.

Note
====
  This routine can only be called after the routine pzgstrs_init(),
  in which the structures of the send and receive buffers are set up.

Arguments
=========

B      (input) doublecomplex*
       The distributed right-hand side matrix of the possibly
       equilibrated system.

m_loc  (input) int (local)
       The local row dimension of matrix B.

nrhs   (input) int (global)
       Number of right-hand sides.

ldb    (input) int (local)
       Leading dimension of matrix B.

fst_row (input) int (global)
       The row number of B's first row in the global matrix.

ilsum  (input) int* (global)
       Starting position of each supernode in a full array.

x      (output) doublecomplex*
       The solution vector. It is valid only on the diagonal processes.

ScalePermstruct (input) zScalePermstruct_t*
       The data structure to store the scaling and permutation vectors
       describing the transformations performed to the original matrix A.

grid   (input) gridinfo_t*
       The 2D process mesh.

SOLVEstruct (input) zSOLVEstruct_t*
       Contains the information for the communication during the
       solution phase.

Return value
============
Here is the call graph for this function:

◆ pzReDistribute_X_to_B()

int_t pzReDistribute_X_to_B ( int_t  n,
doublecomplex B,
int_t  m_loc,
int_t  ldb,
int_t  fst_row,
int_t  nrhs,
doublecomplex x,
int_t ilsum,
zScalePermstruct_t ScalePermstruct,
Glu_persist_t Glu_persist,
gridinfo_t grid,
zSOLVEstruct_t SOLVEstruct 
)
Purpose
=======
  Re-distribute X on the diagonal processes to B distributed on all
  the processes.

Note
====
  This routine can only be called after the routine pzgstrs_init(),
  in which the structures of the send and receive buffers are set up.
Here is the call graph for this function: