This commit is contained in:
Michael Stumpf 2020-04-03 18:32:23 +02:00
parent 2b1d859a7e
commit d485e40ec4
1 changed files with 41 additions and 21 deletions

View File

@ -25,8 +25,9 @@ class ibmppp:
# Declare variables which hold fields
self.field = {'u': None, 'v': None, 'w': None, 'p': None, 's1': None}
# Declare variables which hold statistics
self.spatialMean = {'u': None, 'v': None, 'w': None, 'p': None, 's1': None}
self.spatialMsqr = {'u': None, 'v': None, 'w': None, 'p': None, 's1': None}
self.spatialMean = {'u': None, 'v': None, 'w': None, 'p': None, 's1': None}
self.spatialMsqr = {'u': None, 'v': None, 'w': None, 'p': None, 's1': None}
self.spatialNsample = {'u': None, 'v': None, 'w': None, 'p': None, 's1': None}
# Declare variables for grid (to be set by "loadGrid()")
self.grid = {'u': [None,None,None], 'v':[None,None,None], 'w': [None,None,None], 'p': [None,None,None], 's1': [None,None,None]}
self.domainBounds = (None,None,None,None,None,None)
@ -598,10 +599,13 @@ class ibmppp:
gid = fid.create_group('/'+key)
did = gid.create_dataset('mean',data=self.spatialMean[key].ravel('C'))
did.attrs['F_CONTIGUOUS'] = 0
did.attrs['DIM'] = spatialMean[key].shape
did.attrs['DIM'] = self.spatialMean[key].shape
did = gid.create_dataset('meansquare',data=self.spatialMsqr[key].ravel('C'))
did.attrs['F_CONTIGUOUS'] = 0
did.attrs['DIM'] = spatialMsqr[key].shape
did.attrs['DIM'] = self.spatialMsqr[key].shape
did = gid.create_dataset('nsample',data=self.spatialNsample[key].ravel('C'))
did.attrs['F_CONTIGUOUS'] = 0
did.attrs['DIM'] = self.spatialNsample[key].shape
# Close the file
fid.close()
@ -930,8 +934,11 @@ class ibmppp:
if collocate:
self.shiftField(keyout,self.__getShiftDirection(keyout,'p'))
def spatialStatistics(self,key,homogeneous=(self.__xperiodic,self.__yperiodic,self.__zperiodic)):
def spatialStatistics(self,key,homogeneous=None):
'''Computes spatial mean and meansquare over homogeneous directions.'''
# If no homogeneous directions are provided, use periodicity
if homogeneous is None:
homogeneous = (self.__xperiodic,self.__yperiodic,self.__zperiodic)
# Check if there is any homogeneous direction
if not np.any(homogeneous):
raise ValueError('Cannot compute statistics, since there is no homogeneous direction.')
@ -949,32 +956,43 @@ class ibmppp:
tmpsum = np.nansum(tmpsum,axis=idir,keepdims=True)
tmpssq = np.nansum(tmpssq,axis=idir,keepdims=True)
# Reduce the statistics to processors with col/row/pln=0 in homogeneous directions
proc_layout = np.array(range(0,self.__nproc)).reshape(self.__nxp,self.__nyp,self.__nzp)
for idir in range(0,3):
if homogeneous[idir]:
pos_dst = [self.__icol,self.__jrow,self.__kpln]
pos_dst[idir] = 0
rank_dst = self.__rankFromPosition(pos_dst[0],pos_dst[1],pos_dst[2],self.__nxp,self.__nyp,self.__nzp)
# Create groups for collective communication: get rank of all procs in current direction
proc_slice = [slice(self.__icol,self.__icol+1),slice(self.__jrow,self.__jrow+1),slice(self.__kpln,self.__kpln+1)]
proc_slice[idir] = slice(0,None)
rank_list = proc_layout[proc_slice[0],proc_slice[1],proc_slice[2]].ravel()
comm = self.__comm.Create(MPI.Group.Incl(self.__comm.Get_group(),rank_list))
rank = comm.Get_rank()
# Communicate counter
if self.__rank==rank_dst:
#print('Receiving rank',self.__rank,'is known in this comminucator as',rank,'for idir =',idir)
recvbuf = np.zeros_like(counter)
else:
recvbuf = counter
self.__comm.Reduce(counter,recvbuf,op=MPI.SUM,root=rank_dst)
#print(self.__rank,'Reduce 1: destination',rank_dst,pos_dst)
comm.Reduce(counter,recvbuf,op=MPI.SUM,root=0)
counter = recvbuf
# Communicate tmpsum
if self.__rank==rank_dst:
recvbuf = np.zeros_like(tmpsum)
else:
recvbuf = tmpsum
self.__comm.Reduce(tmpsum,recvbuf,op=MPI.SUM,root=rank_dst)
#print(self.__rank,'Reduce 2')
comm.Reduce(tmpsum,recvbuf,op=MPI.SUM,root=0)
tmpsum = recvbuf
# Communicate tmpssq
if self.__rank==rank_dst:
recvbuf = np.zeros_like(tmpssq)
else:
recvbuf = tmpssq
self.__comm.Reduce(tmpssq,recvbuf,op=MPI.SUM,root=rank_dst)
comm.Reduce(tmpssq,recvbuf,op=MPI.SUM,root=0)
tmpssq = recvbuf
comm.Free()
# All processors with col/row/pln!=0 in homogeneous direction are done.
pos = [self.__icol,self.__jrow,self.__kpln]
for idir in range(0,3):
@ -984,18 +1002,17 @@ class ibmppp:
# the full line or plane. Rank 0 allocates the final array
if self.__rank==0:
# Determine size of final array
nxg = self.__procGrid[key][1][-1]-self.__procGrid[key][0][0]
nyg = self.__procGrid[key][3][-1]-self.__procGrid[key][2][0]
nzg = self.__procGrid[key][5][-1]-self.__procGrid[key][4][0]
nxg = self.__procGrid[key][1][-1]-self.__procGrid[key][0][0]+1
nyg = self.__procGrid[key][3][-1]-self.__procGrid[key][2][0]+1
nzg = self.__procGrid[key][5][-1]-self.__procGrid[key][4][0]+1
dim = [nxg,nyg,nzg]
for idir in range(0,3):
if homogeneous[idir]:
dim[idir] = 1
# Allocate a nice spot of memory
self.spatialMean[key] = np.zeros(dim,dtype=self.__precision)
self.spatialMsqr[key] = np.zeros(dim,dtype=self.__precision)
# Also allocate one helper array
counter_full = np.zeros(dim,dtype=counter.dtype)
self.spatialMean[key] = np.zeros(dim,dtype=self.__precision)
self.spatialMsqr[key] = np.zeros(dim,dtype=self.__precision)
self.spatialNsample[key] = np.zeros(dim,dtype='i')
# Create iterators for receive request
ipend = 0 if homogeneous[0] else self.__nxp-1
jpend = 0 if homogeneous[1] else self.__nyp-1
@ -1017,21 +1034,24 @@ class ibmppp:
if rank_src==0:
recvbuf = counter
else:
recvbuf = self.__comm.Recv(source=rank_src,tag=1)
counter_full[ib:ie+1,jb:je+1,kb:ke+1] = recvbuf
recvbuf = np.empty_like(counter)
self.__comm.Recv(recvbuf,source=rank_src,tag=1)
self.spatialNsample[key][ib:ie+1,jb:je+1,kb:ke+1] = recvbuf
if rank_src==0:
recvbuf = tmpsum
else:
recvbuf = self.__comm.Recv(source=rank_src,tag=2)
recvbuf = np.empty_like(tmpsum)
self.__comm.Recv(recvbuf,source=rank_src,tag=2)
self.spatialMean[key][ib:ie+1,jb:je+1,kb:ke+1] = recvbuf
if rank_src==0:
recvbuf = tmpssq
else:
recvbuf = self.__comm.Recv(source=rank_src,tag=3)
recvbuf = np.empty_like(tmpssq)
self.__comm.Recv(recvbuf,source=rank_src,tag=3)
self.spatialMsqr[key][ib:ie+1,jb:je+1,kb:ke+1] = recvbuf
# Rank 0 got all the data now. Finalize it!
self.spatialMean[key] /= counter_full
self.spatialMsqr[key] /= counter_full
self.spatialMean[key] /= self.spatialNsample[key]
self.spatialMsqr[key] /= self.spatialNsample[key]
else:
# All tanks but rank 0 send their data
self.__comm.Send(counter,dest=0,tag=1)