PmeXPencil Class Reference

Inheritance diagram for PmeXPencil:

PmePencil< CBase_PmeXPencil > CBase_PmeXPencil List of all members.

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, int evirIdx)
void node_process_trans (PmeTransMsg *)
void evir_init ()

Public Attributes

fftw_plan forward_plan
fftw_plan backward_plan
int ny
int nz
int recipEvirPe
PmeKSpacemyKSpace

Detailed Description

Definition at line 4742 of file ComputePme.C.


Constructor & Destructor Documentation

PmeXPencil_SDAG_CODE PmeXPencil::PmeXPencil (  )  [inline]

Definition at line 4745 of file ComputePme.C.

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

04745 { __sdag_init();  myKSpace = 0; setMigratable(false); imsg=imsgb=0; recipEvirPe = -999; }

PmeXPencil::PmeXPencil ( CkMigrateMessage *   )  [inline]

Definition at line 4746 of file ComputePme.C.

04746 { __sdag_init(); }

PmeXPencil::~PmeXPencil (  )  [inline]

Definition at line 4747 of file ComputePme.C.

04747                       {
04748         #ifdef NAMD_FFTW
04749         #ifdef NAMD_FFTW_3
04750                 delete [] forward_plans;
04751                 delete [] backward_plans;
04752         #endif
04753         #endif
04754         }


Member Function Documentation

void PmeXPencil::backward_fft (  ) 

Definition at line 5748 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().

