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 recvNodeAck (PmeAckMsg *)
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)

Detailed Description

Definition at line 4609 of file ComputePme.C.


Constructor & Destructor Documentation

PmeYPencil_SDAG_CODE PmeYPencil::PmeYPencil (  )  [inline]

Definition at line 4612 of file ComputePme.C.

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

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

PmeYPencil::PmeYPencil ( CkMigrateMessage *   )  [inline]

Definition at line 4613 of file ComputePme.C.

04613 { __sdag_init(); }


Member Function Documentation

void PmeYPencil::backward_fft (  ) 

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

05906                               {
05907 #ifdef NAMD_FFTW
05908 #ifdef MANUAL_DEBUG_FFTW3
05909   dumpMatrixFloat3("bw_y_b", data, nx, initdata.grid.K2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05910 #endif
05911 
05912 #ifdef NAMD_FFTW_3
05913 #if     CMK_SMP && USE_CKLOOP
05914   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05915   if(useCkLoop>=CKLOOP_CTRL_PME_BACKWARDFFT
05916      && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
05917           CkLoop_Parallelize(PmeYPencilBackwardFFT, 1, (void *)this, CkMyNodeSize(), 0, nx-1); //sync
05918           return;
05919   }
05920 #endif
05921   //the above is a transformation of the following loop using CkLoop
05922   for ( int i=0; i<nx; ++i ) {
05923 #if CMK_BLUEGENEL
05924         CmiNetworkProgress();
05925 #endif
05926     fftwf_execute_dft(backward_plan,    
05927                                           ((fftwf_complex *) data) + i * nz * initdata.grid.K2,
05928                                           ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
05929   }
05930 #else
05931         for ( int i=0; i<nx; ++i ) {
05932 #if CMK_BLUEGENEL
05933           CmiNetworkProgress();
05934 #endif
05935                 fftw(backward_plan, nz,
05936                 ((fftw_complex *) data) + i * nz * initdata.grid.K2,
05937                 nz, 1, (fftw_complex *) work, 1, 0);
05938         }
05939 #endif
05940 
05941 #ifdef MANUAL_DEBUG_FFTW3
05942   dumpMatrixFloat3("bw_y_a", data, nx, initdata.grid.K2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05943 #endif
05944 
05945 #endif
05946 }

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

Definition at line 5894 of file ComputePme.C.

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

Referenced by PmeYPencilBackwardFFT().

05894                                                            {
05895 #ifdef NAMD_FFTW
05896 #ifdef NAMD_FFTW_3
05897         for(int i=fromIdx; i<=toIdx; i++){
05898                 fftwf_execute_dft(backward_plan,        
05899                                                   ((fftwf_complex *) data) + i * nz * initdata.grid.K2,         
05900                                                   ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
05901         }
05902 #endif
05903 #endif
05904 }

void PmeYPencil::fft_init (  ) 

Definition at line 4839 of file ComputePme.C.

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

04839                           {
04840   CProxy_Node nd(CkpvAccess(BOCclass_group).node);
04841   Node *node = nd.ckLocalBranch();
04842   SimParameters *simParams = node->simParameters;
04843 
04844 #if USE_NODE_PAR_RECEIVE
04845   ((NodePmeMgr *)CkLocalNodeBranch(initdata.pmeNodeProxy))->registerYPencil(thisIndex,this);
04846 #endif
04847 
04848   int K1 = initdata.grid.K1;
04849   int K2 = initdata.grid.K2;
04850   int dim2 = initdata.grid.dim2;
04851   int dim3 = initdata.grid.dim3;
04852   int block1 = initdata.grid.block1;
04853   int block3 = initdata.grid.block3;
04854 
04855   nx = block1;
04856   if ( (thisIndex.x + 1) * block1 > K1 ) nx = K1 - thisIndex.x * block1;
04857   nz = block3;
04858   if ( (thisIndex.z+1)*block3 > dim3/2 ) nz = dim3/2 - thisIndex.z*block3;
04859 
04860 #ifdef NAMD_FFTW
04861   CmiLock(ComputePmeMgr::fftw_plan_lock);
04862 
04863   data = (float *) fftwf_malloc( sizeof(float) * nx*dim2*nz*2);
04864   work = new float[2*K2];
04865 
04866   order_init(initdata.yBlocks);
04867 
04868 #ifdef NAMD_FFTW_3
04869   /* need array of sizes for the dimensions */
04870   /* ideally this should be implementable as a single multidimensional
04871    *  plan, but that has proven tricky to implement, so we maintain the
04872    *  loop of 1d plan executions. */
04873   int sizeLines=nz;
04874   int planLineSizes[1];
04875   planLineSizes[0]=K2;
04876   int fftwFlags = simParams->FFTWPatient ? FFTW_PATIENT  : simParams->FFTWEstimate ? FFTW_ESTIMATE  : FFTW_MEASURE ;
04877   forward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines, 
04878                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04879                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04880                                      FFTW_FORWARD, 
04881                                      fftwFlags);
04882   backward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines, 
04883                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04884                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04885                                      FFTW_BACKWARD, 
04886                                       fftwFlags);
04887   CkAssert(forward_plan != NULL);
04888   CkAssert(backward_plan != NULL);
04889 #else
04890   forward_plan = fftw_create_plan_specific(K2, FFTW_FORWARD,
04891         ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04892         | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
04893         nz, (fftw_complex *) work, 1);
04894   backward_plan = fftw_create_plan_specific(K2, FFTW_BACKWARD,
04895         ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04896         | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
04897         nz, (fftw_complex *) work, 1);
04898 #endif
04899   CmiUnlock(ComputePmeMgr::fftw_plan_lock);
04900 #else
04901   NAMD_die("Sorry, FFTW must be compiled in to use PME.");
04902 #endif
04903 
04904 #if USE_NODE_PAR_RECEIVE
04905   evir = 0;
04906   CmiMemoryWriteFence();
04907 #endif
04908 }

void PmeYPencil::forward_fft (  ) 

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

05396                              {
05397     evir = 0.;
05398 #ifdef NAMD_FFTW
05399 #ifdef MANUAL_DEBUG_FFTW3
05400   dumpMatrixFloat3("fw_y_b", data, nx, initdata.grid.K2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05401 #endif
05402   
05403 #ifdef NAMD_FFTW_3
05404 #if     CMK_SMP && USE_CKLOOP
05405   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05406   if(useCkLoop>=CKLOOP_CTRL_PME_FORWARDFFT
05407      && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
05408           CkLoop_Parallelize(PmeYPencilForwardFFT, 1, (void *)this, CkMyNodeSize(), 0, nx-1); //sync
05409           return;
05410   }
05411 #endif
05412   //the above is a transformation of the following loop using CkLoop
05413   for ( int i=0; i<nx; ++i ) {
05414     fftwf_execute_dft(forward_plan, ((fftwf_complex *) data) + i 
05415                       * nz * initdata.grid.K2,  
05416                       ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
05417   }
05418 #else
05419   for ( int i=0; i<nx; ++i ) {
05420     fftw(forward_plan, nz,
05421         ((fftw_complex *) data) + i * nz * initdata.grid.K2,
05422         nz, 1, (fftw_complex *) work, 1, 0);
05423   }
05424 #endif
05425 #ifdef MANUAL_DEBUG_FFTW3
05426   dumpMatrixFloat3("fw_y_a", data, nx, initdata.grid.dim2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05427 #endif
05428 
05429 #endif
05430 }

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

Definition at line 5384 of file ComputePme.C.

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

Referenced by PmeYPencilForwardFFT().

05384                                                           {
05385 #ifdef NAMD_FFTW
05386 #ifdef NAMD_FFTW_3
05387         for(int i=fromIdx; i<=toIdx; i++){
05388                 fftwf_execute_dft(forward_plan, ((fftwf_complex *) data) + i 
05389                       * nz * initdata.grid.K2,  
05390                       ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
05391         }
05392 #endif
05393 #endif
05394 }

void PmeYPencil::node_process_trans ( PmeTransMsg msg  ) 

Definition at line 4910 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(), PmeTransMsg::sourceNode, and PmePencilInitMsgData::yBlocks.

Referenced by NodePmeMgr::recvYTrans().

04911 {
04912   if ( msg->hasData ) hasData = 1;
04913   needs_reply[msg->sourceNode] = msg->hasData;
04914   recv_trans(msg);
04915   int limsg;
04916   CmiMemoryAtomicFetchAndInc(imsg,limsg);
04917   if(limsg+1 == initdata.yBlocks)
04918     {
04919       if ( hasData ) {
04920         forward_fft();
04921       }
04922       send_trans();
04923       imsg=0;
04924       CmiMemoryWriteFence();
04925     }
04926 }

void PmeYPencil::node_process_untrans ( PmeUntransMsg msg  ) 

Definition at line 4933 of file ComputePme.C.

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

Referenced by recvNodeAck(), and NodePmeMgr::recvYUntrans().

04934 {
04935   if ( msg ) {
04936     if ( ! hasData ) NAMD_bug("PmeYPencil::node_process_untrans non-null msg but not hasData");
04937     recv_untrans(msg);
04938   } else if ( hasData ) NAMD_bug("PmeYPencil::node_process_untrans hasData but null msg");
04939   int limsg;
04940   CmiMemoryAtomicFetchAndInc(imsgb,limsg);
04941   if(limsg+1 == initdata.yBlocks)
04942     {
04943      if ( hasData ) {
04944       backward_fft();
04945      }
04946       hasData=0;
04947       imsgb=0;
04948       CmiMemoryWriteFence();
04949       send_untrans();
04950     }
04951 }

void PmeYPencil::recv_trans ( const PmeTransMsg msg  ) 

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

05343                                                   {
05344   if ( imsg == 0 ) {
05345     lattice = msg->lattice;
05346     sequence = msg->sequence;
05347   }
05348   int block2 = initdata.grid.block2;
05349   int K2 = initdata.grid.K2;
05350   int jb = msg->sourceNode;
05351   int ny = msg->nx;
05352  if ( msg->hasData ) {
05353   const float *md = msg->qgrid;
05354   float *d = data;
05355   for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05356    for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05357     for ( int k=0; k<nz; ++k ) {
05358 #ifdef ZEROCHECK
05359       if ( (*md) == 0. ) CkPrintf("0 in ZY at %d %d %d %d %d %d %d %d %d\n",
05360         thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05361 #endif
05362       d[2*(j*nz+k)] = *(md++);
05363       d[2*(j*nz+k)+1] = *(md++);
05364     }
05365    }
05366   }
05367  } else {
05368   float *d = data;
05369   for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05370    for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05371     for ( int k=0; k<nz; ++k ) {
05372       d[2*(j*nz+k)] = 0;
05373       d[2*(j*nz+k)+1] = 0;
05374     }
05375    }
05376   }
05377  }
05378 }

void PmeYPencil::recv_untrans ( const PmeUntransMsg msg  ) 

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

05865                                                       {
05866   int block2 = initdata.grid.block2;
05867   int K2 = initdata.grid.K2;
05868   int jb = msg->sourceNode;
05869   int ny = msg->ny;
05870   const float *md = msg->qgrid;
05871   float *d = data;
05872   for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05873 #if CMK_BLUEGENEL
05874     CmiNetworkProgress();
05875 #endif   
05876     for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05877       for ( int k=0; k<nz; ++k ) {
05878 #ifdef ZEROCHECK
05879         if ( (*md) == 0. ) CkPrintf("0 in XY at %d %d %d %d %d %d %d %d %d\n",
05880                                     thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05881 #endif
05882         d[2*(j*nz+k)] = *(md++);
05883         d[2*(j*nz+k)+1] = *(md++);
05884       }
05885     }
05886   }
05887 }

void PmeYPencil::recvNodeAck ( PmeAckMsg msg  ) 

Definition at line 4928 of file ComputePme.C.

References node_process_untrans().

04928                                            {
04929   delete msg;
04930   node_process_untrans(0);
04931 }

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

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

05437                                                         {
05438         int yBlocks = initdata.yBlocks;
05439         int block2 = initdata.grid.block2;
05440         int K2 = initdata.grid.K2;
05441     for ( int isend=fromIdx; isend<=toIdx; ++isend ) {
05442           int jb = send_order[isend];
05443           int ny = block2;
05444           if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
05445           int hd = ( hasData ? 1 : 0 );
05446           PmeTransMsg *msg = new (hd*nx*ny*nz*2,PRIORITY_SIZE) PmeTransMsg;
05447           msg->lattice = lattice;
05448           msg->sourceNode = thisIndex.x;
05449           msg->hasData = hasData;
05450           msg->nx = nx;
05451          if ( hasData ) {
05452           float *md = msg->qgrid;
05453           const float *d = data;
05454           for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05455            for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05456                 for ( int k=0; k<nz; ++k ) {
05457                   *(md++) = d[2*(j*nz+k)];
05458                   *(md++) = d[2*(j*nz+k)+1];
05459   #ifdef ZEROCHECK
05460                   if ( *(md-2) == 0. ) CkPrintf("send 0 in YX at %d %d %d %d %d %d %d %d %d\n",
05461           thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05462   #endif
05463                 }
05464            }
05465           }
05466           if ( md != msg->qgrid + nx*ny*nz*2 ) CkPrintf("error in YX at %d %d %d\n",
05467           thisIndex.x, jb, thisIndex.z);
05468          }
05469           msg->sequence = sequence;
05470           SET_PRIORITY(msg,sequence,PME_TRANS2_PRIORITY)
05471       CmiEnableUrgentSend(1);
05472 #if USE_NODE_PAR_RECEIVE
05473       msg->destElem=CkArrayIndex3D(0,jb,thisIndex.z);
05474 #if X_PERSIST 
05475       CmiUsePersistentHandle(&trans_handle[isend], 1);
05476 #endif
05477       initdata.pmeNodeProxy[CmiNodeOf(initdata.xm.ckLocalBranch()->procNum(0,msg->destElem))].recvXTrans(msg);   
05478 #if X_PERSIST 
05479       CmiUsePersistentHandle(NULL, 0);
05480 #endif
05481 #else      
05482 #if X_PERSIST 
05483       CmiUsePersistentHandle(&trans_handle[isend], 1);
05484 #endif
05485       initdata.xPencil(0,jb,thisIndex.z).recvTrans(msg);
05486 #if X_PERSIST 
05487       CmiUsePersistentHandle(NULL, 0);
05488 #endif
05489 #endif
05490       CmiEnableUrgentSend(0);
05491         }
05492 }

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

Definition at line 5953 of file ComputePme.C.

References PmeGrid::block2, PmePencil< CBase_PmeYPencil >::data, PmeUntransMsg::destElem, PmePencilInitMsgData::grid, PmePencil< CBase_PmeYPencil >::initdata, j, PmeGrid::K2, PmePencil< CBase_PmeYPencil >::needs_reply, 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().

05953                                                           {
05954         int yBlocks = initdata.yBlocks;
05955         int block2 = initdata.grid.block2;      
05956         int K2 = initdata.grid.K2;
05957 
05958         for(int isend=fromIdx; isend<=toIdx; isend++) {
05959                 int jb = send_order[isend];
05960                 if ( ! needs_reply[jb] ) {
05961                         PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
05962                         CmiEnableUrgentSend(1);
05963                         SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
05964 #if USE_NODE_PAR_RECEIVE
05965                         initdata.zPencil(thisIndex.x,jb,0).recvNodeAck(msg);
05966 #else
05967                         initdata.zPencil(thisIndex.x,jb,0).recvAck(msg);
05968 #endif
05969                         CmiEnableUrgentSend(0);
05970                         continue;
05971                 }
05972                 int ny = block2;
05973                 if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
05974                 PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
05975                 msg->sourceNode = thisIndex.z;
05976                 msg->ny = nz;
05977                 float *md = msg->qgrid;
05978                 const float *d = data;
05979                 for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05980                         for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05981                                 for ( int k=0; k<nz; ++k ) {
05982                                         *(md++) = d[2*(j*nz+k)];
05983                                         *(md++) = d[2*(j*nz+k)+1];
05984                                 }
05985                         }
05986                 }
05987                 SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
05988             CmiEnableUrgentSend(1);
05989 #if USE_NODE_PAR_RECEIVE
05990         msg->destElem=CkArrayIndex3D( thisIndex.x, jb, 0);
05991         //    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)));
05992 #if Z_PERSIST 
05993         CmiUsePersistentHandle(&untrans_handle[isend], 1);
05994 #endif
05995         initdata.pmeNodeProxy[CmiNodeOf(initdata.zm.ckLocalBranch()->procNum(0,msg->destElem))].recvZUntrans(msg);
05996 #if Z_PERSIST 
05997         CmiUsePersistentHandle(NULL, 0);
05998 #endif
05999 #else
06000 #if Z_PERSIST 
06001         CmiUsePersistentHandle(&untrans_handle[isend], 1);
06002 #endif
06003         initdata.zPencil(thisIndex.x,jb,0).recvUntrans(msg);
06004 #if Z_PERSIST 
06005         CmiUsePersistentHandle(NULL, 0);
06006 #endif
06007 #endif
06008     CmiEnableUrgentSend(0);
06009         }
06010 }

