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 4704 of file ComputePme.C.


Constructor & Destructor Documentation

PmeXPencil_SDAG_CODE PmeXPencil::PmeXPencil (  )  [inline]

Definition at line 4707 of file ComputePme.C.

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

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

PmeXPencil::PmeXPencil ( CkMigrateMessage *   )  [inline]

Definition at line 4708 of file ComputePme.C.

04708 { __sdag_init(); }

PmeXPencil::~PmeXPencil (  )  [inline]

Definition at line 4709 of file ComputePme.C.

04709                       {
04710         #ifdef NAMD_FFTW
04711         #ifdef NAMD_FFTW_3
04712                 delete [] forward_plans;
04713                 delete [] backward_plans;
04714         #endif
04715         #endif
04716         }


Member Function Documentation

void PmeXPencil::backward_fft (  ) 

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

05710                               {
05711 #ifdef NAMD_FFTW
05712 #ifdef MANUAL_DEBUG_FFTW3
05713   dumpMatrixFloat3("bw_x_b", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05714 #endif
05715 
05716 #ifdef NAMD_FFTW_3
05717 #if     CMK_SMP && USE_CKLOOP
05718   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05719   if(useCkLoop>=CKLOOP_CTRL_PME_BACKWARDFFT
05720      && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks) {
05721           //for(int i=0; i<numPlans; i++) fftwf_execute(backward_plans[i]);
05722           //transform the above loop
05723           CkLoop_Parallelize(PmeXZPencilFFT, 1, (void *)backward_plans, CkMyNodeSize(), 0, numPlans-1); //sync
05724           return;
05725   }
05726 #endif
05727   fftwf_execute(backward_plan);
05728 #else
05729   fftw(backward_plan, ny*nz,
05730         ((fftw_complex *) data), ny*nz, 1, (fftw_complex *) work, 1, 0);
05731 #endif
05732 #ifdef MANUAL_DEBUG_FFTW3
05733   dumpMatrixFloat3("bw_x_a", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05734 #endif
05735 #endif
05736 }

void PmeXPencil::evir_init (  ) 

Definition at line 4762 of file ComputePme.C.

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

04762                            {
04763   recipEvirPe = findRecipEvirPe();
04764   initdata.pmeProxy[recipEvirPe].addRecipEvirClient();
04765 }

void PmeXPencil::fft_init (  ) 

Definition at line 5033 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.

05033                           {
05034   CProxy_Node nd(CkpvAccess(BOCclass_group).node);
05035   Node *node = nd.ckLocalBranch();
05036   SimParameters *simParams = node->simParameters;
05037 #if USE_NODE_PAR_RECEIVE
05038   ((NodePmeMgr *)CkLocalNodeBranch(initdata.pmeNodeProxy))->registerXPencil(thisIndex,this);
05039 #endif
05040 
05041   int K1 = initdata.grid.K1;
05042   int K2 = initdata.grid.K2;
05043   int dim3 = initdata.grid.dim3;
05044   int block2 = initdata.grid.block2;
05045   int block3 = initdata.grid.block3;
05046 
05047   ny = block2;
05048   if ( (thisIndex.y + 1) * block2 > K2 ) ny = K2 - thisIndex.y * block2;
05049   nz = block3;
05050   if ( (thisIndex.z+1)*block3 > dim3/2 ) nz = dim3/2 - thisIndex.z*block3;
05051 
05052 #ifdef NAMD_FFTW
05053   CmiLock(ComputePmeMgr::fftw_plan_lock);
05054 
05055   data = (float *) fftwf_malloc( sizeof(float) * K1*ny*nz*2);
05056   work = new float[2*K1];
05057 
05058   order_init(initdata.xBlocks);
05059 
05060 #ifdef NAMD_FFTW_3
05061   /* need array of sizes for the how many */
05062   int fftwFlags = simParams->FFTWPatient ? FFTW_PATIENT  : simParams->FFTWEstimate ? FFTW_ESTIMATE  : FFTW_MEASURE ;
05063   int sizeLines=ny*nz;
05064   int planLineSizes[1];
05065   planLineSizes[0]=K1;
05066   forward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines,
05067                                      (fftwf_complex *) data, NULL, sizeLines, 1,
05068                                      (fftwf_complex *) data, NULL, sizeLines, 1,
05069                                    FFTW_FORWARD,
05070                                      fftwFlags);
05071   backward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines,
05072                                      (fftwf_complex *) data, NULL, sizeLines, 1,
05073                                      (fftwf_complex *) data, NULL, sizeLines, 1,
05074                                           FFTW_BACKWARD,
05075                                       fftwFlags);
05076 
05077 #if     CMK_SMP && USE_CKLOOP
05078   if(simParams->useCkLoop) {
05079           //How many FFT plans to be created? The grain-size issue!!.
05080           //Currently, I am choosing the min(nx, ny) to be coarse-grain
05081           numPlans = (ny<=nz?ny:nz);
05082           // limit attempted parallelism due to false sharing
05083           //if ( numPlans < CkMyNodeSize() ) numPlans = (ny>=nz?ny:nz);
05084           //if ( numPlans < CkMyNodeSize() ) numPlans = sizeLines;
05085           if ( sizeLines/numPlans < 4 ) numPlans = 1;
05086           int howmany = sizeLines/numPlans;
05087           forward_plans = new fftwf_plan[numPlans];
05088           backward_plans = new fftwf_plan[numPlans];
05089           for(int i=0; i<numPlans; i++) {
05090                   int curStride = i*howmany;              
05091                   forward_plans[i] = fftwf_plan_many_dft(1, planLineSizes, howmany,
05092                                                                                                          ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
05093                                                                                                          ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
05094                                                                                                         FFTW_FORWARD,
05095                                                                                                          fftwFlags);
05096 
05097                   backward_plans[i] = fftwf_plan_many_dft(1, planLineSizes, howmany,
05098                                                                                                          ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
05099                                                                                                          ((fftwf_complex *)data)+curStride, NULL, sizeLines, 1,
05100                                                                                                           FFTW_BACKWARD,
05101                                                                                                          fftwFlags);
05102           }
05103   }else
05104 #endif
05105   {
05106           forward_plans = NULL;
05107           backward_plans = NULL;
05108   }
05109 #else
05110   forward_plan = fftw_create_plan_specific(K1, FFTW_FORWARD,
05111         ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
05112         | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
05113         ny*nz, (fftw_complex *) work, 1);
05114   backward_plan = fftw_create_plan_specific(K1, FFTW_BACKWARD,
05115         ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
05116         | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
05117         ny*nz, (fftw_complex *) work, 1);
05118 #endif
05119   CmiUnlock(ComputePmeMgr::fftw_plan_lock);
05120 #else
05121   NAMD_die("Sorry, FFTW must be compiled in to use PME.");
05122 #endif
05123 
05124   myKSpace = new PmeKSpace(initdata.grid,
05125                 thisIndex.y*block2, thisIndex.y*block2 + ny,
05126                 thisIndex.z*block3, thisIndex.z*block3 + nz);
05127 
05128 }