05748                               {
05749 #ifdef NAMD_FFTW
05750 #ifdef MANUAL_DEBUG_FFTW3
05751   dumpMatrixFloat3("bw_x_b", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05752 #endif
05753 
05754 #ifdef NAMD_FFTW_3
05755 #if     CMK_SMP && USE_CKLOOP
05756   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05757   if(useCkLoop>=CKLOOP_CTRL_PME_BACKWARDFFT
05758      && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks) {
05759           //for(int i=0; i<numPlans; i++) fftwf_execute(backward_plans[i]);
05760           //transform the above loop
05761           CkLoop_Parallelize(PmeXZPencilFFT, 1, (void *)backward_plans, CkMyNodeSize(), 0, numPlans-1); //sync
05762           return;
05763   }
05764 #endif
05765   fftwf_execute(backward_plan);
05766 #else
05767   fftw(backward_plan, ny*nz,
05768         ((fftw_complex *) data), ny*nz, 1, (fftw_complex *) work, 1, 0);
05769 #endif
05770 #ifdef MANUAL_DEBUG_FFTW3
05771   dumpMatrixFloat3("bw_x_a", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05772 #endif
05773 #endif
05774 }

void PmeXPencil::evir_init (  ) 

Definition at line 4800 of file ComputePme.C.

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

04800                            {
04801   recipEvirPe = findRecipEvirPe();
04802   initdata.pmeProxy[recipEvirPe].addRecipEvirClient();
04803 }

void PmeXPencil::fft_init (  ) 

Definition at line 5071 of file ComputePme.C.

References backward_plan, PmeGrid::block2, PmeGrid::block3, PmePencil< CBase_PmeXPencil >::data, PmeGrid::dim3, fftwf_malloc, forward_plan, PmePencilInitMsgData::grid, PmePencil< CBase_PmeXPencil >::initdata, PmeGrid::K1, PmeGrid::K2, myKSpace, NAMD_die(), PmePencil< CBase_PmeXPencil >::order_init(), PmePencilInitMsgData::pmeNodeProxy, Node::simParameters, simParams, PmePencil< CBase_PmeXPencil >::work, and PmePencilInitMsgData::xBlocks.

05071                           {
05072   CProxy_Node nd(CkpvAccess(BOCclass_group).node);
05073   Node *node = nd.ckLocalBranch();
05074   SimParameters *simParams = node->simParameters;
05075 #if USE_NODE_PAR_RECEIVE
05076   ((NodePmeMgr *)CkLocalNodeBranch(initdata.pmeNodeProxy))->registerXPencil(thisIndex,this);
05077 #endif
05078 
05079   int K1 = initdata.grid.K1;
05080   int K2 = initdata.grid.K2;
05081   int dim3 = initdata.grid.dim3;
05082   int block2 = initdata.grid.block2;
05083   int block3 = initdata.grid.block3;
05084 
05085   ny = block2;
05086   if ( (thisIndex.y + 1) * block2 > K2 ) ny = K2 - thisIndex.y * block2;
05087   nz = block3;
05088   if ( (thisIndex.z+1)*block3 > dim3/2 ) nz = dim3/2 - thisIndex.z*block3;
05089 
05090 #ifdef NAMD_FFTW
05091   CmiLock(ComputePmeMgr::fftw_plan_lock);
05092 
05093   data = (float *) fftwf_malloc( sizeof(float) * K1*ny*nz*2);
05094   work = new float[2*K1];
05095 
05096   order_init(initdata.xBlocks);
05097 
05098 #ifdef NAMD_FFTW_3
05099   /* need array of sizes for the how many */
05100   int fftwFlags = simParams->FFTWPatient ? FFTW_PATIENT  : simParams->FFTWEstimate ? FFTW_ESTIMATE  : FFTW_MEASURE ;
05101   int sizeLines=ny*nz;
05102   int planLineSizes[1];
05103   planLineSizes[0]=K1;
05104   forward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines,
05105                                      (fftwf_complex *) data, NULL, sizeLines, 1,
05106                                      (fftwf_complex *) data, NULL, sizeLines, 1,
05107                                    FFTW_FORWARD,
05108                                      fftwFlags);
05109   backward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines,
05110                                      (fftwf_complex *) data, NULL, sizeLines, 1,
05111                                      (fftwf_complex *) data, NULL, sizeLines, 1,
05112                                           FFTW_BACKWARD,
05113                                       fftwFlags);
05114 
05115 #if     CMK_SMP && USE_CKLOOP
05116   if(simParams->useCkLoop) {
05117           //How many FFT plans to be created? The grain-size issue!!.
05118           //Currently, I am choosing the min(nx, ny) to be coarse-grain
05119           numPlans = (ny<=nz?ny:nz);
05120           // limit attempted parallelism due to false sharing
05121           //if ( numPlans < CkMyNodeSize() ) numPlans = (ny>=nz?ny:nz);
05122           //if ( numPlans < CkMyNodeSize() ) numPlans = sizeLines;
05123           if ( sizeLines/numPlans < 4 ) numPlans = 1;
05124           int howmany = sizeLines/numPlans;
05125           forward_plans = new fftwf_plan[numPlans];
05126           backward_plans = new fftwf_plan[numPlans];
05127           for(int i=0; i<numPlans; i++) {
05128                   int curStride = i*howmany;              
05129                   forward_plans[i] = fftwf_plan_many_dft(1, planLineSizes, howmany,
05130                                                                                                          ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
05131                                                                                                          ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
05132                                                                                                         FFTW_FORWARD,
05133                                                                                                          fftwFlags);
05134 
05135                   backward_plans[i] = fftwf_plan_many_dft(1, planLineSizes, howmany,
05136                                                                                                          ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
05137                                                                                                          ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
05138                                                                                                           FFTW_BACKWARD,
05139                                                                                                          fftwFlags);
05140           }
05141   }else
05142 #endif
05143   {
05144           forward_plans = NULL;
05145           backward_plans = NULL;
05146   }
05147 #else
05148   forward_plan = fftw_create_plan_specific(K1, FFTW_FORWARD,
05149         ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
05150         | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
05151         ny*nz, (fftw_complex *) work, 1);
05152   backward_plan = fftw_create_plan_specific(K1, FFTW_BACKWARD,
05153         ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
05154         | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
05155         ny*nz, (fftw_complex *) work, 1);
05156 #endif
05157   CmiUnlock(ComputePmeMgr::fftw_plan_lock);
05158 #else
05159   NAMD_die("Sorry, FFTW must be compiled in to use PME.");
05160 #endif
05161 
05162   myKSpace = new PmeKSpace(initdata.grid,
05163                 thisIndex.y*block2, thisIndex.y*block2 + ny,
05164                 thisIndex.z*block3, thisIndex.z*block3 + nz);
05165 
05166 }

