Can some of the VecAXPY() calls in loops in TS be replaced with single VecMAXPY() calls? For example for (j=0; j<r; j++) { ierr = VecAXPY(gl->Z,-shift*u[i*r+j],X[j]);CHKERRQ(ierr); } Just asking, Barry