19#include "fasp_functs.h"
59 INT inroot = -10, nsizei = -10, nsizeall = -10, nlvl = 0;
65 INT *iblock = NULL, *jblock = NULL, *mask = NULL, *maxa = NULL;
73 printf(
"### DEBUG: [-Begin-] %s ...\n", __FUNCTION__);
83 memset(mask, 0,
sizeof(
INT)*n);
84 memset(iblock, 0,
sizeof(
INT)*n);
85 memset(maxa, 0,
sizeof(
INT)*n);
97 for ( i = 0; i < MIS.
row; i++ ) {
99 SWZ_level(inroot,&A,mask,&nlvl,maxa,jblock,maxlev);
105 printf(
"### DEBUG: nsizeall = %d\n", nsizeall);
116 for (i=0;i<MIS.
row;i++) {
118 SWZ_level(inroot,&A,mask,&nlvl,maxa,jb,maxlev);
120 iblock[i+1]=iblock[i]+nsizei;
127 printf(
"### DEBUG: nsizeall = %d, %d\n", nsizeall, iblock[nblk]);
134 memset(mask, 0,
sizeof(
INT)*n);
138 SWZ_block(swzdata, nblk, iblock, jblock, mask);
148 for (i=0; i<nblk; ++i)
149 mumps[i] = fasp_mumps_factorize(&blk[i], NULL, NULL,
PRINT_NONE);
150 swzdata->
mumps = mumps;
162 for (i=0; i<nblk; ++i) {
166 numeric[i] = fasp_umfpack_factorize(&blk[i], 0);
168 swzdata->numeric = numeric;
181 printf(
"### DEBUG: n = %d, #blocks = %d, max block size = %d\n",
182 n, nblk, swzdata->
maxbs);
188 swzdata->
nblk = nblk;
191 swzdata->
mask = mask;
192 swzdata->
maxa = maxa;
196 printf(
"### DEBUG: [--End--] %s ...\n", __FUNCTION__);
221 INT i, j, iblk, ki, kj, kij, is, ibl0, ibl1, nloc, iaa, iab;
242 void **numeric = swzdata->numeric;
249 for (is=0; is<nblk; ++is) {
254 for (i=0; i<nloc; ++i ) {
260 for (i=0; i<nloc; ++i) {
266 for (kij = iaa; kij<iab; ++kij) {
270 rhs.
val[i] -= val[kij]*x->
val[kj];
281 fasp_mumps_solve(&blk[is], &rhs, &u, mumps[is], 0);
289 fasp_umfpack_solve(&blk[is], &rhs, &u, numeric[is], 0);
302 for (i=0; i<nloc; ++i) {
330 INT i, j, iblk, ki, kj, kij, is, ibl0, ibl1, nloc, iaa, iab;
351 void **numeric = swzdata->numeric;
358 for (is=nblk-1; is>=0; --is) {
363 for (i=0; i<nloc; ++i ) {
369 for (i=0; i<nloc; ++i) {
375 for (kij = iaa; kij<iab; ++kij) {
379 rhs.
val[i] -= val[kij]*x->
val[kj];
390 fasp_mumps_solve(&blk[is], &rhs, &u, mumps[is], 0);
398 fasp_umfpack_solve(&blk[is], &rhs, &u, numeric[is], 0);
411 for (i=0; i<nloc; ++i) {
441static void SWZ_level (
const INT inroot,
452 INT i, j, lvl, lbegin, lvlend, nsize, node;
453 INT jstrt, jstop, nbr, lvsize;
456 if (ia[inroot+1]-ia[inroot] <= 1) {
459 jblock[iblock[lvl]] = inroot;
475 while (lvsize > 0 && lvl < maxlev) {
478 iblock[lvl] = lbegin;
480 for(i=lbegin; i<lvlend; ++i) {
483 jstop = ia[node+1]-1;
484 for (j = jstrt; j<jstop; ++j) {
486 if (mask[nbr] == 0) {
493 lvsize = nsize - lvlend;
499 for (i = 0; i< nsize; ++i) {
523static void SWZ_block (
SWZ_data *swzdata,
529 INT i, j, iblk, ki, kj, kij, is, ibl0, ibl1, nloc, iaa, iab;
530 INT maxbs = 0, count, nnz;
540 for (is=0; is<nblk; ++is) {
544 maxbs =
MAX(maxbs, nloc);
547 swzdata->
maxbs = maxbs;
553 for (is=0; is<nblk; ++is) {
558 for (i=0; i<nloc; ++i ) {
571 for (i=0; i<nloc; ++i) {
576 for (kij = iaa; kij<iab; ++kij) {
580 blk[is].
JA[nnz] = j-1;
581 blk[is].
val[nnz] = val[kij];
585 blk[is].
IA[i+1] = nnz;
591 for (i=0; i<nloc; ++i) {
void * fasp_mem_realloc(void *oldmem, const LONGLONG tsize)
Reallocate, initiate, and check memory.
void * fasp_mem_calloc(const unsigned int size, const unsigned int type)
Allocate, initiate, and check memory.
dvector fasp_dvec_create(const INT m)
Create dvector data space of REAL type.
void fasp_dvec_set(INT n, dvector *x, const REAL val)
Initialize dvector x[i]=val for i=0:n-1.
INT fasp_swz_dcsr_setup(SWZ_data *swzdata, SWZ_param *swzparam)
Setup phase for the Schwarz methods.
void fasp_dcsr_swz_forward(SWZ_data *swzdata, SWZ_param *swzparam, dvector *x, dvector *b)
Schwarz smoother: forward sweep.
void fasp_dcsr_swz_backward(SWZ_data *swzdata, SWZ_param *swzparam, dvector *x, dvector *b)
Schwarz smoother: backward sweep.
dCSRmat fasp_dcsr_create(const INT m, const INT n, const INT nnz)
Create CSR sparse matrix data memory space.
void fasp_dcsr_free(dCSRmat *A)
Free CSR sparse matrix data memory space.
void fasp_dcsr_transz(dCSRmat *A, INT *p, dCSRmat *AT)
Generalized transpose of A: (n x m) matrix given in dCSRmat format.
void fasp_dcsr_cp(const dCSRmat *A, dCSRmat *B)
copy a dCSRmat to a new one B=A
ivector fasp_sparse_mis(dCSRmat *A)
Get the maximal independet set of a CSR matrix.
INT fasp_solver_dcsr_pvgmres(dCSRmat *A, dvector *b, dvector *x, precond *pc, const REAL tol, const INT MaxIt, const SHORT restart, const SHORT StopType, const SHORT PrtLvl)
Right preconditioned GMRES method in which the restart parameter can be adaptively modified during it...
Main header file for the FASP project.
#define MAX(a, b)
Definition of max, min, abs.
#define FASP_SUCCESS
Definition of return status and error messages.
#define PRINT_NONE
Print level for all subroutines – not including DEBUG output.
Data for MUMPS interface.
Data for Schwarz methods.
INT maxbs
maximal block size
dCSRmat A
pointer to the original coefficient matrix
dvector rhsloc1
local right hand side
dvector xloc1
local solution
SWZ_param * swzparam
param for Schwarz
dCSRmat * blk_data
matrix for each partition
INT SWZ_type
Schwarz method type.
INT * jblock
column index of blocks
Mumps_data * mumps
param for MUMPS
INT * iblock
row index of blocks
Parameters for Schwarz method.
INT SWZ_maxlvl
maximal level for constructing the blocks
INT SWZ_blksolver
type of Schwarz block solver
SHORT SWZ_type
type for Schwarz method
Sparse matrix of REAL type in CSR format.
REAL * val
nonzero entries of A
INT row
row number of matrix A, m
INT * IA
integer array of row pointers, the size is m+1
INT nnz
number of nonzero entries
INT * JA
integer array of column indexes, the size is nnz
Vector with n entries of REAL type.
REAL * val
actual vector entries
Vector with n entries of INT type.
INT * val
actual vector entries