C++ ParticleSet::rejectMove方法代码示例

  void ThreeDimMomDist::updateDistribution(ParticleSet& PtclSet, TrialWaveFunction& Psi, 
					   IndexType NumCycles) {
    TinyVector<IndexType, 3> Indexes;
    for (IndexType cycle = 0; cycle < NumCycles; cycle++) {

      pcp->NewWalker();
      PosType dr;
      for (Indexes[0] = 0; Indexes[0] < NumPts[0]; Indexes[0]++) {
	dr[1] = 0.0;
	for (Indexes[1] = 0; Indexes[1] < NumPts[1]; Indexes[1]++) {
	  dr[2] = 0.0;
	  for (Indexes[2] = 0; Indexes[2] < NumPts[2]; Indexes[2]++) {
	    IndexType partToDisplace = (*pcp)();

	    PtclSet.makeMove(partToDisplace, dr);

	    if (Indexes[0] == 0 && Indexes[1] == 0 && Indexes[2] == 0) {
	      placeIntsInBin(Indexes, 1.0);
	    } else {
	      placeIntsInBin(Indexes, Psi.ratio(PtclSet, partToDisplace));
	    }
	    totalNumSamples++;
	    PtclSet.rejectMove(partToDisplace);

	    dr[2] += Spacing[2];
	  }      
	  dr[1] += Spacing[1];
	}
	dr[0] += Spacing[0];
      }	

    } 
  }

  NonLocalECPComponent::RealType 
  NonLocalECPComponent::evaluate(ParticleSet& W, TrialWaveFunction& psi,int iat, vector<NonLocalData>& Txy) {
    RealType esum=0.0;

    //int iel=0;
    for(int nn=myTable->M[iat],iel=0; nn<myTable->M[iat+1]; nn++,iel++){

      register RealType r(myTable->r(nn));
      if(r>Rmax) continue;

      register RealType rinv(myTable->rinv(nn));
      register PosType  dr(myTable->dr(nn));

      int txyCounter=Txy.size();
      // Compute ratio of wave functions
      for (int j=0; j < nknot ; j++){ 
        PosType deltar(r*rrotsgrid_m[j]-dr);
        PosType newpos(W.makeMove(iel,deltar)); 
        psiratio[j]=psi.ratio(W,iel)*sgridweight_m[j];
        W.rejectMove(iel);
        //psi.rejectMove(iel);
        //first, add a new NonLocalData with ratio
        Txy.push_back(NonLocalData(iel,psiratio[j],deltar));
      }
      // Compute radial potential
      for(int ip=0;ip< nchannel; ip++){
        vrad[ip]=nlpp_m[ip]->splint(r)*wgt_angpp_m[ip];
      }

      // Compute spherical harmonics on grid
      for (int j=0, jl=0; j<nknot ; j++){ 
        RealType zz=dot(dr,rrotsgrid_m[j])*rinv;
        // Forming the Legendre polynomials
        lpol[0]=1.0;
        RealType lpolprev=0.0;
        for (int l=0 ; l< lmax ; l++){
          //Not a big difference
          //lpol[l+1]=(2*l+1)*zz*lpol[l]-l*lpolprev;
          //lpol[l+1]/=(l+1);
          lpol[l+1]=Lfactor1[l]*zz*lpol[l]-l*lpolprev; 
          lpol[l+1]*=Lfactor2[l]; 
          lpolprev=lpol[l];
        }

        //for(int l=0; l <nchannel; l++,jl++) Amat[jl]=lpol[ angpp_m[l] ]; 
        RealType lsum=0;
        for(int l=0; l <nchannel; l++) lsum += vrad[l]*lpol[ angpp_m[l] ]; 
        esum += Txy[txyCounter++].Weight *= lsum;
      } 
     //BLAS::gemv(nknot, nchannel, &Amat[0], &psiratio[0], &wvec[0]);
     //esum += BLAS::dot(nchannel, &vrad[0], &wvec[0]);
    }   /* end loop over electron */
    return esum;

  }

MomentumEstimator::Return_t MomentumEstimator::evaluate(ParticleSet& P)
{
    const int np=P.getTotalNum();
    nofK=0.0;
    compQ=0.0;

    //will use temp[i].r1 for the Compton profile
    const vector<DistanceTableData::TempDistType>& temp(P.DistTables[0]->Temp);
    Vector<RealType> tmpn_k(nofK);
    for (int s=0; s<M; ++s)
    {
        PosType newpos;

        for (int i=0; i<OHMMS_DIM; ++i) newpos[i]=myRNG();
        //make it cartesian
        newpos=Lattice.toCart(newpos);
        P.makeVirtualMoves(newpos); //updated: temp[i].r1=|newpos-P.R[i]|, temp[i].dr1=newpos-P.R[i]
        refPsi.get_ratios(P,psi_ratios);
//         for (int i=0; i<np; ++i) app_log()<<i<<" "<<psi_ratios[i].real()<<" "<<psi_ratios[i].imag()<<endl;
        P.rejectMove(0); //restore P.R[0] to the orginal position

        for (int ik=0; ik < kPoints.size(); ++ik)
        {
            for (int i=0; i<np; ++i) kdotp[i]=dot(kPoints[ik],temp[i].dr1_nobox);
            eval_e2iphi(np,kdotp.data(),phases.data());
            RealType nofk_here(std::real(BLAS::dot(np,phases.data(),&psi_ratios[0])));//psi_ratios.data())));
            nofK[ik]+= nofk_here;
            tmpn_k[ik]=nofk_here;
        }

        for (int iq=0; iq < compQ.size(); ++iq)
            for (int i=0; i<mappedQtonofK[iq].size(); ++i)
                compQ[iq] += tmpn_k[mappedQtonofK[iq][i]];




    }
    for (int ik=0; ik<nofK.size(); ++ik) nofK[ik] *= norm_nofK;
    for (int iq=0; iq<compQ.size(); ++iq) compQ[iq] *= mappedQnorms[iq];

    if (hdf5_out)
    {
        int j=myIndex;
        for (int ik=0; ik<nofK.size(); ++ik,++j) P.Collectables[j]+= nofK[ik];
        for (int iq=0; iq<compQ.size(); ++iq,++j) P.Collectables[j]+= compQ[iq];
    }
    return 0.0;
}