1*a4963045SJacob Faibussowitsch #pragma once 2421d9b32SPeter Brune 3af0996ceSBarry Smith #include <petsc/private/snesimpl.h> 4af0996ceSBarry Smith #include <petsc/private/linesearchimpl.h> 5539c167fSBarry Smith #include <petscdm.h> 6421d9b32SPeter Brune 7421d9b32SPeter Brune typedef struct { 8421d9b32SPeter Brune /* flags for knowing the global place of this FAS object */ 9421d9b32SPeter Brune PetscInt level; /* level = 0 coarsest level */ 10421d9b32SPeter Brune PetscInt levels; /* if level + 1 = levels; we're the last turtle */ 11421d9b32SPeter Brune 12421d9b32SPeter Brune /* smoothing objects */ 13ab8d36c9SPeter Brune SNES smoothu; /* the SNES for presmoothing */ 14ab8d36c9SPeter Brune SNES smoothd; /* the SNES for postsmoothing */ 1522c1e704SPeter Brune 16421d9b32SPeter Brune /* coarse grid correction objects */ 176273346dSPeter Brune SNES next; /* the SNES instance for the next coarser level in the hierarchy */ 18ab8d36c9SPeter Brune SNES fine; /* the finest SNES instance; used as a reference for prefixes */ 196273346dSPeter Brune SNES previous; /* the SNES instance for the next finer level in the hierarchy */ 20421d9b32SPeter Brune Mat interpolate; /* interpolation */ 21efe1f98aSPeter Brune Mat inject; /* injection operator (unscaled) */ 22421d9b32SPeter Brune Mat restrct; /* restriction operator */ 23421d9b32SPeter Brune Vec rscale; /* the pointwise scaling of the restriction operator */ 24421d9b32SPeter Brune 25ee78dd50SPeter Brune /* method parameters */ 26ee78dd50SPeter Brune PetscInt n_cycles; /* number of cycles on this level */ 2707144faaSPeter Brune SNESFASType fastype; /* FAS type */ 28ee78dd50SPeter Brune PetscInt max_up_it; /* number of pre-smooths */ 29ee78dd50SPeter Brune PetscInt max_down_it; /* number of post-smooth cycles */ 30d1adcc6fSPeter Brune PetscBool usedmfornumberoflevels; /* uses a DM to generate a number of the levels */ 31928e959bSPeter Brune PetscBool full_downsweep; /* smooth on the initial full downsweep */ 326dfbafefSToby Isaac PetscBool full_total; /* use total residual restriction and total solution interpolation on the initial downsweep and upsweep */ 3387f44e3fSPeter Brune PetscBool continuation; /* sets the setup to default to continuation */ 348c94862eSPeter Brune PetscInt full_stage; /* stage of the full cycle -- 0 is the upswing, 1 is the downsweep and final V-cycle */ 356273346dSPeter Brune 366273346dSPeter Brune /* Galerkin FAS state */ 376273346dSPeter Brune PetscBool galerkin; /* use Galerkin formation of the coarse problem */ 386273346dSPeter Brune Vec Xg; /* Galerkin solution projection */ 396273346dSPeter Brune Vec Fg; /* Galerkin function projection */ 406273346dSPeter Brune 410dd27c6cSPeter Brune /* if logging times for each level */ 420dd27c6cSPeter Brune PetscLogEvent eventsmoothsetup; /* level setup */ 430dd27c6cSPeter Brune PetscLogEvent eventsmoothsolve; /* level smoother solves */ 440dd27c6cSPeter Brune PetscLogEvent eventresidual; /* level residual evaluation */ 450dd27c6cSPeter Brune PetscLogEvent eventinterprestrict; /* level interpolation and restriction */ 46421d9b32SPeter Brune } SNES_FAS; 47421d9b32SPeter Brune 481943db53SBarry Smith PETSC_INTERN PetscErrorCode SNESFASCycleCreateSmoother_Private(SNES, SNES *); 49