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


Constructor & Destructor Documentation

PmeYPencil_SDAG_CODE PmeYPencil::PmeYPencil (  )  [inline]

Definition at line 4638 of file ComputePme.C.

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

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

PmeYPencil::PmeYPencil ( CkMigrateMessage *   )  [inline]

Definition at line 4639 of file ComputePme.C.

04639 { __sdag_init(); }


Member Function Documentation

void PmeYPencil::backward_fft (  ) 

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

05932                               {
05933 #ifdef NAMD_FFTW
05934 #ifdef MANUAL_DEBUG_FFTW3
05935   dumpMatrixFloat3("bw_y_b", data, nx, initdata.grid.K2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05936 #endif
05937 
05938 #ifdef NAMD_FFTW_3
05939 #if     CMK_SMP && USE_CKLOOP
05940   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05941   if(useCkLoop>=CKLOOP_CTRL_PME_BACKWARDFFT
05942      && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
05943           CkLoop_Parallelize(PmeYPencilBackwardFFT, 1, (void *)this, CkMyNodeSize(), 0, nx-1); //sync
05944           return;
05945   }
05946 #endif
05947   //the above is a transformation of the following loop using CkLoop
05948   for ( int i=0; i<nx; ++i ) {
05949 #if CMK_BLUEGENEL
05950         CmiNetworkProgress();
05951 #endif
05952     fftwf_execute_dft(backward_plan,    
05953                                           ((fftwf_complex *) data) + i * nz * initdata.grid.K2,
05954                                           ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
05955   }
05956 #else
05957         for ( int i=0; i<nx; ++i ) {
05958 #if CMK_BLUEGENEL
05959           CmiNetworkProgress();
05960 #endif
05961                 fftw(backward_plan, nz,
05962                 ((fftw_complex *) data) + i * nz * initdata.grid.K2,
05963                 nz, 1, (fftw_complex *) work, 1, 0);
05964         }
05965 #endif
05966 
05967 #ifdef MANUAL_DEBUG_FFTW3
05968   dumpMatrixFloat3("bw_y_a", data, nx, initdata.grid.K2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05969 #endif
05970 
05971 #endif
05972 }

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

Definition at line 5920 of file ComputePme.C.

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

Referenced by PmeYPencilBackwardFFT().

05920                                                            {
05921 #ifdef NAMD_FFTW
05922 #ifdef NAMD_FFTW_3
05923         for(int i=fromIdx; i<=toIdx; i++){
05924                 fftwf_execute_dft(backward_plan,        
05925                                                   ((fftwf_complex *) data) + i * nz * initdata.grid.K2,         
05926                                                   ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
05927         }
05928 #endif
05929 #endif
05930 }

void PmeYPencil::fft_init (  ) 

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

04865                           {
04866   CProxy_Node nd(CkpvAccess(BOCclass_group).node);
04867   Node *node = nd.ckLocalBranch();
04868   SimParameters *simParams = node->simParameters;
04869 
04870 #if USE_NODE_PAR_RECEIVE
04871   ((NodePmeMgr *)CkLocalNodeBranch(initdata.pmeNodeProxy))->registerYPencil(thisIndex,this);
04872 #endif
04873 
04874   int K1 = initdata.grid.K1;
04875   int K2 = initdata.grid.K2;
04876   int dim2 = initdata.grid.dim2;
04877   int dim3 = initdata.grid.dim3;
04878   int block1 = initdata.grid.block1;
04879   int block3 = initdata.grid.block3;
04880 
04881   nx = block1;
04882   if ( (thisIndex.x + 1) * block1 > K1 ) nx = K1 - thisIndex.x * block1;
04883   nz = block3;
04884   if ( (thisIndex.z+1)*block3 > dim3/2 ) nz = dim3/2 - thisIndex.z*block3;
04885 
04886 #ifdef NAMD_FFTW
04887   CmiLock(ComputePmeMgr::fftw_plan_lock);
04888 
04889   data = (float *) fftwf_malloc( sizeof(float) * nx*dim2*nz*2);
04890   work = new float[2*K2];
04891 
04892   order_init(initdata.yBlocks);
04893 
04894 #ifdef NAMD_FFTW_3
04895   /* need array of sizes for the dimensions */
04896   /* ideally this should be implementable as a single multidimensional
04897    *  plan, but that has proven tricky to implement, so we maintain the
04898    *  loop of 1d plan executions. */
04899   int sizeLines=nz;
04900   int planLineSizes[1];
04901   planLineSizes[0]=K2;
04902   int fftwFlags = simParams->FFTWPatient ? FFTW_PATIENT  : simParams->FFTWEstimate ? FFTW_ESTIMATE  : FFTW_MEASURE ;
04903   forward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines, 
04904                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04905                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04906                                      FFTW_FORWARD, 
04907                                      fftwFlags);
04908   backward_plan = fftwf_plan_many_dft(1, planLineSizes, sizeLines, 
04909                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04910                                      (fftwf_complex *) data, NULL, sizeLines, 1,
04911                                      FFTW_BACKWARD, 
04912                                       fftwFlags);
04913   CkAssert(forward_plan != NULL);
04914   CkAssert(backward_plan != NULL);
04915 #else
04916   forward_plan = fftw_create_plan_specific(K2, FFTW_FORWARD,
04917         ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04918         | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
04919         nz, (fftw_complex *) work, 1);
04920   backward_plan = fftw_create_plan_specific(K2, FFTW_BACKWARD,
04921         ( simParams->FFTWEstimate ? FFTW_ESTIMATE : FFTW_MEASURE )
04922         | FFTW_IN_PLACE | FFTW_USE_WISDOM, (fftw_complex *) data,
04923         nz, (fftw_complex *) work, 1);
04924 #endif
04925   CmiUnlock(ComputePmeMgr::fftw_plan_lock);
04926 #else
04927   NAMD_die("Sorry, FFTW must be compiled in to use PME.");
04928 #endif
04929 
04930 #if USE_NODE_PAR_RECEIVE
04931   evir = 0;
04932   CmiMemoryWriteFence();
04933 #endif
04934 }

void PmeYPencil::forward_fft (  ) 

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

05422                              {
05423     evir = 0.;
05424 #ifdef NAMD_FFTW
05425 #ifdef MANUAL_DEBUG_FFTW3
05426   dumpMatrixFloat3("fw_y_b", data, nx, initdata.grid.K2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05427 #endif
05428   
05429 #ifdef NAMD_FFTW_3
05430 #if     CMK_SMP && USE_CKLOOP
05431   int useCkLoop = Node::Object()->simParameters->useCkLoop;
05432   if(useCkLoop>=CKLOOP_CTRL_PME_FORWARDFFT
05433      && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
05434           CkLoop_Parallelize(PmeYPencilForwardFFT, 1, (void *)this, CkMyNodeSize(), 0, nx-1); //sync
05435           return;
05436   }
05437 #endif
05438   //the above is a transformation of the following loop using CkLoop
05439   for ( int i=0; i<nx; ++i ) {
05440     fftwf_execute_dft(forward_plan, ((fftwf_complex *) data) + i 
05441                       * nz * initdata.grid.K2,  
05442                       ((fftwf_complex *) data) + i * nz * initdata.grid.K2);
05443   }
05444 #else
05445   for ( int i=0; i<nx; ++i ) {
05446     fftw(forward_plan, nz,
05447         ((fftw_complex *) data) + i * nz * initdata.grid.K2,
05448         nz, 1, (fftw_complex *) work, 1, 0);
05449   }
05450 #endif
05451 #ifdef MANUAL_DEBUG_FFTW3
05452   dumpMatrixFloat3("fw_y_a", data, nx, initdata.grid.dim2, nz, thisIndex.x, thisIndex.y, thisIndex.z);
05453 #endif
05454 
05455 #endif
05456 }

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

Definition at line 5410 of file ComputePme.C.

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

Referenced by PmeYPencilForwardFFT().

05410                                                           {
05411 #ifdef NAMD_FFTW
05412 #ifdef NAMD_FFTW_3
05413         for(int i=fromIdx; i<=toIdx; 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 #endif
05419 #endif
05420 }

void PmeYPencil::node_process_trans ( PmeTransMsg msg  ) 

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

04937 {
04938   if ( msg->hasData ) hasData = 1;
04939   needs_reply[msg->sourceNode] = msg->hasData;
04940   recv_trans(msg);
04941   int limsg;
04942   CmiMemoryAtomicFetchAndInc(imsg,limsg);
04943   if(limsg+1 == initdata.yBlocks)
04944     {
04945       if ( hasData ) {
04946         forward_fft();
04947       }
04948       send_trans();
04949       imsg=0;
04950       CmiMemoryWriteFence();
04951     }
04952 }

void PmeYPencil::node_process_untrans ( PmeUntransMsg msg  ) 

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

04960 {
04961   if ( msg ) {
04962     if ( ! hasData ) NAMD_bug("PmeYPencil::node_process_untrans non-null msg but not hasData");
04963     recv_untrans(msg);
04964   } else if ( hasData ) NAMD_bug("PmeYPencil::node_process_untrans hasData but null msg");
04965   int limsg;
04966   CmiMemoryAtomicFetchAndInc(imsgb,limsg);
04967   if(limsg+1 == initdata.yBlocks)
04968     {
04969      if ( hasData ) {
04970       backward_fft();
04971      }
04972       hasData=0;
04973       imsgb=0;
04974       CmiMemoryWriteFence();
04975       send_untrans();
04976     }
04977 }

void PmeYPencil::recv_trans ( const PmeTransMsg msg  ) 

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

05369                                                   {
05370   if ( imsg == 0 ) {
05371     lattice = msg->lattice;
05372     sequence = msg->sequence;
05373   }
05374   int block2 = initdata.grid.block2;
05375   int K2 = initdata.grid.K2;
05376   int jb = msg->sourceNode;
05377   int ny = msg->nx;
05378  if ( msg->hasData ) {
05379   const float *md = msg->qgrid;
05380   float *d = data;
05381   for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05382    for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05383     for ( int k=0; k<nz; ++k ) {
05384 #ifdef ZEROCHECK
05385       if ( (*md) == 0. ) CkPrintf("0 in ZY at %d %d %d %d %d %d %d %d %d\n",
05386         thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05387 #endif
05388       d[2*(j*nz+k)] = *(md++);
05389       d[2*(j*nz+k)+1] = *(md++);
05390     }
05391    }
05392   }
05393  } else {
05394   float *d = data;
05395   for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05396    for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05397     for ( int k=0; k<nz; ++k ) {
05398       d[2*(j*nz+k)] = 0;
05399       d[2*(j*nz+k)+1] = 0;
05400     }
05401    }
05402   }
05403  }
05404 }

void PmeYPencil::recv_untrans ( const PmeUntransMsg msg  ) 

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

05891                                                       {
05892   int block2 = initdata.grid.block2;
05893   int K2 = initdata.grid.K2;
05894   int jb = msg->sourceNode;
05895   int ny = msg->ny;
05896   const float *md = msg->qgrid;
05897   float *d = data;
05898   for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05899 #if CMK_BLUEGENEL
05900     CmiNetworkProgress();
05901 #endif   
05902     for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05903       for ( int k=0; k<nz; ++k ) {
05904 #ifdef ZEROCHECK
05905         if ( (*md) == 0. ) CkPrintf("0 in XY at %d %d %d %d %d %d %d %d %d\n",
05906                                     thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05907 #endif
05908         d[2*(j*nz+k)] = *(md++);
05909         d[2*(j*nz+k)+1] = *(md++);
05910       }
05911     }
05912   }
05913 }

void PmeYPencil::recvNodeAck ( PmeAckMsg msg  ) 

Definition at line 4954 of file ComputePme.C.

References node_process_untrans().

04954                                            {
04955   delete msg;
04956   node_process_untrans(0);
04957 }

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

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

05463                                                         {
05464         int yBlocks = initdata.yBlocks;
05465         int block2 = initdata.grid.block2;
05466         int K2 = initdata.grid.K2;
05467     for ( int isend=fromIdx; isend<=toIdx; ++isend ) {
05468           int jb = send_order[isend];
05469           int ny = block2;
05470           if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
05471           int hd = ( hasData ? 1 : 0 );
05472           PmeTransMsg *msg = new (hd*nx*ny*nz*2,PRIORITY_SIZE) PmeTransMsg;
05473           msg->lattice = lattice;
05474           msg->sourceNode = thisIndex.x;
05475           msg->hasData = hasData;
05476           msg->nx = nx;
05477          if ( hasData ) {
05478           float *md = msg->qgrid;
05479           const float *d = data;
05480           for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05481            for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05482                 for ( int k=0; k<nz; ++k ) {
05483                   *(md++) = d[2*(j*nz+k)];
05484                   *(md++) = d[2*(j*nz+k)+1];
05485   #ifdef ZEROCHECK
05486                   if ( *(md-2) == 0. ) CkPrintf("send 0 in YX at %d %d %d %d %d %d %d %d %d\n",
05487           thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05488   #endif
05489                 }
05490            }
05491           }
05492           if ( md != msg->qgrid + nx*ny*nz*2 ) CkPrintf("error in YX at %d %d %d\n",
05493           thisIndex.x, jb, thisIndex.z);
05494          }
05495           msg->sequence = sequence;
05496           SET_PRIORITY(msg,sequence,PME_TRANS2_PRIORITY)
05497       CmiEnableUrgentSend(1);
05498 #if USE_NODE_PAR_RECEIVE
05499       msg->destElem=CkArrayIndex3D(0,jb,thisIndex.z);
05500 #if X_PERSIST 
05501       CmiUsePersistentHandle(&trans_handle[isend], 1);
05502 #endif
05503       initdata.pmeNodeProxy[CmiNodeOf(initdata.xm.ckLocalBranch()->procNum(0,msg->destElem))].recvXTrans(msg);   
05504 #if X_PERSIST 
05505       CmiUsePersistentHandle(NULL, 0);
05506 #endif
05507 #else      
05508 #if X_PERSIST 
05509       CmiUsePersistentHandle(&trans_handle[isend], 1);
05510 #endif
05511       initdata.xPencil(0,jb,thisIndex.z).recvTrans(msg);
05512 #if X_PERSIST 
05513       CmiUsePersistentHandle(NULL, 0);
05514 #endif
05515 #endif
05516       CmiEnableUrgentSend(0);
05517         }
05518 }

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

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

