关于negSumDiag()的一个问题



  • 在看关于negSumDiag()的代码的时候遇到一个问题。diag上面的项是由该“行”上面的的off-diag的项相加后取负值得到。为什么下面的代码感觉像是把该“列”的数相加了呢?

    template<class Type, class DType, class LUType>
    void Foam::LduMatrix<Type, DType, LUType>::negSumDiag()
    {
    const Field<LUType>& Lower = const_cast<const LduMatrix&>(*this).lower();
    const Field<LUType>& Upper = const_cast<const LduMatrix&>(*this).upper();
    Field<DType>& Diag = diag();

    const unallocLabelList& l = lduAddr().lowerAddr();
    const unallocLabelList& u = lduAddr().upperAddr();
    
    for (label face=0; face<l.size(); face++)
    {
        Diag[l[face]] -= Lower[face];
        Diag[u[face]] -= Upper[face];
    }
    

    }

    0_1475689903884_Capture.JPG
    以这个网格为例
    0_1475690047966_Capture.JPG
    lowerAdd 为 (0 1 2 1 2 3 4 3 )
    upperAdd 为 (2 2 3 6 4 5 6 5 )

    矩阵的形式为

    0_1475689800365_Capture2.JPG

    那么Diag[0] 应该是由 Upper[1]得到。这个和上面的代码不符。如果矩阵是对称的,这个计算没有问题。而如果是 非对称的,是不是就有问题了?

    谢谢大家!



  • 为了进一步测试,建立一个3X3的网格,并创建一个fvMatrix。

        fvVectorMatrix UEqn_lap
        (
          fvm::laplacian(nu, U)
        );
    

    网格编号如下:

    |6|7|8|
    |3|4|5|
    |0|1|2|

    输出 upper 和lower 的address,没有问题

        Info<< "UEqn_lap upper add before sum   " << UEqn_lap.lduAddr().upperAddr() << nl << endl;
        Info<< "UEqn_lap lower add before sum   " << UEqn_lap.lduAddr().lowerAddr() << nl << endl;
    

    UEqn_lap upper add before sum
    12
    (
    1
    3
    2
    4
    5
    4
    6
    5
    7
    8
    7
    8
    )

    UEqn_lap lower add before sum
    12
    (
    0
    0
    1
    1
    2
    3
    3
    4
    4
    5
    6
    7
    )

    下面把upper lower 和diag重新赋值 并输出

        UEqn_lap.upper() = 100;
        UEqn_lap.diag() = 0;
        UEqn_lap.lower() = 1;
    
        Info<< "UEqn_lap diag after change   " << UEqn_lap.diag() << nl << endl;
        Info<< "UEqn_lap upper after change   " << UEqn_lap.upper() << nl << endl;
        Info<< "UEqn_lap lower after change   " << UEqn_lap.lower() << nl << endl;
    

    输出为:

    UEqn_lap diag after change 9{0}

    UEqn_lap upper after change 12{100}

    UEqn_lap lower after change 12{1}

    调用

        UEqn_lap.negSumDiag();
    
        Info<< "UEqn_lap diag after negSumDiag   " << UEqn_lap.diag() << nl << endl;
    

    输出结果为:

    UEqn_lap diag after negSumDiag 9(-2 -102 -101 -102 -202 -201 -101 -201 -200)

    我特意把upper里面的数值设的很大,从输出的结果看,diag 是把lower里面的数加和了,也就是把列中的数加和了。

    这个结果让我有点儿弄不明白。


  • OpenFOAM教授

    @mengweilm425 这个问题cfd-online上也有人讨论过,不过没有结果。我个人觉得这确实是OpenFOAM的bug,不过由于fvm离散后的系数矩阵的特殊性质,列元素相加和行元素相加得到的结果是一样的,因此这个bug不影响计算。有机会细说。



  • @wwzhao 谢谢您的回复!

    我找到了这个帖子

    http://www.cfd-online.com/Forums/openfoam-programming-development/72456-matrix-storage-diagonal-coeffients-calculation.html

    但是也没说清楚。不过系数矩阵是对称的话,确实不影响结果。


Log in to reply