[SirIsaac90]

Hello Foamers!

I have again a problem. I think its easy to solve for you guys, but my code can't be correct.
The implicit Euler method shall be implemented. Therefore, I have to solve a linear system Ax = b.

'b' is in this case the velocity vector of the actual time step 'n'. 'Ax' contains everything with respect the new 'U' of time step n+1 that I want to solve this system for.
This part of my code looks as follows:

        fvVectorMatrix UEqn
        (
            U - (-runTime.deltaT()*fvm::div(phi, U)+
            runTime.deltaT()*nu * fvm::laplacian(U))
        );

        solve (UEqn == U);

Sorry, but I thought nothing else has to be done. In 'fvSolution' I chose the PBiCG with the DILU preconditioner.

Can you tell me how to solve that system to get the new velocity field 'U'.

Thank you and best regards!

[mkraposhin]

Hi,

why you don't want to write your matrix as in standard OpenFOAM solver:

    fvVectorMatrix UEqn
        (
            fvm::ddt(U)
             + fvm::div(phi, U) -
             fvm::laplacian(nu, U))
        );

        UEqn.solve();

But if you want to implement Euler time scheme as you suggested, you need to insert U at 'n+1' slice as implicit source at lhs:

    volVectorField Uold = U.oldTime()

        fvVectorMatrix UEqn
        (
            fvm::Sp(geometricOneField(),U)
            ==
            runTime.deltaT()*
            (
               -fvm::div(phi, U)
                +
                nu * fvm::laplacian(U))
            )
            +
            Uold
         );
        
         UEqn.solve();

[SirIsaac90]

Thank you very much! That helped me a lot. I want to do this, because of a Benchmark of methods based on the projection method in OpenFOAM.

Instead of 'geometricOneField()' I just wrote '1', because an error occured in the compiling process. But in my opinion, it works fine now.

Now, I want to implement a DIRK method, so I have to compute velocity slopes, let's say 'ki'.

Is it possible to use 'UEqn' to solve the linear system for a volVectorField 'ki' somehow? Or can I only do this always with 'U'? Is there a general method to do that with any volVectorFields without writing a completely new code for that? The 'ki's should be stored in the PtrList<volVectorField> ki

It would be nice to define the method for the linear system in fvSolution, too.

Best regards!

[SirIsaac90]

Okay, maybe the question is too stupid I am a beginner and it is quite difficult to get the knowledge from other sources..

Maybe someone can help me with this:
('lauf3, lauf4' are just running numbers.., 'k' is a volVectorField,
Uold is derived from volVectorField 'U')
I want to solve the system for k

        fvVectorMatrix kEqn
        (
                fvm::Sp(1, k)
                ==
                - fvm::div
                (phi, Uold + runTime.deltaT()*A[lauf3][lauf4] * k)
                +
                nu * fvm::laplacian
                (Uold + runTime.deltaT()*A[lauf3][lauf4] * k)
        );
        kEqn.solve()

Can someone tell me where is my mistake? Actually, it should be the same as for 'UEqn' from my first post, but the compiling does not work.. The error has to do something with k.

Thanks and best regards!

[chriss85]

You need to only have the field you want to solve for in the second argument of the fvm:: operators, not others. You can put the other fields in the first argument.