/* Serially Checked ! */ void Dtm_psi(spinor * const P, spinor * const Q){ if(P==Q){ printf("Error in Dtm_psi (D_psi.c):\n"); printf("Arguments must be differen spinor fields\n"); printf("Program aborted\n"); exit(1); } #ifdef _GAUGE_COPY2 if(g_update_gauge_copy) { update_backward_gauge(g_gauge_field); } #endif # if defined TM_USE_MPI xchange_lexicfield(Q); # endif #ifdef TM_USE_OMP #pragma omp parallel { #endif int ix,iy,iz; su3 *up,*um; spinor *s,*sp,*sm,*rn; _Complex double fact1, fact2; spinor rs __attribute__ ((aligned (16))); fact1 = 1. + g_mu * I; fact2 = conj(fact1); #ifndef TM_USE_OMP iy=g_iup[0][0]; sp=(spinor *) Q + iy; up=&g_gauge_field[0][0]; #endif /************************ loop over all lattice sites *************************/ #ifdef TM_USE_OMP #pragma omp for #endif for (ix=0;ix<VOLUME;ix++){ #ifdef TM_USE_OMP iy=g_iup[ix][0]; up=&g_gauge_field[ix][0]; sp=(spinor *) Q + iy; #endif s=(spinor *) Q + ix; _prefetch_spinor(s); /******************************* direction +0 *********************************/ iy=g_idn[ix][0]; sm = (spinor *) Q + iy; _prefetch_spinor(sm); _sse_load(sp->s0); _sse_load_up(sp->s2); _sse_vector_add(); _sse_su3_multiply((*up)); _sse_vector_cmplx_mul(phase_0); _sse_store_up(rs.s2); // the diagonal bit _sse_load_up(s->s0); _sse_vector_cmplx_mul(fact1); _sse_load(rs.s2); _sse_vector_add(); _sse_store(rs.s0); // g5 in the twisted term _sse_load_up(s->s2); _sse_vector_cmplx_mul(fact2); _sse_load(rs.s2); _sse_vector_add(); _sse_store(rs.s2); um=&g_gauge_field[iy][0]; _prefetch_su3(um); _sse_load(sp->s1); _sse_load_up(sp->s3); _sse_vector_add(); _sse_su3_multiply((*up)); _sse_vector_cmplx_mul(phase_0); _sse_store_up(rs.s3); // the diagonal bit _sse_load_up(s->s1); _sse_vector_cmplx_mul(fact1); _sse_load(rs.s3); _sse_vector_add(); _sse_store(rs.s1); // g5 in the twisted term _sse_load_up(s->s3); _sse_vector_cmplx_mul(fact2); _sse_load(rs.s3); _sse_vector_add(); _sse_store(rs.s3); /******************************* direction -0 *********************************/ iy=g_iup[ix][1]; sp = (spinor *) Q + iy; _prefetch_spinor(sp); _sse_load(sm->s0); _sse_load_up(sm->s2); _sse_vector_sub(); _sse_su3_inverse_multiply((*um)); _sse_vector_cmplxcg_mul(phase_0); _sse_load(rs.s0); _sse_vector_add(); _sse_store(rs.s0); _sse_load(rs.s2); _sse_vector_sub(); _sse_store(rs.s2); up+=1; _prefetch_su3(up); _sse_load(sm->s1); _sse_load_up(sm->s3); _sse_vector_sub(); _sse_su3_inverse_multiply((*um)); _sse_vector_cmplxcg_mul(phase_0); _sse_load(rs.s1); _sse_vector_add(); _sse_store(rs.s1); _sse_load(rs.s3); _sse_vector_sub(); _sse_store(rs.s3); /******************************* direction +1 *********************************/ iy=g_idn[ix][1]; sm = (spinor *) Q + iy; _prefetch_spinor(sm); _sse_load(sp->s0); _sse_load_up(sp->s3); _sse_vector_i_mul(); _sse_vector_add(); _sse_su3_multiply((*up)); _sse_vector_cmplx_mul(phase_1); _sse_load(rs.s0); _sse_vector_add(); _sse_store(rs.s0); _sse_load(rs.s3); _sse_vector_i_mul(); _sse_vector_sub(); _sse_store(rs.s3); um=&g_gauge_field[iy][1]; _prefetch_su3(um); _sse_load(sp->s1); _sse_load_up(sp->s2); _sse_vector_i_mul(); _sse_vector_add(); _sse_su3_multiply((*up)); _sse_vector_cmplx_mul(phase_1); _sse_load(rs.s1); _sse_vector_add(); _sse_store(rs.s1); _sse_load(rs.s2); _sse_vector_i_mul(); _sse_vector_sub(); _sse_store(rs.s2); /******************************* direction -1 *********************************/ iy=g_iup[ix][2]; sp = (spinor *) Q + iy; _prefetch_spinor(sp); _sse_load(sm->s0); _sse_load_up(sm->s3); _sse_vector_i_mul(); _sse_vector_sub(); _sse_su3_inverse_multiply((*um)); _sse_vector_cmplxcg_mul(phase_1); _sse_load(rs.s0); _sse_vector_add(); _sse_store(rs.s0); _sse_load(rs.s3); _sse_vector_i_mul(); _sse_vector_add(); _sse_store(rs.s3); up+=1; _prefetch_su3(up); _sse_load(sm->s1); _sse_load_up(sm->s2); _sse_vector_i_mul(); _sse_vector_sub(); _sse_su3_inverse_multiply((*um)); _sse_vector_cmplxcg_mul(phase_1); _sse_load(rs.s1); _sse_vector_add(); _sse_store(rs.s1); _sse_load(rs.s2); _sse_vector_i_mul(); _sse_vector_add(); _sse_store(rs.s2); /******************************* direction +2 *********************************/ iy=g_idn[ix][2]; sm = (spinor *) Q + iy; _prefetch_spinor(sm); _sse_load(sp->s0); _sse_load_up(sp->s3); _sse_vector_add(); _sse_su3_multiply((*up)); _sse_vector_cmplx_mul(phase_2); _sse_load(rs.s0); _sse_vector_add(); _sse_store(rs.s0); _sse_load(rs.s3); _sse_vector_add(); _sse_store(rs.s3); um=&g_gauge_field[iy][2]; _prefetch_su3(um); _sse_load(sp->s1); _sse_load_up(sp->s2); _sse_vector_sub(); _sse_su3_multiply((*up)); _sse_vector_cmplx_mul(phase_2); _sse_load(rs.s1); _sse_vector_add(); _sse_store(rs.s1); _sse_load(rs.s2); _sse_vector_sub(); _sse_store(rs.s2); /******************************* direction -2 *********************************/ iy=g_iup[ix][3]; sp = (spinor *) Q + iy; _prefetch_spinor(sp); _sse_load(sm->s0); _sse_load_up(sm->s3); _sse_vector_sub(); _sse_su3_inverse_multiply((*um)); _sse_vector_cmplxcg_mul(phase_2); _sse_load(rs.s0); _sse_vector_add(); _sse_store(rs.s0); _sse_load(rs.s3); _sse_vector_sub(); _sse_store(rs.s3); up+=1; _prefetch_su3(up); _sse_load(sm->s1); _sse_load_up(sm->s2); _sse_vector_add(); _sse_su3_inverse_multiply((*um)); _sse_vector_cmplxcg_mul(phase_2); _sse_load(rs.s1); _sse_vector_add(); _sse_store(rs.s1); _sse_load(rs.s2); _sse_vector_add(); _sse_store(rs.s2); /******************************* direction +3 *********************************/ iy=g_idn[ix][3]; sm = (spinor *) Q + iy; _prefetch_spinor(sm); _sse_load(sp->s0); _sse_load_up(sp->s2); _sse_vector_i_mul(); _sse_vector_add(); _sse_su3_multiply((*up)); _sse_vector_cmplx_mul(phase_3); _sse_load(rs.s0); _sse_vector_add(); _sse_store(rs.s0); _sse_load(rs.s2); _sse_vector_i_mul(); _sse_vector_sub(); _sse_store(rs.s2); um=&g_gauge_field[iy][3]; _prefetch_su3(um); _sse_load(sp->s1); _sse_load_up(sp->s3); _sse_vector_i_mul(); _sse_vector_sub(); _sse_su3_multiply((*up)); _sse_vector_cmplx_mul(phase_3); _sse_load(rs.s1); _sse_vector_add(); _sse_store(rs.s1); _sse_load(rs.s3); _sse_vector_i_mul(); _sse_vector_add(); _sse_store(rs.s3); /******************************* direction -3 *********************************/ iz=(ix+1+VOLUME)%VOLUME; iy=g_iup[iz][0]; sp = (spinor *) Q + iy; _prefetch_spinor(sp); _sse_load(sm->s0); _sse_load_up(sm->s2); _sse_vector_i_mul(); _sse_vector_sub(); _sse_su3_inverse_multiply((*um)); _sse_vector_cmplxcg_mul(phase_3); rn = (spinor *) P + ix; _sse_load(rs.s0); _sse_vector_add(); _sse_store_nt(rn->s0); _sse_load(rs.s2); _sse_vector_i_mul(); _sse_vector_add(); _sse_store_nt(rn->s2); up=&g_gauge_field[iz][0]; _prefetch_su3(up); _sse_load(sm->s1); _sse_load_up(sm->s3); _sse_vector_i_mul(); _sse_vector_add(); _sse_su3_inverse_multiply((*um)); _sse_vector_cmplxcg_mul(phase_3); _sse_load(rs.s1); _sse_vector_add(); _sse_store_nt(rn->s1); _sse_load(rs.s3); _sse_vector_i_mul(); _sse_vector_sub(); _sse_store_nt(rn->s3); /******************************** end of loop *********************************/ } #ifdef TM_USE_OMP } /* OpenMP closing brace */ #endif }
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 }