1 #include <petsc-private/snesimpl.h> /*I "petscsnes.h" I*/ 2 3 typedef struct { 4 PetscInt sweeps; 5 PetscBool norms; 6 } SNES_GS; 7 8 #undef __FUNCT__ 9 #define __FUNCT__ "SNESReset_GS" 10 PetscErrorCode SNESReset_GS(SNES snes) 11 { 12 PetscFunctionBegin; 13 PetscFunctionReturn(0); 14 } 15 16 #undef __FUNCT__ 17 #define __FUNCT__ "SNESDestroy_GS" 18 PetscErrorCode SNESDestroy_GS(SNES snes) 19 { 20 PetscErrorCode ierr; 21 22 PetscFunctionBegin; 23 ierr = SNESReset_GS(snes);CHKERRQ(ierr); 24 ierr = PetscFree(snes->data); 25 PetscFunctionReturn(0); 26 } 27 28 #undef __FUNCT__ 29 #define __FUNCT__ "SNESSetUp_GS" 30 PetscErrorCode SNESSetUp_GS(SNES snes) 31 { 32 PetscFunctionBegin; 33 PetscFunctionReturn(0); 34 } 35 36 #undef __FUNCT__ 37 #define __FUNCT__ "SNESSetFromOptions_GS" 38 PetscErrorCode SNESSetFromOptions_GS(SNES snes) 39 { 40 PetscErrorCode ierr; 41 SNES_GS *gs = (SNES_GS *)snes->data; 42 43 PetscFunctionBegin; 44 ierr = PetscOptionsHead("SNES GS options");CHKERRQ(ierr); 45 ierr = PetscOptionsBool("-snes_gs_norms","Compute norms for monitoring","none",gs->norms,&gs->norms,PETSC_NULL);CHKERRQ(ierr); 46 ierr = PetscOptionsTail();CHKERRQ(ierr); 47 PetscFunctionReturn(0); 48 } 49 50 #undef __FUNCT__ 51 #define __FUNCT__ "SNESView_GS" 52 PetscErrorCode SNESView_GS(SNES snes, PetscViewer viewer) 53 { 54 PetscFunctionBegin; 55 PetscFunctionReturn(0); 56 } 57 58 #undef __FUNCT__ 59 #define __FUNCT__ "SNESSolve_GS" 60 PetscErrorCode SNESSolve_GS(SNES snes) 61 { 62 Vec F; 63 Vec X; 64 Vec B; 65 PetscInt i; 66 PetscReal fnorm; 67 PetscErrorCode ierr; 68 SNES_GS *gs = (SNES_GS *)snes->data; 69 70 PetscFunctionBegin; 71 72 X = snes->vec_sol; 73 F = snes->vec_func; 74 B = snes->vec_rhs; 75 76 fnorm = 1.; 77 78 snes->reason = SNES_CONVERGED_ITERATING; 79 /* compute the initial function and preconditioned update delX */ 80 if (!snes->vec_func_init_set) { 81 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 82 if (snes->domainerror) { 83 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 84 PetscFunctionReturn(0); 85 } 86 } else { 87 snes->vec_func_init_set = PETSC_FALSE; 88 } 89 90 /* convergence test */ 91 if (!snes->norm_init_set) { 92 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); /* fnorm <- ||F|| */ 93 if (PetscIsInfOrNanReal(fnorm)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"Infinite or not-a-number generated in norm"); 94 } else { 95 fnorm = snes->norm_init; 96 snes->norm_init_set = PETSC_FALSE; 97 } 98 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 99 snes->iter = 0; 100 snes->norm = fnorm; 101 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 102 SNESLogConvHistory(snes,snes->norm,0); 103 ierr = SNESMonitor(snes,0,snes->norm);CHKERRQ(ierr); 104 105 /* set parameter for default relative tolerance convergence test */ 106 snes->ttol = fnorm*snes->rtol; 107 /* test convergence */ 108 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 109 if (snes->reason) PetscFunctionReturn(0); 110 111 /* Call general purpose update function */ 112 if (snes->ops->update) { 113 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 114 } 115 for (i = 0; i < snes->max_its; i++) { 116 ierr = SNESComputeGS(snes, B, X);CHKERRQ(ierr); 117 /* only compute norms if requested or about to exit due to maximum iterations */ 118 if (gs->norms || i == snes->max_its - 1) { 119 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 120 if (snes->domainerror) { 121 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 122 PetscFunctionReturn(0); 123 } 124 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); /* fnorm <- ||F|| */ 125 if (PetscIsInfOrNanReal(fnorm)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"Infinite or not-a-number generated in norm"); 126 /* Monitor convergence */ 127 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 128 snes->iter = i+1; 129 snes->norm = fnorm; 130 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 131 } 132 SNESLogConvHistory(snes,snes->norm,0); 133 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 134 /* Test for convergence */ 135 ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 136 if (snes->reason) PetscFunctionReturn(0); 137 138 /* Call general purpose update function */ 139 if (snes->ops->update) { 140 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 141 } 142 } 143 if (!snes->reason) snes->reason = SNES_CONVERGED_ITS; /* GS is meant to be used as a preconditioner */ 144 PetscFunctionReturn(0); 145 } 146 147 /*MC 148 SNESGS - Just calls the user-provided solution routine provided with SNESSetGS() 149 150 Level: advanced 151 152 Notes: 153 the Gauss-Seidel smoother is inherited through composition. If a solver has been created with SNESGetPC(), it will have 154 its parent's Gauss-Seidel routine associated with it. 155 156 .seealso: SNESCreate(), SNES, SNESSetType(), SNESSetGS(), SNESType (for list of available types) 157 M*/ 158 159 EXTERN_C_BEGIN 160 #undef __FUNCT__ 161 #define __FUNCT__ "SNESCreate_GS" 162 PetscErrorCode SNESCreate_GS(SNES snes) 163 { 164 SNES_GS *gs; 165 PetscErrorCode ierr; 166 167 PetscFunctionBegin; 168 snes->ops->destroy = SNESDestroy_GS; 169 snes->ops->setup = SNESSetUp_GS; 170 snes->ops->setfromoptions = SNESSetFromOptions_GS; 171 snes->ops->view = SNESView_GS; 172 snes->ops->solve = SNESSolve_GS; 173 snes->ops->reset = SNESReset_GS; 174 175 snes->usesksp = PETSC_FALSE; 176 snes->usespc = PETSC_FALSE; 177 178 if (!snes->tolerancesset) { 179 snes->max_its = 10000; 180 snes->max_funcs = 10000; 181 } 182 183 ierr = PetscNewLog(snes, SNES_GS, &gs);CHKERRQ(ierr); 184 snes->data = (void*) gs; 185 gs->norms = PETSC_FALSE; 186 PetscFunctionReturn(0); 187 } 188 EXTERN_C_END 189