New automatic determination of Jacobi weights. Also different weights

are used on different levels (pfmg).
This commit is contained in:
barrylee 2007-03-30 17:22:42 +00:00
parent 87f895f913
commit 956fcf613c
2 changed files with 51 additions and 51 deletions

View File

@ -42,7 +42,6 @@ typedef struct
void *relax_data;
void *rb_relax_data;
int relax_type;
int usr_jacobi_weight;
double jacobi_weight;
} hypre_PFMGRelaxData;
@ -60,7 +59,6 @@ hypre_PFMGRelaxCreate( MPI_Comm comm )
(pfmg_relax_data -> relax_data) = hypre_PointRelaxCreate(comm);
(pfmg_relax_data -> rb_relax_data) = hypre_RedBlackGSCreate(comm);
(pfmg_relax_data -> relax_type) = 0; /* Weighted Jacobi */
(pfmg_relax_data -> usr_jacobi_weight) = 0;
(pfmg_relax_data -> jacobi_weight) = 0.0;
return (void *) pfmg_relax_data;
@ -138,9 +136,7 @@ hypre_PFMGRelaxSetup( void *pfmg_relax_vdata,
{
hypre_PFMGRelaxData *pfmg_relax_data = pfmg_relax_vdata;
int relax_type = (pfmg_relax_data -> relax_type);
int usr_jacobi_weight= (pfmg_relax_data -> usr_jacobi_weight);
double jacobi_weight = (pfmg_relax_data -> jacobi_weight);
int ndim = hypre_StructGridDim(hypre_StructVectorGrid(x));
int ierr;
switch(relax_type)
@ -159,33 +155,9 @@ hypre_PFMGRelaxSetup( void *pfmg_relax_vdata,
if (relax_type==1)
{
if (usr_jacobi_weight)
{
hypre_PointRelaxSetWeight(pfmg_relax_data -> relax_data, jacobi_weight);
}
else /* weights dimensionally dependent */
{
switch(ndim)
{
case 1: /* Weighted Jacobi (weight = 2/3)- already set */
{
break;
}
case 2: /* Weighted Jacobi (weight = 0.80) */
{
hypre_PointRelaxSetWeight(pfmg_relax_data -> relax_data, 0.80);
break;
}
case 3: /* Weighted Jacobi (weight = 6/7) */
{
hypre_PointRelaxSetWeight(pfmg_relax_data -> relax_data, 0.857142857);
break;
}
}
}
hypre_PointRelaxSetWeight(pfmg_relax_data -> relax_data, jacobi_weight);
}
return ierr;
}
@ -201,29 +173,17 @@ hypre_PFMGRelaxSetType( void *pfmg_relax_vdata,
hypre_PFMGRelaxData *pfmg_relax_data = pfmg_relax_vdata;
void *relax_data = (pfmg_relax_data -> relax_data);
int ierr = 0;
double weight = (pfmg_relax_data->jacobi_weight);
int usr_jacobi_weight = (pfmg_relax_data -> usr_jacobi_weight);
(pfmg_relax_data -> relax_type) = relax_type;
hypre_PointRelaxSetWeight(relax_data, 1.0);
switch(relax_type)
{
case 1: /* Weighted Jacobi (weight = 2/3) */
{
if ( !usr_jacobi_weight )
{
weight = 0.6666667;
(pfmg_relax_data->jacobi_weight) = weight;
}
hypre_PointRelaxSetWeight(relax_data, weight);
}
case 0: /* Jacobi */
{
hypre_Index stride;
hypre_Index indices[1];
hypre_PointRelaxSetWeight(relax_data, 1.0);
hypre_PointRelaxSetNumPointsets(relax_data, 1);
hypre_SetIndex(stride, 1, 1, 1);
@ -253,7 +213,6 @@ hypre_PFMGRelaxSetJacobiWeight(void *pfmg_relax_vdata,
hypre_PFMGRelaxData *pfmg_relax_data = pfmg_relax_vdata;
(pfmg_relax_data -> jacobi_weight) = weight;
(pfmg_relax_data -> usr_jacobi_weight)= 1;
return hypre_error_flag;
}

View File

@ -94,6 +94,9 @@ hypre_PFMGSetup( void *pfmg_vdata,
double *data;
int data_size = 0;
double *relax_weights;
double alpha;
hypre_StructMatrix **A_l;
hypre_StructMatrix **P_l;
hypre_StructMatrix **RT_l;
@ -159,12 +162,14 @@ hypre_PFMGSetup( void *pfmg_vdata,
P_grid_l[0] = NULL;
cdir_l = hypre_TAlloc(int, max_levels);
active_l = hypre_TAlloc(int, max_levels);
relax_weights= hypre_CTAlloc(double, max_levels);
hypre_SetIndex(coarsen, 1, 1, 1); /* forces relaxation on finest grid */
for (l = 0; ; l++)
{
/* determine cdir */
min_dxyz = dxyz[0] + dxyz[1] + dxyz[2] + 1;
cdir = -1;
alpha= 0.0;
for (d = 0; d < dim; d++)
{
if ((hypre_BoxIMaxD(cbox, d) > hypre_BoxIMinD(cbox, d)) &&
@ -173,6 +178,29 @@ hypre_PFMGSetup( void *pfmg_vdata,
min_dxyz = dxyz[d];
cdir = d;
}
alpha+= 1.0/(dxyz[d]*dxyz[d]);
}
if (cdir != -1)
{
for (d = 0; d < dim; d++)
{
if (d != cdir)
{
relax_weights[l]+= 1.0/(dxyz[d]*dxyz[d]);
}
}
alpha= relax_weights[l]/alpha;
/* determine level Jacobi weights */
if (dim > 1)
{
relax_weights[l]= 2.0/(3.0 - alpha);
}
else
{
relax_weights[l]= 2.0/3.0; /* always 2/3 for 1-d */
}
}
/* if cannot coarsen in any direction, stop */
@ -231,10 +259,10 @@ hypre_PFMGSetup( void *pfmg_vdata,
}
}
(pfmg_data -> num_levels) = num_levels;
(pfmg_data -> cdir_l) = cdir_l;
(pfmg_data -> grid_l) = grid_l;
(pfmg_data -> P_grid_l) = P_grid_l;
(pfmg_data -> num_levels) = num_levels;
(pfmg_data -> cdir_l) = cdir_l;
(pfmg_data -> grid_l) = grid_l;
(pfmg_data -> P_grid_l) = P_grid_l;
/*-----------------------------------------------------
* Set up matrix and vector structures
@ -424,6 +452,10 @@ hypre_PFMGSetup( void *pfmg_vdata,
{
hypre_PFMGRelaxSetJacobiWeight(relax_data_l[0], jacobi_weight);
}
else
{
hypre_PFMGRelaxSetJacobiWeight(relax_data_l[0], relax_weights[0]);
}
hypre_PFMGRelaxSetType(relax_data_l[0], relax_type);
hypre_PFMGRelaxSetTempVec(relax_data_l[0], tx_l[0]);
hypre_PFMGRelaxSetup(relax_data_l[0], A_l[0], b_l[0], x_l[0]);
@ -440,6 +472,10 @@ hypre_PFMGSetup( void *pfmg_vdata,
{
hypre_PFMGRelaxSetJacobiWeight(relax_data_l[l], jacobi_weight);
}
else
{
hypre_PFMGRelaxSetJacobiWeight(relax_data_l[l], relax_weights[l]);
}
hypre_PFMGRelaxSetType(relax_data_l[l], relax_type);
hypre_PFMGRelaxSetTempVec(relax_data_l[l], tx_l[l]);
hypre_PFMGRelaxSetup(relax_data_l[l], A_l[l], b_l[l], x_l[l]);
@ -455,11 +491,16 @@ hypre_PFMGSetup( void *pfmg_vdata,
{
hypre_PFMGRelaxSetJacobiWeight(relax_data_l[l], jacobi_weight);
}
else
{
hypre_PFMGRelaxSetJacobiWeight(relax_data_l[l], relax_weights[l]);
}
hypre_PFMGRelaxSetType(relax_data_l[l], 0);
hypre_PFMGRelaxSetTempVec(relax_data_l[l], tx_l[l]);
hypre_PFMGRelaxSetup(relax_data_l[l], A_l[l], b_l[l], x_l[l]);
}
}
hypre_TFree(relax_weights);
for (l = 0; l < num_levels; l++)
{
@ -754,9 +795,9 @@ hypre_PFMGComputeDxyz( hypre_StructMatrix *A,
}
}
cx += tcxyz[0];
cy += tcxyz[1];
cz += tcxyz[2];
cx += tcxyz[0]/2.0;
cy += tcxyz[1]/2.0;
cz += tcxyz[2]/2.0;
}
hypre_BoxLoop1End(Ai);
}