PmeYPencil Class Reference

Inheritance diagram for PmeYPencil:

PmePencil< CBase_PmeYPencil > CBase_PmeYPencil List of all members.

Public Member Functions

PmeYPencil_SDAG_CODE PmeYPencil ()
 PmeYPencil (CkMigrateMessage *)
void fft_init ()
void recv_trans (const PmeTransMsg *)
void forward_fft ()
void forward_subset_fft (int fromIdx, int toIdx)
void send_trans ()
void send_subset_trans (int fromIdx, int toIdx)
void recv_untrans (const PmeUntransMsg *)
void node_process_trans (PmeTransMsg *)
void node_process_untrans (PmeUntransMsg *)
void backward_fft ()
void backward_subset_fft (int fromIdx, int toIdx)
void send_untrans ()
void send_subset_untrans (int fromIdx, int toIdx, int evirIdx)

Detailed Description

Definition at line 4679 of file ComputePme.C.


Constructor & Destructor Documentation

PmeYPencil_SDAG_CODE PmeYPencil::PmeYPencil (  )  [inline]

Definition at line 4682 of file ComputePme.C.

References PmePencil< CBase_PmeYPencil >::imsg, and PmePencil< CBase_PmeYPencil >::imsgb.

04682 { __sdag_init(); setMigratable(false); imsg=imsgb=0;}

PmeYPencil::PmeYPencil ( CkMigrateMessage *   )  [inline]

Definition at line 4683 of file ComputePme.C.

04683 { __sdag_init(); }


Member Function Documentation

void PmeYPencil::backward_fft (  ) 

Definition at line 6030 of file ComputePme.C.

References CKLOOP_CTRL_PME_BACKWARDFFT, PmePencil< CBase_PmeYPencil >::data, PmePencilInitMsgData::grid, PmePencil< CBase_PmeYPencil >::initdata, PmeGrid::K2, Node::Object(), PmeYPencilBackwardFFT(), Node::simParameters, SimParameters::useCkLoop, PmePencil< CBase_PmeYPencil >::work, PmePencilInitMsgData::xBlocks, and PmePencilInitMsgData::zBlocks.

Referenced by node_process_untrans().