05979                                                           {
05980         int yBlocks = initdata.yBlocks;
05981         int block2 = initdata.grid.block2;      
05982         int K2 = initdata.grid.K2;
05983 
05984         for(int isend=fromIdx; isend<=toIdx; isend++) {
05985                 int jb = send_order[isend];
05986                 if ( ! needs_reply[jb] ) {
05987                         PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
05988                         CmiEnableUrgentSend(1);
05989                         SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
05990 #if USE_NODE_PAR_RECEIVE
05991                         initdata.zPencil(thisIndex.x,jb,0).recvNodeAck(msg);
05992 #else
05993                         initdata.zPencil(thisIndex.x,jb,0).recvAck(msg);
05994 #endif
05995                         CmiEnableUrgentSend(0);
05996                         continue;
05997                 }
05998                 int ny = block2;
05999                 if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
06000                 PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
06001                 msg->sourceNode = thisIndex.z;
06002                 msg->ny = nz;
06003                 float *md = msg->qgrid;
06004                 const float *d = data;
06005                 for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
06006                         for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
06007                                 for ( int k=0; k<nz; ++k ) {
06008                                         *(md++) = d[2*(j*nz+k)];
06009                                         *(md++) = d[2*(j*nz+k)+1];
06010                                 }
06011                         }
06012                 }
06013                 SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
06014             CmiEnableUrgentSend(1);
06015 #if USE_NODE_PAR_RECEIVE
06016         msg->destElem=CkArrayIndex3D( thisIndex.x, jb, 0);
06017         //    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)));
06018 #if Z_PERSIST 
06019         CmiUsePersistentHandle(&untrans_handle[isend], 1);
06020 #endif
06021         initdata.pmeNodeProxy[CmiNodeOf(initdata.zm.ckLocalBranch()->procNum(0,msg->destElem))].recvZUntrans(msg);
06022 #if Z_PERSIST 
06023         CmiUsePersistentHandle(NULL, 0);
06024 #endif
06025 #else
06026 #if Z_PERSIST 
06027         CmiUsePersistentHandle(&untrans_handle[isend], 1);
06028 #endif
06029         initdata.zPencil(thisIndex.x,jb,0).recvUntrans(msg);
06030 #if Z_PERSIST 
06031         CmiUsePersistentHandle(NULL, 0);
06032 #endif
06033 #endif
06034     CmiEnableUrgentSend(0);
06035         }
06036 }

