NAMD
Public Member Functions | Public Attributes | List of all members
PmeXPencil Class Reference
Inheritance diagram for PmeXPencil:
PmePencil< CBase_PmeXPencil >

Public Member Functions

PmeXPencil_SDAG_CODE PmeXPencil ()
 
 PmeXPencil (CkMigrateMessage *)
 
 ~PmeXPencil ()
 
void fft_init ()
 
void recv_trans (const PmeTransMsg *)
 
void forward_fft ()
 
void pme_kspace ()
 
void backward_fft ()
 
void send_untrans ()
 
void send_subset_untrans (int fromIdx, int toIdx)
 
void node_process_trans (PmeTransMsg *)
 
void evir_init ()
 
- Public Member Functions inherited from PmePencil< CBase_PmeXPencil >
 PmePencil ()
 
 ~PmePencil ()
 
void base_init (PmePencilInitMsg *msg)
 
void order_init (int nBlocks)
 

Public Attributes

fftw_plan forward_plan
 
fftw_plan backward_plan
 
int ny
 
int nz
 
int recipEvirPe
 
PmeKSpacemyKSpace
 
- Public Attributes inherited from PmePencil< CBase_PmeXPencil >
PmePencilInitMsgData initdata
 
Lattice lattice
 
PmeReduction evir
 
int sequence
 
AtomicInt imsg
 
AtomicInt imsgb
 
int hasData
 
int offload
 
float * data
 
float * work
 
int * send_order
 
int * needs_reply
 

Additional Inherited Members

- Public Types inherited from PmePencil< CBase_PmeXPencil >
typedef int AtomicInt
 

Detailed Description

Definition at line 4699 of file ComputePme.C.

Constructor & Destructor Documentation

PmeXPencil_SDAG_CODE PmeXPencil::PmeXPencil ( )
inline

Definition at line 4702 of file ComputePme.C.

References PmePencil< CBase_PmeXPencil >::imsg, PmePencil< CBase_PmeXPencil >::imsgb, myKSpace, and recipEvirPe.

4702 { __sdag_init(); myKSpace = 0; setMigratable(false); imsg=imsgb=0; recipEvirPe = -999; }
PmeKSpace * myKSpace
Definition: ComputePme.C:4733
int recipEvirPe
Definition: ComputePme.C:4731
PmeXPencil::PmeXPencil ( CkMigrateMessage *  )
inline

Definition at line 4703 of file ComputePme.C.

4703 { __sdag_init(); }
PmeXPencil::~PmeXPencil ( )
inline

Definition at line 4704 of file ComputePme.C.

4704  {
4705  #ifdef NAMD_FFTW
4706  #ifdef NAMD_FFTW_3
4707  delete [] forward_plans;
4708  delete [] backward_plans;
4709  #endif
4710  #endif
4711  }

Member Function Documentation

void PmeXPencil::backward_fft ( )

Definition at line 5712 of file ComputePme.C.

References backward_plan, CKLOOP_CTRL_PME_BACKWARDFFT, PmePencil< CBase_PmeXPencil >::data, PmePencilInitMsgData::grid, PmePencil< CBase_PmeXPencil >::initdata, PmeGrid::K1, Node::Object(), PmeXZPencilFFT(), Node::simParameters, SimParameters::useCkLoop, PmePencil< CBase_PmeXPencil >::work, PmePencilInitMsgData::yBlocks, and PmePencilInitMsgData::zBlocks.

Referenced by node_process_trans().