void PmeXPencil::forward_fft (  ) 

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

05651                              {
05652 #ifdef NAMD_FFTW
05653 
05654 #ifdef MANUAL_DEBUG_FFTW3
05655   dumpMatrixFloat3("fw_x_b", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05656 #endif
05657 
05658 #ifdef NAMD_FFTW_3
05659 #if     CMK_SMP && USE_CKLOOP
05660   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05661   if(useCkLoop>=CKLOOP_CTRL_PME_FORWARDFFT
05662      && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks) {
05663           //for(int i=0; i<numPlans; i++) fftwf_execute(forward_plans[i]);
05664           //transform the above loop
05665           CkLoop_Parallelize(PmeXZPencilFFT, 1, (void *)forward_plans, CkMyNodeSize(), 0, numPlans-1); //sync
05666           return;
05667   }
05668 #endif
05669   fftwf_execute(forward_plan);
05670 #else
05671   fftw(forward_plan, ny*nz,
05672         ((fftw_complex *) data), ny*nz, 1, (fftw_complex *) work, 1, 0);
05673 #endif
05674 #ifdef MANUAL_DEBUG_FFTW3
05675   dumpMatrixFloat3("fw_x_a", data, initdata.grid.K1, ny, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05676 #endif
05677 
05678 #endif
05679 }

void PmeXPencil::node_process_trans ( PmeTransMsg  ) 

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

05595 {
05596   if(msg->hasData) hasData=1;
05597   needs_reply[msg->sourceNode] = msg->hasData;
05598   recv_trans(msg);
05599   int limsg;
05600   CmiMemoryAtomicFetchAndInc(imsg,limsg);
05601   if(limsg+1 == initdata.xBlocks)
05602     {
05603       if(hasData){
05604         forward_fft();
05605         pme_kspace();
05606         backward_fft();
05607       }
05608       send_untrans();
05609       imsg=0;
05610       CmiMemoryWriteFence();
05611     }
05612 }

void PmeXPencil::pme_kspace (  ) 

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

05681                             {
05682 
05683   evir = 0.;
05684 
05685 #ifdef FFTCHECK
05686   return;
05687 #endif
05688 
05689   BigReal ewaldcof = ComputeNonbondedUtil::ewaldcof;
05690 
05691   int useCkLoop = 0;
05692 #if CMK_SMP && USE_CKLOOP
05693   if ( Node::Object()->simParameters->useCkLoop >= CKLOOP_CTRL_PME_KSPACE
05694        && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks ) {
05695     useCkLoop = 1;
05696   }
05697 #endif
05698 
05699   int numGrids = 1;
05700   for ( int g=0; g<numGrids; ++g ) {
05701     evir[0] = myKSpace->compute_energy(data+0*g,
05702                 lattice, ewaldcof, &(evir[1]), useCkLoop);
05703   }
05704   
05705 #if USE_NODE_PAR_RECEIVE
05706     CmiMemoryWriteFence();
05707 #endif
05708 }

void PmeXPencil::recv_trans ( const PmeTransMsg  ) 

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

05614                                                   {
05615   if ( imsg == 0 ) {
05616     lattice = msg->lattice;
05617     sequence = msg->sequence;
05618   }
05619   int block1 = initdata.grid.block1;
05620   int K1 = initdata.grid.K1;
05621   int ib = msg->sourceNode;
05622   int nx = msg->nx;
05623  if ( msg->hasData ) {
05624   const float *md = msg->qgrid;
05625   for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
05626    float *d = data + i*ny*nz*2;
05627    for ( int j=0; j<ny; ++j, d += nz*2 ) {
05628     for ( int k=0; k<nz; ++k ) {
05629 #ifdef ZEROCHECK
05630       if ( (*md) == 0. ) CkPrintf("0 in YX at %d %d %d %d %d %d %d %d %d\n",
05631         ib, thisIndex.y, thisIndex.z, i, j, k, nx, ny, nz);
05632 #endif
05633       d[2*k] = *(md++);
05634       d[2*k+1] = *(md++);
05635     }
05636    }
05637   }
05638  } else {
05639   for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
05640    float *d = data + i*ny*nz*2;
05641    for ( int j=0; j<ny; ++j, d += nz*2 ) {
05642     for ( int k=0; k<nz; ++k ) {
05643       d[2*k] = 0;
05644       d[2*k+1] = 0;
05645     }
05646    }
05647   }
05648  }
05649 }

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

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

