C++ Site::End方法代码示例

/*!
  Non-Zero Momentum nucleon  from wall sources
*/
void Nuc2pt::calcWallMomNucleon(QPropW& quark, QPropWMomSrc& mom_quark)
{
  char *fname = "Nuc2pt()";
  ThreeMom sink_mom(mom_quark.Mom()) ;
  sink_mom.conj();
  //The momentum on the sink has to be complex conjugated!!
  // p --> - p 
  
  srcdesc[0]=srcdesc[1]=srcdesc[2]=srcdesc[3]=-1 ;
  srcdesc[4]=srcdesc[5]=-1 ; // TIZB
  switch (quark.SrcType())
    {
    case WALL:
      srcdesc[0] = quark.SourceTime() ;
      source="WALL" ;
      break ;
    default:
      ERR.General(cname,fname,"Bad source type\n") ;
    }
  
  quark_mass = quark.Mass() ;
    
  WilsonMatrix Q1;  
  WilsonMatrix Q2;
  WilsonMatrix Q3;  
  
  SnkType=POINT;
  sprintf(mom,"%i %i %i", mom_quark.Mom().cmp(0), 
	  mom_quark.Mom().cmp(1), mom_quark.Mom().cmp(2) );
  Site s ;
  for(s.Begin();s.End();s.nextSite())
    {
      int t = s.physT() ;
      Q1 =  quark[s.Index()] ;  
      if(HalfFerm) Q1.PParProjectSink() ;
      Q2 = Q1 ;
      // here is the sink so I need exp(ipx) 
      // I did conjugate the momentum at the construction
      Q3 = sink_mom.Fact(s)*mom_quark[s.Index()]  ; 
      if(HalfFerm) Q3.PParProjectSink() ;
      calc_nuc_(Q1, Q2, Q3,t) ;	  
    }
  plus_parity.GlobalSum() ;
  minus_parity.GlobalSum() ;
}

/*! 
  Computes
  \f[
   Trace(seqQ \, 
     \gamma_{\mu_1} \, \gamma_{\mu_2} \, \gamma_{\mu_3}\, 
     e^{-i mom . x}\, Quark)
  \f]
  where \f$\mu\f$ is the member int gamma and ThreeMom mom is
  the momentum injected at the operator. 
**/
void Nuc3ptGamma::InsertOp(CorrFunc& tmp,QPropW& seqQ, QPropW& Quark)
{
  char *fname = "InsertOp()";

  VRB.Func(cname, fname);

  Site s ;
  for(s.Begin();s.End();s.nextSite())
    {
      int t(s.physT()) ;
      WilsonMatrix sq(seqQ[s.Index()]) ;
      for(int mu(0);mu<G.N();mu++)
	sq.gr(G[mu]) ; // Multiply by gamma_mu 
      Complex cc(mom.Fact(s)) ; // The momentum factor exp(-ipx)
      cc *= Trace(sq,Quark[s.Index()]) ;
      tmp[t] += cc ;
    }
}

/*!
  Calculate zero momentum nucleon two point functions
*/
void Nuc2pt::calcNucleon(QPropW& quark, int dohalf){
  char *cname = "Nuc2pt" ;
  char *fname = "calcNucleon(QPropW&)" ;

  strcpy(mom,"0 0 0");  

  HalfFerm=dohalf;

  srcdesc[0]=srcdesc[1]=srcdesc[2]=srcdesc[3]=-1 ;
  srcdesc[4]=srcdesc[5]=-1 ; // TIZB
  switch (quark.SrcType())
    {
    case WALL:
      srcdesc[0] = quark.SourceTime() ;
      source="WALL" ;
      break ;
    case BOX:
      srcdesc[0] = quark.BoxSrcStart() ; 
      srcdesc[1] = quark.BoxSrcEnd() ;
      srcdesc[2] = quark.SourceTime() ;
      if( quark. BoxSrcUseXYZOffset() ){ // use QPropWArg.{x,y,z} as offset
	srcdesc[3] = quark.PointSrcX() ; 
	srcdesc[4] = quark.PointSrcX() ; 
	srcdesc[5] = quark.PointSrcX() ; 
      }
      source="BOX" ;
      break ;
    case POINT:
      srcdesc[0] = quark.PointSrcX() ; 
      srcdesc[1] = quark.PointSrcY() ;
      srcdesc[2] = quark.PointSrcZ() ;
      srcdesc[3] = quark.SourceTime() ;
      source = "POINT" ;
      break ;
    case GAUSS_GAUGE_INV:
      srcdesc[0] = quark.Gauss_N() ; 
      srcdesc[1] = int(100*quark.Gauss_W()) ;
      srcdesc[2] = quark.SourceTime() ;
      source = "GAUSS" ;
      break ;
    default:
      ERR.General(cname,fname,"Bad source type\n") ;
    }
  quark_mass = quark.Mass() ;
    
  WilsonMatrix Q1;  
  WilsonMatrix Q2;
  WilsonMatrix Q3;
  WilsonMatrix Q4;
  WilsonMatrix Q5;
  
    switch (SnkType)
      {
      case POINT:
      case GAUSS_GAUGE_INV:	
	{
	  if(SnkType==GAUSS_GAUGE_INV)
	    sink="GAUSS" ;
	  else
	    sink="POINT" ;
	  Site s ;
	  for(s.Begin();s.End();s.nextSite())
	    {
	      int t = s.physT() ;
	      Q1 = quark[s.Index()] ;
	      if(HalfFerm) Q1.PParProjectSink() ;
	      Q2 = Q1 ; 
	      Q3 = Q1 ; 
	      Q4 = Q1 ; 
	      Q5 = Q1 ; 
	      calc_nuc_(Q1, Q2, Q3, Q4, Q5, t, dohalf) ;
	    }
	  plus_parity.GlobalSum() ;
	  minus_parity.GlobalSum() ;
	}
	break ;
      case WALL:
	{
	  sink="WALL" ;

	  if(!quark.GFixedSnk())
	    ERR.General(cname,fname,"Wall sink needs gauge fixing\n") ;

	  for(int t = 0 ; t<plus_parity.TimeSize();t++)
	    {
	      Q1 =  quark.WallSinkProp(t) ;
	      if(HalfFerm) Q1.PParProjectSink() ;
	      Q2 = Q1 ;
	      calc_nuc_(Q1, Q2, Q1, t) ;
	    }
	}
	break ;
      default:
	ERR.General(cname,fname,"Unknown Sink Type\n") ;
      } 
}