5712  {
5713 #ifdef NAMD_FFTW
5714 #ifdef MANUAL_DEBUG_FFTW3
5715  dumpMatrixFloat3("bw_x_b", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
5716 #endif
5717 
5718 #ifdef NAMD_FFTW_3
5719 #if CMK_SMP && USE_CKLOOP
5720  int useCkLoop = Node::Object()->simParameters->useCkLoop;
5721  if(useCkLoop>=CKLOOP_CTRL_PME_BACKWARDFFT
5722  && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks) {
5723  //for(int i=0; i<numPlans; i++) fftwf_execute(backward_plans[i]);
5724  //transform the above loop
5725  CkLoop_Parallelize(PmeXZPencilFFT, 1, (void *)backward_plans, CkMyNodeSize(), 0, numPlans-1); //sync
5726  return;
5727  }
5728 #endif
5729  fftwf_execute(backward_plan);
5730 #else
5731  fftw(backward_plan, ny*nz,
5732  ((fftw_complex *) data), ny*nz, 1, (fftw_complex *) work, 1, 0);
5733 #endif
5734 #ifdef MANUAL_DEBUG_FFTW3
5735  dumpMatrixFloat3("bw_x_a", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
5736 #endif
5737 #endif
5738 }
static Node * Object()
Definition: Node.h:86
PmePencilInitMsgData initdata
Definition: ComputePme.C:4535
SimParameters * simParameters
Definition: Node.h:178
int K1
Definition: PmeBase.h:18
static void PmeXZPencilFFT(int first, int last, void *result, int paraNum, void *param)
Definition: ComputePme.C:5177
fftw_plan backward_plan
Definition: ComputePme.C:4727
#define CKLOOP_CTRL_PME_BACKWARDFFT
Definition: SimParameters.h:97
void PmeXPencil::evir_init ( )

Definition at line 4757 of file ComputePme.C.

References findRecipEvirPe(), PmePencil< CBase_PmeXPencil >::initdata, PmePencilInitMsgData::pmeProxy, and recipEvirPe.

4757  {
4759  initdata.pmeProxy[recipEvirPe].addRecipEvirClient();
4760 }
PmePencilInitMsgData initdata
Definition: ComputePme.C:4535
CProxy_ComputePmeMgr pmeProxy
Definition: ComputePme.C:217
static int findRecipEvirPe()
Definition: ComputePme.C:241
int recipEvirPe
Definition: ComputePme.C:4731
void PmeXPencil::fft_init ( )

Definition at line 5035 of file ComputePme.C.

References backward_plan, PmeGrid::block2, PmeGrid::block3, PmePencil< CBase_PmeXPencil >::data, PmeGrid::dim3, ComputePmeMgr::fftw_plan_lock, SimParameters::FFTWEstimate, fftwf_malloc, SimParameters::FFTWPatient, forward_plan, PmePencilInitMsgData::grid, if(), PmePencil< CBase_PmeXPencil >::initdata, PmeGrid::K1, PmeGrid::K2, myKSpace, NAMD_die(), nz, PmePencil< CBase_PmeXPencil >::order_init(), PmePencilInitMsgData::pmeNodeProxy, Node::simParameters, SimParameters::useCkLoop, PmePencil< CBase_PmeXPencil >::work, and PmePencilInitMsgData::xBlocks.

5035  {
5036  CProxy_Node nd(CkpvAccess(BOCclass_group).node);
5037  Node *node = nd.ckLocalBranch();
5039 #if USE_NODE_PAR_RECEIVE
5040  ((NodePmeMgr *)CkLocalNodeBranch(initdata.pmeNodeProxy))->registerXPencil(thisIndex,this);
5041 #endif
5042 
5043  int K1 = initdata.grid.K1;
5044  int K2 = initdata.grid.K2;
5045  int dim3 = initdata.grid.dim3;
5046  int block2 = initdata.grid.block2;
5047  int block3 = initdata.grid.block3;
5048 
5049  ny = block2;
5050  if ( (thisIndex.y + 1) * block2 > K2 ) ny = K2 - thisIndex.y * block2;
5051  nz = block3;
5052  if ( (thisIndex.z+1)*block3 > dim3/2 ) nz = dim3/2 - thisIndex.z*block3;
5053 
5054 #ifdef NAMD_FFTW
5056 
5057  data = (float *) fftwf_malloc( sizeof(float) * K1*ny*nz*2);
5058  work = new float[2*K1];
5059 
5061 
5062 #ifdef NAMD_FFTW_3
5063  /* need array of sizes for the how many */
5064  int fftwFlags = simParams->FFTWPatient ? FFTW_PATIENT : simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE ;
5065  int sizeLines=ny*nz;
5066  int planLineSizes[1];
5067  planLineSizes[0]=K1;
5068  forward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines,
5069  (fftwf_complex *) data, NULL, sizeLines, 1,
5070  (fftwf_complex *) data, NULL, sizeLines, 1,
5071  FFTW_FORWARD,
5072  fftwFlags);
5073  backward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines,
5074  (fftwf_complex *) data, NULL, sizeLines, 1,
5075  (fftwf_complex *) data, NULL, sizeLines, 1,
5076  FFTW_BACKWARD,
5077  fftwFlags);
5078 
5079 #if CMK_SMP && USE_CKLOOP
5080  if(simParams->useCkLoop) {
5081  //How many FFT plans to be created? The grain-size issue!!.
5082  //Currently, I am choosing the min(nx, ny) to be coarse-grain
5083  numPlans = (ny<=nz?ny:nz);
5084  // limit attempted parallelism due to false sharing
5085  //if ( numPlans < CkMyNodeSize() ) numPlans = (ny>=nz?ny:nz);
5086  //if ( numPlans < CkMyNodeSize() ) numPlans = sizeLines;
5087  if ( sizeLines/numPlans < 4 ) numPlans = 1;
5088  int howmany = sizeLines/numPlans;
5089  forward_plans = new fftwf_plan[numPlans];
5090  backward_plans = new fftwf_plan[numPlans];
5091  for(int i=0; i<numPlans; i++) {
5092  int curStride = i*howmany;
5093  forward_plans[i] = fftwf_plan_many_dft(1, planLineSizes, howmany,
5094  ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
5095  ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
5096  FFTW_FORWARD,
5097  fftwFlags);
5098 
5099  backward_plans[i] = fftwf_plan_many_dft(1, planLineSizes, howmany,
5100  ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
5101  ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
5102  FFTW_BACKWARD,
5103  fftwFlags);
5104  }
5105  }else
5106 #endif
5107  {
5108  forward_plans = NULL;
5109  backward_plans = NULL;
5110  }
5111 #else
5112  forward_plan = fftw_create_plan_specific(K1, FFTW_FORWARD,
5113  ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
5114  | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
5115  ny*nz, (fftw_complex *) work, 1);
5116  backward_plan = fftw_create_plan_specific(K1, FFTW_BACKWARD,
5117  ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
5118  | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
5119  ny*nz, (fftw_complex *) work, 1);
5120 #endif
5121  CmiUnlock(ComputePmeMgr::fftw_plan_lock);
5122 #else
5123  NAMD_die("Sorry, FFTW must be compiled in to use PME.");
5124 #endif
5125 
5127  thisIndex.y*block2, thisIndex.y*block2 + ny,
5128  thisIndex.z*block3, thisIndex.z*block3 + nz);
5129 
5130 }
static CmiNodeLock fftw_plan_lock
Definition: ComputePme.C:414
PmePencilInitMsgData initdata
Definition: ComputePme.C:4535
int dim3
Definition: PmeBase.h:19
Definition: Node.h:78
void order_init(int nBlocks)
Definition: ComputePme.C:4523
int K2
Definition: PmeBase.h:18
SimParameters * simParameters
Definition: Node.h:178
int K1
Definition: PmeBase.h:18
if(ComputeNonbondedUtil::goMethod==2)
CProxy_NodePmeMgr pmeNodeProxy
Definition: ComputePme.C:218
int block2
Definition: PmeBase.h:21
PmeKSpace * myKSpace
Definition: ComputePme.C:4733
int block3
Definition: PmeBase.h:21
void NAMD_die(const char *err_msg)
Definition: common.C:83
#define simParams
Definition: Output.C:127
fftw_plan backward_plan
Definition: ComputePme.C:4727
fftw_plan forward_plan
Definition: ComputePme.C:4727
#define fftwf_malloc
Definition: ComputePme.C:13
void PmeXPencil::forward_fft ( )

