I want to implement an iterative method for linear systems with the GS rule in MATLAB: at each step we pick the equation whose corresponding component of the residual is the largest in absolute value. For special cases like discretisation of Poisson equation one direction is to use a max-heap for tracking the residuals. I'm not sure if I understand it correctly:
for an initial guess ($x_{0}=0$) calculate the residual and construct a max-heap with its components. For each iteration, find the largest residual, pop it from the heap, update that component of vector x and calculate components of residual which are affected by changing the component we picked and replace the old ones for those in heap (e.g. for Poisson in one dimension I change 2-3 components).
Also, do I have to keep indices of the equation in the max-heap or can I use functions like find in Matlab to find the equation whose residual component is equal to the maximum in the heap?