1 #include <petsc-private/snesimpl.h> /*I "petscsnes.h" I*/ 2 #include <petscdm.h> 3 4 typedef struct { 5 PetscInt n; /* local subdomains */ 6 SNES *subsnes; /* nonlinear solvers for each subdomain */ 7 8 Vec *x; /* solution vectors */ 9 Vec *xl; /* solution local vectors */ 10 Vec *y; /* step vectors */ 11 Vec *b; /* rhs vectors */ 12 13 VecScatter *oscatter; /* scatter from global space to the subdomain global space */ 14 VecScatter *iscatter; /* scatter from global space to the nonoverlapping subdomain space */ 15 VecScatter *gscatter; /* scatter from global space to the subdomain local space */ 16 17 PCASMType type; /* ASM type */ 18 19 PetscBool usesdm; /* use the DM for setting up the subproblems */ 20 } SNES_NASM; 21 22 #undef __FUNCT__ 23 #define __FUNCT__ "SNESReset_NASM" 24 PetscErrorCode SNESReset_NASM(SNES snes) 25 { 26 SNES_NASM *nasm = (SNES_NASM *)snes->data; 27 PetscErrorCode ierr; 28 PetscInt i; 29 30 PetscFunctionBegin; 31 for (i=0;i<nasm->n;i++) { 32 if (nasm->xl) { ierr = VecDestroy(&nasm->xl[i]);CHKERRQ(ierr); } 33 if (nasm->x) { ierr = VecDestroy(&nasm->x[i]);CHKERRQ(ierr); } 34 if (nasm->y) { ierr = VecDestroy(&nasm->y[i]);CHKERRQ(ierr); } 35 if (nasm->b) { ierr = VecDestroy(&nasm->b[i]);CHKERRQ(ierr); } 36 37 if (nasm->subsnes) { ierr = SNESDestroy(&nasm->subsnes[i]);CHKERRQ(ierr); } 38 if (nasm->oscatter) { ierr = VecScatterDestroy(&nasm->oscatter[i]);CHKERRQ(ierr); } 39 if (nasm->iscatter) { ierr = VecScatterDestroy(&nasm->iscatter[i]);CHKERRQ(ierr); } 40 if (nasm->gscatter) { ierr = VecScatterDestroy(&nasm->gscatter[i]);CHKERRQ(ierr); } 41 } 42 43 if (nasm->x) {ierr = PetscFree(nasm->x);CHKERRQ(ierr);} 44 if (nasm->xl) {ierr = PetscFree(nasm->xl);CHKERRQ(ierr);} 45 if (nasm->y) {ierr = PetscFree(nasm->y);CHKERRQ(ierr);} 46 if (nasm->b) {ierr = PetscFree(nasm->b);CHKERRQ(ierr);} 47 48 if (nasm->subsnes) {ierr = PetscFree(nasm->subsnes);CHKERRQ(ierr);} 49 if (nasm->oscatter) {ierr = PetscFree(nasm->oscatter);CHKERRQ(ierr);} 50 if (nasm->iscatter) {ierr = PetscFree(nasm->iscatter);CHKERRQ(ierr);} 51 if (nasm->gscatter) {ierr = PetscFree(nasm->gscatter);CHKERRQ(ierr);} 52 PetscFunctionReturn(0); 53 } 54 55 #undef __FUNCT__ 56 #define __FUNCT__ "SNESDestroy_NASM" 57 PetscErrorCode SNESDestroy_NASM(SNES snes) 58 { 59 PetscErrorCode ierr; 60 61 PetscFunctionBegin; 62 ierr = SNESReset_NASM(snes);CHKERRQ(ierr); 63 ierr = PetscFree(snes->data);CHKERRQ(ierr); 64 PetscFunctionReturn(0); 65 } 66 67 #undef __FUNCT__ 68 #define __FUNCT__ "DMGlobalToLocalSubDomainDirichletHook_Private" 69 PetscErrorCode DMGlobalToLocalSubDomainDirichletHook_Private(DM dm,Vec g,InsertMode mode,Vec l,void *ctx) 70 { 71 PetscErrorCode ierr; 72 Vec bcs = (Vec)ctx; 73 74 PetscFunctionBegin; 75 ierr = VecCopy(bcs,l);CHKERRQ(ierr); 76 PetscFunctionReturn(0); 77 } 78 79 #undef __FUNCT__ 80 #define __FUNCT__ "SNESSetUp_NASM" 81 PetscErrorCode SNESSetUp_NASM(SNES snes) 82 { 83 SNES_NASM *nasm = (SNES_NASM *)snes->data; 84 PetscErrorCode ierr; 85 DM dm,ddm; 86 DM *subdms; 87 PetscInt i; 88 const char *optionsprefix; 89 Vec F; 90 91 PetscFunctionBegin; 92 if (!nasm->subsnes) { 93 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 94 if (dm) { 95 nasm->usesdm = PETSC_TRUE; 96 ierr = DMCreateDomainDecomposition(dm,&nasm->n,PETSC_NULL,PETSC_NULL,PETSC_NULL,&subdms);CHKERRQ(ierr); 97 if (!subdms) { 98 ierr = DMCreateDomainDecompositionDM(dm,"default",&ddm);CHKERRQ(ierr); 99 if (!ddm) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE,"DM has no default decomposition defined. Set subsolves manually with SNESNASMSetSubdomains()."); 100 ierr = SNESSetDM(snes,ddm);CHKERRQ(ierr); 101 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 102 ierr = DMCreateDomainDecomposition(dm,&nasm->n,PETSC_NULL,PETSC_NULL,PETSC_NULL,&subdms);CHKERRQ(ierr); 103 } 104 if (!subdms) SETERRQ(((PetscObject)dm)->comm,PETSC_ERR_ARG_WRONGSTATE,"DM has no default decomposition defined. Set subsolves manually with SNESNASMSetSubdomains()."); 105 ierr = DMCreateDomainDecompositionScatters(dm,nasm->n,subdms,&nasm->iscatter,&nasm->oscatter,&nasm->gscatter);CHKERRQ(ierr); 106 107 ierr = SNESGetOptionsPrefix(snes, &optionsprefix);CHKERRQ(ierr); 108 ierr = PetscMalloc(nasm->n*sizeof(SNES),&nasm->subsnes);CHKERRQ(ierr); 109 110 for (i=0;i<nasm->n;i++) { 111 ierr = SNESCreate(PETSC_COMM_SELF,&nasm->subsnes[i]);CHKERRQ(ierr); 112 ierr = SNESAppendOptionsPrefix(nasm->subsnes[i],optionsprefix);CHKERRQ(ierr); 113 ierr = SNESAppendOptionsPrefix(nasm->subsnes[i],"sub_");CHKERRQ(ierr); 114 ierr = SNESSetDM(nasm->subsnes[i],subdms[i]);CHKERRQ(ierr); 115 ierr = SNESSetFromOptions(nasm->subsnes[i]);CHKERRQ(ierr); 116 ierr = DMDestroy(&subdms[i]);CHKERRQ(ierr); 117 } 118 ierr = PetscFree(subdms);CHKERRQ(ierr); 119 } else SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_ARG_WRONGSTATE,"Cannot construct local problems automatically without a DM!"); 120 } else SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_ARG_WRONGSTATE,"Must set subproblems manually if there is no DM!"); 121 /* allocate the global vectors */ 122 ierr = PetscMalloc(nasm->n*sizeof(Vec),&nasm->x);CHKERRQ(ierr); 123 ierr = PetscMalloc(nasm->n*sizeof(Vec),&nasm->xl);CHKERRQ(ierr); 124 ierr = PetscMalloc(nasm->n*sizeof(Vec),&nasm->y);CHKERRQ(ierr); 125 ierr = PetscMalloc(nasm->n*sizeof(Vec),&nasm->b);CHKERRQ(ierr); 126 127 for (i=0;i<nasm->n;i++) { 128 DM subdm; 129 ierr = SNESGetFunction(nasm->subsnes[i],&F,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr); 130 ierr = VecDuplicate(F,&nasm->x[i]);CHKERRQ(ierr); 131 ierr = VecDuplicate(F,&nasm->y[i]);CHKERRQ(ierr); 132 ierr = VecDuplicate(F,&nasm->b[i]);CHKERRQ(ierr); 133 ierr = SNESGetDM(nasm->subsnes[i],&subdm);CHKERRQ(ierr); 134 ierr = DMCreateLocalVector(subdm,&nasm->xl[i]);CHKERRQ(ierr); 135 ierr = DMGlobalToLocalHookAdd(subdm,DMGlobalToLocalSubDomainDirichletHook_Private,PETSC_NULL,nasm->xl[i]);CHKERRQ(ierr); 136 } 137 PetscFunctionReturn(0); 138 } 139 140 #undef __FUNCT__ 141 #define __FUNCT__ "SNESSetFromOptions_NASM" 142 PetscErrorCode SNESSetFromOptions_NASM(SNES snes) 143 { 144 PetscErrorCode ierr; 145 DM dm,ddm; 146 char ddm_name[1024]; 147 PCASMType asmtype; 148 const char *const SNESNASMTypes[] = {"NONE","RESTRICT","INTERPOLATE","BASIC","PCASMType","PC_ASM_",0}; 149 PetscBool flg; 150 SNES_NASM *nasm = (SNES_NASM *)snes->data; 151 152 PetscFunctionBegin; 153 ierr = PetscOptionsHead("Nonlinear Additive Schwartz options");CHKERRQ(ierr); 154 ierr = PetscOptionsEnum("-snes_nasm_type","Type of restriction/extension","",SNESNASMTypes,(PetscEnum)nasm->type,(PetscEnum*)&asmtype,&flg);CHKERRQ(ierr); 155 if (flg) {nasm->type = asmtype;} 156 ierr = PetscOptionsString("-snes_nasm_decomposition", "Name of the DM defining the composition", "SNESSetDM", ddm_name, ddm_name,1024,&flg);CHKERRQ(ierr); 157 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 158 if (flg) { 159 if (dm) { 160 ierr = DMCreateDomainDecompositionDM(dm, ddm_name, &ddm);CHKERRQ(ierr); 161 if (!ddm) SETERRQ1(((PetscObject)snes)->comm, PETSC_ERR_ARG_WRONGSTATE, "Unknown DM decomposition name %s", ddm_name); 162 ierr = PetscInfo(snes,"Using domain decomposition DM defined using options database\n");CHKERRQ(ierr); 163 ierr = SNESSetDM(snes,ddm);CHKERRQ(ierr); 164 } else SETERRQ(((PetscObject)snes)->comm, PETSC_ERR_ARG_WRONGSTATE, "No DM to decompose"); 165 } 166 ierr = PetscOptionsTail();CHKERRQ(ierr); 167 PetscFunctionReturn(0); 168 } 169 170 #undef __FUNCT__ 171 #define __FUNCT__ "SNESView_NASM" 172 PetscErrorCode SNESView_NASM(SNES snes, PetscViewer viewer) 173 { 174 PetscFunctionBegin; 175 PetscFunctionReturn(0); 176 } 177 178 #undef __FUNCT__ 179 #define __FUNCT__ "SNESNASMSetSubdomains" 180 PetscErrorCode SNESNASMSetSubdomains(SNES snes,PetscInt n,SNES subsnes[],VecScatter iscatter[],VecScatter oscatter[],VecScatter gscatter[]) 181 { 182 PetscErrorCode ierr; 183 PetscErrorCode (*f)(SNES,PetscInt,SNES*,VecScatter*,VecScatter*,VecScatter*); 184 185 PetscFunctionBegin; 186 ierr = PetscObjectQueryFunction((PetscObject)snes,"SNESNASMSetSubdomains_C",(void (**)(void))&f);CHKERRQ(ierr); 187 ierr = (f)(snes,n,subsnes,iscatter,oscatter,gscatter);CHKERRQ(ierr); 188 PetscFunctionReturn(0); 189 } 190 191 EXTERN_C_BEGIN 192 #undef __FUNCT__ 193 #define __FUNCT__ "SNESNASMSetSubdomains_NASM" 194 PetscErrorCode SNESNASMSetSubdomains_NASM(SNES snes,PetscInt n,SNES subsnes[],VecScatter iscatter[],VecScatter oscatter[],VecScatter gscatter[]) 195 { 196 PetscInt i; 197 PetscErrorCode ierr; 198 SNES_NASM *nasm = (SNES_NASM *)snes->data; 199 200 PetscFunctionBegin; 201 if (snes->setupcalled)SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_ARG_WRONGSTATE,"SNESNASMSetSubdomains() should be called before calling SNESSetUp()."); 202 203 /* tear down the previously set things */ 204 ierr = SNESReset(snes);CHKERRQ(ierr); 205 206 nasm->n = n; 207 if (oscatter) { 208 for (i=0; i<n; i++) {ierr = PetscObjectReference((PetscObject)oscatter[i]);CHKERRQ(ierr);} 209 } 210 if (iscatter) { 211 for (i=0; i<n; i++) {ierr = PetscObjectReference((PetscObject)iscatter[i]);CHKERRQ(ierr);} 212 } 213 if (gscatter) { 214 for (i=0; i<n; i++) {ierr = PetscObjectReference((PetscObject)gscatter[i]);CHKERRQ(ierr);} 215 } 216 if (oscatter) { 217 ierr = PetscMalloc(n*sizeof(IS),&nasm->oscatter);CHKERRQ(ierr); 218 for (i=0; i<n; i++) { 219 nasm->oscatter[i] = oscatter[i]; 220 } 221 } 222 if (iscatter) { 223 ierr = PetscMalloc(n*sizeof(IS),&nasm->iscatter);CHKERRQ(ierr); 224 for (i=0; i<n; i++) { 225 nasm->iscatter[i] = iscatter[i]; 226 } 227 } 228 if (gscatter) { 229 ierr = PetscMalloc(n*sizeof(IS),&nasm->gscatter);CHKERRQ(ierr); 230 for (i=0; i<n; i++) { 231 nasm->gscatter[i] = gscatter[i]; 232 } 233 } 234 235 if (subsnes) { 236 ierr = PetscMalloc(n*sizeof(SNES),&nasm->subsnes);CHKERRQ(ierr); 237 for (i=0; i<n; i++) { 238 nasm->subsnes[i] = subsnes[i]; 239 } 240 } 241 PetscFunctionReturn(0); 242 } 243 EXTERN_C_END 244 245 #undef __FUNCT__ 246 #define __FUNCT__ "SNESNASMSolveLocal_Private" 247 PetscErrorCode SNESNASMSolveLocal_Private(SNES snes,Vec B,Vec Y,Vec X) 248 { 249 SNES_NASM *nasm = (SNES_NASM *)snes->data; 250 SNES subsnes; 251 PetscInt i; 252 PetscErrorCode ierr; 253 Vec Xlloc,Xl,Bl,Yl; 254 VecScatter iscat,oscat,gscat; 255 DM dm,subdm; 256 257 PetscFunctionBegin; 258 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 259 ierr = VecSet(Y,0);CHKERRQ(ierr); 260 for (i=0;i<nasm->n;i++) { 261 subsnes = nasm->subsnes[i]; 262 ierr = SNESGetDM(subsnes,&subdm);CHKERRQ(ierr); 263 iscat = nasm->iscatter[i]; 264 oscat = nasm->oscatter[i]; 265 gscat = nasm->gscatter[i]; 266 ierr = DMSubDomainRestrict(dm,oscat,gscat,subdm);CHKERRQ(ierr); 267 268 /* scatter the solution to the local solution */ 269 Xl = nasm->x[i]; 270 Xlloc = nasm->xl[i]; 271 Yl = nasm->y[i]; 272 273 ierr = VecScatterBegin(oscat,X,Xl,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 274 ierr = VecScatterEnd(oscat,X,Xl,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 275 276 ierr = VecScatterBegin(gscat,X,Xlloc,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 277 ierr = VecScatterEnd(gscat,X,Xlloc,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 278 279 ierr = VecCopy(Xl,Yl);CHKERRQ(ierr); 280 281 if (B) { 282 /* scatter the RHS to the local RHS */ 283 Bl = nasm->b[i]; 284 ierr = VecScatterBegin(oscat,B,Bl,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 285 ierr = VecScatterEnd(oscat,B,Bl,INSERT_VALUES,SCATTER_FORWARD);CHKERRQ(ierr); 286 } else { 287 Bl = PETSC_NULL; 288 } 289 ierr = SNESSolve(subsnes,Bl,Yl);CHKERRQ(ierr); 290 291 ierr = VecAXPY(Yl,-1.0,Xl);CHKERRQ(ierr); 292 293 if (nasm->type == PC_ASM_BASIC) { 294 ierr = VecScatterBegin(oscat,Yl,Y,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 295 ierr = VecScatterEnd(oscat,Yl,Y,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 296 } else if (nasm->type == PC_ASM_RESTRICT) { 297 ierr = VecScatterBegin(iscat,Yl,Y,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 298 ierr = VecScatterEnd(iscat,Yl,Y,ADD_VALUES,SCATTER_REVERSE);CHKERRQ(ierr); 299 } else SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_ARG_WRONGSTATE,"Only basic and interpolate types are supported for SNESNASM"); 300 } 301 302 ierr = VecAssemblyBegin(Y);CHKERRQ(ierr); 303 ierr = VecAssemblyEnd(Y);CHKERRQ(ierr); 304 305 ierr = VecAXPY(X,1.0,Y);CHKERRQ(ierr); 306 PetscFunctionReturn(0); 307 } 308 309 #undef __FUNCT__ 310 #define __FUNCT__ "SNESSolve_NASM" 311 PetscErrorCode SNESSolve_NASM(SNES snes) 312 { 313 Vec F; 314 Vec X; 315 Vec B; 316 Vec Y; 317 PetscInt i; 318 PetscReal fnorm; 319 PetscErrorCode ierr; 320 SNESNormType normtype; 321 322 PetscFunctionBegin; 323 X = snes->vec_sol; 324 Y = snes->vec_sol_update; 325 F = snes->vec_func; 326 B = snes->vec_rhs; 327 328 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 329 snes->iter = 0; 330 snes->norm = 0.; 331 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 332 snes->reason = SNES_CONVERGED_ITERATING; 333 ierr = SNESGetNormType(snes, &normtype);CHKERRQ(ierr); 334 if (normtype == SNES_NORM_FUNCTION || normtype == SNES_NORM_INITIAL_ONLY || normtype == SNES_NORM_INITIAL_FINAL_ONLY) { 335 /* compute the initial function and preconditioned update delX */ 336 if (!snes->vec_func_init_set) { 337 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 338 if (snes->domainerror) { 339 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 340 PetscFunctionReturn(0); 341 } 342 } else { 343 snes->vec_func_init_set = PETSC_FALSE; 344 } 345 346 /* convergence test */ 347 if (!snes->norm_init_set) { 348 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); /* fnorm <- ||F|| */ 349 if (PetscIsInfOrNanReal(fnorm)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"Infinite or not-a-number generated in norm"); 350 } else { 351 fnorm = snes->norm_init; 352 snes->norm_init_set = PETSC_FALSE; 353 } 354 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 355 snes->iter = 0; 356 snes->norm = fnorm; 357 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 358 SNESLogConvHistory(snes,snes->norm,0); 359 ierr = SNESMonitor(snes,0,snes->norm);CHKERRQ(ierr); 360 361 /* set parameter for default relative tolerance convergence test */ 362 snes->ttol = fnorm*snes->rtol; 363 364 /* test convergence */ 365 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 366 if (snes->reason) PetscFunctionReturn(0); 367 } else { 368 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 369 SNESLogConvHistory(snes,snes->norm,0); 370 ierr = SNESMonitor(snes,0,snes->norm);CHKERRQ(ierr); 371 } 372 373 /* Call general purpose update function */ 374 if (snes->ops->update) { 375 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 376 } 377 378 for (i = 0; i < snes->max_its; i++) { 379 ierr = SNESNASMSolveLocal_Private(snes,B,Y,X);CHKERRQ(ierr); 380 if (normtype == SNES_NORM_FUNCTION || ((i == snes->max_its - 1) && (normtype == SNES_NORM_INITIAL_FINAL_ONLY || normtype == SNES_NORM_FINAL_ONLY))) { 381 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 382 if (snes->domainerror) { 383 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 384 PetscFunctionReturn(0); 385 } 386 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); /* fnorm <- ||F|| */ 387 if (PetscIsInfOrNanReal(fnorm)) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_FP,"Infinite or not-a-number generated in norm"); 388 } 389 /* Monitor convergence */ 390 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 391 snes->iter = i+1; 392 snes->norm = fnorm; 393 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 394 SNESLogConvHistory(snes,snes->norm,0); 395 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 396 /* Test for convergence */ 397 if (normtype == SNES_NORM_FUNCTION) ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 398 if (snes->reason) PetscFunctionReturn(0); 399 /* Call general purpose update function */ 400 if (snes->ops->update) { 401 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 402 } 403 } 404 if (normtype == SNES_NORM_FUNCTION) { 405 if (i == snes->max_its) { 406 ierr = PetscInfo1(snes,"Maximum number of iterations has been reached: %D\n",snes->max_its);CHKERRQ(ierr); 407 if (!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT; 408 } 409 } else { 410 if (!snes->reason) snes->reason = SNES_CONVERGED_ITS; /* NASM is meant to be used as a preconditioner */ 411 } 412 PetscFunctionReturn(0); 413 } 414 415 /*MC 416 SNESNASM - Nonlinear Additive Schwartz 417 418 Level: advanced 419 420 .seealso: SNESCreate(), SNES, SNESSetType(), SNESType (for list of available types) 421 M*/ 422 423 EXTERN_C_BEGIN 424 #undef __FUNCT__ 425 #define __FUNCT__ "SNESCreate_NASM" 426 PetscErrorCode SNESCreate_NASM(SNES snes) 427 { 428 SNES_NASM *nasm; 429 PetscErrorCode ierr; 430 431 PetscFunctionBegin; 432 ierr = PetscNewLog(snes, SNES_NASM, &nasm);CHKERRQ(ierr); 433 snes->data = (void*)nasm; 434 435 nasm->n = PETSC_DECIDE; 436 nasm->subsnes = 0; 437 nasm->x = 0; 438 nasm->xl = 0; 439 nasm->y = 0; 440 nasm->b = 0; 441 nasm->oscatter = 0; 442 nasm->iscatter = 0; 443 nasm->gscatter = 0; 444 445 nasm->type = PC_ASM_BASIC; 446 447 snes->ops->destroy = SNESDestroy_NASM; 448 snes->ops->setup = SNESSetUp_NASM; 449 snes->ops->setfromoptions = SNESSetFromOptions_NASM; 450 snes->ops->view = SNESView_NASM; 451 snes->ops->solve = SNESSolve_NASM; 452 snes->ops->reset = SNESReset_NASM; 453 454 snes->usesksp = PETSC_FALSE; 455 snes->usespc = PETSC_FALSE; 456 457 if (!snes->tolerancesset) { 458 snes->max_its = 10000; 459 snes->max_funcs = 10000; 460 } 461 462 ierr = PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESNASMSetSubdomains_C","SNESNASMSetSubdomains_NASM", 463 SNESNASMSetSubdomains_NASM);CHKERRQ(ierr); 464 PetscFunctionReturn(0); 465 } 466 EXTERN_C_END 467