Definition at line 5653 of file ComputePme.C.

References CKLOOP_CTRL_PME_FORWARDFFT, PmePencil< CBase_PmeXPencil >::data, forward_plan, PmePencilInitMsgData::grid, PmePencil< CBase_PmeXPencil >::initdata, PmeGrid::K1, Node::Object(), PmeXZPencilFFT(), Node::simParameters, SimParameters::useCkLoop, PmePencil< CBase_PmeXPencil >::work, PmePencilInitMsgData::yBlocks, and PmePencilInitMsgData::zBlocks.

Referenced by node_process_trans().

5653  {
5654 #ifdef NAMD_FFTW
5655 
5656 #ifdef MANUAL_DEBUG_FFTW3
5657  dumpMatrixFloat3("fw_x_b", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
5658 #endif
5659 
5660 #ifdef NAMD_FFTW_3
5661 #if CMK_SMP && USE_CKLOOP
5662  int useCkLoop = Node::Object()->simParameters->useCkLoop;
5663  if(useCkLoop>=CKLOOP_CTRL_PME_FORWARDFFT
5664  && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks) {
5665  //for(int i=0; i<numPlans; i++) fftwf_execute(forward_plans[i]);
5666  //transform the above loop
5667  CkLoop_Parallelize(PmeXZPencilFFT, 1, (void *)forward_plans, CkMyNodeSize(), 0, numPlans-1); //sync
5668  return;
5669  }
5670 #endif
5671  fftwf_execute(forward_plan);
5672 #else
5673  fftw(forward_plan, ny*nz,
5674  ((fftw_complex *) data), ny*nz, 1, (fftw_complex *) work, 1, 0);
5675 #endif
5676 #ifdef MANUAL_DEBUG_FFTW3
5677  dumpMatrixFloat3("fw_x_a", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
5678 #endif
5679 
5680 #endif
5681 }
static Node * Object()
Definition: Node.h:86
PmePencilInitMsgData initdata
Definition: ComputePme.C:4535
SimParameters * simParameters
Definition: Node.h:178
int K1
Definition: PmeBase.h:18
static void PmeXZPencilFFT(int first, int last, void *result, int paraNum, void *param)
Definition: ComputePme.C:5177
#define CKLOOP_CTRL_PME_FORWARDFFT
Definition: SimParameters.h:94
fftw_plan forward_plan
Definition: ComputePme.C:4727
void PmeXPencil::node_process_trans ( PmeTransMsg msg)

Definition at line 5596 of file ComputePme.C.

References backward_fft(), forward_fft(), PmeTransMsg::hasData, PmePencil< CBase_PmeXPencil >::hasData, PmePencil< CBase_PmeXPencil >::imsg, PmePencil< CBase_PmeXPencil >::initdata, PmePencil< CBase_PmeXPencil >::needs_reply, pme_kspace(), recv_trans(), send_untrans(), PmeTransMsg::sourceNode, and PmePencilInitMsgData::xBlocks.

Referenced by NodePmeMgr::recvXTrans().

5597 {
5598  if(msg->hasData) hasData=1;
5599  needs_reply[msg->sourceNode] = msg->hasData;
5600  recv_trans(msg);
5601  int limsg;
5602  CmiMemoryAtomicFetchAndInc(imsg,limsg);
5603  if(limsg+1 == initdata.xBlocks)
5604  {
5605  if(hasData){
5606  forward_fft();
5607  pme_kspace();
5608  backward_fft();
5609  }
5610  send_untrans();
5611  imsg=0;
5612  CmiMemoryWriteFence();
5613  }
5614 }
PmePencilInitMsgData initdata
Definition: ComputePme.C:4535
void recv_trans(const PmeTransMsg *)
Definition: ComputePme.C:5616
void pme_kspace()
Definition: ComputePme.C:5683
void forward_fft()
Definition: ComputePme.C:5653
void send_untrans()
Definition: ComputePme.C:5803
int sourceNode
Definition: ComputePme.C:131
void backward_fft()
Definition: ComputePme.C:5712
void PmeXPencil::pme_kspace ( )

Definition at line 5683 of file ComputePme.C.

References CKLOOP_CTRL_PME_KSPACE, PmeKSpace::compute_energy(), PmePencil< CBase_PmeXPencil >::data, PmePencil< CBase_PmeXPencil >::evir, ComputeNonbondedUtil::ewaldcof, PmePencil< CBase_PmeXPencil >::initdata, PmePencil< CBase_PmeXPencil >::lattice, myKSpace, Node::Object(), PmePencilInitMsgData::yBlocks, and PmePencilInitMsgData::zBlocks.

Referenced by node_process_trans().

5683  {
5684 
5685  evir = 0.;
5686 
5687 #ifdef FFTCHECK
5688  return;
5689 #endif
5690 
5692 
5693  int useCkLoop = 0;
5694 #if CMK_SMP && USE_CKLOOP
5695  if ( Node::Object()->simParameters->useCkLoop >= CKLOOP_CTRL_PME_KSPACE
5696  && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks ) {
5697  useCkLoop = 1;
5698  }
5699 #endif
5700 
5701  int numGrids = 1;
5702  for ( int g=0; g<numGrids; ++g ) {
5703  evir[0] = myKSpace->compute_energy(data+0*g,
5704  lattice, ewaldcof, &(evir[1]), useCkLoop);
5705  }
5706 
5707 #if USE_NODE_PAR_RECEIVE
5708  CmiMemoryWriteFence();
5709 #endif
5710 }
static Node * Object()
Definition: Node.h:86
PmePencilInitMsgData initdata
Definition: ComputePme.C:4535
PmeKSpace * myKSpace
Definition: ComputePme.C:4733
double compute_energy(float q_arr[], const Lattice &lattice, double ewald, double virial[], int useCkLoop)
Definition: PmeKSpace.C:321
#define CKLOOP_CTRL_PME_KSPACE
Definition: SimParameters.h:96
double BigReal
Definition: common.h:114
void PmeXPencil::recv_trans ( const PmeTransMsg msg)

Definition at line 5616 of file ComputePme.C.

References PmeGrid::block1, PmePencil< CBase_PmeXPencil >::data, PmePencilInitMsgData::grid, PmeTransMsg::hasData, PmePencil< CBase_PmeXPencil >::imsg, PmePencil< CBase_PmeXPencil >::initdata, PmeGrid::K1, PmeTransMsg::lattice, PmePencil< CBase_PmeXPencil >::lattice, PmeTransMsg::nx, ny, nz, PmeTransMsg::qgrid, PmeTransMsg::sequence, PmePencil< CBase_PmeXPencil >::sequence, and PmeTransMsg::sourceNode.

Referenced by node_process_trans().

5616  {
5617  if ( imsg == 0 ) {
5618  lattice = msg->lattice;
5619  sequence = msg->sequence;
5620  }
5621  int block1 = initdata.grid.block1;
5622  int K1 = initdata.grid.K1;
5623  int ib = msg->sourceNode;
5624  int nx = msg->nx;
5625  if ( msg->hasData ) {
5626  const float *md = msg->qgrid;
5627  for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
5628  float *d = data + i*ny*nz*2;
5629  for ( int j=0; j<ny; ++j, d += nz*2 ) {
5630  for ( int k=0; k<nz; ++k ) {
5631 #ifdef ZEROCHECK
5632  if ( (*md) == 0. ) CkPrintf("0 in YX at %d %d %d %d %d %d %d %d %d\n",
5633  ib, thisIndex.y, thisIndex.z, i, j, k, nx, ny, nz);
5634 #endif
5635  d[2*k] = *(md++);
5636  d[2*k+1] = *(md++);
5637  }
5638  }
5639  }
5640  } else {
5641  for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
5642  float *d = data + i*ny*nz*2;
5643  for ( int j=0; j<ny; ++j, d += nz*2 ) {
5644  for ( int k=0; k<nz; ++k ) {
5645  d[2*k] = 0;
5646  d[2*k+1] = 0;
5647  }
5648  }
5649  }
5650  }
5651 }
PmePencilInitMsgData initdata
Definition: ComputePme.C:4535
int K1
Definition: PmeBase.h:18
int block1
Definition: PmeBase.h:21
float * qgrid
Definition: ComputePme.C:137
int sourceNode
Definition: ComputePme.C:131
Lattice lattice
Definition: ComputePme.C:134
void PmeXPencil::send_subset_untrans ( int  fromIdx,
int  toIdx 
)

