Actual source code: tmap.c

petsc-3.3-p5 2012-12-01
  1: #include <petsc-private/vecimpl.h>  /*I   "petscvec.h"  I*/

  5: /*
  6:      PetscThreadsLayoutCreate - Allocates PetsThreadscLayout space and sets the map contents to the default.


  9:    Input Parameters:
 10: .    map - pointer to the map

 12:    Level: developer

 14: .seealso: PetscThreadsLayoutSetLocalSizes(), PetscThreadsLayoutGetLocalSizes(), PetscThreadsLayout, 
 15:           PetscThreadsLayoutDestroy(), PetscThreadsLayoutSetUp()
 16: */
 17: PetscErrorCode PetscThreadsLayoutCreate(PetscThreadsLayout *tmap)
 18: {

 22:   PetscNew(struct _n_PetscThreadsLayout,tmap);
 23:   (*tmap)->nthreads = -1;
 24:   (*tmap)->N        = -1;
 25:   (*tmap)->trstarts =  0;
 26:   (*tmap)->affinity =  0;
 27:   return(0);
 28: }

 32: /*
 33:      PetscThreadsLayoutDestroy - Frees a map object and frees its range if that exists.

 35:    Input Parameters:
 36: .    map - the PetscThreadsLayout

 38:    Level: developer

 40:       The PetscThreadsLayout object and methods are intended to be used in the PETSc threaded Vec and Mat implementions; it is 
 41:       recommended they not be used in user codes unless you really gain something in their use.

 43: .seealso: PetscThreadsLayoutSetLocalSizes(), PetscThreadsLayoutGetLocalSizes(), PetscThreadsLayout, 
 44:           PetscThreadsLayoutSetUp()
 45: */
 46: PetscErrorCode PetscThreadsLayoutDestroy(PetscThreadsLayout *tmap)
 47: {
 49: 
 51:   if(!*tmap) return(0);
 52:   PetscFree((*tmap)->trstarts);
 53:   PetscFree((*tmap)->affinity);
 54:   PetscFree((*tmap));
 55:   *tmap = PETSC_NULL;
 56:   return(0);
 57: }

 61: /*
 62:      PetscThreadsLayoutSetUp - given a map where you have set the thread count, either global size or
 63:            local sizes sets up the map so that it may be used.

 65:    Input Parameters:
 66: .    map - pointer to the map

 68:    Level: developer

 70:    Notes: Typical calling sequence
 71:       PetscThreadsLayoutCreate(PetscThreadsLayout *);
 72:       PetscThreadsLayoutSetNThreads(PetscThreadsLayout,nthreads);
 73:       PetscThreadsLayoutSetSize(PetscThreadsLayout,N) or PetscThreadsLayoutSetLocalSizes(PetscThreadsLayout, *n); or both
 74:       PetscThreadsLayoutSetUp(PetscThreadsLayout);

 76:        If the local sizes, global size are already set and row offset exists then this does nothing.

 78: .seealso: PetscThreadsLayoutSetLocalSizes(), PetscThreadsLayoutGetLocalSizes(), PetscThreadsLayout, 
 79:           PetscThreadsLayoutDestroy()
 80: */
 81: PetscErrorCode PetscThreadsLayoutSetUp(PetscThreadsLayout tmap)
 82: {
 83:   PetscErrorCode     ierr;
 84:   PetscInt           t,rstart=0,n,Q,R;
 85:   PetscBool          S;
 86: 
 88:   if(!tmap->nthreads) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Number of threads not set yet");
 89:   if((tmap->N >= 0) && (tmap->trstarts)) return(0);
 90:   PetscMalloc((tmap->nthreads+1)*sizeof(PetscInt),&tmap->trstarts);

 92:   Q = tmap->N/tmap->nthreads;
 93:   R = tmap->N - Q*tmap->nthreads;
 94:   for(t=0;t < tmap->nthreads;t++) {
 95:     tmap->trstarts[t] = rstart;
 96:     S               = (PetscBool)(t<R);
 97:     n               = S?Q+1:Q;
 98:     rstart         += n;
 99:   }
100:   tmap->trstarts[tmap->nthreads] = rstart;
101:   return(0);
102: }
103: 
106: /*

108:     PetscThreadsLayoutDuplicate - creates a new PetscThreadsLayout with the same information as a given one. If the PetscThreadsLayout already exists it is destroyed first.

110:      Collective on PetscThreadsLayout

112:     Input Parameter:
113: .     in - input PetscThreadsLayout to be copied

115:     Output Parameter:
116: .     out - the copy

118:    Level: developer

120:     Notes: PetscThreadsLayoutSetUp() does not need to be called on the resulting PetscThreadsLayout

122: .seealso: PetscThreadsLayoutCreate(), PetscThreadsLayoutDestroy(), PetscThreadsLayoutSetUp()
123: */
124: PetscErrorCode PetscThreadsLayoutDuplicate(PetscThreadsLayout in,PetscThreadsLayout *out)
125: {

129:   PetscThreadsLayoutDestroy(out);
130:   PetscThreadsLayoutCreate(out);
131:   PetscMemcpy(*out,in,sizeof(struct _n_PetscThreadsLayout));

133:   PetscMalloc(in->nthreads*sizeof(PetscInt),&(*out)->trstarts);
134:   PetscMemcpy((*out)->trstarts,in->trstarts,in->nthreads*sizeof(PetscInt));
135: 
136:   PetscMalloc(in->nthreads*sizeof(PetscInt),&(*out)->affinity);
137:   PetscMemcpy((*out)->affinity,in->affinity,in->nthreads*sizeof(PetscInt));

139:   return(0);
140: }

