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


Constructor & Destructor Documentation

PmeYPencil_SDAG_CODE PmeYPencil::PmeYPencil (  )  [inline]

Definition at line 4675 of file ComputePme.C.

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

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

PmeYPencil::PmeYPencil ( CkMigrateMessage *   )  [inline]

Definition at line 4676 of file ComputePme.C.

04676 { __sdag_init(); }


Member Function Documentation

void PmeYPencil::backward_fft (  ) 

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

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

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

Definition at line 6011 of file ComputePme.C.

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

Referenced by PmeYPencilBackwardFFT().

06011                                                            {
06012 #ifdef NAMD_FFTW
06013 #ifdef NAMD_FFTW_3
06014         for(int i=fromIdx; i<=toIdx; i++){
06015                 fftwf_execute_dft(backward_plan,        
06016                                                   ((fftwf_complex *) data) + i * nz * initdata.grid.K2,         
06017                                                   ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
06018         }
06019 #endif
06020 #endif
06021 }

void PmeYPencil::fft_init (  ) 

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

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

void PmeYPencil::forward_fft (  ) 

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

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

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

Definition at line 5439 of file ComputePme.C.

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

Referenced by PmeYPencilForwardFFT().

05439                                                           {
05440 #ifdef NAMD_FFTW
05441 #ifdef NAMD_FFTW_3
05442         for(int i=fromIdx; i<=toIdx; i++){
05443                 fftwf_execute_dft(forward_plan, ((fftwf_complex *) data) + i 
05444                       * nz * initdata.grid.K2,  
05445                       ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
05446         }
05447 #endif
05448 #endif
05449 }

void PmeYPencil::node_process_trans ( PmeTransMsg  ) 

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

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

void PmeYPencil::node_process_untrans ( PmeUntransMsg  ) 

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

04995 {
04996   recv_untrans(msg);
04997   int limsg;
04998   CmiMemoryAtomicFetchAndInc(imsgb,limsg);
04999   if(limsg+1 == initdata.yBlocks)
05000     {
05001       backward_fft();
05002       send_untrans();
05003       imsgb=0;
05004       CmiMemoryWriteFence();
05005     }
05006 }

void PmeYPencil::recv_trans ( const PmeTransMsg  ) 

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

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

void PmeYPencil::recv_untrans ( const PmeUntransMsg  ) 

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

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

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

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

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

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

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

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

void PmeYPencil::send_trans (  ) 

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

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

void PmeYPencil::send_untrans (  ) 

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

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


The documentation for this class was generated from the following file:
Generated on Sat Sep 23 01:17:21 2017 for NAMD by  doxygen 1.4.7