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