144: /*
145:      PetscThreadsLayoutSetLocalSizes - Sets the local size for each thread 

147:    Input Parameters:
148: +    map - pointer to the map
149: -    n - local sizes

151:    Level: developer

153: .seealso: PetscThreadsLayoutCreate(), PetscThreadsLayoutDestroy(), PetscThreadsLayoutGetLocalSizes()

155: */
156: PetscErrorCode PetscThreadsLayoutSetLocalSizes(PetscThreadsLayout tmap,PetscInt n[])
157: {
159:   PetscInt       i;

162:   if(!tmap->nthreads) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Number of threads not set yet");
163:   if (tmap->trstarts) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_WRONGSTATE,"Already set local sizes");
164:   PetscMalloc((tmap->nthreads+1)*sizeof(PetscInt),&tmap->trstarts);
165:   tmap->trstarts[0] = 0;
166:   for(i=1;i < tmap->nthreads+1;i++) tmap->trstarts[i] += n[i-1];
167:   return(0);
168: }

172: /*
173:      PetscThreadsLayoutGetLocalSizes - Gets the local size for each thread 

175:    Input Parameters:
176: .    map - pointer to the map

178:    Output Parameters:
179: .    n - array to hold the local sizes (must be allocated)

181:    Level: developer

183: .seealso: PetscThreadsLayoutCreate(), PetscThreadsLayoutDestroy(), PetscThreadsLayoutSetLocalSizes()
184: */
185: PetscErrorCode PetscThreadsLayoutGetLocalSizes(PetscThreadsLayout tmap,PetscInt *n[])
186: {
187:   PetscInt i;
188:   PetscInt *tn=*n;
190:   for(i=0;i < tmap->nthreads;i++) tn[i] = tmap->trstarts[i+1] - tmap->trstarts[i];
191:   return(0);
192: }
193: 
196: /*
197:      PetscThreadsLayoutSetSize - Sets the global size for PetscThreadsLayout object

199:    Input Parameters:
200: +    map - pointer to the map
201: -    n -   global size

203:    Level: developer

205: .seealso: PetscThreadsLayoutCreate(), PetscThreadsLayoutSetLocalSizes(), PetscThreadsLayoutGetSize()
206: */
207: PetscErrorCode PetscThreadsLayoutSetSize(PetscThreadsLayout tmap,PetscInt N)
208: {
210:   tmap->N = N;
211:   return(0);
212: }

216: /*
217:      PetscThreadsLayoutGetSize - Gets the global size 

219:    Input Parameters:
220: .    map - pointer to the map

222:    Output Parameters:
223: .    n - global size

225:    Level: developer

227: .seealso: PetscThreadsLayoutCreate(), PetscThreadsLayoutSetSize(), PetscThreadsLayoutGetLocalSizes()
228: */
229: PetscErrorCode PetscThreadsLayoutGetSize(PetscThreadsLayout tmap,PetscInt *N)
230: {
232:   *N = tmap->N;
233:   return(0);
234: }

238: /*
239:      PetscThreadsLayoutSetNThreads - Sets the thread count for PetscThreadsLayout object

241:    Input Parameters:
242: +    map - pointer to the map
243: -    nthreads -   number of threads to be used with the map

245:    Level: developer

247: .seealso: PetscThreadsLayoutCreate(), PetscThreadsLayoutSetLocalSizes(), PetscThreadsLayoutGetSize()
248: */
249: PetscErrorCode PetscThreadsLayoutSetNThreads(PetscThreadsLayout tmap,PetscInt nthreads)
250: {
252:   tmap->nthreads = nthreads;
253:   return(0);
254: }

258: /*
259:      PetscThreadsLayoutSetLocalSizes - Sets the core affinities for PetscThreadsLayout object

261:    Input Parameters:
262: +    map - pointer to the map
263: -    affinities - core affinities for PetscThreadsLayout 

265:    Level: developer

267: .seealso: PetscThreadsLayoutGetThreadAffinities()

269: */
270: PetscErrorCode PetscThreadsLayoutSetThreadAffinities(PetscThreadsLayout tmap, PetscInt affinities[])
271: {

275:   PetscMalloc(tmap->nthreads*sizeof(PetscInt),&tmap->affinity);
276:   PetscMemcpy(tmap->affinity,affinities,tmap->nthreads*sizeof(PetscInt));
277:   return(0);
278: }

282: /*
283:      PetscThreadsLayoutGetThreadAffinities - Gets the core affinities of threads

285:    Input Parameters:
286: .    map - pointer to the map

288:    Output Parameters:
289: .    affinity - core affinities of threads

291:    Level: developer

293: .seealso: PetscThreadsLayoutSetThreadAffinities()
294: */
295: PetscErrorCode PetscThreadsLayoutGetThreadAffinities(PetscThreadsLayout tmap,const PetscInt *affinity[])
296: {
298:   *affinity = tmap->affinity;
299:   return(0);
300: }