Definition at line 5745 of file ComputePme.C.

References PmeGrid::block1, PmePencil< CBase_PmeXPencil >::data, PmeUntransMsg::destElem, PmePencilInitMsgData::grid, PmePencil< CBase_PmeXPencil >::initdata, PmeGrid::K1, PmePencil< CBase_PmeXPencil >::needs_reply, PmeUntransMsg::ny, ny, nz, PME_UNTRANS_PRIORITY, PmePencilInitMsgData::pmeNodeProxy, PRIORITY_SIZE, PmeUntransMsg::qgrid, PmePencil< CBase_PmeXPencil >::send_order, PmePencil< CBase_PmeXPencil >::sequence, SET_PRIORITY, PmeUntransMsg::sourceNode, PmePencilInitMsgData::xBlocks, PmePencilInitMsgData::ym, and PmePencilInitMsgData::yPencil.

Referenced by PmeXPencilSendUntrans().

5745  {
5746  int xBlocks = initdata.xBlocks;
5747  int block1 = initdata.grid.block1;
5748  int K1 = initdata.grid.K1;
5749 
5750  for(int isend=fromIdx; isend<=toIdx; isend++) {
5751  int ib = send_order[isend];
5752  if ( ! needs_reply[ib] ) {
5753  PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
5754  CmiEnableUrgentSend(1);
5756 #if USE_NODE_PAR_RECEIVE
5757  initdata.yPencil(ib,0,thisIndex.z).recvNodeAck(msg);
5758 #else
5759  initdata.yPencil(ib,0,thisIndex.z).recvAck(msg);
5760 #endif
5761  CmiEnableUrgentSend(0);
5762  continue;
5763  }
5764  int nx = block1;
5765  if ( (ib+1)*block1 > K1 ) nx = K1 - ib*block1;
5766  PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
5767  msg->sourceNode = thisIndex.y;
5768  msg->ny = ny;
5769  float *md = msg->qgrid;
5770  for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
5771  float *d = data + i*ny*nz*2;
5772  for ( int j=0; j<ny; ++j, d += nz*2 ) {
5773  for ( int k=0; k<nz; ++k ) {
5774  *(md++) = d[2*k];
5775  *(md++) = d[2*k+1];
5776  }
5777  }
5778  }
5780  CmiEnableUrgentSend(1);
5781 #if USE_NODE_PAR_RECEIVE
5782  msg->destElem=CkArrayIndex3D(ib,0, thisIndex.z);
5783 #if Y_PERSIST
5784  CmiUsePersistentHandle(&untrans_handle[isend], 1);
5785 #endif
5786  initdata.pmeNodeProxy[CmiNodeOf(initdata.ym.ckLocalBranch()->procNum(0,msg->destElem))].recvYUntrans(msg);
5787 #if Y_PERSIST
5788  CmiUsePersistentHandle(NULL, 0);
5789 #endif
5790 #else
5791 #if Y_PERSIST
5792  // CmiUsePersistentHandle(&untrans_handle[isend], 1);
5793 #endif
5794  initdata.yPencil(ib,0,thisIndex.z).recvUntrans(msg);
5795 #if Y_PERSIST
5796  // CmiUsePersistentHandle(NULL, 0);
5797 #endif
5798 #endif
5799  CmiEnableUrgentSend(0);
5800  }
5801 }
float * qgrid
Definition: ComputePme.C:154
PmePencilInitMsgData initdata
Definition: ComputePme.C:4535
int K1
Definition: PmeBase.h:18
int block1
Definition: PmeBase.h:21
CProxy_PmeYPencil yPencil
Definition: ComputePme.C:215
CProxy_NodePmeMgr pmeNodeProxy
Definition: ComputePme.C:218
#define PRIORITY_SIZE
Definition: Priorities.h:13
#define PME_UNTRANS_PRIORITY
Definition: Priorities.h:33
CProxy_PmePencilMap ym
Definition: ComputePme.C:220
#define SET_PRIORITY(MSG, SEQ, PRIO)
Definition: Priorities.h:18
CkArrayIndex3D destElem
Definition: ComputePme.C:155
void PmeXPencil::send_untrans ( )

