Home Forums Preconditioned Steepest Descent Method (PSDM)

Viewing 1 post (of 1 total)
  • Author
  • #1828

    When I was implementing an iterative solver using Steepest Descent Method (SDM), I tried to search on the web how to implement a preconditioned version of the algorithm.

    However, I only found posts in some forums asking the same question.

    Therefore, I started reading this famous document on Conjugate Gradient (CG) – An Introduction to the Conjugate Gradient Method Without the Agonizing Pain by Jonathan Richard Shewchuk (http://www.cs.cmu.edu/~quake-papers/painless-conjugate-gradient.pdf).

    Then, I derived PSDM by myself.

    The original SDM looks like this:

    Based on a symmetric positive definite preconditioner M,
    M = EE[sup:20nmn0md]T[/sup:20nmn0md] by Cholesky factorization
    the preconditioned coefficient matrix should be E[sup:20nmn0md]-1[/sup:20nmn0md]AE[sup:20nmn0md]-T[/sup:20nmn0md] because it is also symmetric positive definite
    the original system Ax = b is transformed to
    E[sup:20nmn0md]-1[/sup:20nmn0md]AE[sup:20nmn0md]-T[/sup:20nmn0md]y = E[sup:20nmn0md]-1[/sup:20nmn0md]b
    where y = E[sup:20nmn0md]T[/sup:20nmn0md]x

    Thus, given z = M[sup:20nmn0md]-1[/sup:20nmn0md]r, we can update the symbols in our original algorithm:
    A -> A’ = E[sup:20nmn0md]-1[/sup:20nmn0md]AE[sup:20nmn0md]-T[/sup:20nmn0md]
    b -> b’ = E[sup:20nmn0md]-1[/sup:20nmn0md]b
    r -> r’ = E[sup:20nmn0md]-1[/sup:20nmn0md]r
    q -> q’ = A’r’ = E[sup:20nmn0md]-1[/sup:20nmn0md]Az
    -> r'[sup:20nmn0md]T[/sup:20nmn0md]r’ = r[sup:20nmn0md]T[/sup:20nmn0md]z
    -> r'[sup:20nmn0md]T[/sup:20nmn0md]q’ = (E[sup:20nmn0md]-1[/sup:20nmn0md]r)[sup:20nmn0md]T[/sup:20nmn0md] E[sup:20nmn0md]-1[/sup:20nmn0md]Az = r[sup:20nmn0md]T[/sup:20nmn0md]E[sup:20nmn0md]-T[/sup:20nmn0md]E[sup:20nmn0md]-1[/sup:20nmn0md]Az = z[sup:20nmn0md]T[/sup:20nmn0md]Az

    The update for x becomes
    x = x + alpha * r -> x’ = x’ + alpha’ * r’ -> E[sup:20nmn0md]T[/sup:20nmn0md]x = E[sup:20nmn0md]T[/sup:20nmn0md]x + alpha’ * E[sup:20nmn0md]-1[/sup:20nmn0md]r
    multiply on both sides by E[sup:20nmn0md]-T[/sup:20nmn0md], becomes
    x = x + alpha’ z   (since E[sup:20nmn0md]-T[/sup:20nmn0md]E[sup:20nmn0md]-1[/sup:20nmn0md] = M[sup:20nmn0md]-1[/sup:20nmn0md])

    Similarly, the update for r becomes
    r = r – alpha * q -> r’ = r’ – alpha’ * q’ -> E[sup:20nmn0md]-1[/sup:20nmn0md]r = E[sup:20nmn0md]-1[/sup:20nmn0md]r – alpha’ * E[sup:20nmn0md]-1[/sup:20nmn0md]Az
    multiply on both sides by E, becomes
    r = r – alpha’ * q’   (since Az = q’)

    And, PSDM looks like this:

    I hope this would help others who are also implementing the same method. Correct me if I miss anything.

Viewing 1 post (of 1 total)
  • You must be logged in to reply to this topic.