06030                               {
06031 #ifdef NAMD_FFTW
06032 #ifdef MANUAL_DEBUG_FFTW3
06033   dumpMatrixFloat3("bw_y_b", data, nx, initdata.grid.K2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
06034 #endif
06035 
06036 #ifdef NAMD_FFTW_3
06037 #if     CMK_SMP && USE_CKLOOP
06038   int useCkLoop = Node::Object()->simParameters->useCkLoop;
06039   if(useCkLoop>=CKLOOP_CTRL_PME_BACKWARDFFT
06040      && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
06041           CkLoop_Parallelize(PmeYPencilBackwardFFT, 1, (void *)this, CkMyNodeSize(), 0, nx-1); //sync
06042           return;
06043   }
06044 #endif
06045   //the above is a transformation of the following loop using CkLoop
06046   for ( int i=0; i<nx; ++i ) {
06047 #if CMK_BLUEGENEL
06048         CmiNetworkProgress();
06049 #endif
06050     fftwf_execute_dft(backward_plan,    
06051                                           ((fftwf_complex *) data) + i * nz * initdata.grid.K2,
06052                                           ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
06053   }
06054 #else
06055         for ( int i=0; i<nx; ++i ) {
06056 #if CMK_BLUEGENEL
06057           CmiNetworkProgress();
06058 #endif
06059                 fftw(backward_plan, nz,
06060                 ((fftw_complex *) data) + i * nz * initdata.grid.K2,
06061                 nz, 1, (fftw_complex *) work, 1, 0);
06062         }
06063 #endif
06064 
06065 #ifdef MANUAL_DEBUG_FFTW3
06066   dumpMatrixFloat3("bw_y_a", data, nx, initdata.grid.K2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
06067 #endif
06068 
06069 #endif
06070 }

void PmeYPencil::backward_subset_fft ( int  fromIdx,
int  toIdx 
)

Definition at line 6018 of file ComputePme.C.

References PmePencil< CBase_PmeYPencil >::data, PmePencilInitMsgData::grid, PmePencil< CBase_PmeYPencil >::initdata, and PmeGrid::K2.

Referenced by PmeYPencilBackwardFFT().

06018                                                            {
06019 #ifdef NAMD_FFTW
06020 #ifdef NAMD_FFTW_3
06021         for(int i=fromIdx; i<=toIdx; i++){
06022                 fftwf_execute_dft(backward_plan,        
06023                                                   ((fftwf_complex *) data) + i * nz * initdata.grid.K2,         
06024                                                   ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
06025         }
06026 #endif
06027 #endif
06028 }

void PmeYPencil::fft_init (  ) 

Definition at line 4908 of file ComputePme.C.

References PmeGrid::block1, PmeGrid::block3, PmePencil< CBase_PmeYPencil >::data, PmeGrid::dim2, PmeGrid::dim3, PmePencil< CBase_PmeYPencil >::evir, fftwf_malloc, PmePencilInitMsgData::grid, PmePencil< CBase_PmeYPencil >::initdata, PmeGrid::K1, PmeGrid::K2, NAMD_die(), PmePencil< CBase_PmeYPencil >::order_init(), PmePencilInitMsgData::pmeNodeProxy, Node::simParameters, simParams, PmePencil< CBase_PmeYPencil >::work, and PmePencilInitMsgData::yBlocks.

04908                           {
04909   CProxy_Node nd(CkpvAccess(BOCclass_group).node);
04910   Node *node = nd.ckLocalBranch();
04911   SimParameters *simParams = node->simParameters;
04912 
04913 #if USE_NODE_PAR_RECEIVE
04914   ((NodePmeMgr *)CkLocalNodeBranch(initdata.pmeNodeProxy))->registerYPencil(thisIndex,this);
04915 #endif
04916 
04917   int K1 = initdata.grid.K1;
04918   int K2 = initdata.grid.K2;
04919   int dim2 = initdata.grid.dim2;
04920   int dim3 = initdata.grid.dim3;
04921   int block1 = initdata.grid.block1;
04922   int block3 = initdata.grid.block3;
04923 
04924   nx = block1;
04925   if ( (thisIndex.x + 1) * block1 > K1 ) nx = K1 - thisIndex.x * block1;
04926   nz = block3;
04927   if ( (thisIndex.z+1)*block3 > dim3/2 ) nz = dim3/2 - thisIndex.z*block3;
04928 
04929 #ifdef NAMD_FFTW
04930   CmiLock(ComputePmeMgr::fftw_plan_lock);
04931 
04932   data = (float *) fftwf_malloc( sizeof(float) * nx*dim2*nz*2);
04933   work = new float[2*K2];
04934 
04935   order_init(initdata.yBlocks);
04936 
04937 #ifdef NAMD_FFTW_3
04938   /* need array of sizes for the dimensions */
04939   /* ideally this should be implementable as a single multidimensional
04940    *  plan, but that has proven tricky to implement, so we maintain the
04941    *  loop of 1d plan executions. */
04942   int sizeLines=nz;
04943   int planLineSizes[1];
04944   planLineSizes[0]=K2;
04945   int fftwFlags = simParams->FFTWPatient ? FFTW_PATIENT  : simParams->FFTWEstimate ? FFTW_ESTIMATE  : FFTW_MEASURE ;
04946   forward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines, 
04947                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04948                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04949                                      FFTW_FORWARD, 
04950                                      fftwFlags);
04951   backward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines, 
04952                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04953                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04954                                      FFTW_BACKWARD, 
04955                                       fftwFlags);
04956   CkAssert(forward_plan != NULL);
04957   CkAssert(backward_plan != NULL);
04958 #else
04959   forward_plan = fftw_create_plan_specific(K2, FFTW_FORWARD,
04960         ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04961         | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
04962         nz, (fftw_complex *) work, 1);
04963   backward_plan = fftw_create_plan_specific(K2, FFTW_BACKWARD,
04964         ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04965         | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
04966         nz, (fftw_complex *) work, 1);
04967 #endif
04968   CmiUnlock(ComputePmeMgr::fftw_plan_lock);
04969 #else
04970   NAMD_die("Sorry, FFTW must be compiled in to use PME.");
04971 #endif
04972 
04973 #if USE_NODE_PAR_RECEIVE
04974   evir = 0;
04975   CmiMemoryWriteFence();
04976 #endif
04977 }

void PmeYPencil::forward_fft (  ) 

Definition at line 5458 of file ComputePme.C.

References CKLOOP_CTRL_PME_FORWARDFFT, PmePencil< CBase_PmeYPencil >::data, PmeGrid::dim2, PmePencil< CBase_PmeYPencil >::evir, PmePencilInitMsgData::grid, PmePencil< CBase_PmeYPencil >::initdata, PmeGrid::K2, Node::Object(), PmeYPencilForwardFFT(), Node::simParameters, SimParameters::useCkLoop, PmePencil< CBase_PmeYPencil >::work, PmePencilInitMsgData::xBlocks, and PmePencilInitMsgData::zBlocks.

Referenced by node_process_trans().

05458                              {
05459     evir = 0.;
05460 #ifdef NAMD_FFTW
05461 #ifdef MANUAL_DEBUG_FFTW3
05462   dumpMatrixFloat3("fw_y_b", data, nx, initdata.grid.K2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05463 #endif
05464   
05465 #ifdef NAMD_FFTW_3
05466 #if     CMK_SMP && USE_CKLOOP
05467   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05468   if(useCkLoop>=CKLOOP_CTRL_PME_FORWARDFFT
05469      && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
05470           CkLoop_Parallelize(PmeYPencilForwardFFT, 1, (void *)this, CkMyNodeSize(), 0, nx-1); //sync
05471           return;
05472   }
05473 #endif
05474   //the above is a transformation of the following loop using CkLoop
05475   for ( int i=0; i<nx; ++i ) {
05476     fftwf_execute_dft(forward_plan, ((fftwf_complex *) data) + i 
05477                       * nz * initdata.grid.K2,  
05478                       ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
05479   }
05480 #else
05481   for ( int i=0; i<nx; ++i ) {
05482     fftw(forward_plan, nz,
05483         ((fftw_complex *) data) + i * nz * initdata.grid.K2,
05484         nz, 1, (fftw_complex *) work, 1, 0);
05485   }
05486 #endif
05487 #ifdef MANUAL_DEBUG_FFTW3
05488   dumpMatrixFloat3("fw_y_a", data, nx, initdata.grid.dim2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05489 #endif
05490 
05491 #endif
05492 }

void PmeYPencil::forward_subset_fft ( int  fromIdx,
int  toIdx 
)

Definition at line 5446 of file ComputePme.C.

References PmePencil< CBase_PmeYPencil >::data, PmePencilInitMsgData::grid, PmePencil< CBase_PmeYPencil >::initdata, and PmeGrid::K2.

Referenced by PmeYPencilForwardFFT().

05446                                                           {
05447 #ifdef NAMD_FFTW
05448 #ifdef NAMD_FFTW_3
05449         for(int i=fromIdx; i<=toIdx; i++){
05450                 fftwf_execute_dft(forward_plan, ((fftwf_complex *) data) + i 
05451                       * nz * initdata.grid.K2,  
05452                       ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
05453         }
05454 #endif
05455 #endif
05456 }

void PmeYPencil::node_process_trans ( PmeTransMsg  ) 

Definition at line 4979 of file ComputePme.C.

References forward_fft(), PmePencil< CBase_PmeYPencil >::hasData, PmeTransMsg::hasData, PmePencil< CBase_PmeYPencil >::imsg, PmePencil< CBase_PmeYPencil >::initdata, PmePencil< CBase_PmeYPencil >::needs_reply, recv_trans(), send_trans(), send_untrans(), PmeTransMsg::sourceNode, and PmePencilInitMsgData::yBlocks.

Referenced by NodePmeMgr::recvYTrans().

04980 {
04981   if ( msg->hasData ) hasData = 1;
04982   needs_reply[msg->sourceNode] = msg->hasData;
04983   recv_trans(msg);
04984   int limsg;
04985   CmiMemoryAtomicFetchAndInc(imsg,limsg);
04986   if(limsg+1 == initdata.yBlocks)
04987     {
04988       if ( hasData ) {
04989         forward_fft();
04990       }
04991       send_trans();
04992       if( ! hasData)
04993         {
04994           send_untrans(); //todo, what is up with the recvAck in SDAG version?
04995         }
04996       imsg=0;
04997       CmiMemoryWriteFence();
04998     }
04999 }

void PmeYPencil::node_process_untrans ( PmeUntransMsg  ) 

Definition at line 5001 of file ComputePme.C.

References backward_fft(), PmePencil< CBase_PmeYPencil >::imsgb, PmePencil< CBase_PmeYPencil >::initdata, recv_untrans(), send_untrans(), and PmePencilInitMsgData::yBlocks.

Referenced by NodePmeMgr::recvYUntrans().

05002 {
05003   recv_untrans(msg);
05004   int limsg;
05005   CmiMemoryAtomicFetchAndInc(imsgb,limsg);
05006   if(limsg+1 == initdata.yBlocks)
05007     {
05008       backward_fft();
05009       send_untrans();
05010       imsgb=0;
05011       CmiMemoryWriteFence();
05012     }
05013 }

void PmeYPencil::recv_trans ( const PmeTransMsg  ) 

Definition at line 5405 of file ComputePme.C.

References PmeGrid::block2, PmePencil< CBase_PmeYPencil >::data, PmePencilInitMsgData::grid, PmeTransMsg::hasData, PmePencil< CBase_PmeYPencil >::imsg, PmePencil< CBase_PmeYPencil >::initdata, j, PmeGrid::K2, PmeTransMsg::lattice, PmePencil< CBase_PmeYPencil >::lattice, PmeTransMsg::nx, PmeTransMsg::qgrid, PmeTransMsg::sequence, PmePencil< CBase_PmeYPencil >::sequence, and PmeTransMsg::sourceNode.

Referenced by node_process_trans().

05405                                                   {
05406   if ( imsg == 0 ) {
05407     lattice = msg->lattice;
05408     sequence = msg->sequence;
05409   }
05410   int block2 = initdata.grid.block2;
05411   int K2 = initdata.grid.K2;
05412   int jb = msg->sourceNode;
05413   int ny = msg->nx;
05414  if ( msg->hasData ) {
05415   const float *md = msg->qgrid;
05416   float *d = data;
05417   for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05418    for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05419     for ( int k=0; k<nz; ++k ) {
05420 #ifdef ZEROCHECK
05421       if ( (*md) == 0. ) CkPrintf("0 in ZY at %d %d %d %d %d %d %d %d %d\n",
05422         thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05423 #endif
05424       d[2*(j*nz+k)] = *(md++);
05425       d[2*(j*nz+k)+1] = *(md++);
05426     }
05427    }
05428   }
05429  } else {
05430   float *d = data;
05431   for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05432    for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05433     for ( int k=0; k<nz; ++k ) {
05434       d[2*(j*nz+k)] = 0;
05435       d[2*(j*nz+k)+1] = 0;
05436     }
05437    }
05438   }
05439  }
05440 }

void PmeYPencil::recv_untrans ( const PmeUntransMsg  ) 

Definition at line 5989 of file ComputePme.C.

References PmeGrid::block2, PmePencil< CBase_PmeYPencil >::data, PmePencilInitMsgData::grid, PmePencil< CBase_PmeYPencil >::initdata, j, PmeGrid::K2, PmeUntransMsg::ny, PmeUntransMsg::qgrid, and PmeUntransMsg::sourceNode.

Referenced by node_process_untrans().

05989                                                       {
05990   int block2 = initdata.grid.block2;
05991   int K2 = initdata.grid.K2;
05992   int jb = msg->sourceNode;
05993   int ny = msg->ny;
05994   const float *md = msg->qgrid;
05995   float *d = data;
05996   for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05997 #if CMK_BLUEGENEL
05998     CmiNetworkProgress();
05999 #endif   
06000     for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
06001       for ( int k=0; k<nz; ++k ) {
06002 #ifdef ZEROCHECK
06003         if ( (*md) == 0. ) CkPrintf("0 in XY at %d %d %d %d %d %d %d %d %d\n",
06004                                     thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
06005 #endif
06006         d[2*(j*nz+k)] = *(md++);
06007         d[2*(j*nz+k)+1] = *(md++);
06008       }
06009     }
06010   }
06011 }

void PmeYPencil::send_subset_trans ( int  fromIdx,
int  toIdx 
)

Definition at line 5499 of file ComputePme.C.

References PmeGrid::block2, PmePencil< CBase_PmeYPencil >::data, PmeTransMsg::destElem, PmePencilInitMsgData::grid, PmeTransMsg::hasData, PmePencil< CBase_PmeYPencil >::hasData, PmePencil< CBase_PmeYPencil >::initdata, j, PmeGrid::K2, PmePencil< CBase_PmeYPencil >::lattice, PmeTransMsg::lattice, PmeTransMsg::nx, PME_TRANS2_PRIORITY, PmePencilInitMsgData::pmeNodeProxy, PRIORITY_SIZE, PmeTransMsg::qgrid, PmePencil< CBase_PmeYPencil >::send_order, PmePencil< CBase_PmeYPencil >::sequence, PmeTransMsg::sequence, SET_PRIORITY, PmeTransMsg::sourceNode, PmePencilInitMsgData::xm, PmePencilInitMsgData::xPencil, and PmePencilInitMsgData::yBlocks.

Referenced by PmeYPencilSendTrans().

05499                                                         {
05500         int yBlocks = initdata.yBlocks;
05501         int block2 = initdata.grid.block2;
05502         int K2 = initdata.grid.K2;
05503     for ( int isend=fromIdx; isend<=toIdx; ++isend ) {
05504           int jb = send_order[isend];
05505           int ny = block2;
05506           if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
05507           int hd = ( hasData ? 1 : 0 );
05508           PmeTransMsg *msg = new (hd*nx*ny*nz*2,PRIORITY_SIZE) PmeTransMsg;
05509           msg->lattice = lattice;
05510           msg->sourceNode = thisIndex.x;
05511           msg->hasData = hasData;
05512           msg->nx = nx;
05513          if ( hasData ) {
05514           float *md = msg->qgrid;
05515           const float *d = data;
05516           for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05517            for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05518                 for ( int k=0; k<nz; ++k ) {
05519                   *(md++) = d[2*(j*nz+k)];
05520                   *(md++) = d[2*(j*nz+k)+1];
05521   #ifdef ZEROCHECK
05522                   if ( *(md-2) == 0. ) CkPrintf("send 0 in YX at %d %d %d %d %d %d %d %d %d\n",
05523           thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05524   #endif
05525                 }
05526            }
05527           }
05528           if ( md != msg->qgrid + nx*ny*nz*2 ) CkPrintf("error in YX at %d %d %d\n",
05529           thisIndex.x, jb, thisIndex.z);
05530          }
05531           msg->sequence = sequence;
05532           SET_PRIORITY(msg,sequence,PME_TRANS2_PRIORITY)
05533       CmiEnableUrgentSend(1);
05534 #if USE_NODE_PAR_RECEIVE
05535       msg->destElem=CkArrayIndex3D(0,jb,thisIndex.z);
05536 #if X_PERSIST 
05537       CmiUsePersistentHandle(&trans_handle[isend], 1);
05538 #endif
05539       initdata.pmeNodeProxy[CmiNodeOf(initdata.xm.ckLocalBranch()->procNum(0,msg->destElem))].recvXTrans(msg);   
05540 #if X_PERSIST 
05541       CmiUsePersistentHandle(NULL, 0);
05542 #endif
05543 #else      
05544 #if X_PERSIST 
05545       CmiUsePersistentHandle(&trans_handle[isend], 1);
05546 #endif
05547       initdata.xPencil(0,jb,thisIndex.z).recvTrans(msg);
05548 #if X_PERSIST 
05549       CmiUsePersistentHandle(NULL, 0);
05550 #endif
05551 #endif
05552       CmiEnableUrgentSend(0);
05553         }
05554 }

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

Definition at line 6078 of file ComputePme.C.

References PmeGrid::block2, PmePencil< CBase_PmeYPencil >::data, PmeUntransMsg::destElem, PmePencilInitMsgData::grid, PmePencil< CBase_PmeYPencil >::initdata, j, PmeGrid::K2, PmeUntransMsg::ny, PME_UNTRANS2_PRIORITY, PmePencilInitMsgData::pmeNodeProxy, PRIORITY_SIZE, PmeUntransMsg::qgrid, PmePencil< CBase_PmeYPencil >::send_order, PmePencil< CBase_PmeYPencil >::sequence, SET_PRIORITY, PmeUntransMsg::sourceNode, PmePencilInitMsgData::yBlocks, PmePencilInitMsgData::zm, and PmePencilInitMsgData::zPencil.

Referenced by PmeYPencilSendUntrans().

06078                                                                        {
06079         int yBlocks = initdata.yBlocks;
06080         int block2 = initdata.grid.block2;      
06081         int K2 = initdata.grid.K2;
06082 
06083         int ackL=0, ackH=-1;
06084         int unL=0, unH=-1;
06085         int send_evir=0;
06086         if(fromIdx >= evirIdx+1) {
06087                 //send PmeUntransMsg with has_evir=0
06088                 unL = fromIdx;
06089                 unH = toIdx;            
06090         } else if(toIdx <= evirIdx-1) {
06091                 //send PmeAckMsg
06092                 ackL=fromIdx;
06093                 ackH=toIdx;             
06094         } else {
06095                 //partially send PmeAckMsg and partially send PmeUntransMsg
06096                 ackL=fromIdx;
06097                 ackH=evirIdx-1;
06098                 send_evir=1;
06099                 unL=evirIdx+1;
06100                 unH=toIdx;
06101         }
06102 
06103         for(int isend=ackL; isend<=ackH; isend++) {
06104                 //send PmeAckMsg
06105         CmiEnableUrgentSend(1);
06106                 int jb = send_order[isend];
06107                 PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
06108                 SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
06109                 initdata.zPencil(thisIndex.x,jb,0).recvAck(msg);
06110         CmiEnableUrgentSend(0);
06111         }
06112 
06113     CmiEnableUrgentSend(1);
06114         //send PmeUntransMsg with has_evir=1
06115         if(send_evir) {
06116                 int jb = send_order[evirIdx];
06117                 int ny = block2;
06118                 if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
06119                 PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;              
06120                 msg->sourceNode = thisIndex.z;
06121                 msg->ny = nz;
06122                 float *md = msg->qgrid;
06123                 const float *d = data;
06124                 for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
06125                         for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
06126                                 for ( int k=0; k<nz; ++k ) {
06127                                         *(md++) = d[2*(j*nz+k)];
06128                                         *(md++) = d[2*(j*nz+k)+1];
06129                                 }
06130                         }
06131                 }
06132                 SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
06133 #if USE_NODE_PAR_RECEIVE
06134         msg->destElem=CkArrayIndex3D( thisIndex.x, jb, 0);
06135     //    CkPrintf("[%d] sending to %d %d %d recvZUntrans on node %d\n", CkMyPe(), thisIndex.x, jb, 0, CmiNodeOf(initdata.zm.ckLocalBranch()->procNum(0,msg->destElem)));
06136         initdata.pmeNodeProxy[CmiNodeOf(initdata.zm.ckLocalBranch()->procNum(0,msg->destElem))].recvZUntrans(msg);
06137 #else
06138         initdata.zPencil(thisIndex.x,jb,0).recvUntrans(msg);
06139 #endif
06140         }
06141 
06142     CmiEnableUrgentSend(0);
06143         //send PmeUntransMsg with has_evir=0
06144         for(int isend=unL; isend<=unH; isend++) {
06145                 int jb = send_order[isend];
06146                 int ny = block2;
06147                 if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
06148                 PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
06149                 msg->sourceNode = thisIndex.z;
06150                 msg->ny = nz;
06151                 float *md = msg->qgrid;
06152                 const float *d = data;
06153                 for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
06154                         for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
06155                                 for ( int k=0; k<nz; ++k ) {
06156                                         *(md++) = d[2*(j*nz+k)];
06157                                         *(md++) = d[2*(j*nz+k)+1];
06158                                 }
06159                         }
06160                 }
06161                 SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
06162             CmiEnableUrgentSend(1);
06163 #if USE_NODE_PAR_RECEIVE
06164         msg->destElem=CkArrayIndex3D( thisIndex.x, jb, 0);
06165         //    CkPrintf("[%d] sending to %d %d %d recvZUntrans on node %d\n", CkMyPe(), thisIndex.x, jb, 0, CmiNodeOf(initdata.zm.ckLocalBranch()->procNum(0,msg->destElem)));
06166 #if Z_PERSIST 
06167         CmiUsePersistentHandle(&untrans_handle[isend], 1);
06168 #endif
06169         initdata.pmeNodeProxy[CmiNodeOf(initdata.zm.ckLocalBranch()->procNum(0,msg->destElem))].recvZUntrans(msg);
06170 #if Z_PERSIST 
06171         CmiUsePersistentHandle(NULL, 0);
06172 #endif
06173 #else
06174 #if Z_PERSIST 
06175         CmiUsePersistentHandle(&untrans_handle[isend], 1);
06176 #endif
06177         initdata.zPencil(thisIndex.x,jb,0).recvUntrans(msg);
06178 #if Z_PERSIST 
06179         CmiUsePersistentHandle(NULL, 0);
06180 #endif
06181 #endif
06182     CmiEnableUrgentSend(0);
06183         }
06184 }