void PmeYPencil::send_trans (  ) 

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

05494                             {
05495 #if USE_PERSISTENT
05496     if (trans_handle == NULL) setup_persistent();
05497 #endif
05498 #if     CMK_SMP && USE_CKLOOP
05499         int useCkLoop = Node::Object()->simParameters->useCkLoop;
05500         if(useCkLoop>=CKLOOP_CTRL_PME_SENDTRANS
05501            && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
05508                 //send_subset_trans(0, initdata.yBlocks-1);
05509                 CkLoop_Parallelize(PmeYPencilSendTrans, 1, (void *)this, CkMyNodeSize(), 0, initdata.yBlocks-1, 1); //not sync
05510                 return;
05511         }
05512 #endif
05513   int yBlocks = initdata.yBlocks;
05514   int block2 = initdata.grid.block2;
05515   int K2 = initdata.grid.K2;
05516   for ( int isend=0; isend<yBlocks; ++isend ) {
05517     int jb = send_order[isend];
05518     int ny = block2;
05519     if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
05520     int hd = ( hasData ? 1 : 0 );
05521     PmeTransMsg *msg = new (hd*nx*ny*nz*2,PRIORITY_SIZE) PmeTransMsg;
05522     msg->lattice = lattice;
05523     msg->sourceNode = thisIndex.x;
05524     msg->hasData = hasData;
05525     msg->nx = nx;
05526    if ( hasData ) {
05527     float *md = msg->qgrid;
05528     const float *d = data;
05529     for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05530      for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05531       for ( int k=0; k<nz; ++k ) {
05532         *(md++) = d[2*(j*nz+k)];
05533         *(md++) = d[2*(j*nz+k)+1];
05534 #ifdef ZEROCHECK
05535         if ( *(md-2) == 0. ) CkPrintf("send 0 in YX at %d %d %d %d %d %d %d %d %d\n",
05536         thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05537 #endif
05538       }
05539      }
05540     }
05541     if ( md != msg->qgrid + nx*ny*nz*2 ) CkPrintf("error in YX at %d %d %d\n",
05542         thisIndex.x, jb, thisIndex.z);
05543    }
05544     msg->sequence = sequence;
05545     SET_PRIORITY(msg,sequence,PME_TRANS2_PRIORITY)
05546     CmiEnableUrgentSend(1);
05547 #if USE_NODE_PAR_RECEIVE
05548     msg->destElem=CkArrayIndex3D(0,jb,thisIndex.z);
05549 #if X_PERSIST 
05550         CmiUsePersistentHandle(&trans_handle[isend], 1);
05551 #endif
05552     initdata.pmeNodeProxy[CmiNodeOf(initdata.xm.ckLocalBranch()->procNum(0,msg->destElem))].recvXTrans(msg);   
05553 #if X_PERSIST 
05554         CmiUsePersistentHandle(NULL, 0);
05555 #endif
05556 #else
05557 #if X_PERSIST 
05558         CmiUsePersistentHandle(&trans_handle[isend], 1);
05559 #endif
05560     initdata.xPencil(0,jb,thisIndex.z).recvTrans(msg);
05561 #if X_PERSIST 
05562         CmiUsePersistentHandle(NULL, 0);
05563 #endif
05564     
05565 #endif
05566     CmiEnableUrgentSend(0);
05567   }
05568 }