void PmeXPencil::forward_fft (  ) 

Definition at line 5689 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().

05689                              {
05690 #ifdef NAMD_FFTW
05691 
05692 #ifdef MANUAL_DEBUG_FFTW3
05693   dumpMatrixFloat3("fw_x_b", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05694 #endif
05695 
05696 #ifdef NAMD_FFTW_3
05697 #if     CMK_SMP && USE_CKLOOP
05698   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05699   if(useCkLoop>=CKLOOP_CTRL_PME_FORWARDFFT
05700      && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks) {
05701           //for(int i=0; i<numPlans; i++) fftwf_execute(forward_plans[i]);
05702           //transform the above loop
05703           CkLoop_Parallelize(PmeXZPencilFFT, 1, (void *)forward_plans, CkMyNodeSize(), 0, numPlans-1); //sync
05704           return;
05705   }
05706 #endif
05707   fftwf_execute(forward_plan);
05708 #else
05709   fftw(forward_plan, ny*nz,
05710         ((fftw_complex *) data), ny*nz, 1, (fftw_complex *) work, 1, 0);
05711 #endif
05712 #ifdef MANUAL_DEBUG_FFTW3
05713   dumpMatrixFloat3("fw_x_a", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05714 #endif
05715 
05716 #endif
05717 }

void PmeXPencil::node_process_trans ( PmeTransMsg  ) 

Definition at line 5632 of file ComputePme.C.

References backward_fft(), forward_fft(), PmePencil< CBase_PmeXPencil >::hasData, PmeTransMsg::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().

05633 {
05634   if(msg->hasData) hasData=1;
05635   needs_reply[msg->sourceNode] = msg->hasData;
05636   recv_trans(msg);
05637   int limsg;
05638   CmiMemoryAtomicFetchAndInc(imsg,limsg);
05639   if(limsg+1 == initdata.xBlocks)
05640     {
05641       if(hasData){
05642         forward_fft();
05643         pme_kspace();
05644         backward_fft();
05645       }
05646       send_untrans();
05647       imsg=0;
05648       CmiMemoryWriteFence();
05649     }
05650 }

void PmeXPencil::pme_kspace (  ) 

Definition at line 5719 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().

05719                             {
05720 
05721   evir = 0.;
05722 
05723 #ifdef FFTCHECK
05724   return;
05725 #endif
05726 
05727   BigReal ewaldcof = ComputeNonbondedUtil::ewaldcof;
05728 
05729   int useCkLoop = 0;
05730 #if CMK_SMP && USE_CKLOOP
05731   if ( Node::Object()->simParameters->useCkLoop >= CKLOOP_CTRL_PME_KSPACE
05732        && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks ) {
05733     useCkLoop = 1;
05734   }
05735 #endif
05736 
05737   int numGrids = 1;
05738   for ( int g=0; g<numGrids; ++g ) {
05739     evir[0] = myKSpace->compute_energy(data+0*g,
05740                 lattice, ewaldcof, &(evir[1]), useCkLoop);
05741   }
05742   
05743 #if USE_NODE_PAR_RECEIVE
05744     CmiMemoryWriteFence();
05745 #endif
05746 }

void PmeXPencil::recv_trans ( const PmeTransMsg  ) 

Definition at line 5652 of file ComputePme.C.

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

Referenced by node_process_trans().

05652                                                   {
05653   if ( imsg == 0 ) {
05654     lattice = msg->lattice;
05655     sequence = msg->sequence;
05656   }
05657   int block1 = initdata.grid.block1;
05658   int K1 = initdata.grid.K1;
05659   int ib = msg->sourceNode;
05660   int nx = msg->nx;
05661  if ( msg->hasData ) {
05662   const float *md = msg->qgrid;
05663   for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
05664    float *d = data + i*ny*nz*2;
05665    for ( int j=0; j<ny; ++j, d += nz*2 ) {
05666     for ( int k=0; k<nz; ++k ) {
05667 #ifdef ZEROCHECK
05668       if ( (*md) == 0. ) CkPrintf("0 in YX at %d %d %d %d %d %d %d %d %d\n",
05669         ib, thisIndex.y, thisIndex.z, i, j, k, nx, ny, nz);
05670 #endif
05671       d[2*k] = *(md++);
05672       d[2*k+1] = *(md++);
05673     }
05674    }
05675   }
05676  } else {
05677   for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
05678    float *d = data + i*ny*nz*2;
05679    for ( int j=0; j<ny; ++j, d += nz*2 ) {
05680     for ( int k=0; k<nz; ++k ) {
05681       d[2*k] = 0;
05682       d[2*k+1] = 0;
05683     }
05684    }
05685   }
05686  }
05687 }