void PmeYPencil::send_trans (  ) 

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

05520                             {
05521 #if USE_PERSISTENT
05522     if (trans_handle == NULL) setup_persistent();
05523 #endif
05524 #if     CMK_SMP && USE_CKLOOP
05525         int useCkLoop = Node::Object()->simParameters->useCkLoop;
05526         if(useCkLoop>=CKLOOP_CTRL_PME_SENDTRANS
05527            && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
05534                 //send_subset_trans(0, initdata.yBlocks-1);
05535                 CkLoop_Parallelize(PmeYPencilSendTrans, 1, (void *)this, CkMyNodeSize(), 0, initdata.yBlocks-1, 1); //not sync
05536                 return;
05537         }
05538 #endif
05539   int yBlocks = initdata.yBlocks;
05540   int block2 = initdata.grid.block2;
05541   int K2 = initdata.grid.K2;
05542   for ( int isend=0; isend<yBlocks; ++isend ) {
05543     int jb = send_order[isend];
05544     int ny = block2;
05545     if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
05546     int hd = ( hasData ? 1 : 0 );
05547     PmeTransMsg *msg = new (hd*nx*ny*nz*2,PRIORITY_SIZE) PmeTransMsg;
05548     msg->lattice = lattice;
05549     msg->sourceNode = thisIndex.x;
05550     msg->hasData = hasData;
05551     msg->nx = nx;
05552    if ( hasData ) {
05553     float *md = msg->qgrid;
05554     const float *d = data;
05555     for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
05556      for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
05557       for ( int k=0; k<nz; ++k ) {
05558         *(md++) = d[2*(j*nz+k)];
05559         *(md++) = d[2*(j*nz+k)+1];
05560 #ifdef ZEROCHECK
05561         if ( *(md-2) == 0. ) CkPrintf("send 0 in YX at %d %d %d %d %d %d %d %d %d\n",
05562         thisIndex.x, jb, thisIndex.z, i, j, k, nx, ny, nz);
05563 #endif
05564       }
05565      }
05566     }
05567     if ( md != msg->qgrid + nx*ny*nz*2 ) CkPrintf("error in YX at %d %d %d\n",
05568         thisIndex.x, jb, thisIndex.z);
05569    }
05570     msg->sequence = sequence;
05571     SET_PRIORITY(msg,sequence,PME_TRANS2_PRIORITY)
05572     CmiEnableUrgentSend(1);
05573 #if USE_NODE_PAR_RECEIVE
05574     msg->destElem=CkArrayIndex3D(0,jb,thisIndex.z);
05575 #if X_PERSIST 
05576         CmiUsePersistentHandle(&trans_handle[isend], 1);
05577 #endif
05578     initdata.pmeNodeProxy[CmiNodeOf(initdata.xm.ckLocalBranch()->procNum(0,msg->destElem))].recvXTrans(msg);   
05579 #if X_PERSIST 
05580         CmiUsePersistentHandle(NULL, 0);
05581 #endif
05582 #else
05583 #if X_PERSIST 
05584         CmiUsePersistentHandle(&trans_handle[isend], 1);
05585 #endif
05586     initdata.xPencil(0,jb,thisIndex.z).recvTrans(msg);
05587 #if X_PERSIST 
05588         CmiUsePersistentHandle(NULL, 0);
05589 #endif
05590     
05591 #endif
05592     CmiEnableUrgentSend(0);
05593   }
05594 }

