xref: /petsc/src/ksp/pc/impls/bddc/bddcscalingbasic.c (revision ab5c6b0cbfbda615248e56bcc3c683a1052cebab)
1*ab5c6b0cSJed Brown #include <../src/ksp/pc/impls/bddc/bddc.h>
2*ab5c6b0cSJed Brown #include <../src/ksp/pc/impls/bddc/bddcprivate.h>
3674ae819SStefano Zampini 
4674ae819SStefano Zampini /* prototypes for deluxe public functions */
5681e7c04SStefano Zampini #if 0
6674ae819SStefano Zampini extern PetscErrorCode PCBDDCScalingCreateDeluxe(PC);
7674ae819SStefano Zampini extern PetscErrorCode PCBDDCScalingDestroyDeluxe(PC);
8674ae819SStefano Zampini extern PetscErrorCode PCBDDCScalingSetUpDeluxe(PC);
9681e7c04SStefano Zampini #endif
10674ae819SStefano Zampini 
11674ae819SStefano Zampini #undef __FUNCT__
12674ae819SStefano Zampini #define __FUNCT__ "PCBDDCScalingExtension_Basic"
13674ae819SStefano Zampini static PetscErrorCode PCBDDCScalingExtension_Basic(PC pc, Vec local_interface_vector, Vec global_vector)
14674ae819SStefano Zampini {
15674ae819SStefano Zampini   PC_IS* pcis = (PC_IS*)pc->data;
16674ae819SStefano Zampini   PC_BDDC* pcbddc = (PC_BDDC*)pc->data;
17674ae819SStefano Zampini   PetscErrorCode ierr;
18674ae819SStefano Zampini 
19674ae819SStefano Zampini   PetscFunctionBegin;
20674ae819SStefano Zampini   /* Apply partition of unity */
21674ae819SStefano Zampini   ierr = VecPointwiseMult(pcbddc->work_scaling,pcis->D,local_interface_vector);CHKERRQ(ierr);
22674ae819SStefano Zampini   ierr = VecSet(global_vector,0.0);CHKERRQ(ierr);
23674ae819SStefano Zampini   ierr = VecScatterBegin(pcis->global_to_B,pcbddc->work_scaling,global_vector,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
24674ae819SStefano Zampini   ierr = VecScatterEnd(pcis->global_to_B,pcbddc->work_scaling,global_vector,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
25674ae819SStefano Zampini   PetscFunctionReturn(0);
26674ae819SStefano Zampini }
27674ae819SStefano Zampini 
28681e7c04SStefano Zampini #if 0
29674ae819SStefano Zampini #undef __FUNCT__
30674ae819SStefano Zampini #define __FUNCT__ "PCBDDCScalingExtension_Deluxe"
31674ae819SStefano Zampini static PetscErrorCode PCBDDCScalingExtension_Deluxe(PC pc, Vec x, Vec y)
32674ae819SStefano Zampini {
33674ae819SStefano Zampini   PC_IS*              pcis=(PC_IS*)pc->data;
34674ae819SStefano Zampini   PC_BDDC*            pcbddc=(PC_BDDC*)pc->data;
35674ae819SStefano Zampini   PCBDDCDeluxeScaling deluxe_ctx = pcbddc->deluxe_ctx;
36674ae819SStefano Zampini   PetscScalar         *array_x,*array_D,*array;
37674ae819SStefano Zampini   PetscScalar         zero=0.0;
38674ae819SStefano Zampini   PetscInt            i;
39674ae819SStefano Zampini   PetscMPIInt         color_rank;
40674ae819SStefano Zampini   PetscErrorCode      ierr;
41674ae819SStefano Zampini 
42674ae819SStefano Zampini   /* TODO CHECK STUFF RELATED WITH FAKE WORK */
43674ae819SStefano Zampini   PetscFunctionBegin;
44674ae819SStefano Zampini   ierr = VecSet(pcbddc->work_scaling,zero);CHKERRQ(ierr); /* needed by the fake work below */
45674ae819SStefano Zampini   /* scale deluxe vertices using diagonal scaling */
46674ae819SStefano Zampini   ierr = VecGetArray(x,&array_x);CHKERRQ(ierr);
47674ae819SStefano Zampini   ierr = VecGetArray(pcis->D,&array_D);CHKERRQ(ierr);
48674ae819SStefano Zampini   ierr = VecGetArray(pcbddc->work_scaling,&array);CHKERRQ(ierr);
49674ae819SStefano Zampini   for (i=0;i<deluxe_ctx->n_simple;i++) {
50674ae819SStefano Zampini     array[deluxe_ctx->idx_simple_B[i]] = array_x[deluxe_ctx->idx_simple_B[i]]*array_D[deluxe_ctx->idx_simple_B[i]];
51674ae819SStefano Zampini   }
52674ae819SStefano Zampini   ierr = VecRestoreArray(pcbddc->work_scaling,&array);CHKERRQ(ierr);
53674ae819SStefano Zampini   ierr = VecRestoreArray(pcis->D,&array_D);CHKERRQ(ierr);
54674ae819SStefano Zampini   ierr = VecRestoreArray(x,&array_x);CHKERRQ(ierr);
55674ae819SStefano Zampini   /* sequential part : all problems and Schur applications collapsed into seq_mat at setup phase */
56674ae819SStefano Zampini   if (deluxe_ctx->seq_mat) {
57674ae819SStefano Zampini     ierr = VecScatterBegin(deluxe_ctx->seq_scctx,x,deluxe_ctx->seq_work1,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
58674ae819SStefano Zampini     ierr = VecScatterEnd(deluxe_ctx->seq_scctx,x,deluxe_ctx->seq_work1,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
59674ae819SStefano Zampini     ierr = MatMult(deluxe_ctx->seq_mat,deluxe_ctx->seq_work1,deluxe_ctx->seq_work2);CHKERRQ(ierr);
60674ae819SStefano Zampini     ierr = KSPSolve(deluxe_ctx->seq_ksp,deluxe_ctx->seq_work2,deluxe_ctx->seq_work1);CHKERRQ(ierr);
61674ae819SStefano Zampini     /* fake work due to final ADD VALUES and vertices scaling needed? TODO: check it */
62674ae819SStefano Zampini     ierr = VecScatterBegin(deluxe_ctx->seq_scctx,deluxe_ctx->seq_work1,pcbddc->work_scaling,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
63674ae819SStefano Zampini     ierr = VecScatterEnd(deluxe_ctx->seq_scctx,deluxe_ctx->seq_work1,pcbddc->work_scaling,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
64674ae819SStefano Zampini   }
65674ae819SStefano Zampini   /* parallel part */
66674ae819SStefano Zampini   for (i=0;i<deluxe_ctx->par_colors;i++) {
67674ae819SStefano Zampini     if (deluxe_ctx->par_ksp[i]) {
68674ae819SStefano Zampini       ierr = MPI_Comm_rank(deluxe_ctx->par_subcomm[i]->comm,&color_rank);CHKERRQ(ierr);
69674ae819SStefano Zampini       ierr = VecSet(deluxe_ctx->work1_B,zero);CHKERRQ(ierr);
70674ae819SStefano Zampini       ierr = VecScatterBegin(deluxe_ctx->par_scctx_s[i],x,deluxe_ctx->work1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
71674ae819SStefano Zampini       ierr = VecScatterEnd(deluxe_ctx->par_scctx_s[i],x,deluxe_ctx->work1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
72674ae819SStefano Zampini       /* apply local schur on subset S_j^-1 */
73674ae819SStefano Zampini       ierr = PCBDDCApplySchur(pc,deluxe_ctx->work1_B,deluxe_ctx->work2_B,(Vec)0,deluxe_ctx->work1_D,deluxe_ctx->work2_D);CHKERRQ(ierr);
74674ae819SStefano Zampini       /* parallel transpose solve (\sum_j S_j)^-1 */
75674ae819SStefano Zampini       ierr = VecSet(deluxe_ctx->par_vec[i],zero);CHKERRQ(ierr);
76674ae819SStefano Zampini       ierr = VecScatterBegin(deluxe_ctx->par_scctx_p[i],deluxe_ctx->work2_B,deluxe_ctx->par_vec[i],ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
77674ae819SStefano Zampini       ierr = VecScatterEnd(deluxe_ctx->par_scctx_p[i],deluxe_ctx->work2_B,deluxe_ctx->par_vec[i],ADD_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
78674ae819SStefano Zampini       ierr = KSPSolve(deluxe_ctx->par_ksp[i],deluxe_ctx->par_vec[i],deluxe_ctx->par_vec[i]);CHKERRQ(ierr);
79674ae819SStefano Zampini       if (!color_rank) {
80674ae819SStefano Zampini         ierr = VecScatterBegin(deluxe_ctx->par_scctx_p[i],deluxe_ctx->par_vec[i],pcbddc->work_scaling,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
81674ae819SStefano Zampini         ierr = VecScatterEnd(deluxe_ctx->par_scctx_p[i],deluxe_ctx->par_vec[i],pcbddc->work_scaling,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
82674ae819SStefano Zampini       } else { /* fake work due to final ADD VALUES and vertices scaling */
83674ae819SStefano Zampini         ierr = VecScatterBegin(deluxe_ctx->par_scctx_p[i],deluxe_ctx->par_vec[i],deluxe_ctx->work1_B,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
84674ae819SStefano Zampini         ierr = VecScatterEnd(deluxe_ctx->par_scctx_p[i],deluxe_ctx->par_vec[i],deluxe_ctx->work1_B,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
85674ae819SStefano Zampini       }
86674ae819SStefano Zampini     }
87674ae819SStefano Zampini   }
88674ae819SStefano Zampini   /* put local boundary part in global vector */
89674ae819SStefano Zampini   ierr = VecSet(y,zero);CHKERRQ(ierr);
90674ae819SStefano Zampini   ierr = VecScatterBegin(pcis->global_to_B,pcbddc->work_scaling,y,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
91674ae819SStefano Zampini   ierr = VecScatterEnd  (pcis->global_to_B,pcbddc->work_scaling,y,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
92674ae819SStefano Zampini   PetscFunctionReturn(0);
93674ae819SStefano Zampini }
94681e7c04SStefano Zampini #endif
95674ae819SStefano Zampini 
96674ae819SStefano Zampini #undef __FUNCT__
97674ae819SStefano Zampini #define __FUNCT__ "PCBDDCScalingExtension"
98674ae819SStefano Zampini PetscErrorCode PCBDDCScalingExtension(PC pc, Vec local_interface_vector, Vec global_vector)
99674ae819SStefano Zampini {
100674ae819SStefano Zampini   PC_BDDC *pcbddc=(PC_BDDC*)pc->data;
101674ae819SStefano Zampini   PetscErrorCode ierr;
102674ae819SStefano Zampini 
103674ae819SStefano Zampini   PetscFunctionBegin;
104674ae819SStefano Zampini   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
105674ae819SStefano Zampini   PetscValidHeaderSpecific(local_interface_vector,VEC_CLASSID,2);
106674ae819SStefano Zampini   PetscValidHeaderSpecific(global_vector,VEC_CLASSID,3);
107674ae819SStefano Zampini   if (local_interface_vector == pcbddc->work_scaling) {
108674ae819SStefano Zampini     SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Local vector cannot be pcbddc->work_scaling!\n");
109674ae819SStefano Zampini   }
110674ae819SStefano Zampini   ierr = PetscTryMethod(pc,"PCBDDCScalingExtension_C",(PC,Vec,Vec),(pc,local_interface_vector,global_vector));CHKERRQ(ierr);
111674ae819SStefano Zampini   PetscFunctionReturn(0);
112674ae819SStefano Zampini }
113674ae819SStefano Zampini 
114674ae819SStefano Zampini #undef __FUNCT__
115674ae819SStefano Zampini #define __FUNCT__ "PCBDDCScalingRestriction_Basic"
116674ae819SStefano Zampini static PetscErrorCode PCBDDCScalingRestriction_Basic(PC pc, Vec global_vector, Vec local_interface_vector)
117674ae819SStefano Zampini {
118674ae819SStefano Zampini   PetscErrorCode ierr;
119674ae819SStefano Zampini   PC_IS* pcis = (PC_IS*)pc->data;
120674ae819SStefano Zampini 
121674ae819SStefano Zampini   PetscFunctionBegin;
122674ae819SStefano Zampini   ierr = VecScatterBegin(pcis->global_to_B,global_vector,local_interface_vector,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
123674ae819SStefano Zampini   ierr = VecScatterEnd(pcis->global_to_B,global_vector,local_interface_vector,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
124674ae819SStefano Zampini   /* Apply partition of unity */
125674ae819SStefano Zampini   ierr = VecPointwiseMult(local_interface_vector,pcis->D,local_interface_vector);CHKERRQ(ierr);
126674ae819SStefano Zampini   PetscFunctionReturn(0);
127674ae819SStefano Zampini }
128674ae819SStefano Zampini 
129681e7c04SStefano Zampini #if 0
130674ae819SStefano Zampini #undef __FUNCT__
131674ae819SStefano Zampini #define __FUNCT__ "PCBDDCScalingRestriction_Deluxe"
132674ae819SStefano Zampini static PetscErrorCode PCBDDCScalingRestriction_Deluxe(PC pc, Vec x, Vec y)
133674ae819SStefano Zampini {
134674ae819SStefano Zampini   PC_IS*              pcis=(PC_IS*)pc->data;
135674ae819SStefano Zampini   PC_BDDC*            pcbddc=(PC_BDDC*)pc->data;
136674ae819SStefano Zampini   PCBDDCDeluxeScaling deluxe_ctx = pcbddc->deluxe_ctx;
137674ae819SStefano Zampini   PetscScalar         *array_y,*array_D,zero=0.0;
138674ae819SStefano Zampini   PetscInt            i;
139674ae819SStefano Zampini   PetscErrorCode      ierr;
140674ae819SStefano Zampini 
141674ae819SStefano Zampini   PetscFunctionBegin;
142674ae819SStefano Zampini   /* get local boundary part of global vector */
143674ae819SStefano Zampini   ierr = VecScatterBegin(pcis->global_to_B,x,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
144674ae819SStefano Zampini   ierr = VecScatterEnd  (pcis->global_to_B,x,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
145674ae819SStefano Zampini   /* scale vertices using diagonal scaling -> every scaling perform the same */
146674ae819SStefano Zampini   ierr = VecGetArray(y,&array_y);CHKERRQ(ierr);
147674ae819SStefano Zampini   ierr = VecGetArray(pcis->D,&array_D);CHKERRQ(ierr);
148674ae819SStefano Zampini   for (i=0;i<deluxe_ctx->n_simple;i++) {
149674ae819SStefano Zampini     array_y[deluxe_ctx->idx_simple_B[i]] *= array_D[deluxe_ctx->idx_simple_B[i]];
150674ae819SStefano Zampini   }
151674ae819SStefano Zampini   ierr = VecRestoreArray(pcis->D,&array_D);CHKERRQ(ierr);
152674ae819SStefano Zampini   ierr = VecRestoreArray(y,&array_y);CHKERRQ(ierr);
153674ae819SStefano Zampini   /* sequential part : all problems and Schur applications collapsed into seq_mat at setup phase */
154674ae819SStefano Zampini   if (deluxe_ctx->seq_mat) {
155674ae819SStefano Zampini     ierr = VecScatterBegin(deluxe_ctx->seq_scctx,y,deluxe_ctx->seq_work1,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
156674ae819SStefano Zampini     ierr = VecScatterEnd(deluxe_ctx->seq_scctx,y,deluxe_ctx->seq_work1,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
157674ae819SStefano Zampini     ierr = KSPSolveTranspose(deluxe_ctx->seq_ksp,deluxe_ctx->seq_work1,deluxe_ctx->seq_work2);CHKERRQ(ierr);
158674ae819SStefano Zampini     ierr = MatMultTranspose(deluxe_ctx->seq_mat,deluxe_ctx->seq_work2,deluxe_ctx->seq_work1);CHKERRQ(ierr);
159674ae819SStefano Zampini     ierr = VecScatterBegin(deluxe_ctx->seq_scctx,deluxe_ctx->seq_work1,y,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
160674ae819SStefano Zampini     ierr = VecScatterEnd(deluxe_ctx->seq_scctx,deluxe_ctx->seq_work1,y,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
161674ae819SStefano Zampini   }
162674ae819SStefano Zampini   /* parallel part */
163674ae819SStefano Zampini   for (i=0;i<deluxe_ctx->par_colors;i++) {
164674ae819SStefano Zampini     if (deluxe_ctx->par_ksp[i]) {
165674ae819SStefano Zampini       /* parallel solve (\sum_j S_j)^-1 */
166674ae819SStefano Zampini       ierr = VecSet(deluxe_ctx->par_vec[i],zero);CHKERRQ(ierr);
167674ae819SStefano Zampini       ierr = VecScatterBegin(deluxe_ctx->par_scctx_p[i],y,deluxe_ctx->par_vec[i],INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
168674ae819SStefano Zampini       ierr = VecScatterEnd(deluxe_ctx->par_scctx_p[i],y,deluxe_ctx->par_vec[i],INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
169674ae819SStefano Zampini       ierr = KSPSolveTranspose(deluxe_ctx->par_ksp[i],deluxe_ctx->par_vec[i],deluxe_ctx->par_vec[i]);CHKERRQ(ierr);
170674ae819SStefano Zampini       /* apply local schur S_j^-1 */
171674ae819SStefano Zampini       ierr = VecSet(deluxe_ctx->work1_B,zero);CHKERRQ(ierr);
172674ae819SStefano Zampini       ierr = VecScatterBegin(deluxe_ctx->par_scctx_p[i],deluxe_ctx->par_vec[i],deluxe_ctx->work1_B,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
173674ae819SStefano Zampini       ierr = VecScatterEnd(deluxe_ctx->par_scctx_p[i],deluxe_ctx->par_vec[i],deluxe_ctx->work1_B,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
174674ae819SStefano Zampini       ierr = PCBDDCApplySchurTranspose(pc,deluxe_ctx->work1_B,deluxe_ctx->work2_B,(Vec)0,deluxe_ctx->work1_D,deluxe_ctx->work2_D);CHKERRQ(ierr);
175674ae819SStefano Zampini       ierr = VecScatterBegin(deluxe_ctx->par_scctx_s[i],deluxe_ctx->work2_B,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
176674ae819SStefano Zampini       ierr = VecScatterEnd(deluxe_ctx->par_scctx_s[i],deluxe_ctx->work2_B,y,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
177674ae819SStefano Zampini     }
178674ae819SStefano Zampini   }
179674ae819SStefano Zampini   PetscFunctionReturn(0);
180674ae819SStefano Zampini }
181681e7c04SStefano Zampini #endif
182674ae819SStefano Zampini 
183674ae819SStefano Zampini #undef __FUNCT__
184674ae819SStefano Zampini #define __FUNCT__ "PCBDDCScalingRestriction"
185674ae819SStefano Zampini PetscErrorCode PCBDDCScalingRestriction(PC pc, Vec global_vector, Vec local_interface_vector)
186674ae819SStefano Zampini {
187674ae819SStefano Zampini   PC_BDDC        *pcbddc=(PC_BDDC*)pc->data;
188674ae819SStefano Zampini   PetscErrorCode ierr;
189674ae819SStefano Zampini 
190674ae819SStefano Zampini   PetscFunctionBegin;
191674ae819SStefano Zampini   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
192674ae819SStefano Zampini   PetscValidHeaderSpecific(global_vector,VEC_CLASSID,2);
193674ae819SStefano Zampini   PetscValidHeaderSpecific(local_interface_vector,VEC_CLASSID,3);
194674ae819SStefano Zampini   if (local_interface_vector == pcbddc->work_scaling) {
195674ae819SStefano Zampini     SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Local vector should cannot be pcbddc->work_scaling!\n");
196674ae819SStefano Zampini   }
197674ae819SStefano Zampini   ierr = PetscTryMethod(pc,"PCBDDCScalingRestriction_C",(PC,Vec,Vec),(pc,global_vector,local_interface_vector));CHKERRQ(ierr);
198674ae819SStefano Zampini   PetscFunctionReturn(0);
199674ae819SStefano Zampini }
200674ae819SStefano Zampini 
201674ae819SStefano Zampini #undef __FUNCT__
202674ae819SStefano Zampini #define __FUNCT__ "PCBDDCScalingSetUp"
203674ae819SStefano Zampini PetscErrorCode PCBDDCScalingSetUp(PC pc)
204674ae819SStefano Zampini {
205674ae819SStefano Zampini   PC_IS* pcis=(PC_IS*)pc->data;
206674ae819SStefano Zampini   PC_BDDC* pcbddc=(PC_BDDC*)pc->data;
207674ae819SStefano Zampini   PetscErrorCode ierr;
208674ae819SStefano Zampini 
209674ae819SStefano Zampini   PetscFunctionBegin;
210674ae819SStefano Zampini   PetscValidHeaderSpecific(pc,PC_CLASSID,1);
211674ae819SStefano Zampini   /* create work vector for the operator */
212674ae819SStefano Zampini   ierr = VecDuplicate(pcis->vec1_B,&pcbddc->work_scaling);CHKERRQ(ierr);
21328d874f6SStefano Zampini   /* rebuild pcis->D if stiffness scaling has been requested */
21428d874f6SStefano Zampini   if (pcis->use_stiffness_scaling) {
215674ae819SStefano Zampini     ierr = MatGetDiagonal(pcbddc->local_mat,pcis->vec1_N);CHKERRQ(ierr);
216674ae819SStefano Zampini     ierr = VecScatterBegin(pcis->N_to_B,pcis->vec1_N,pcis->D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
217674ae819SStefano Zampini     ierr = VecScatterEnd(pcis->N_to_B,pcis->vec1_N,pcis->D,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
218674ae819SStefano Zampini   }
219674ae819SStefano Zampini   ierr = VecCopy(pcis->D,pcis->vec1_B);CHKERRQ(ierr);
220674ae819SStefano Zampini   ierr = VecSet(pcis->vec1_global,0.0);CHKERRQ(ierr);
221674ae819SStefano Zampini   ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_B,pcis->vec1_global,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
222674ae819SStefano Zampini   ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_B,pcis->vec1_global,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
223674ae819SStefano Zampini   ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_global,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
224674ae819SStefano Zampini   ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_global,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
225674ae819SStefano Zampini   ierr = VecPointwiseDivide(pcis->D,pcis->D,pcis->vec1_B);CHKERRQ(ierr);
226674ae819SStefano Zampini   /* now setup */
227681e7c04SStefano Zampini #if 0
228681e7c04SStefano Zampini   if (pcbddc->use_deluxe_scaling) {
229674ae819SStefano Zampini     ierr = PCBDDCScalingCreateDeluxe(pc);CHKERRQ(ierr);
230674ae819SStefano Zampini     ierr = PCBDDCScalingSetUpDeluxe(pc);CHKERRQ(ierr);
231674ae819SStefano Zampini     ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCScalingRestriction_C",PCBDDCScalingRestriction_Deluxe);CHKERRQ(ierr);
232674ae819SStefano Zampini     ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCScalingExtension_C",PCBDDCScalingExtension_Deluxe);CHKERRQ(ierr);
233674ae819SStefano Zampini   } else {
234681e7c04SStefano Zampini #endif
235674ae819SStefano Zampini     ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCScalingRestriction_C",PCBDDCScalingRestriction_Basic);CHKERRQ(ierr);
236674ae819SStefano Zampini     ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCScalingExtension_C",PCBDDCScalingExtension_Basic);CHKERRQ(ierr);
237681e7c04SStefano Zampini #if 0
238674ae819SStefano Zampini   }
239681e7c04SStefano Zampini #endif
240674ae819SStefano Zampini   /* test */
241674ae819SStefano Zampini   if (pcbddc->dbg_flag) {
242674ae819SStefano Zampini     PetscViewer viewer=pcbddc->dbg_viewer;
243674ae819SStefano Zampini     PetscReal   error,gerror;
244674ae819SStefano Zampini     MPI_Comm    test_comm;
245674ae819SStefano Zampini 
246674ae819SStefano Zampini     /* extension -> from local to parallel */
247674ae819SStefano Zampini     ierr = PetscObjectGetComm((PetscObject)pc,&test_comm);CHKERRQ(ierr);
248674ae819SStefano Zampini     ierr = VecSetRandom(pcis->vec1_global,NULL);CHKERRQ(ierr);
249674ae819SStefano Zampini     ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_global,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
250674ae819SStefano Zampini     ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_global,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
251674ae819SStefano Zampini     ierr = PCBDDCScalingExtension(pc,pcis->vec1_B,pcis->vec1_global);CHKERRQ(ierr);
252674ae819SStefano Zampini     ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_global,pcbddc->work_scaling,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
253674ae819SStefano Zampini     ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_global,pcbddc->work_scaling,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
254674ae819SStefano Zampini     ierr = VecAXPY(pcis->vec1_B,-1.0,pcbddc->work_scaling);CHKERRQ(ierr);
255674ae819SStefano Zampini     ierr = VecNorm(pcis->vec1_B,NORM_INFINITY,&error);CHKERRQ(ierr);
256674ae819SStefano Zampini     ierr = MPI_Allreduce(&error,&gerror,1,MPIU_REAL,MPI_MAX,test_comm);CHKERRQ(ierr);
257674ae819SStefano Zampini     ierr = PetscViewerASCIIPrintf(viewer,"Error scaling extension %1.14e\n",error);CHKERRQ(ierr);
258674ae819SStefano Zampini     if (PetscAbsReal(gerror)>1.e-8) {
259674ae819SStefano Zampini       ierr = VecSet(pcis->vec1_global,0.0);CHKERRQ(ierr);
260674ae819SStefano Zampini       ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_B,pcis->vec1_global,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
261674ae819SStefano Zampini       ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_B,pcis->vec1_global,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
262674ae819SStefano Zampini       ierr = VecView(pcis->vec1_global,viewer);CHKERRQ(ierr);
263674ae819SStefano Zampini     }
264674ae819SStefano Zampini     /* restriction -> from parallel to local */
265674ae819SStefano Zampini     ierr = VecSetRandom(pcis->vec1_global,NULL);CHKERRQ(ierr);
266674ae819SStefano Zampini     ierr = PCBDDCScalingRestriction(pc,pcis->vec1_global,pcis->vec1_B);CHKERRQ(ierr);
267674ae819SStefano Zampini     ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_global,pcbddc->work_scaling,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
268674ae819SStefano Zampini     ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_global,pcbddc->work_scaling,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
269674ae819SStefano Zampini     ierr = VecSet(pcis->vec1_global,0.0);CHKERRQ(ierr);
270674ae819SStefano Zampini     ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_B,pcis->vec1_global,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
271674ae819SStefano Zampini     ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_B,pcis->vec1_global,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
272674ae819SStefano Zampini     ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_global,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
273674ae819SStefano Zampini     ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_global,pcis->vec1_B,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr);
274674ae819SStefano Zampini     ierr = VecAXPY(pcis->vec1_B,-1.0,pcbddc->work_scaling);CHKERRQ(ierr);
275674ae819SStefano Zampini     ierr = VecNorm(pcis->vec1_B,NORM_INFINITY,&error);CHKERRQ(ierr);
276674ae819SStefano Zampini     ierr = MPI_Allreduce(&error,&gerror,1,MPIU_REAL,MPI_MAX,test_comm);CHKERRQ(ierr);
277674ae819SStefano Zampini     ierr = PetscViewerASCIIPrintf(viewer,"Error scaling restriction %1.14e\n",gerror);CHKERRQ(ierr);
278674ae819SStefano Zampini     if (PetscAbsReal(gerror)>1.e-8) {
279674ae819SStefano Zampini       ierr = VecSet(pcis->vec1_global,0.0);CHKERRQ(ierr);
280674ae819SStefano Zampini       ierr = VecScatterBegin(pcis->global_to_B,pcis->vec1_B,pcis->vec1_global,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
281674ae819SStefano Zampini       ierr = VecScatterEnd(pcis->global_to_B,pcis->vec1_B,pcis->vec1_global,INSERT_VALUES,SCATTER_REVERSE);CHKERRQ(ierr);
282674ae819SStefano Zampini       ierr = VecView(pcis->vec1_global,viewer);CHKERRQ(ierr);
283674ae819SStefano Zampini     }
284674ae819SStefano Zampini   }
285674ae819SStefano Zampini   PetscFunctionReturn(0);
286674ae819SStefano Zampini }
287674ae819SStefano Zampini 
288674ae819SStefano Zampini #undef __FUNCT__
289674ae819SStefano Zampini #define __FUNCT__ "PCBDDCScalingDestroy"
290674ae819SStefano Zampini PetscErrorCode PCBDDCScalingDestroy(PC pc)
291674ae819SStefano Zampini {
292674ae819SStefano Zampini   PC_BDDC* pcbddc=(PC_BDDC*)pc->data;
293674ae819SStefano Zampini   PetscErrorCode ierr;
294674ae819SStefano Zampini 
295674ae819SStefano Zampini   PetscFunctionBegin;
296681e7c04SStefano Zampini #if 0
297681e7c04SStefano Zampini   if (pcbddc->use_deluxe_scaling) {
298674ae819SStefano Zampini     ierr = PCBDDCScalingDestroyDeluxe(pc);CHKERRQ(ierr);
299674ae819SStefano Zampini   }
300681e7c04SStefano Zampini #endif
301674ae819SStefano Zampini   ierr = VecDestroy(&pcbddc->work_scaling);CHKERRQ(ierr);
302674ae819SStefano Zampini   /* remove functions */
303674ae819SStefano Zampini   ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCScalingRestriction_C",NULL);CHKERRQ(ierr);
304674ae819SStefano Zampini   ierr = PetscObjectComposeFunction((PetscObject)pc,"PCBDDCScalingExtension_C",NULL);CHKERRQ(ierr);
305674ae819SStefano Zampini   PetscFunctionReturn(0);
306674ae819SStefano Zampini }
307674ae819SStefano Zampini 
308