void PmeXPencil::send_subset_untrans ( int  fromIdx,
int  toIdx,
int  evirIdx 
)

Definition at line 5782 of file ComputePme.C.

References PmeGrid::block1, PmePencil< CBase_PmeXPencil >::data, PmeUntransMsg::destElem, PmePencilInitMsgData::grid, PmePencil< CBase_PmeXPencil >::initdata, j, PmeGrid::K1, PmeUntransMsg::ny, 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().

05782                                                                        {
05783         int xBlocks = initdata.xBlocks;
05784         int block1 = initdata.grid.block1;      
05785         int K1 = initdata.grid.K1;
05786 
05787         int ackL=0, ackH=-1;
05788         int unL=0, unH=-1;
05789         int send_evir=0;
05790         if(fromIdx >= evirIdx+1) {
05791                 //send PmeUntransMsg with has_evir=0
05792                 unL = fromIdx;
05793                 unH = toIdx;            
05794         } else if(toIdx <= evirIdx-1) {
05795                 //send PmeAckMsg
05796                 ackL=fromIdx;
05797                 ackH=toIdx;             
05798         } else {
05799                 //partially send PmeAckMsg and partially send PmeUntransMsg
05800                 ackL=fromIdx;
05801                 ackH=evirIdx-1;
05802                 send_evir=1;
05803                 unL=evirIdx+1;
05804                 unH=toIdx;
05805         }
05806 
05807         for(int isend=ackL; isend<=ackH; isend++) {
05808                 //send PmeAckMsg
05809         CmiEnableUrgentSend(1);
05810                 int ib = send_order[isend];
05811                 PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
05812                 SET_PRIORITY(msg,sequence,PME_UNTRANS_PRIORITY)
05813                 initdata.yPencil(ib,0,thisIndex.z).recvAck(msg);
05814         CmiEnableUrgentSend(0);
05815     }
05816 
05817     CmiEnableUrgentSend(1);
05818         //send PmeUntransMsg with has_evir=1
05819         if(send_evir) {
05820                 int ib = send_order[evirIdx];
05821                 int nx = block1;
05822                 if ( (ib+1)*block1 > K1 ) nx = K1 - ib*block1;
05823                 PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;              
05824                 msg->sourceNode = thisIndex.y;
05825                 msg->ny = ny;
05826                 float *md = msg->qgrid;
05827                 for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
05828                         float *d = data + i*ny*nz*2;
05829                         for ( int j=0; j<ny; ++j, d += nz*2 ) {
05830                                 for ( int k=0; k<nz; ++k ) {
05831                                         *(md++) = d[2*k];
05832                                         *(md++) = d[2*k+1];
05833                                 }
05834                         }
05835                 }
05836                 SET_PRIORITY(msg,sequence,PME_UNTRANS_PRIORITY)
05837 #if USE_NODE_PAR_RECEIVE
05838         msg->destElem=CkArrayIndex3D(ib,0, thisIndex.z);
05839         initdata.pmeNodeProxy[CmiNodeOf(initdata.ym.ckLocalBranch()->procNum(0,msg->destElem))].recvYUntrans(msg);
05840 #else
05841         initdata.yPencil(ib,0,thisIndex.z).recvUntrans(msg);
05842 #endif
05843          }
05844     CmiEnableUrgentSend(0);
05845         
05846         //send PmeUntransMsg with has_evir=0
05847         for(int isend=unL; isend<=unH; isend++) {
05848                 int ib = send_order[isend];
05849                 int nx = block1;
05850                 if ( (ib+1)*block1 > K1 ) nx = K1 - ib*block1;
05851                 PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
05852                 msg->sourceNode = thisIndex.y;
05853                 msg->ny = ny;
05854                 float *md = msg->qgrid;
05855                 for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
05856                         float *d = data + i*ny*nz*2;
05857                         for ( int j=0; j<ny; ++j, d += nz*2 ) {
05858                                 for ( int k=0; k<nz; ++k ) {
05859                                         *(md++) = d[2*k];
05860                                         *(md++) = d[2*k+1];
05861                                 }
05862                         }
05863                 }
05864                 SET_PRIORITY(msg,sequence,PME_UNTRANS_PRIORITY)
05865         CmiEnableUrgentSend(1);
05866 #if USE_NODE_PAR_RECEIVE
05867         msg->destElem=CkArrayIndex3D(ib,0, thisIndex.z);
05868 #if Y_PERSIST 
05869         CmiUsePersistentHandle(&untrans_handle[isend], 1);
05870 #endif
05871         initdata.pmeNodeProxy[CmiNodeOf(initdata.ym.ckLocalBranch()->procNum(0,msg->destElem))].recvYUntrans(msg);
05872 #if Y_PERSIST 
05873         CmiUsePersistentHandle(NULL, 0);
05874 #endif
05875 #else
05876 #if Y_PERSIST 
05877   //      CmiUsePersistentHandle(&untrans_handle[isend], 1);
05878 #endif
05879         initdata.yPencil(ib,0,thisIndex.z).recvUntrans(msg);
05880 #if Y_PERSIST 
05881    //     CmiUsePersistentHandle(NULL, 0);
05882 #endif
05883 #endif
05884         CmiEnableUrgentSend(0);
05885         }
05886 }