05744                                                                        {
05745         int xBlocks = initdata.xBlocks;
05746         int block1 = initdata.grid.block1;      
05747         int K1 = initdata.grid.K1;
05748 
05749         int ackL=0, ackH=-1;
05750         int unL=0, unH=-1;
05751         int send_evir=0;
05752         if(fromIdx >= evirIdx+1) {
05753                 //send PmeUntransMsg with has_evir=0
05754                 unL = fromIdx;
05755                 unH = toIdx;            
05756         } else if(toIdx <= evirIdx-1) {
05757                 //send PmeAckMsg
05758                 ackL=fromIdx;
05759                 ackH=toIdx;             
05760         } else {
05761                 //partially send PmeAckMsg and partially send PmeUntransMsg
05762                 ackL=fromIdx;
05763                 ackH=evirIdx-1;
05764                 send_evir=1;
05765                 unL=evirIdx+1;
05766                 unH=toIdx;
05767         }
05768 
05769         for(int isend=ackL; isend<=ackH; isend++) {
05770                 //send PmeAckMsg
05771         CmiEnableUrgentSend(1);
05772                 int ib = send_order[isend];
05773                 PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
05774                 SET_PRIORITY(msg,sequence,PME_UNTRANS_PRIORITY)
05775                 initdata.yPencil(ib,0,thisIndex.z).recvAck(msg);
05776         CmiEnableUrgentSend(0);
05777     }
05778 
05779     CmiEnableUrgentSend(1);
05780         //send PmeUntransMsg with has_evir=1
05781         if(send_evir) {
05782                 int ib = send_order[evirIdx];
05783                 int nx = block1;
05784                 if ( (ib+1)*block1 > K1 ) nx = K1 - ib*block1;
05785                 PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;              
05786                 msg->sourceNode = thisIndex.y;
05787                 msg->ny = ny;
05788                 float *md = msg->qgrid;
05789                 for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
05790                         float *d = data + i*ny*nz*2;
05791                         for ( int j=0; j<ny; ++j, d += nz*2 ) {
05792                                 for ( int k=0; k<nz; ++k ) {
05793                                         *(md++) = d[2*k];
05794                                         *(md++) = d[2*k+1];
05795                                 }
05796                         }
05797                 }
05798                 SET_PRIORITY(msg,sequence,PME_UNTRANS_PRIORITY)
05799 #if USE_NODE_PAR_RECEIVE
05800         msg->destElem=CkArrayIndex3D(ib,0, thisIndex.z);
05801         initdata.pmeNodeProxy[CmiNodeOf(initdata.ym.ckLocalBranch()->procNum(0,msg->destElem))].recvYUntrans(msg);
05802 #else
05803         initdata.yPencil(ib,0,thisIndex.z).recvUntrans(msg);
05804 #endif
05805          }
05806     CmiEnableUrgentSend(0);
05807         
05808         //send PmeUntransMsg with has_evir=0
05809         for(int isend=unL; isend<=unH; isend++) {
05810                 int ib = send_order[isend];
05811                 int nx = block1;
05812                 if ( (ib+1)*block1 > K1 ) nx = K1 - ib*block1;
05813                 PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
05814                 msg->sourceNode = thisIndex.y;
05815                 msg->ny = ny;
05816                 float *md = msg->qgrid;
05817                 for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
05818                         float *d = data + i*ny*nz*2;
05819                         for ( int j=0; j<ny; ++j, d += nz*2 ) {
05820                                 for ( int k=0; k<nz; ++k ) {
05821                                         *(md++) = d[2*k];
05822                                         *(md++) = d[2*k+1];
05823                                 }
05824                         }
05825                 }
05826                 SET_PRIORITY(msg,sequence,PME_UNTRANS_PRIORITY)
05827         CmiEnableUrgentSend(1);
05828 #if USE_NODE_PAR_RECEIVE
05829         msg->destElem=CkArrayIndex3D(ib,0, thisIndex.z);
05830 #if Y_PERSIST 
05831         CmiUsePersistentHandle(&untrans_handle[isend], 1);
05832 #endif
05833         initdata.pmeNodeProxy[CmiNodeOf(initdata.ym.ckLocalBranch()->procNum(0,msg->destElem))].recvYUntrans(msg);
05834 #if Y_PERSIST 
05835         CmiUsePersistentHandle(NULL, 0);
05836 #endif
05837 #else
05838 #if Y_PERSIST 
05839   //      CmiUsePersistentHandle(&untrans_handle[isend], 1);
05840 #endif
05841         initdata.yPencil(ib,0,thisIndex.z).recvUntrans(msg);
05842 #if Y_PERSIST 
05843    //     CmiUsePersistentHandle(NULL, 0);
05844 #endif
05845 #endif
05846         CmiEnableUrgentSend(0);
05847         }
05848 }

