Actual source code: dacreate.c

  1: #define PETSCDM_DLL
 2:  #include private/daimpl.h

  4: #undef  __FUNCT__
  6: /*@
  7:   DAViewFromOptions - This function visualizes the DA based upon user options.

  9:   Collective on DA

 11:   Input Parameters:
 12: + da   - The DA
 13: - title - The title (currently ignored)

 15:   Level: intermediate

 17: .keywords: DA, view, options, database
 18: .seealso: DASetFromOptions(), DAView()
 19: @*/
 20: PetscErrorCode  DAViewFromOptions(DA da, const char title[])
 21: {

 25:   DAView_Private(da);
 26:   return(0);
 27: }

 31: /*
 32:   DASetTypeFromOptions_Private - Sets the type of DA from user options. Defaults to a 1D DA.

 34:   Collective on Vec

 36:   Input Parameter:
 37: . da - The DA

 39:   Level: intermediate

 41: .keywords: DA, set, options, database, type
 42: .seealso: DASetFromOptions(), DASetType()
 43: */
 44: static PetscErrorCode DASetTypeFromOptions_Private(DA da)
 45: {
 46:   const DAType   defaultType = DA1D;
 47:   char           typeName[256];
 48:   PetscTruth     opt;

 52:   switch (da->dim) {
 53:     case 1: defaultType = DA1D; break;
 54:     case 2: defaultType = DA2D; break;
 55:     case 3: defaultType = DA3D; break;
 56:   }
 57:   if (((PetscObject)da)->type_name) {
 58:     defaultType = ((PetscObject)da)->type_name;
 59:   }
 60:   if (!DARegisterAllCalled) {DARegisterAll(PETSC_NULL);}
 61:   PetscOptionsList("-da_type","DA type","DASetType",DAList,defaultType,typeName,256,&opt);
 62:   if (opt) {
 63:     DASetType(da, typeName);
 64:   } else {
 65:     DASetType(da, defaultType);
 66:   }
 67:   return(0);
 68: }

 72: /*@
 73:   DASetFromOptions - Configures the vector from the options database.

 75:   Collective on DA

 77:   Input Parameter:
 78: . da - The DA

 80:   Notes:  To see all options, run your program with the -help option, or consult the users manual.
 81:           Must be called after DACreate() but before the DA is used.

 83:   Level: beginner

 85:   Concepts: DA^setting options
 86:   Concepts: DA^setting type

 88: .keywords: DA, set, options, database
 89: .seealso: DACreate(), DASetOptionsPrefix()
 90: @*/
 91: PetscErrorCode  DASetFromOptions(DA da)
 92: {
 94:   PetscInt       dim;
 95:   PetscTruth     flg;


100:   PetscOptionsBegin(((PetscObject)da)->comm,((PetscObject)da)->prefix,"DA Options","DA");
101:     /* Handle DA dimensions */
102:     PetscOptionsInt("-da_dim","Number of dimensions","DASetDim",da->dim,&dim,&flg);
103:     if (flg) {
104:       DASetDim(da,dim);
105:     }
106:     /* Handle DA grid sizes */
107:     if (da->M < 0) {
108:       PetscInt newM = -da->M;
109:       PetscOptionsInt("-da_grid_x","Number of grid points in x direction","DASetSizes",newM,&newM,PETSC_NULL);
110:       da->M = newM;
111:     }
112:     if (da->N < 0) {
113:       PetscInt newN = -da->N;
114:       PetscOptionsInt("-da_grid_y","Number of grid points in y direction","DASetSizes",newN,&newN,PETSC_NULL);
115:       da->N = newN;
116:     }
117:     if (da->P < 0) {
118:       PetscInt newP = -da->P;
119:       PetscOptionsInt("-da_grid_z","Number of grid points in z direction","DASetSizes",newP,&newP,PETSC_NULL);
120:       da->P = newP;
121:     }
122:     /* Handle DA parallel distibution */
123:     PetscOptionsInt("-da_processors_x","Number of processors in x direction","DASetNumProcs",da->m,&da->m,PETSC_NULL);
124:     PetscOptionsInt("-da_processors_y","Number of processors in y direction","DASetNumProcs",da->n,&da->n,PETSC_NULL);
125:     PetscOptionsInt("-da_processors_z","Number of processors in z direction","DASetNumProcs",da->p,&da->p,PETSC_NULL);
126:     /* Handle DA refinement */
127:     PetscOptionsInt("-da_refine_x","Refinement ratio in x direction","DASetRefinementFactor",da->refine_x,&da->refine_x,PETSC_NULL);
128:     PetscOptionsInt("-da_refine_y","Refinement ratio in y direction","DASetRefinementFactor",da->refine_y,&da->refine_y,PETSC_NULL);
129:     PetscOptionsInt("-da_refine_z","Refinement ratio in z direction","DASetRefinementFactor",da->refine_z,&da->refine_z,PETSC_NULL);
130:     /* Handle DA type options */
131:     DASetTypeFromOptions_Private(da);
132:     /* Handle specific DA options */
133:     if (da->ops->setfromoptions) {
134:       (*da->ops->setfromoptions)(da);
135:     }
136:   PetscOptionsEnd();

138:   DAViewFromOptions(da, ((PetscObject)da)->name);
139:   return(0);
140: }

144: /*@
145:   DACreate - Creates an empty DA object. The type can then be set with DASetType(),
146:   or DASetFromOptions().

148:    If you never  call DASetType() or DASetFromOptions() it will generate an 
149:    error when you try to use the DA.

151:   Collective on MPI_Comm

153:   Input Parameter:
154: . comm - The communicator for the DA object

156:   Output Parameter:
157: . da  - The DA object

159:   Level: beginner

161: .keywords: DA, create
162: .seealso: DASetType(), DASetSizes(), DADuplicate()
163: @*/
164: PetscErrorCode  DACreate(MPI_Comm comm, DA *da)
165: {
166:   DA             d;

171:   *da = PETSC_NULL;
172: #ifndef PETSC_USE_DYNAMIC_LIBRARIES
173:   DMInitializePackage(PETSC_NULL);
174: #endif

176:   PetscHeaderCreate(d, _p_DA, struct _DAOps, DM_COOKIE, 0, "DM", comm, DADestroy, DAView);
177:   PetscMemzero(d->ops, sizeof(struct _DAOps));

179:   d->dim        = -1;
180:   d->interptype = DA_Q1;
181:   d->refine_x   = 2;
182:   d->refine_y   = 2;
183:   d->refine_z   = 2;
184:   d->fieldname  = PETSC_NULL;
185:   d->nlocal     = -1;
186:   d->Nlocal     = -1;
187:   d->M          = -1;
188:   d->N          = -1;
189:   d->P          = -1;
190:   d->m          = -1;
191:   d->n          = -1;
192:   d->p          = -1;
193:   d->w          = -1;
194:   d->s          = -1;
195:   d->xs = -1; d->xe = -1; d->ys = -1; d->ye = -1; d->zs = -1; d->ze = -1;
196:   d->Xs = -1; d->Xe = -1; d->Ys = -1; d->Ye = -1; d->Zs = -1; d->Ze = -1;

198:   d->gtol         = PETSC_NULL;
199:   d->ltog         = PETSC_NULL;
200:   d->ltol         = PETSC_NULL;
201:   d->ltogmap      = PETSC_NULL;
202:   d->ltogmapb     = PETSC_NULL;
203:   d->ao           = PETSC_NULL;
204:   d->base         = -1;
205:   d->wrap         = DA_NONPERIODIC;
206:   d->stencil_type = DA_STENCIL_BOX;
207:   d->interptype   = DA_Q1;
208:   d->idx          = PETSC_NULL;
209:   d->Nl           = -1;
210:   d->lx           = PETSC_NULL;
211:   d->ly           = PETSC_NULL;
212:   d->lz           = PETSC_NULL;

214:   d->ops->globaltolocalbegin = DAGlobalToLocalBegin;
215:   d->ops->globaltolocalend   = DAGlobalToLocalEnd;
216:   d->ops->localtoglobal      = DALocalToGlobal;
217:   d->ops->createglobalvector = DACreateGlobalVector;
218:   d->ops->createlocalvector  = DACreateLocalVector;
219:   d->ops->getinterpolation   = DAGetInterpolation;
220:   d->ops->getcoloring        = DAGetColoring;
221:   d->ops->getmatrix          = DAGetMatrix;
222:   d->ops->refine             = DARefine;
223:   d->ops->coarsen            = DACoarsen;
224:   d->ops->refinehierarchy    = DARefineHierarchy;
225:   d->ops->coarsenhierarchy   = DACoarsenHierarchy;
226:   d->ops->getinjection       = DAGetInjection;
227:   d->ops->getaggregates      = DAGetAggregates;

229:   *da = d;
230:   return(0);
231: }