[mahdi_kh8]

Hello my friends,

I have a simple question but I have a difficulty solving that and need some help and advice.

I am using finite volume method to solve N-S equations.

Assume I have the coefficients: AP, AW, AWW, AE, AEE
(For 1D case)

How I can construct coefficient matrix using these coefficients?
For better understanding please apply your answer in the sample code.

Sample code:

// CONSTRUCT MATRIX OF COEFFICIENT ( X DIRECTION )
// U: ( 3 --> NIM1 , 2 --> NJM1 )
double[,] CoefAx = new double[NI, NI];

for (int JJ = 2; JJ <= NJM1; JJ++)
{
// matrix of coefficient
    for (int I = 3; I <= NIM1; I++)
{
    {
        for (int J = 3; J <= NIM1; J++)
        {
            // Inside matrix
            if (I == J && I >= 5 && I <= NIM1 - 2)
            {
                CoefAx[I, J] = AP[I, JJ];
                CoefAx[I, J + 1] = -AE[I, JJ];
                CoefAx[I, J - 1] = -AW[I, JJ];
                CoefAx[I, J + 2] = -AEE[I, JJ];
                CoefAx[I, J - 2] = -AWW[I, JJ];
            }
            // last raw of matrix
            if (I == J && I == NIM1)
            {
                CoefAx[I, J] = AP[I, JJ];
                CoefAx[I, J - 1] = -AW[I, JJ];
                CoefAx[I, J - 2] = -AWW[I, JJ];
            }
            if (I == J && I == NIM1 - 1)
            {
                CoefAx[I, J] = AP[I, JJ];
                CoefAx[I, J - 1] = -AW[I, JJ];
                CoefAx[I, J - 2] = -AWW[I, JJ];
                CoefAx[I, J + 1] = -AE[I, JJ];
            }
            // first raw of matrix
            if (I == J && I == 3)
            {
                CoefAx[I, J] = G.AP[I, JJ];
                CoefAx[I, J + 1] = -G.AE[I, JJ];
                CoefAx[I, J + 2] = -G.AEE[I, JJ];
            }
            if (I == J && I == 4)
            {
                CoefAx[I, J] = AP[I, JJ];
                CoefAx[I, J + 1] = -AE[I, JJ];
                CoefAx[I, J + 2] = -AEE[I, JJ];
                CoefAx[I, J - 1] = -AW[I, JJ];
            }
        }
    double[] Bx = new double[G.NI];
    }
    // RHS of equation
    for (int I = 3; I <= NIM1; I++)
    {
        Bx[I] = SU[I, JJ] + ANN[I, JJ] * U[I, JJ + 2] + ASS[I, JJ] * U[I, JJ - 2];
    }

    Matrix A = new Matrix(NI, NI);
    Matrix b = new Matrix(NI, 1);
    double temm;

    for (int i = 3; i <= NIM1; i++)
    {
        for (int j = 3; j <= NIM1; j++)
        {
            A[i, j] = CoefAx[i, j];
            temm = Bx[i];

            b[i, 0] = temm;       
        }
    }

     // Solving equation
    Matrix xx = (A).SolveWith(b);
    for (int i = 0; i <= NIM1; i++)    
        { 
            U[i, JJ] = xx[i, 0];    
        } 
}

[FMDenaro]

In 1D is very simple, AP is in the main diagonal, AE, AEE are in the next upper diagonals and AW, AWW in the next lower one. I see you already did that.

[mahdi_kh8]

Originally Posted by FMDenaro

In 1D is very simple, AP is in the main diagonal, AE, AEE are in the next upper diagonals and AW, AWW in the next lower one. I see you already did that.

So my code is correct ?
What about boundary condition ?
I am a little bit confused about constructing matrix and applying BC.

[FMDenaro]

Originally Posted by mahdi_kh8

So my code is correct ?
What about boundary condition ?
I am a little bit confused about constructing matrix and applying BC.

you have to modify the corresponding term in the diagonal and, depending on either Neumann or Dirichlet BC.s add a term to the known terms vector

[mahdi_kh8]

Originally Posted by FMDenaro

you have to modify the corresponding term in the diagonal and, depending on either Neumann or Dirichlet BC.s add a term to the known terms vector

If the variable was Pressure and BC was as follow:

Left: dP/dx = 0
Right: P=0

Could you please apply these BC in the code?

Thank you