void PmeYPencil::send_trans (  ) 

Definition at line 5556 of file ComputePme.C.

References PmeGrid::block2, CKLOOP_CTRL_PME_SENDTRANS, PmePencil< CBase_PmeYPencil >::data, PmeTransMsg::destElem, PmePencilInitMsgData::grid, PmeTransMsg::hasData, PmePencil< CBase_PmeYPencil >::hasData, PmePencil< CBase_PmeYPencil >::initdata, j, PmeGrid::K2, PmePencil< CBase_PmeYPencil >::lattice, PmeTransMsg::lattice, PmeTransMsg::nx, Node::Object(), PME_TRANS2_PRIORITY, PmePencilInitMsgData::pmeNodeProxy, PmeYPencilSendTrans(), PRIORITY_SIZE, PmeTransMsg::qgrid, PmePencil< CBase_PmeYPencil >::send_order, PmePencil< CBase_PmeYPencil >::sequence, PmeTransMsg::sequence, SET_PRIORITY, Node::simParameters, PmeTransMsg::sourceNode, SimParameters::useCkLoop, PmePencilInitMsgData::xBlocks, PmePencilInitMsgData::xm, PmePencilInitMsgData::xPencil, PmePencilInitMsgData::yBlocks, and PmePencilInitMsgData::zBlocks.

