CFD Online Discussion Forums - Adding U*U as a source term

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)

- OpenFOAM Programming & Development (https://www.cfd-online.com/Forums/openfoam-programming-development/)

- - Adding U*U as a source term (https://www.cfd-online.com/Forums/openfoam-programming-development/79666-adding-u-u-source-term.html)

Adding U*U as a source term

Hi Foamers,

I am trying to solve the following equation:

dU/dt = grad(p) - 0.5*(f/L)*(U*U)

here, U(velocity vector), p (pressure), and L is length, (f is dimensionless scalar).

In fact, I made an attempt and it seems working for 1D .. However, I can not be very sure about the approach that I use (especially for 2D and 3D problems), and I need some comments from more experienced people.

I constructed my equations as follows:

Code:

        ...



        vector identity(1,1,1);



        ...

    

        volVectorField gradp = fvc::grad(p);



        ...



        fvVectorMatrix UEqn

        (

                 fvm::ddt(U) + fvm::Sp((U & identity)/2*(1/L+f/D), U)

        );



        ...



        solve(UEqn == -fvc::grad(p));



        ...

The reason that I am not sure is about the term "U*U" in the original equation. Since U is a vector, U*U will produce a matrix, and therefore it was not possible to write "U*U" as it is to the code as a source term.

Instead, assuming that only the i component of the vectors will be used while solving Ux, I wrote "U & identity" to obtain the value of the x component of U and put it in to the first parameter place of Sp function.

This works for 1D but I cannot be sure about my assumption about the behavior of U & identity expression in 2D and 3D problems..

In short, say that U(1, 2, 3) then will "U & identity" produce 1+2+3=6 in ALL Ux, Uy and Uz solutions OR will it produce 1 in Ux solution, 2 in Uy solution and 3 in Uz solution??

Kind regards,

Ozgur

Hi,

PS: The above is not a realistic CFD problem, but I am just playing a little with the code so to get some understanding...

Anyway, I think I got the answer of my question by constructing a 2D case. And, I saw that unfortunately (in fact fortunately) my assumption was wrong, i.e. use of "U & identity" produced Ux+Uy+Uz.

So, in this problem, what is needed is something like scale(U,U) function, which gives the product of corresponding components of vectors U and U; so it will produce another vector (not a matrix which only the main diagonal is useful).

I couldn't still make scale(U,U) work. But asked about it under another title.

Ozgur

scale(U,U) doesn't work and unfortunately it has not been updated in programmers guide. You should use cmptMultiply(U,U) instead.