Definition at line 5803 of file ComputePme.C.

References PmeGrid::block1, CKLOOP_CTRL_PME_SENDUNTRANS, PmePencil< CBase_PmeXPencil >::data, PmeUntransMsg::destElem, PmeEvirMsg::evir, PmePencil< CBase_PmeXPencil >::evir, PmePencilInitMsgData::grid, PmePencil< CBase_PmeXPencil >::initdata, PmeGrid::K1, PmePencil< CBase_PmeXPencil >::needs_reply, PmeUntransMsg::ny, ny, nz, Node::Object(), PME_UNGRID_PRIORITY, PME_UNTRANS_PRIORITY, PmePencilInitMsgData::pmeNodeProxy, PmePencilInitMsgData::pmeProxy, PmeXPencilSendUntrans(), PRIORITY_SIZE, PmeUntransMsg::qgrid, recipEvirPe, PmePencil< CBase_PmeXPencil >::send_order, PmePencil< CBase_PmeXPencil >::sequence, SET_PRIORITY, Node::simParameters, PmeUntransMsg::sourceNode, SimParameters::useCkLoop, PmePencilInitMsgData::xBlocks, PmePencilInitMsgData::yBlocks, PmePencilInitMsgData::ym, PmePencilInitMsgData::yPencil, and PmePencilInitMsgData::zBlocks.