Referenced by node_process_trans().

05556                             {
05557 #if USE_PERSISTENT
05558     if (trans_handle == NULL) setup_persistent();
05559 #endif
05560 #if     CMK_SMP && USE_CKLOOP
05561         int useCkLoop = Node::Object()->simParameters->useCkLoop;
05562         if(useCkLoop>=CKLOOP_CTRL_PME_SENDTRANS
05563            && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
05570                 //send_subset_trans(0, initdata.yBlocks-1);
05571                 CkLoop_Parallelize(PmeYPencilSendTrans, 1, (void *)this, CkMyNodeSize(), 0, initdata.yBlocks-1, 1); //not sync
05572                 return;
05573         }
05574 #endif
05575   int yBlocks = initdata.yBlocks;
05576   int block2 = initdata.grid.block2;
05577   int K2 = initdata.grid.K2;
05578   for ( int isend=0; isend<yBlocks; ++isend ) {
05579     int jb = send_order[isend];
05580     int ny = block2;
05581     if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
05582     int hd = ( hasData ? 1 : 0 );
05583     PmeTransMsg *msg = new (hd*nx*ny*nz*2,PRIORITY_SIZE) PmeTransMsg;
05584     msg->lattice = lattice;
05585     msg->sourceNode = thisIndex.x;
05586     msg->hasData = hasData;
05587     msg->nx = nx;
05588    if ( hasData ) {
05589     float *md = msg->qgrid;
05590     const float *d = data;
05591     for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05592      for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05593       for ( int k=0; k<nz; ++k ) {
05594         *(md++) = d[2*(j*nz+k)];
05595         *(md++) = d[2*(j*nz+k)+1];
05596 #ifdef ZEROCHECK
05597         if ( *(md-2) == 0. ) CkPrintf("send 0 in YX at %d %d %d %d %d %d %d %d %d\n",
05598         thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05599 #endif
05600       }
05601      }
05602     }
05603     if ( md != msg->qgrid + nx*ny*nz*2 ) CkPrintf("error in YX at %d %d %d\n",
05604         thisIndex.x, jb, thisIndex.z);
05605    }
05606     msg->sequence = sequence;
05607     SET_PRIORITY(msg,sequence,PME_TRANS2_PRIORITY)
05608     CmiEnableUrgentSend(1);
05609 #if USE_NODE_PAR_RECEIVE
05610     msg->destElem=CkArrayIndex3D(0,jb,thisIndex.z);
05611 #if X_PERSIST 
05612         CmiUsePersistentHandle(&trans_handle[isend], 1);
05613 #endif
05614     initdata.pmeNodeProxy[CmiNodeOf(initdata.xm.ckLocalBranch()->procNum(0,msg->destElem))].recvXTrans(msg);   
05615 #if X_PERSIST 
05616         CmiUsePersistentHandle(NULL, 0);
05617 #endif
05618 #else
05619 #if X_PERSIST 
05620         CmiUsePersistentHandle(&trans_handle[isend], 1);
05621 #endif
05622     initdata.xPencil(0,jb,thisIndex.z).recvTrans(msg);
05623 #if X_PERSIST 
05624         CmiUsePersistentHandle(NULL, 0);
05625 #endif
05626     
05627 #endif
05628     CmiEnableUrgentSend(0);
05629   }
05630 }