void PmeXPencil::send_untrans (  ) 

Definition at line 5888 of file ComputePme.C.

References PmeGrid::block1, CKLOOP_CTRL_PME_SENDUNTRANS, PmePencil< CBase_PmeXPencil >::data, PmeUntransMsg::destElem, PmePencil< CBase_PmeXPencil >::evir, PmeEvirMsg::evir, PmePencilInitMsgData::grid, PmePencil< CBase_PmeXPencil >::initdata, j, PmeGrid::K1, PmePencil< CBase_PmeXPencil >::needs_reply, PmeUntransMsg::ny, 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().

05888                               {
05889 
05890   { // send energy and virial
05891     int numGrids = 1;
05892     PmeEvirMsg *newmsg = new (numGrids, PRIORITY_SIZE) PmeEvirMsg;
05893     newmsg->evir[0] = evir;
05894     SET_PRIORITY(newmsg,sequence,PME_UNGRID_PRIORITY)
05895     CmiEnableUrgentSend(1);
05896     initdata.pmeProxy[recipEvirPe].recvRecipEvir(newmsg);
05897     CmiEnableUrgentSend(0);
05898   }
05899 
05900 #if USE_PERSISTENT
05901   if (untrans_handle == NULL) setup_persistent();
05902 #endif
05903 #if     CMK_SMP && USE_CKLOOP
05904   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05905   if(useCkLoop>=CKLOOP_CTRL_PME_SENDUNTRANS
05906      && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks) {
05907                 int xBlocks = initdata.xBlocks;
05908                 int evirIdx = 0;
05909                 for ( int isend=0; isend<xBlocks; ++isend ) {
05910                         int ib = send_order[isend];
05911                         if (needs_reply[ib]) {
05912                                 evirIdx = isend;
05913                                 break;
05914                         }
05915                 }
05916 
05917                 //basically: 
05918                 //[0,evirIdx-1]->send PmeAckMsg
05919                 //evirIdx->send PmeUntransMsg with has_evir=1
05920                 //[evirIdx+1, xBlocks-1]->send PmeUntransMsg with has_evir=0
05921                 //send_subset_untrans(0, xBlocks-1, evirIdx);
05922 #if USE_NODE_PAR_RECEIVE
05923                 //CkLoop_Parallelize(PmeXPencilSendUntrans, evirIdx, (void *)this, CkMyNodeSize(), 0, xBlocks-1, 1); //has to sync
05924                 CkLoop_Parallelize(PmeXPencilSendUntrans, evirIdx, (void *)this, xBlocks, 0, xBlocks-1, 1); //has to sync
05925 #else
05926         //CkLoop_Parallelize(PmeXPencilSendUntrans, evirIdx, (void *)this, CkMyNodeSize(), 0, xBlocks-1, 0); //not sync
05927                 CkLoop_Parallelize(PmeXPencilSendUntrans, evirIdx, (void *)this, xBlocks, 0, xBlocks-1, 0); //not sync
05928 #endif        
05929                 return;
05930   }
05931 #endif
05932   int xBlocks = initdata.xBlocks;
05933   int block1 = initdata.grid.block1;
05934   int K1 = initdata.grid.K1;
05935   int send_evir = 1;
05936   for ( int isend=0; isend<xBlocks; ++isend ) {
05937     int ib = send_order[isend];
05938     if ( ! needs_reply[ib] ) {
05939       PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
05940       CmiEnableUrgentSend(1);
05941       SET_PRIORITY(msg,sequence,PME_UNTRANS_PRIORITY)
05942       initdata.yPencil(ib,0,thisIndex.z).recvAck(msg);
05943       CmiEnableUrgentSend(0);
05944       continue;
05945     }
05946     int nx = block1;
05947     if ( (ib+1)*block1 > K1 ) nx = K1 - ib*block1;
05948     PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
05949     if ( send_evir ) {
05950       send_evir = 0;
05951     }
05952     msg->sourceNode = thisIndex.y;
05953     msg->ny = ny;
05954     float *md = msg->qgrid;
05955     for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
05956      float *d = data + i*ny*nz*2;
05957      for ( int j=0; j<ny; ++j, d += nz*2 ) {
05958       for ( int k=0; k<nz; ++k ) {
05959         *(md++) = d[2*k];
05960         *(md++) = d[2*k+1];
05961       }
05962      }
05963     }
05964     SET_PRIORITY(msg,sequence,PME_UNTRANS_PRIORITY)
05965 
05966     CmiEnableUrgentSend(1);
05967 #if USE_NODE_PAR_RECEIVE
05968     msg->destElem=CkArrayIndex3D(ib,0, thisIndex.z);
05969 #if Y_PERSIST 
05970     CmiUsePersistentHandle(&untrans_handle[isend], 1);
05971 #endif
05972     initdata.pmeNodeProxy[CmiNodeOf(initdata.ym.ckLocalBranch()->procNum(0,msg->destElem))].recvYUntrans(msg);
05973 #if Y_PERSIST 
05974     CmiUsePersistentHandle(NULL, 0);
05975 #endif
05976 #else
05977 #if Y_PERSIST 
05978     CmiUsePersistentHandle(&untrans_handle[isend], 1);
05979 #endif
05980     initdata.yPencil(ib,0,thisIndex.z).recvUntrans(msg);
05981 #if Y_PERSIST 
05982     CmiUsePersistentHandle(NULL, 0);
05983 #endif
05984 #endif
05985     CmiEnableUrgentSend(0);
05986   }
05987 }


Member Data Documentation

fftw_plan PmeXPencil::backward_plan

Definition at line 4770 of file ComputePme.C.

Referenced by backward_fft(), and fft_init().

fftw_plan PmeXPencil::forward_plan

Definition at line 4770 of file ComputePme.C.

Referenced by fft_init(), and forward_fft().

PmeKSpace* PmeXPencil::myKSpace

Definition at line 4776 of file ComputePme.C.

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

int PmeXPencil::ny

Definition at line 4773 of file ComputePme.C.

int PmeXPencil::nz

Definition at line 4773 of file ComputePme.C.

int PmeXPencil::recipEvirPe

Definition at line 4774 of file ComputePme.C.

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


The documentation for this class was generated from the following file:
Generated on Thu Nov 23 01:17:19 2017 for NAMD by  doxygen 1.4.7