void PmeYPencil::send_untrans (  ) 

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

06012                               {
06013 #if USE_PERSISTENT
06014   if (untrans_handle == NULL) setup_persistent();
06015 #endif
06016 #if     CMK_SMP && USE_CKLOOP
06017   int useCkLoop = Node::Object()->simParameters->useCkLoop;
06018   if(useCkLoop>=CKLOOP_CTRL_PME_SENDUNTRANS
06019      && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
06020           int yBlocks = initdata.yBlocks;
06021 
06022 #if USE_NODE_PAR_RECEIVE      
06023           //CkLoop_Parallelize(PmeYPencilSendUntrans, 1, (void *)this, CkMyNodeSize(), 0, yBlocks-1, 1); //sync
06024           CkLoop_Parallelize(PmeYPencilSendUntrans, 1, (void *)this, yBlocks, 0, yBlocks-1, 1);
06025 #else
06026       //CkLoop_Parallelize(PmeYPencilSendUntrans, 1, (void *)this, CkMyNodeSize(), 0, yBlocks-1, 0); //not sync
06027           CkLoop_Parallelize(PmeYPencilSendUntrans, 1, (void *)this, yBlocks, 0, yBlocks-1, 0); //not sync
06028 #endif
06029           return;
06030   }
06031 #endif
06032   int yBlocks = initdata.yBlocks;
06033   int block2 = initdata.grid.block2;
06034   int K2 = initdata.grid.K2;
06035   for ( int isend=0; isend<yBlocks; ++isend ) {
06036     int jb = send_order[isend];
06037     if ( ! needs_reply[jb] ) {
06038       PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
06039       CmiEnableUrgentSend(1);
06040       SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
06041 #if USE_NODE_PAR_RECEIVE
06042       initdata.zPencil(thisIndex.x,jb,0).recvNodeAck(msg);
06043 #else
06044       initdata.zPencil(thisIndex.x,jb,0).recvAck(msg);
06045 #endif
06046       CmiEnableUrgentSend(0);
06047       continue;
06048     }
06049     int ny = block2;
06050     if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
06051     PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
06052     msg->sourceNode = thisIndex.z;
06053     msg->ny = nz;
06054     float *md = msg->qgrid;
06055     const float *d = data;
06056     for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
06057      for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
06058       for ( int k=0; k<nz; ++k ) {
06059         *(md++) = d[2*(j*nz+k)];
06060         *(md++) = d[2*(j*nz+k)+1];
06061       }
06062      }
06063     }
06064     SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
06065 
06066     CmiEnableUrgentSend(1);
06067 #if USE_NODE_PAR_RECEIVE
06068     msg->destElem=CkArrayIndex3D( thisIndex.x, jb, 0);
06069     //    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)));
06070 #if Z_PERSIST 
06071     CmiUsePersistentHandle(&untrans_handle[isend], 1);
06072 #endif
06073     initdata.pmeNodeProxy[CmiNodeOf(initdata.zm.ckLocalBranch()->procNum(0,msg->destElem))].recvZUntrans(msg);
06074 #if Z_PERSIST
06075     CmiUsePersistentHandle(NULL, 0);
06076 #endif
06077 #else
06078 #if Z_PERSIST 
06079     CmiUsePersistentHandle(&untrans_handle[isend], 1);
06080 #endif
06081     initdata.zPencil(thisIndex.x,jb,0).recvUntrans(msg);
06082 #if Z_PERSIST 
06083     CmiUsePersistentHandle(NULL, 0);
06084 #endif
06085 #endif    
06086     CmiEnableUrgentSend(0);
06087   }
06088   
06089 #if USE_NODE_PAR_RECEIVE
06090   evir = 0.;
06091   CmiMemoryWriteFence();
06092 #endif
06093 }


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

Generated on 8 Dec 2019 for NAMD by  doxygen 1.6.1