void PmeYPencil::send_untrans (  ) 

Definition at line 6186 of file ComputePme.C.

References PmeGrid::block2, CKLOOP_CTRL_PME_SENDUNTRANS, PmePencil< CBase_PmeYPencil >::data, PmeUntransMsg::destElem, PmePencil< CBase_PmeYPencil >::evir, PmePencilInitMsgData::grid, PmePencil< CBase_PmeYPencil >::initdata, j, PmeGrid::K2, PmePencil< CBase_PmeYPencil >::needs_reply, PmeUntransMsg::ny, Node::Object(), PME_UNTRANS2_PRIORITY, PmePencilInitMsgData::pmeNodeProxy, PmeYPencilSendUntrans(), PRIORITY_SIZE, PmeUntransMsg::qgrid, PmePencil< CBase_PmeYPencil >::send_order, PmePencil< CBase_PmeYPencil >::sequence, SET_PRIORITY, Node::simParameters, PmeUntransMsg::sourceNode, SimParameters::useCkLoop, PmePencilInitMsgData::xBlocks, PmePencilInitMsgData::yBlocks, PmePencilInitMsgData::zBlocks, PmePencilInitMsgData::zm, and PmePencilInitMsgData::zPencil.

Referenced by node_process_trans(), and node_process_untrans().

