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