Referenced by node_process_trans().

5803  {
5804 
5805  { // send energy and virial
5806  int numGrids = 1;
5807  PmeEvirMsg *newmsg = new (numGrids, PRIORITY_SIZE) PmeEvirMsg;
5808  newmsg->evir[0] = evir;
5810  CmiEnableUrgentSend(1);
5811  initdata.pmeProxy[recipEvirPe].recvRecipEvir(newmsg);
5812  CmiEnableUrgentSend(0);
5813  }
5814 
5815 #if USE_PERSISTENT
5816  if (untrans_handle == NULL) setup_persistent();
5817 #endif
5818 #if CMK_SMP && USE_CKLOOP
5819  int useCkLoop = Node::Object()->simParameters->useCkLoop;
5820  if(useCkLoop>=CKLOOP_CTRL_PME_SENDUNTRANS
5821  && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks) {
5822  int xBlocks = initdata.xBlocks;
5823 
5824 #if USE_NODE_PAR_RECEIVE
5825  //CkLoop_Parallelize(PmeXPencilSendUntrans, 1, (void *)this, CkMyNodeSize(), 0, xBlocks-1, 1); //has to sync
5826  CkLoop_Parallelize(PmeXPencilSendUntrans, 1, (void *)this, xBlocks, 0, xBlocks-1, 1); //has to sync
5827 #else
5828  //CkLoop_Parallelize(PmeXPencilSendUntrans, 1, (void *)this, CkMyNodeSize(), 0, xBlocks-1, 0); //not sync
5829  CkLoop_Parallelize(PmeXPencilSendUntrans, 1, (void *)this, xBlocks, 0, xBlocks-1, 0); //not sync
5830 #endif
5831  return;
5832  }
5833 #endif
5834  int xBlocks = initdata.xBlocks;
5835  int block1 = initdata.grid.block1;
5836  int K1 = initdata.grid.K1;
5837  for ( int isend=0; isend<xBlocks; ++isend ) {
5838  int ib = send_order[isend];
5839  if ( ! needs_reply[ib] ) {
5840  PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
5841  CmiEnableUrgentSend(1);
5843 #if USE_NODE_PAR_RECEIVE
5844  initdata.yPencil(ib,0,thisIndex.z).recvNodeAck(msg);
5845 #else
5846  initdata.yPencil(ib,0,thisIndex.z).recvAck(msg);
5847 #endif
5848  CmiEnableUrgentSend(0);
5849  continue;
5850  }
5851  int nx = block1;
5852  if ( (ib+1)*block1 > K1 ) nx = K1 - ib*block1;
5853  PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
5854  msg->sourceNode = thisIndex.y;
5855  msg->ny = ny;
5856  float *md = msg->qgrid;
5857  for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
5858  float *d = data + i*ny*nz*2;
5859  for ( int j=0; j<ny; ++j, d += nz*2 ) {
5860  for ( int k=0; k<nz; ++k ) {
5861  *(md++) = d[2*k];
5862  *(md++) = d[2*k+1];
5863  }
5864  }
5865  }
5867 
5868  CmiEnableUrgentSend(1);
5869 #if USE_NODE_PAR_RECEIVE
5870  msg->destElem=CkArrayIndex3D(ib,0, thisIndex.z);
5871 #if Y_PERSIST
5872  CmiUsePersistentHandle(&untrans_handle[isend], 1);
5873 #endif
5874  initdata.pmeNodeProxy[CmiNodeOf(initdata.ym.ckLocalBranch()->procNum(0,msg->destElem))].recvYUntrans(msg);
5875 #if Y_PERSIST
5876  CmiUsePersistentHandle(NULL, 0);
5877 #endif
5878 #else
5879 #if Y_PERSIST
5880  CmiUsePersistentHandle(&untrans_handle[isend], 1);
5881 #endif
5882  initdata.yPencil(ib,0,thisIndex.z).recvUntrans(msg);
5883 #if Y_PERSIST
5884  CmiUsePersistentHandle(NULL, 0);
5885 #endif
5886 #endif
5887  CmiEnableUrgentSend(0);
5888  }
5889 }
static Node * Object()
Definition: Node.h:86
float * qgrid
Definition: ComputePme.C:154
PmePencilInitMsgData initdata
Definition: ComputePme.C:4535
#define PME_UNGRID_PRIORITY
Definition: Priorities.h:74
SimParameters * simParameters
Definition: Node.h:178
int K1
Definition: PmeBase.h:18
int block1
Definition: PmeBase.h:21
CProxy_PmeYPencil yPencil
Definition: ComputePme.C:215
CProxy_NodePmeMgr pmeNodeProxy
Definition: ComputePme.C:218
PmeReduction * evir
Definition: ComputePme.C:167
#define CKLOOP_CTRL_PME_SENDUNTRANS
Definition: SimParameters.h:98
#define PRIORITY_SIZE
Definition: Priorities.h:13
static void PmeXPencilSendUntrans(int first, int last, void *result, int paraNum, void *param)
Definition: ComputePme.C:5740
int recipEvirPe
Definition: ComputePme.C:4731
#define PME_UNTRANS_PRIORITY
Definition: Priorities.h:33
CProxy_PmePencilMap ym
Definition: ComputePme.C:220
#define SET_PRIORITY(MSG, SEQ, PRIO)
Definition: Priorities.h:18
CkArrayIndex3D destElem
Definition: ComputePme.C:155

Member Data Documentation

fftw_plan PmeXPencil::backward_plan

Definition at line 4727 of file ComputePme.C.

Referenced by backward_fft(), and fft_init().

fftw_plan PmeXPencil::forward_plan

Definition at line 4727 of file ComputePme.C.

Referenced by fft_init(), and forward_fft().

PmeKSpace* PmeXPencil::myKSpace

Definition at line 4733 of file ComputePme.C.

Referenced by fft_init(), pme_kspace(), and PmeXPencil().

int PmeXPencil::ny

Definition at line 4730 of file ComputePme.C.

Referenced by recv_trans(), send_subset_untrans(), and send_untrans().

int PmeXPencil::nz

Definition at line 4730 of file ComputePme.C.

Referenced by fft_init(), recv_trans(), send_subset_untrans(), and send_untrans().

int PmeXPencil::recipEvirPe

Definition at line 4731 of file ComputePme.C.

Referenced by evir_init(), PmeXPencil(), and send_untrans().


The documentation for this class was generated from the following file: