Actual source code: dadist.c
1: /*$Id: dadist.c,v 1.29 2001/03/23 23:25:00 balay Exp $*/
2:
3: /*
4: Code for manipulating distributed regular arrays in parallel.
5: */
7: #include src/dm/da/daimpl.h
9: /*@C
10: DACreateGlobalVector - Creates a parallel PETSc vector that
11: may be used with the DAXXX routines.
13: Collective on DA
15: Input Parameter:
16: . da - the distributed array
18: Output Parameter:
19: . g - the distributed global vector
21: Level: beginner
23: Note:
24: The output parameter, g, is a regular PETSc vector that should be destroyed
25: with a call to VecDestroy() when usage is finished.
27: .keywords: distributed array, create, global, distributed, vector
29: .seealso: DACreateLocalVector(), VecDuplicate(), VecDuplicateVecs(),
30: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
31: DAGlobalToLocalEnd(), DALocalToGlobal(), DACreateNaturalVector()
32: @*/
33: int DACreateGlobalVector(DA da,Vec* g)
34: {
39: VecDuplicate(da->global,g);
40: PetscObjectCompose((PetscObject)*g,"DA",(PetscObject)da);
41: return(0);
42: }
44: /*@C
45: DACreateNaturalVector - Creates a parallel PETSc vector that
46: will hold vector values in the natural numbering, rather than in
47: the PETSc parallel numbering associated with the DA.
49: Collective on DA
51: Input Parameter:
52: . da - the distributed array
54: Output Parameter:
55: . g - the distributed global vector
57: Level: developer
59: Note:
60: The output parameter, g, is a regular PETSc vector that should be destroyed
61: with a call to VecDestroy() when usage is finished.
63: .keywords: distributed array, create, global, distributed, vector
65: .seealso: DACreateLocalVector(), VecDuplicate(), VecDuplicateVecs(),
66: DACreate1d(), DACreate2d(), DACreate3d(), DAGlobalToLocalBegin(),
67: DAGlobalToLocalEnd(), DALocalToGlobal()
68: @*/
69: int DACreateNaturalVector(DA da,Vec* g)
70: {
71: int cnt,m,ierr;
75: if (da->natural) {
76: PetscObjectGetReference((PetscObject)da->natural,&cnt);
77: if (cnt == 1) { /* object is not currently used by anyone */
78: PetscObjectReference((PetscObject)da->natural);
79: *g = da->natural;
80: } else {
81: VecDuplicate(da->natural,g);
82: }
83: } else { /* create the first version of this guy */
84: VecGetLocalSize(da->global,&m);
85: VecCreateMPI(da->comm,m,PETSC_DETERMINE,g);
86: PetscObjectReference((PetscObject)*g);
87: da->natural = *g;
88: }
89: return(0);
90: }