/*! 
  Computes
  \f[
   Trace(seqQ \, 
     \gamma_{\mu_1} \, \gamma_{\mu_2} \, \gamma_{\mu_3}\, 
     e^{-i mom . x}\, Quark)
  \f]
  where \f$\mu\f$ is the member int gamma and ThreeMom mom is
  the momentum injected at the operator. 
**/
void Nuc3ptCons::InsertOp(CorrFunc* tmp,QPropW& seqQ, QPropW& Quark, int Nmom, ThreeMom* momo)
{
  char *fname = "InsertOp()";

  VRB.Func(cname, fname);

  WilsonMatrix v0, v1, v0_next, v1_next;
  WilsonMatrix temp_v0D, temp_v0D_next;
  WilsonMatrix temp_v0_1pD, temp_v0_1mD_next;
  WilsonMatrix temp_v1_u;
  WilsonMatrix temp_v1_next_u;


  int shift_t=GJP.TnodeCoor()*GJP.TnodeSites();
  int vol= (GJP.VolNodeSites()/GJP.TnodeSites());

  int dir = G[G.N()-1];


  const int LORENTZs(4); 
  int lcl_sites[LORENTZs]; 
  lcl_sites[0] = GJP.XnodeSites();
  lcl_sites[1] = GJP.YnodeSites();
  lcl_sites[2] = GJP.ZnodeSites();
  lcl_sites[3] = GJP.TnodeSites();

  int lcl_node[LORENTZs];
  lcl_node[0] = GJP.XnodeCoor();
  lcl_node[1] = GJP.YnodeCoor();
  lcl_node[2] = GJP.ZnodeCoor();
  lcl_node[3] = GJP.TnodeCoor();

  int lcl2glb_offset[LORENTZs];
  for (int i = 0; i < LORENTZs; ++i) 
    lcl2glb_offset[i] = lcl_sites[i] * lcl_node[i];

 
  Lattice& lat = seqQ.AlgLattice();
  Matrix *gauge_field = lat.GaugeField();  

  int ls_glb = GJP.SnodeSites() * GJP.Snodes();
  int ls_ini = 0;
  int ls_end = ls_glb;
  if(GJP.Snodes()==2) {
    ls_ini = GJP.SnodeCoor() * GJP.SnodeSites();
    ls_end = ( GJP.SnodeCoor() + 1 ) * GJP.SnodeSites();
  }

  char dummy[50];
  Site site;

  for(int s=ls_ini;s<ls_end;s++){ 
    //for(int s=0;s<ls_glb;s++){ 
    Quark.RestoreQPropLs(dummy, s);
    if(GJP.Snodes()==2) {
    seqQ.RestoreQPropLs(dummy, ls_end-1-s+ls_ini);
    } else {
    seqQ.RestoreQPropLs(dummy, ls_end-1-s);
    }
    //seqQ.RestoreQPropLs(dummy, ls_glb-1-s);

    // This operations are included in QPropWSeqProtU(D)Src.
    // When we restore 5d prop the operations are nessesary.
    //Multiply by gamma5 and take the dagger to make it in to quark.
    for (site.Begin();site.End();site.nextSite()) {
      seqQ[site.Index()].gl(-5);
      seqQ[site.Index()].hconj();
    }


    int lcl[LORENTZs];
    int lcl_next[LORENTZs]; // Next site along propagation direction


    VRB.Clock(cname,fname,"RO\n");
    for(lcl[3]=0; lcl[3]<=GJP.TnodeSites()-1; lcl[3]++) {
    int t=lcl[3]+shift_t;

    for(lcl[0]=0; lcl[0]<=GJP.XnodeSites()-1; lcl[0]++) 
    for(lcl[1]=0; lcl[1]<=GJP.YnodeSites()-1; lcl[1]++) 
    for(lcl[2]=0; lcl[2]<=GJP.ZnodeSites()-1; lcl[2]++) {

      Float coeff = 1.0;
      if( G.N() == 2 ) {
      coeff = -1.0;
      if(GJP.Snodes()==2) {
	if(GJP.SnodeCoor()==1) coeff = -coeff;
      } else {
	if( s > ls_glb/2 - 1 ) coeff = -coeff;
      }
//.........这里部分代码省略.........

/*! 
  Computes
  \f[
   \mbox{Trace}(\mbox{seqQ} \, 
     \gamma_{\mu_1} \, \gamma_{\mu_2} \, \gamma_{\mu_3}\, 
     \stackrel{\displaystyle \leftrightarrow}{D}_\nu 
     \stackrel{\displaystyle \leftrightarrow}{D}_\rho ...
     e^{-i mom . x}\, \mbox{Quark})
  \f]
  where \f$\mu\f$ is the member int gamma and ThreeMom mom is
  the momentum injected at the operator.

  The string of derivatives gets broken up in strings of forward
  and backward hops. For each term the quark and anti-quark are sitting
  on different sites and they are connected with the gauge connection U.
  The gauge connection U is calculated using the PathOrdPlus routine
  given the location of the anti-quark (starting point) and a list
  of hops. The Derivative class is used not only  to do the break-up
  of the string of derivatives in strings of hops, but also to 
  compute the locations of the quark and the anti-quark and 
  provide the list of links needed to compute the gauge connection U
  for each term.

  In order to speed up the code first we compute for all sites
  all derivative terms with the derivatives applied on the quark i.e
  the anti-quark is not displaced. 
  These terms are not summed up but they are stored individually in the
  DTerms Terms class. 
  Then I we go again through the all sites and construct the derivatives.
  If for a given term the anti-quark is displaced we just have to communicate
  a single Complex number from that site. The DTerms class handles all that.
  The routine Derivative :: DTermIndx() is the mapping of the derivative
  terms to a particular storage order. The DTerm class does not have
  to know the storage order. Only the Derivative class knows about it.
 
 NOTE: The quark is at the sink of the Quark propagator
        while the anti-quark is at to source of the sequential propagator.
	
**/
void Nuc3ptStru::InsertOp(CorrFunc& tmp,QPropW& seqQ, QPropW& Quark)
{
  char *fname = "InsertOp()";

  VRB.Func(cname, fname);
  
  if(seqQ.GFixedSnk() || Quark.GFixedSnk())
    ERR.General(cname,fname,"Gauge Fixed sinks!\n") ;

  Lattice& lat=seqQ.AlgLattice() ; // The lattice

#ifdef DEBUG  
printf("%s:: NDTerms= %i\n",cname,D.NDTerms());
#endif

  DTerms Terms(D.NDTerms()); // The derivative terms
  Matrix U ;
  WilsonMatrix aQ ; // the anti-quark
  WilsonMatrix  Q ; // quark
  WilsonMatrix aQU ; // anti-quark times Gauge path
  WilsonMatrix buff ; // Buffer needed for communication 
  int aq_vec[4]; // Anti-quark location
  int q_vec[4] ; // quark location
  int *dirs = new int[D.N()];
  if(dirs == 0) ERR.Pointer(cname, fname, "dirs");

  Site s ;
  //first compute all derivatives terms on quarks
  for(s.Begin();s.End();s.nextSite()) // Loop over sites
    {
      aQ=seqQ[s.Index()] ;
      // Apply the gamma matrices on the sequential  propagator
      // We multiply the source indices of the sequential propagator
      // which are the anti-quark
      for(int mu(0);mu<G.N();mu++)
	aQ.gr(G[mu]) ; // Multiply by gamma_mu 

      // Loop over all derivative terms on the Quark
      for(D.Start();D.NotEnd();D.NextQuark()) 
	{
	  // Calculate the endpoints of the shift
	  D.CalcEndPoints(aq_vec,q_vec,s);
	  // The sink indices of the quark propagator are the quark
	  Q=Quark.GetMatrix(q_vec,buff) ;
	  
	 
	  D.CalcListOfLinks(dirs);

#ifdef DEBUG
printf("\n%s:: dirs ",cname);
for(int i(0);i<D.N();i++) printf("%i ",dirs[i]);
printf(" aq_vec: ");
for(int i(0);i<4;i++) printf("%i ",aq_vec[i]);
printf(" q_vec: ") ; 
for(int i(0);i<4;i++) printf("%i ",q_vec[i]);
#endif

	  // Start from the anti-quark and move along the path described
	  // by the list of dirs multiplying all the encountered links
	  // at the end of the day you should be at the position of the quark
          if(D.N()>1)  
//.........这里部分代码省略.........