void PmeYPencil::send_untrans (  ) 

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

06038                               {
06039 #if USE_PERSISTENT
06040   if (untrans_handle == NULL) setup_persistent();
06041 #endif
06042 #if     CMK_SMP && USE_CKLOOP
06043   int useCkLoop = Node::Object()->simParameters->useCkLoop;
06044   if(useCkLoop>=CKLOOP_CTRL_PME_SENDUNTRANS
06045      && CkNumPes() >= 2 * initdata.xBlocks * initdata.zBlocks) {
06046           int yBlocks = initdata.yBlocks;
06047 
06048 #if USE_NODE_PAR_RECEIVE      
06049           //CkLoop_Parallelize(PmeYPencilSendUntrans, 1, (void *)this, CkMyNodeSize(), 0, yBlocks-1, 1); //sync
06050           CkLoop_Parallelize(PmeYPencilSendUntrans, 1, (void *)this, yBlocks, 0, yBlocks-1, 1);
06051 #else
06052       //CkLoop_Parallelize(PmeYPencilSendUntrans, 1, (void *)this, CkMyNodeSize(), 0, yBlocks-1, 0); //not sync
06053           CkLoop_Parallelize(PmeYPencilSendUntrans, 1, (void *)this, yBlocks, 0, yBlocks-1, 0); //not sync
06054 #endif
06055           return;
06056   }
06057 #endif
06058   int yBlocks = initdata.yBlocks;
06059   int block2 = initdata.grid.block2;
06060   int K2 = initdata.grid.K2;
06061   for ( int isend=0; isend<yBlocks; ++isend ) {
06062     int jb = send_order[isend];
06063     if ( ! needs_reply[jb] ) {
06064       PmeAckMsg *msg = new (PRIORITY_SIZE) PmeAckMsg;
06065       CmiEnableUrgentSend(1);
06066       SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
06067 #if USE_NODE_PAR_RECEIVE
06068       initdata.zPencil(thisIndex.x,jb,0).recvNodeAck(msg);
06069 #else
06070       initdata.zPencil(thisIndex.x,jb,0).recvAck(msg);
06071 #endif
06072       CmiEnableUrgentSend(0);
06073       continue;
06074     }
06075     int ny = block2;
06076     if ( (jb+1)*block2 > K2 ) ny = K2 - jb*block2;
06077     PmeUntransMsg *msg = new (nx*ny*nz*2,PRIORITY_SIZE) PmeUntransMsg;
06078     msg->sourceNode = thisIndex.z;
06079     msg->ny = nz;
06080     float *md = msg->qgrid;
06081     const float *d = data;
06082     for ( int i=0; i<nx; ++i, d += K2*nz*2 ) {
06083      for ( int j=jb*block2; j<(jb*block2+ny); ++j ) {
06084       for ( int k=0; k<nz; ++k ) {
06085         *(md++) = d[2*(j*nz+k)];
06086         *(md++) = d[2*(j*nz+k)+1];
06087       }
06088      }
06089     }
06090     SET_PRIORITY(msg,sequence,PME_UNTRANS2_PRIORITY)
06091 
06092     CmiEnableUrgentSend(1);
06093 #if USE_NODE_PAR_RECEIVE
06094     msg->destElem=CkArrayIndex3D( thisIndex.x, jb, 0);
06095     //    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)));
06096 #if Z_PERSIST 
06097     CmiUsePersistentHandle(&untrans_handle[isend], 1);
06098 #endif
06099     initdata.pmeNodeProxy[CmiNodeOf(initdata.zm.ckLocalBranch()->procNum(0,msg->destElem))].recvZUntrans(msg);
06100 #if Z_PERSIST
06101     CmiUsePersistentHandle(NULL, 0);
06102 #endif
06103 #else
06104 #if Z_PERSIST 
06105     CmiUsePersistentHandle(&untrans_handle[isend], 1);
06106 #endif
06107     initdata.zPencil(thisIndex.x,jb,0).recvUntrans(msg);
06108 #if Z_PERSIST 
06109     CmiUsePersistentHandle(NULL, 0);
06110 #endif
06111 #endif    
06112     CmiEnableUrgentSend(0);
06113   }
06114   
06115 #if USE_NODE_PAR_RECEIVE
06116   evir = 0.;
06117   CmiMemoryWriteFence();
06118 #endif
06119 }


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

Generated on 1 Jun 2020 for NAMD by  doxygen 1.6.1