1 / 27

A Dimension Abstraction Approach to Vectorization in Matlab

A Dimension Abstraction Approach to Vectorization in Matlab. Neil Birkbeck Jonathan Levesque Jose Nelson Amaral. Computing Science University of Alberta Edmonton, Alberta, Canada. Problem. Problem Statement:

Download Presentation

A Dimension Abstraction Approach to Vectorization in Matlab

An Image/Link below is provided (as is) to download presentation Download Policy: Content on the Website is provided to you AS IS for your information and personal use and may not be sold / licensed / shared on other websites without getting consent from its author. Content is provided to you AS IS for your information and personal use only. Download presentation by click this link. While downloading, if for some reason you are not able to download a presentation, the publisher may have deleted the file from their server. During download, if you can't get a presentation, the file might be deleted by the publisher.

E N D

Presentation Transcript


  1. A Dimension Abstraction Approach to Vectorization in Matlab Neil Birkbeck Jonathan Levesque Jose Nelson Amaral Computing ScienceUniversity of AlbertaEdmonton, Alberta, Canada

  2. Problem • Problem Statement: • Generate equivalent, error-free vectorized source code for Matlab source while utilizing higher level matrix operations when possible to improve efficiency.

  3. Motivation • Loop-based code is slower than vector code in Matlab. • Why? • interpretive overhead (type/shape checking,…) • resizing of arrays in loops n=1000; for i=1:n, A(i)=B(i)+C(i); end 5x faster! n=1000; A(1:n)=B(1:n)+C(1:n); • Vectorization also useful for compiled Matlab code, where optimized vector routines could be substituted.

  4. Related Work • Data dependence vectorization • Allen & Kennedy’s Codegen algorithm • Build data dependence graph • Topological visit strongly connected components • Abstract Matrix Form (AMF) [Menon & Pingali] • axioms used to transform array code • take advantage of matrix multiplication • Not clear if it is easily extensible or allows for vectorization of irregular access (e.g., access to the diagonal)

  5. If this is not true the vectorized code will introduce an error! Incorrect Vectorization • Example 1: Pull out of loop. Index variable substitution (i1:n) for i=1:n, a(i)=b(i)+c(i); end a(1:n)=b(1:n)+c(1:n) • Vectorization correct if a,b, and c are row vectors or column vectors

  6. Incorrect Vectorization for i=1:n, x(i)=y(i,h)*z(h,i); end • Example 2: • Matlab is untyped • Vectorization depends on whether h is a vector or scalar. • If h is a scalar: • Otherwise: x(1:n)=y(1:n,h).*z(h,1:n)’; x(1:n)=sum(y(1:n,h).*z(h,1:n)’,2);

  7. Overview of Solution Data dependence-based vectorizer Vectorizable statement Knowledge of Shape of variables Propagate dimensionality up parse tree Dimensions Agree? Yes Perform Transformations No Leave statement in loop Output Vector statement

  8. More Specifically Examples: • Represent dimensionality of expressions as list of symbols 1 or “*” (>1) • Assume known for variables. • Propagate up parse tree according to Matlab rules • Compatibility: • dim(A)≈dim(B) when the lists are equivalent (after removal of redundant 1’s)

  9. Vectorized Dimensionality • Vectorized dimensionality: • representation of dimensions after vectorization of a loop • denoted dimi for loop with index variable i • Introduce new symbol ri for index variable i for i=1:n, a(i)=10+i; end

  10. Vectorized Dimensionality • Expressions with incompatible vectorized dimensionality should not be vectorized. • When do dimensionalities agree? • Assignment expressions: elhs=erhs • dimi(elhs)≈dimi(erhs) || erhs≈(1) • Element-wise binary operators: e=elhsΘerhs • dimi(elhs) ≈(1)||dimi(erhs)≈(1)||dimi(elhs)≈dim(erhs) Θ in {+,-,.*,…} Dimensionalities are compatible or erhs is scalar Dimensionalities are compatible or either is scalar

  11. Vectorized Dimensionality • Rules very restrictive: • Assume dim(A)=dim(B)=dim(C)=(*,*) for i=1:100, for j=1:100 A(i,j)=B(j,i)+C(i,j); end end dimi,j(B)=(rj,ri) dimi,j(C)=(ri,rj) Vectorization fails because (ri,rj) is not compatible with (rj,ri)

  12. Transpose Transformation • Extension to utilize transpose when necessary is straightforward: • For assignment: • if dimi(A)≈reverse(dimi(B)) then A=BT is allowable for i=1:m, for j=1:n A(i,j)=B(j,i); end end dimi,j(A)=reverse(dimi,j(B))=(ri,rj) A(1:m,1:n)=(B(1:n,1:m))’

  13. Transpose Transformation • Extension to utilize transpose when necessary is straightforward: • Similar for pointwise operations: • if dimi(A)≈reverse(dimi(B)) then AΘBT is allowable, propagate dimi(AΘBT)=dimi(A) • if dimi(reverse(A))≈dimi(A) then ATΘB is allowable, propagate dimi(ATΘB)=dimi(B)

  14. Pattern Database • Dimensionality disagreement at binary operators inhibits vectorization. • Recognizing patterns (consisting of operator type and operand dimensionalities) can be used to identify a transformation enabling vectorization. • lhs operation rhs output • (ri, rj) Θ(ri,1) (ri, rj) Pattern: for i=1:m, for j=1:n, A(i,j)=B(i,j)+C(i); end end Transformed Result B(i,j)+C(i); B(1:m,1:n)+repmat(C(1:m),1,n);

  15. Pattern Database • Diagonal access pattern: • lhs operation rhs output • (ri, ri) (index) nil (1, ri) Pattern: for i=1:n, a(i)=A(i,i)*b(i); end a(1:n)=A((1:n)+size(A,1)*((1:n)-1)).*b(1:n); Column major indexing of A

  16. for i1=…, for i2=…, … for ik=… A(J)=A(J)+E; … end end end Loop nest variables I={i1,i2,…,ik} J is a subset of E for i=1:m, for j=1:n, a(i)=a(i)+B(i,j); end end I={i,j} J={i} Additive Reduction Statements • Additive-reduction statements use a loop variable to perform an accumulation. • Not all loop nest index variables appear in output dimensionality

  17. for i=1:m a=a+10; end I={i},J={} I-J={i} ρ(10)={} ri not in dimi(10) Reduce: 10m*10, ρ(m*10)={ri} Vectorize: a=a+m*10; for i=1:m a=a+b(i); end I={i},J={} I-J={i} ρ(b(i))={} ri in dimi(b(i))=(ri,1) Reduce: b(i)sum(b(i),1); Vectorize: a=a+sum(b(1:m)); Additive Reduction (Solution) • Maintain/propagate dimensionality and reduced variables for an expression. • ρ(E) denotes the reduced variables for expression E • When checking statement A(J)=A(J)+E • ensure dimi1,i2,…,ikA(J)≈dimi1,i2,…,ik(E) and ρ(E)=I-J • any variable ri in I-J but not in ρ(E) must be reduced

  18. for i=1:m for j=1:n a(i)=a(i)+B(i,j)*x(j); end end • j is used for reduction • dimi,j(B(i,j))=(ri,rj) • dimi,j (x(j))=(rj) a(1:m)=a(1:m)+… B(1:m,1:n)*x(1:n); Additive Reduction via Matrix Multiplication • Matrix multiplication can be used to perform reductions on e=elhs*erhs , provided: • dimi1,…,ik(elhs)=(Sl,rk) • dimi1,…,ik(erhs)=(rk,Sr) • rk is a reduction variable. • Implies: • dimi1,…,ik(e)=(Sl,Sr) • ρ(e)=union(ρ(elhs), ρ(erhs),{rk})

  19. Use matrix multiplication to reduce rj ρ(a(i,j)*b(j))={rj}, dimi,j(a(i,j)*b(j))=(ri) Additive Reduction Example • Additive reduction example: ρ(a(i,j))={}, dimi,j(a(i,j))=(ri,rj) ρ(b(j))={}, dimi,j(b(j))=(rj) rj is reduction variable for i=1:m, for j=1:n, d(i)=d(i)+a(i,j)*b(j)+c(i,j) end end ρ(c(i,j))={}, dimi,j(c(i,j))={ri,rj} Need to reduce rj: c(i,j)sum(c(i,j),2); ρ(a(i,j)*b(j)+sum(c(i,j),2))={rj}, dimi,j(a(i,j)*b(j)+sum(c(i,j),2)=(ri,rj) Dimensionality and reduced variables agree, now replace index variables: d(1:m)=d(1:m)+a(1:m,1:n)*b(1:n)+sum(c(1:m,1:n),2);

  20. Implementation Prototype Vectorized Loop Original Loop • Pattern database and corresponding transformations are specified in modular end-user extensible manner. Vectorizer Embedded Control Statements Code Generator Octave Parser yes Dimension Check Success no no Create DDG Vectorize Statement yes

  21. Results • Source-to-source transformation • Timing results averaged over 100 runs: • Platform: • Matlab 7.2.0.283 • 3.0 GHz Pentium D Processor

  22. Results • Histogram Equalization: Input source Vectorized Result h=hist(im(:),[0:255]);%histogram heq=255*cumsum(h(:))/sum(h(:)); for i=1:size(im,1), for j=1:size(im,2), im2(i,j)=heq(im(i,j)+1); end end h=hist(im(:),[(0:255)]); heq=255*cumsum(h(:))/sum(h(:)); im2(1:size(im,1),1:size(im,2))=... heq(im(1:size(im,1),1:size(im,2))+1); • For monochrome 8-bit 800x600 image: • original/vectorized: • Entire routine: 0.178s/0.114s (speedup: 1.56) • Loop Portion only: 0.0814s/0.0176s (speedup: 4.6)

  23. Results (Menon & Pingali Examples) for k=1:p, for j=1:(i-1), X(i,k)=X(i,k)-L(i,j)*X(j,k); end end X(i,1:p)=X(i,1:p)-L(i,1:i-1)*X(1:i-1,1:p); for i=1:N,for j=1:N phi(k)=phi(k)+a(i,j)*x_se(i)*f(j); end end phi(k)=phi(k)+sum(a(1:N,1:N)’* x_se(1:N).*f(1:N),1); for i=1:n,for j=1:n, for k=1:n,for l=1:n y(i)=y(i)+x(j)*A(i,k)* B(l,k)*C(l,j); end end end end y(1:n)=y(1:n)+x(1:n)’*... (A(1:n,1:n)*B(1:n,1:n)’*C(1:n,1:n))’;

  24. Remaining Issues/Future Work • Each pattern transformation is local; no optimization over entire statement. • e.g., we do not optimize and distribute transposes • Control flow within loop • Function calls • functions are treated as pointwise operators (correct for many predefined arithmetic functions) • Incorporate our analysis directly with shape analysis

  25. Summary • Contributions: • A simple method to prevent incorrect vectorization in Matlab • A user extensible operator/dimensionality pattern database can be used to improve vectorization • These patterns can make use of higher level semantics (e.g., matrix multiplication) or diagonal accesses in vectorization.

  26. Acknowledgements • Funding provided by NSERC • Grateful for reviewers comments and suggestions

  27. Thank You Questions?

More Related