Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Conversion from Centered Diff Scheme to Ax = b

  1. Jun 18, 2012 #1
    Hey all, I've been working with an in-house code for a while and have decided to pursue a different method of solving the pde inside. Currently a spatial centered difference scheme is used to give us our equation to solve with the Jacobi iterative method. I want to investigate the conjugate gradient method as it pertains to the problem, but have been struggling. I need a new set of eyes.
    So first off, here is the equation inside the Jacobi: (Let i denote x direction, j the y direction)

    [itex]\Psi_{i,j} = \kappa_A \left(\left(\Delta x\right)^2 * \left(x_{i,j}^2 + y_{i,j}^2\right) * \omega_{i,j} + \Psi_{i+1,j} + \Psi_{i-1,j} + \kappa_2*\left(\Psi_{i,j+1} + \Psi_{i,j-1}\right)\right)[/itex]

    Where [itex]\Delta x[/itex], [itex]\kappa_A[/itex], and [itex]\kappa_2[/itex] are scalar constants.
    [itex]\Psi[/itex] is what I'm solving for, [itex]\omega[/itex] is the RHS.

    So for conjugate gradient method (CGM) we need to rearrange this to the form Ax = b, where A is an NxN symmetric positive definite matrix. So this is where I need someone to double check my work. I started with the above equation and rearranged it to:

    [itex]\Psi_{i,j}/\kappa_A - \Psi_{i+1,j} - \Psi_{i-1,j} - \kappa_2*\left(\Psi_{i,j+1} + \Psi_{i,j-1}\right) = \left(\Delta x\right)^2\left(x_{i,j}^2 + y_{i,j}^2\right)*\omega_{i,j}[/itex]

    Assuming that is correct, then the matrix A should consist of [itex]1/\kappa_A[/itex] along the diagonal, -1's on both sides of the diagonal and then [itex]-\kappa_2[/itex] at each side spaced out from the diagonal based on the size of [itex]\Psi[/itex].
    Might look something like this:

    The vector, x, would then be [itex] \left[ \Psi_{1:end,1};\Psi_{1:end,2};... \Psi_{1:end,end}\right] [/itex]

    And vector, b, would be [itex]\left(\Delta x\right)^2 * \left[\left(x_{1:end,1}^2 + y_{1:end,1}^2\right)*\omega_{1:end,1};... \left(x_{1:end,end}^2 + y_{1:end,end}^2\right)*\omega_{1:end,end}\right] [/itex]

    So how have I done? Am I good up to here?
  2. jcsd
Share this great discussion with others via Reddit, Google+, Twitter, or Facebook

Can you offer guidance or do you also need help?
Draft saved Draft deleted