void sf_gauge_derivative(const int id, hamiltonian_field_t * const hf) { int i, mu; static su3 v, w; su3 *z; su3adj *xm; monomial * mnl = &monomial_list[id]; double factor = -1. * g_beta/3.0; if(mnl->use_rectangles) { factor = -mnl->c0 * g_beta/3.0; } for(i = 0; i < VOLUME; i++) { for(mu=0;mu<4;mu++) { z=&hf->gaugefield[i][mu]; xm=&hf->derivative[i][mu]; get_staples(&v,i,mu, (const su3**) hf->gaugefield); _su3_times_su3d(w,*z,v); _trace_lambda_mul_add_assign((*xm), factor, w); if(mnl->use_rectangles) { get_rectangle_staples(&v, i, mu); _su3_times_su3d(w, *z, v); _trace_lambda_mul_add_assign((*xm), factor*mnl->c1/mnl->c0, w); } } } return; }
/* this function calculates the derivative of the momenta: equation 13 of Gottlieb */ void gauge_derivative(const int id, hamiltonian_field_t * const hf) { monomial * mnl = &monomial_list[id]; double factor = -1. * g_beta/3.0; if(mnl->use_rectangles) { mnl->forcefactor = 1.; factor = -mnl->c0 * g_beta/3.0; } double atime, etime; atime = gettime(); #ifdef OMP #pragma omp parallel { #endif su3 ALIGN v, w; int i, mu; su3 *z; su3adj *xm; #ifdef OMP #pragma omp for #endif for(i = 0; i < VOLUME; i++) { for(mu=0;mu<4;mu++) { z=&hf->gaugefield[i][mu]; xm=&hf->derivative[i][mu]; get_staples(&v,i,mu, (const su3**) hf->gaugefield); _su3_times_su3d(w,*z,v); _trace_lambda_mul_add_assign((*xm), factor, w); if(mnl->use_rectangles) { get_rectangle_staples(&v, i, mu); _su3_times_su3d(w, *z, v); _trace_lambda_mul_add_assign((*xm), factor*mnl->c1/mnl->c0, w); } } } #ifdef OMP } /* OpenMP closing brace */ #endif etime = gettime(); if(g_debug_level > 1 && g_proc_id == 0) { printf("# Time for %s monomial derivative: %e s\n", mnl->name, etime-atime); } return; }
void sw_all(hamiltonian_field_t * const hf, const double kappa, const double c_sw) { int k,l; int x,xpk,xpl,xmk,xml,xpkml,xplmk,xmkml; su3 *w1,*w2,*w3,*w4; double ka_csw_8 = kappa*c_sw/8.; static su3 v1,v2,vv1,vv2,plaq; static su3 vis[4][4]; for(x = 0; x < VOLUME; x++) { _minus_itimes_su3_plus_su3(vis[0][1],swm[x][1],swm[x][3]); _su3_minus_su3(vis[0][2],swm[x][1],swm[x][3]); _itimes_su3_minus_su3(vis[0][3],swm[x][2],swm[x][0]); _minus_itimes_su3_plus_su3(vis[2][3],swp[x][1],swp[x][3]); _su3_minus_su3(vis[1][3],swp[x][3],swp[x][1]); _itimes_su3_minus_su3(vis[1][2],swp[x][2],swp[x][0]); // project to the traceless anti-hermitian part _su3_dagger(v1,vis[0][1]); _su3_minus_su3(vis[0][1],vis[0][1],v1); _su3_dagger(v1,vis[0][2]); _su3_minus_su3(vis[0][2],vis[0][2],v1); _su3_dagger(v1,vis[0][3]); _su3_minus_su3(vis[0][3],vis[0][3],v1); _su3_dagger(v1,vis[2][3]); _su3_minus_su3(vis[2][3],vis[2][3],v1); _su3_dagger(v1,vis[1][3]); _su3_minus_su3(vis[1][3],vis[1][3],v1); _su3_dagger(v1,vis[1][2]); _su3_minus_su3(vis[1][2],vis[1][2],v1); for(k = 0; k < 4; k++) { for(l = k+1; l < 4; l++) { xpk=g_iup[x][k]; xpl=g_iup[x][l]; xmk=g_idn[x][k]; xml=g_idn[x][l]; xpkml=g_idn[xpk][l]; xplmk=g_idn[xpl][k]; xmkml=g_idn[xml][k]; w1=&hf->gaugefield[x][k]; w2=&hf->gaugefield[xpk][l]; w3=&hf->gaugefield[xpl][k]; /*dag*/ w4=&hf->gaugefield[x][l]; /*dag*/ _su3_times_su3(v1,*w1,*w2); _su3_times_su3(v2,*w4,*w3); _su3_times_su3d(plaq,v1,v2); _su3_times_su3(vv1,plaq,vis[k][l]); _trace_lambda_mul_add_assign(hf->derivative[x][k], -2.*ka_csw_8, vv1); _su3d_times_su3(vv2,*w1,vv1); _su3_times_su3(vv1,vv2,*w1); _trace_lambda_mul_add_assign(hf->derivative[xpk][l], -2.*ka_csw_8, vv1); _su3_times_su3(vv2,vis[k][l],plaq); _su3_dagger(vv1,vv2); _trace_lambda_mul_add_assign(hf->derivative[x][l], -2.*ka_csw_8, vv1); _su3d_times_su3(vv2,*w4,vv1); _su3_times_su3(vv1,vv2,*w4); _trace_lambda_mul_add_assign(hf->derivative[xpl][k], -2.*ka_csw_8, vv1); w1=&hf->gaugefield[x][l]; w2=&hf->gaugefield[xplmk][k]; /*dag*/ w3=&hf->gaugefield[xmk][l]; /*dag*/ w4=&hf->gaugefield[xmk][k]; _su3_times_su3d(v1,*w1,*w2); _su3d_times_su3(v2,*w3,*w4); _su3_times_su3(plaq,v1,v2); _su3_times_su3(vv1,plaq,vis[k][l]); _trace_lambda_mul_add_assign(hf->derivative[x][l], -2.*ka_csw_8, vv1); _su3_dagger(vv1,v1); _su3_times_su3d(vv2,vv1,vis[k][l]); _su3_times_su3d(vv1,vv2,v2); _trace_lambda_mul_add_assign(hf->derivative[xplmk][k], -2.*ka_csw_8, vv1); _su3_times_su3(vv2,*w3,vv1); _su3_times_su3d(vv1,vv2,*w3); _trace_lambda_mul_add_assign(hf->derivative[xmk][l], -2.*ka_csw_8, vv1); _su3_dagger(vv2,vv1); _trace_lambda_mul_add_assign(hf->derivative[xmk][k], -2.*ka_csw_8, vv2); w1=&hf->gaugefield[xmk][k]; /*dag*/ w2=&hf->gaugefield[xmkml][l]; /*dag*/ w3=&hf->gaugefield[xmkml][k]; w4=&hf->gaugefield[xml][l]; _su3_times_su3(v1,*w2,*w1); _su3_times_su3(v2,*w3,*w4); _su3_times_su3d(vv1,*w1,vis[k][l]); _su3_times_su3d(vv2,vv1,v2); _su3_times_su3(vv1,vv2,*w2); _trace_lambda_mul_add_assign(hf->derivative[xmk][k], -2.*ka_csw_8, vv1); _su3_times_su3(vv2,*w2,vv1); _su3_times_su3d(vv1,vv2,*w2); _trace_lambda_mul_add_assign(hf->derivative[xmkml][l], -2.*ka_csw_8, vv1); _su3_dagger(vv2,vv1); _trace_lambda_mul_add_assign(hf->derivative[xmkml][k], -2.*ka_csw_8, vv2); _su3d_times_su3(vv1,*w3,vv2); _su3_times_su3(vv2,vv1,*w3); _trace_lambda_mul_add_assign(hf->derivative[xml][l], -2.*ka_csw_8, vv2); w1=&hf->gaugefield[xml][l]; /*dag*/ w2=&hf->gaugefield[xml][k]; w3=&hf->gaugefield[xpkml][l]; w4=&hf->gaugefield[x][k]; /*dag*/ _su3d_times_su3(v1,*w1,*w2); _su3_times_su3d(v2,*w3,*w4); _su3_times_su3d(vv1,*w1,vis[k][l]); _su3_times_su3d(vv2,vv1,v2); _su3_times_su3d(vv1,vv2,*w2); _trace_lambda_mul_add_assign(hf->derivative[xml][l], -2.*ka_csw_8, vv1); _su3_dagger(vv2,vv1); _trace_lambda_mul_add_assign(hf->derivative[xml][k], -2.*ka_csw_8, vv2); _su3d_times_su3(vv1,*w2,vv2); _su3_times_su3(vv2,vv1,*w2); _trace_lambda_mul_add_assign(hf->derivative[xpkml][l], -2.*ka_csw_8, vv2); _su3_dagger(vv2,v2); _su3_times_su3d(vv1,vv2,v1); _su3_times_su3d(vv2,vv1,vis[k][l]); _trace_lambda_mul_add_assign(hf->derivative[x][k], -2.*ka_csw_8, vv2); } } } return; }
void deriv_Sb_D_psi(spinor * const l, spinor * const k, hamiltonian_field_t * const hf, const double factor) { int ix,iy; su3 * restrict up ALIGN; su3 * restrict um ALIGN; static su3 v1,v2; static su3_vector psia,psib,phia,phib; static spinor rr; /* spinor * restrict r ALIGN; */ spinor * restrict sp ALIGN; spinor * restrict sm ALIGN; #ifdef _KOJAK_INST #pragma pomp inst begin(derivSb) #endif #ifdef XLC #pragma disjoint(*sp, *sm, *up, *um) #endif #ifdef BGL __alignx(16, l); __alignx(16, k); #endif /* for parallelization */ #ifdef MPI xchange_lexicfield(k); xchange_lexicfield(l); #endif /************** loop over all lattice sites ****************/ for(ix = 0; ix < (VOLUME); ix++){ rr = (*(l + ix)); /* rr=g_spinor_field[l][icx-ioff]; */ /*multiply the left vector with gamma5*/ _vector_minus_assign(rr.s2, rr.s2); _vector_minus_assign(rr.s3, rr.s3); /*********************** direction +0 ********************/ iy=g_iup[ix][0]; sp = k + iy; up=&hf->gaugefield[ix][0]; _vector_add(psia,(*sp).s0,(*sp).s2); _vector_add(psib,(*sp).s1,(*sp).s3); _vector_add(phia,rr.s0,rr.s2); _vector_add(phib,rr.s1,rr.s3); _vector_tensor_vector_add(v1, phia, psia, phib, psib); _su3_times_su3d(v2,*up,v1); _complex_times_su3(v1,ka0,v2); _trace_lambda_mul_add_assign(hf->derivative[ix][0], 2.*factor, v1); /************** direction -0 ****************************/ iy=g_idn[ix][0]; sm = k + iy; um=&hf->gaugefield[iy][0]; _vector_sub(psia,(*sm).s0,(*sm).s2); _vector_sub(psib,(*sm).s1,(*sm).s3); _vector_sub(phia,rr.s0,rr.s2); _vector_sub(phib,rr.s1,rr.s3); _vector_tensor_vector_add(v1, psia, phia, psib, phib); _su3_times_su3d(v2,*um,v1); _complex_times_su3(v1,ka0,v2); _trace_lambda_mul_add_assign(hf->derivative[iy][0], 2.*factor, v1); /*************** direction +1 **************************/ iy=g_iup[ix][1]; sp = k + iy; up=&hf->gaugefield[ix][1]; _vector_i_add(psia,(*sp).s0,(*sp).s3); _vector_i_add(psib,(*sp).s1,(*sp).s2); _vector_i_add(phia,rr.s0,rr.s3); _vector_i_add(phib,rr.s1,rr.s2); _vector_tensor_vector_add(v1, phia, psia, phib, psib); _su3_times_su3d(v2,*up,v1); _complex_times_su3(v1,ka1,v2); _trace_lambda_mul_add_assign(hf->derivative[ix][1], 2.*factor, v1); /**************** direction -1 *************************/ iy=g_idn[ix][1]; sm = k + iy; um=&hf->gaugefield[iy][1]; _vector_i_sub(psia,(*sm).s0,(*sm).s3); _vector_i_sub(psib,(*sm).s1,(*sm).s2); _vector_i_sub(phia,rr.s0,rr.s3); _vector_i_sub(phib,rr.s1,rr.s2); _vector_tensor_vector_add(v1, psia, phia, psib, phib); _su3_times_su3d(v2,*um,v1); _complex_times_su3(v1,ka1,v2); _trace_lambda_mul_add_assign(hf->derivative[iy][1], 2.*factor, v1); /*************** direction +2 **************************/ iy=g_iup[ix][2]; sp = k + iy; up=&hf->gaugefield[ix][2]; _vector_add(psia,(*sp).s0,(*sp).s3); _vector_sub(psib,(*sp).s1,(*sp).s2); _vector_add(phia,rr.s0,rr.s3); _vector_sub(phib,rr.s1,rr.s2); _vector_tensor_vector_add(v1, phia, psia, phib, psib); _su3_times_su3d(v2,*up,v1); _complex_times_su3(v1,ka2,v2); _trace_lambda_mul_add_assign(hf->derivative[ix][2], 2.*factor, v1); /***************** direction -2 ************************/ iy=g_idn[ix][2]; sm = k + iy; um=&hf->gaugefield[iy][2]; _vector_sub(psia,(*sm).s0,(*sm).s3); _vector_add(psib,(*sm).s1,(*sm).s2); _vector_sub(phia,rr.s0,rr.s3); _vector_add(phib,rr.s1,rr.s2); _vector_tensor_vector_add(v1, psia, phia, psib, phib); _su3_times_su3d(v2,*um,v1); _complex_times_su3(v1,ka2,v2); _trace_lambda_mul_add_assign(hf->derivative[iy][2], 2.*factor, v1); /****************** direction +3 ***********************/ iy=g_iup[ix][3]; sp = k + iy; up=&hf->gaugefield[ix][3]; _vector_i_add(psia,(*sp).s0,(*sp).s2); _vector_i_sub(psib,(*sp).s1,(*sp).s3); _vector_i_add(phia,rr.s0,rr.s2); _vector_i_sub(phib,rr.s1,rr.s3); _vector_tensor_vector_add(v1, phia, psia, phib, psib); _su3_times_su3d(v2,*up,v1); _complex_times_su3(v1,ka3,v2); _trace_lambda_mul_add_assign(hf->derivative[ix][3], 2.*factor, v1); /***************** direction -3 ************************/ iy=g_idn[ix][3]; sm = k + iy; um=&hf->gaugefield[iy][3]; _vector_i_sub(psia,(*sm).s0,(*sm).s2); _vector_i_add(psib,(*sm).s1,(*sm).s3); _vector_i_sub(phia,rr.s0,rr.s2); _vector_i_add(phib,rr.s1,rr.s3); _vector_tensor_vector_add(v1, psia, phia, psib, phib); _su3_times_su3d(v2,*um,v1); _complex_times_su3(v1,ka3,v2); _trace_lambda_mul_add_assign(hf->derivative[iy][3], 2.*factor, v1); /****************** end of loop ************************/ } #ifdef _KOJAK_INST #pragma pomp inst end(derivSb) #endif }