06186                               {
06187 #if USE_PERSISTENT
06188   if (untrans_handle == NULL) setup_persistent();
06189 #endif
06190 #if     CMK_SMP && USE_CKLOOP
06191   int useCkLoop = Node::Object()->simParameters->useCkLoop;
06192   if(useCkLoop>=CKLOOP_CTRL_PME_SENDUNTRANS
06193      && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
06194           int yBlocks = initdata.yBlocks;
06195           int evirIdx = 0;
06196           for ( int isend=0; isend<yBlocks; ++isend ) {
06197                   int jb = send_order[isend];
06198                   if (needs_reply[jb]) {
06199                           evirIdx = isend;
06200                           break;
06201                   }
06202           }
06203 
06204           //basically: 
06205           //[0,evirIdx-1]->send PmeAckMsg
06206           //evirIdx->send PmeUntransMsg with has_evir=1
06207           //[evirIdx+1, yBlocks-1]->send PmeUntransMsg with has_evir=0
06208           //send_subset_untrans(0, yBlocks-1, evirIdx);
06209 #if USE_NODE_PAR_RECEIVE      
06210           //CkLoop_Parallelize(PmeYPencilSendUntrans, evirIdx, (void *)this, CkMyNodeSize(), 0, yBlocks-1, 1); //sync
06211           CkLoop_Parallelize(PmeYPencilSendUntrans, evirIdx, (void *)this, yBlocks, 0, yBlocks-1, 1);
06212       evir = 0.;
06213       CmiMemoryWriteFence();
06214 #else
06215       //CkLoop_Parallelize(PmeYPencilSendUntrans, evirIdx, (void *)this, CkMyNodeSize(), 0, yBlocks-1, 0); //not sync
06216           CkLoop_Parallelize(PmeYPencilSendUntrans, evirIdx, (void *)this, yBlocks, 0, yBlocks-1, 0); //not sync
06217 #endif
06218           return;
06219   }
06220 #endif
06221   int yBlocks = initdata.yBlocks;
06222   int block2 = initdata.grid.block2;
06223   int K2 = initdata.grid.K2;
06224   int send_evir = 1;
06225   for ( int isend=0; isend<yBlocks; ++isend ) {
06226     int jb = send_order[isend];
06227     if ( ! needs_reply[jb] ) {
06228       PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
06229       CmiEnableUrgentSend(1);
06230       SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
06231       initdata.zPencil(thisIndex.x,jb,0).recvAck(msg);
06232       CmiEnableUrgentSend(0);
06233       continue;
06234     }
06235     int ny = block2;
06236     if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
06237     PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
06238     if ( send_evir ) {
06239       send_evir = 0;
06240     }
06241     msg->sourceNode = thisIndex.z;
06242     msg->ny = nz;
06243     float *md = msg->qgrid;
06244     const float *d = data;
06245     for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
06246      for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
06247       for ( int k=0; k<nz; ++k ) {
06248         *(md++) = d[2*(j*nz+k)];
06249         *(md++) = d[2*(j*nz+k)+1];
06250       }
06251      }
06252     }
06253     SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
06254 
06255     CmiEnableUrgentSend(1);
06256 #if USE_NODE_PAR_RECEIVE
06257     msg->destElem=CkArrayIndex3D( thisIndex.x, jb, 0);
06258     //    CkPrintf("[%d] sending to %d %d %d recvZUntrans on node %d\n", CkMyPe(), thisIndex.x, jb, 0, CmiNodeOf(initdata.zm.ckLocalBranch()->procNum(0,msg->destElem)));
06259 #if Z_PERSIST 
06260     CmiUsePersistentHandle(&untrans_handle[isend], 1);
06261 #endif
06262     initdata.pmeNodeProxy[CmiNodeOf(initdata.zm.ckLocalBranch()->procNum(0,msg->destElem))].recvZUntrans(msg);
06263 #if Z_PERSIST
06264     CmiUsePersistentHandle(NULL, 0);
06265 #endif
06266 #else
06267 #if Z_PERSIST 
06268     CmiUsePersistentHandle(&untrans_handle[isend], 1);
06269 #endif
06270     initdata.zPencil(thisIndex.x,jb,0).recvUntrans(msg);
06271 #if Z_PERSIST 
06272     CmiUsePersistentHandle(NULL, 0);
06273 #endif
06274 #endif    
06275     CmiEnableUrgentSend(0);
06276   }
06277   
06278 #if USE_NODE_PAR_RECEIVE
06279   evir = 0.;
06280   CmiMemoryWriteFence();
06281 #endif
06282 }


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