static #endif void exponent_from_coefficients(su3 *out, _Complex double f0, _Complex double f1, _Complex double f2, su3 const *in) { su3 ALIGN tmp; _complex_times_su3(tmp, f2, *in); _su3_add_equals_complex_identity(tmp, f1); _su3_times_su3(*out, tmp, *in); _su3_add_equals_complex_identity(*out, f0); }
void deriv_Sb_D_psi(spinor * const l, spinor * const k, hamiltonian_field_t * const hf, const double factor) { #ifdef BGL __alignx(16, l); __alignx(16, k); #endif /* for parallelization */ #ifdef MPI xchange_lexicfield(k); xchange_lexicfield(l); #endif #ifdef OMP #define static #pragma omp parallel { #endif 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 OMP #undef static #endif #ifdef _KOJAK_INST #pragma pomp inst begin(derivSb) #endif #ifdef XLC #pragma disjoint(*sp, *sm, *up, *um) #endif /************** loop over all lattice sites ****************/ #ifdef OMP #pragma omp for #endif 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_nonlocal(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_nonlocal(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_nonlocal(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_nonlocal(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_nonlocal(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_nonlocal(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_nonlocal(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_nonlocal(hf->derivative[iy][3], 2.*factor, v1); /****************** end of loop ************************/ } #ifdef _KOJAK_INST #pragma pomp inst end(derivSb) #endif #ifdef OMP } /*OpenMP closing brace */ #endif }
void deriv_Sb_D_psi(spinor * const l, spinor * const k) { /* const int l, const int k){ */ int ix,iy; su3 * restrict up ALIGN; su3 * restrict um ALIGN; /* su3adj * restrict ddd; */ /* static su3adj der; */ 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; /* sp=&g_spinor_field[k][icy]; */ up=&g_gauge_field[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); /* _vector_tensor_vector(v1,phia,psia); */ /* _vector_tensor_vector(v2,phib,psib); */ /* _su3_plus_su3(v1,v1,v2); */ _su3_times_su3d(v2,*up,v1); _complex_times_su3(v1,ka0,v2); _trace_lambda_add_assign(df0[ix][0], v1); /* _trace_lambda(der,v1); */ /* ddd=&df0[ix][0]; */ /* _add_su3adj(*ddd,der); */ /************** direction -0 ****************************/ iy=g_idn[ix][0]; sm = k + iy; /* sm=&g_spinor_field[k][icy]; */ um=&g_gauge_field[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); /* _vector_tensor_vector(v1,psia,phia); */ /* _vector_tensor_vector(v2,psib,phib); */ /* _su3_plus_su3(v1,v1,v2); */ _su3_times_su3d(v2,*um,v1); _complex_times_su3(v1,ka0,v2); _trace_lambda_add_assign(df0[iy][0], v1); /* _trace_lambda(der,v1); */ /* ddd=&df0[iy][0]; */ /* _add_su3adj(*ddd,der); */ /*************** direction +1 **************************/ iy=g_iup[ix][1]; sp = k + iy; /* sp=&g_spinor_field[k][icy]; */ up=&g_gauge_field[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); /* _vector_tensor_vector(v1,phia,psia); */ /* _vector_tensor_vector(v2,phib,psib); */ /* _su3_plus_su3(v1,v1,v2); */ _su3_times_su3d(v2,*up,v1); _complex_times_su3(v1,ka1,v2); _trace_lambda_add_assign(df0[ix][1], v1); /* _trace_lambda(der,v1); */ /* ddd=&df0[ix][1]; */ /* _add_su3adj(*ddd,der); */ /**************** direction -1 *************************/ iy=g_idn[ix][1]; sm = k + iy; /* sm=&g_spinor_field[k][icy]; */ um=&g_gauge_field[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); /* _vector_tensor_vector(v1,psia,phia); */ /* _vector_tensor_vector(v2,psib,phib); */ /* _su3_plus_su3(v1,v1,v2); */ _su3_times_su3d(v2,*um,v1); _complex_times_su3(v1,ka1,v2); _trace_lambda_add_assign(df0[iy][1], v1); /* _trace_lambda(der,v1); */ /* ddd=&df0[iy][1]; */ /* _add_su3adj(*ddd,der); */ /*************** direction +2 **************************/ iy=g_iup[ix][2]; sp = k + iy; /* sp=&g_spinor_field[k][icy]; */ up=&g_gauge_field[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); /* _vector_tensor_vector(v1,phia,psia); */ /* _vector_tensor_vector(v2,phib,psib); */ /* _su3_plus_su3(v1,v1,v2); */ _su3_times_su3d(v2,*up,v1); _complex_times_su3(v1,ka2,v2); _trace_lambda_add_assign(df0[ix][2], v1); /* _trace_lambda(der,v1); */ /* ddd=&df0[ix][2]; */ /* _add_su3adj(*ddd,der); */ /***************** direction -2 ************************/ iy=g_idn[ix][2]; sm = k + iy; /* sm=&g_spinor_field[k][icy]; */ um=&g_gauge_field[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); /* _vector_tensor_vector(v1,psia,phia); */ /* _vector_tensor_vector(v2,psib,phib); */ /* _su3_plus_su3(v1,v1,v2); */ _su3_times_su3d(v2,*um,v1); _complex_times_su3(v1,ka2,v2); _trace_lambda_add_assign(df0[iy][2], v1); /* _trace_lambda(der,v1); */ /* ddd=&df0[iy][2]; */ /* _add_su3adj(*ddd,der); */ /****************** direction +3 ***********************/ iy=g_iup[ix][3]; sp = k + iy; /* sp=&g_spinor_field[k][icy]; */ up=&g_gauge_field[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); /* _vector_tensor_vector(v1,phia,psia); */ /* _vector_tensor_vector(v2,phib,psib); */ /* _su3_plus_su3(v1,v1,v2); */ _su3_times_su3d(v2,*up,v1); _complex_times_su3(v1,ka3,v2); _trace_lambda_add_assign(df0[ix][3], v1); /* _trace_lambda(der,v1); */ /* ddd=&df0[ix][3]; */ /* _add_su3adj(*ddd,der); */ /***************** direction -3 ************************/ iy=g_idn[ix][3]; sm = k + iy; /* sm=&g_spinor_field[k][icy]; */ um=&g_gauge_field[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); /* _vector_tensor_vector(v1,psia,phia); */ /* _vector_tensor_vector(v2,psib,phib); */ /* _su3_plus_su3(v1,v1,v2); */ _su3_times_su3d(v2,*um,v1); _complex_times_su3(v1,ka3,v2); _trace_lambda_add_assign(df0[iy][3], v1); /* _trace_lambda(der,v1); */ /* ddd=&df0[iy][3]; */ /* _add_su3adj(*ddd,der); */ /****************** end of loop ************************/ } #ifdef _KOJAK_INST #pragma pomp inst end(derivSb) #endif }