void PmeXPencil::send_untrans (  ) 

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

05850                               {
05851 
05852   { // send energy and virial
05853     int numGrids = 1;
05854     PmeEvirMsg *newmsg = new (numGrids, PRIORITY_SIZE) PmeEvirMsg;
05855     newmsg->evir[0] = evir;
05856     SET_PRIORITY(newmsg,sequence,PME_UNGRID_PRIORITY)
05857     CmiEnableUrgentSend(1);
05858     initdata.pmeProxy[recipEvirPe].recvRecipEvir(newmsg);
05859     CmiEnableUrgentSend(0);
05860   }
05861 
05862 #if USE_PERSISTENT
05863   if (untrans_handle == NULL) setup_persistent();
05864 #endif
05865 #if     CMK_SMP && USE_CKLOOP
05866   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05867   if(useCkLoop>=CKLOOP_CTRL_PME_SENDUNTRANS
05868      && CkNumPes() >= 2 * initdata.yBlocks * initdata.zBlocks) {
05869                 int xBlocks = initdata.xBlocks;
05870                 int evirIdx = 0;
05871                 for ( int isend=0; isend<xBlocks; ++isend ) {
05872                         int ib = send_order[isend];
05873                         if (needs_reply[ib]) {
05874                                 evirIdx = isend;
05875                                 break;
05876                         }
05877                 }
05878 
05879                 //basically: 
05880                 //[0,evirIdx-1]->send PmeAckMsg
05881                 //evirIdx->send PmeUntransMsg with has_evir=1
05882                 //[evirIdx+1, xBlocks-1]->send PmeUntransMsg with has_evir=0
05883                 //send_subset_untrans(0, xBlocks-1, evirIdx);
05884 #if USE_NODE_PAR_RECEIVE
05885                 //CkLoop_Parallelize(PmeXPencilSendUntrans, evirIdx, (void *)this, CkMyNodeSize(), 0, xBlocks-1, 1); //has to sync
05886                 CkLoop_Parallelize(PmeXPencilSendUntrans, evirIdx, (void *)this, xBlocks, 0, xBlocks-1, 1); //has to sync
05887 #else
05888         //CkLoop_Parallelize(PmeXPencilSendUntrans, evirIdx, (void *)this, CkMyNodeSize(), 0, xBlocks-1, 0); //not sync
05889                 CkLoop_Parallelize(PmeXPencilSendUntrans, evirIdx, (void *)this, xBlocks, 0, xBlocks-1, 0); //not sync
05890 #endif        
05891                 return;
05892   }
05893 #endif
05894   int xBlocks = initdata.xBlocks;
05895   int block1 = initdata.grid.block1;
05896   int K1 = initdata.grid.K1;
05897   int send_evir = 1;
05898   for ( int isend=0; isend<xBlocks; ++isend ) {
05899     int ib = send_order[isend];
05900     if ( ! needs_reply[ib] ) {
05901       PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
05902       CmiEnableUrgentSend(1);
05903       SET_PRIORITY(msg,sequence,PME_UNTRANS_PRIORITY)
05904       initdata.yPencil(ib,0,thisIndex.z).recvAck(msg);
05905       CmiEnableUrgentSend(0);
05906       continue;
05907     }
05908     int nx = block1;
05909     if ( (ib+1)*block1 > K1 ) nx = K1 - ib*block1;
05910     PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
05911     if ( send_evir ) {
05912       send_evir = 0;
05913     }
05914     msg->sourceNode = thisIndex.y;
05915     msg->ny = ny;
05916     float *md = msg->qgrid;
05917     for ( int i=ib*block1; i<(ib*block1+nx); ++i ) {
05918      float *d = data + i*ny*nz*2;
05919      for ( int j=0; j<ny; ++j, d += nz*2 ) {
05920       for ( int k=0; k<nz; ++k ) {
05921         *(md++) = d[2*k];
05922         *(md++) = d[2*k+1];
05923       }
05924      }
05925     }
05926     SET_PRIORITY(msg,sequence,PME_UNTRANS_PRIORITY)
05927 
05928     CmiEnableUrgentSend(1);
05929 #if USE_NODE_PAR_RECEIVE
05930     msg->destElem=CkArrayIndex3D(ib,0, thisIndex.z);
05931 #if Y_PERSIST 
05932     CmiUsePersistentHandle(&untrans_handle[isend], 1);
05933 #endif
05934     initdata.pmeNodeProxy[CmiNodeOf(initdata.ym.ckLocalBranch()->procNum(0,msg->destElem))].recvYUntrans(msg);
05935 #if Y_PERSIST 
05936     CmiUsePersistentHandle(NULL, 0);
05937 #endif
05938 #else
05939 #if Y_PERSIST 
05940     CmiUsePersistentHandle(&untrans_handle[isend], 1);
05941 #endif
05942     initdata.yPencil(ib,0,thisIndex.z).recvUntrans(msg);
05943 #if Y_PERSIST 
05944     CmiUsePersistentHandle(NULL, 0);
05945 #endif
05946 #endif
05947     CmiEnableUrgentSend(0);
05948   }
05949 }


Member Data Documentation

fftw_plan PmeXPencil::backward_plan

Definition at line 4732 of file ComputePme.C.

Referenced by backward_fft(), and fft_init().

fftw_plan PmeXPencil::forward_plan

Definition at line 4732 of file ComputePme.C.

Referenced by fft_init(), and forward_fft().

PmeKSpace* PmeXPencil::myKSpace

Definition at line 4738 of file ComputePme.C.

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

int PmeXPencil::ny

Definition at line 4735 of file ComputePme.C.

int PmeXPencil::nz

Definition at line 4735 of file ComputePme.C.

int PmeXPencil::recipEvirPe

Definition at line 4736 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 Mon Jun 18 01:17:20 2018 for NAMD by  doxygen 1.4.7