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