[KTG]

Hi all,


I am trying to solve poisson's equation by adding a forcing term to laplacianFoam and removing the time derivative. My code runs, but the solution around the forcing is way off - I get large numbers that depend on the size of the mesh, but the solutions at the Dirichlet boundaries make sense. This lead me to believe that when I write:

           

 fvScalarMatrix psixEqn
            (
                fvm::laplacian(DT, psix) + f1
            );

openfoam does not multiply by the cell volume but rather inserts the forcing term f1 directly into the matrix problem? After a lot of internet searching, I modified the forcing this way in my createFields.H:

Info<< "Creating forcing term F1"<<endl;
volScalarField F1 (
        IOobject (
            "F1",
            runTime.timeName(),
            mesh,
            IOobject::NO_READ
       ),
       mesh,
       dimensionedScalar("F1",dimensionSet(0,-3,0,0,0,0,0), 0.0),
       calculatedFvPatchScalarField::typeName
);

F1.ref() = f1/mesh.V();

which gave results closer to what I expected, though it is still off depending on the mesh size. I also feel like I should be multiplying by the cell volume rather than dividing, but that gave really wild results. Clearly there is something I am not getting here - there is probably an easier way to do this that will make it work properly...help.



Thanks

[KTG]

This is such a basic thing, its very strange that I can't find a casefile anywhere that shows how to do it! Its seems like I should be using fvc::sp(f1) or fvc::su(f1), but when I use these, I get an error:

error: no matching function for call to ‘Sp(Foam::volScalarField&)’
                 fvm::laplacian(DT, psix) == fvc::Sp(f1)

where the ^ points to the last ) after f1


I have defined f1 this way:

   volScalarField f1
    (
        IOobject
        (
            "f1",
            runTime.timeName(),
            mesh,
            IOobject::MUST_READ,
            IOobject::AUTO_WRITE
        ),
        mesh
    );

Maybe volScalarField is the problem? I don't know what else to put there, and am not finding any clues in fvcSup.H (if that is even the right file).

[clapointe]

What version of OpenFOAM are you using? I quick search of the code (version 4) returns many uses of fvm::Sp, but none of fvc::Sp, which I think is the reason you are getting the error.

Caelan

[KTG]

Yes, I think that you are correct for the reason for the error - maybe it is outdated, but I read that fvc:: should explicitly add the source term to the right hand side. My original question stands though. I used this code:

            fvScalarMatrix psixEqn
            (
                fvm::laplacian(DT, psix) == fvm::Su(f1,psix)
            );

and am getting the same results as before - the forcing terms seems to vary with the size of the mesh, and is negative when I put it on the right hand side. When I set f1=1, it comes out to about 1800 in the solution for psix. Do you know if openfoam is just setting the source term equal to f1, or actually multiplying it by the cell volume?

[clapointe]

It looks like the strength is scaled by cell volume -- take a look at (for Openfoam 4x): https://github.com/OpenFOAM/OpenFOAM...e/fvm/fvmSup.C.

Caelan

[KTG]

Yes, I think you are right, which makes sense in terms of integrating over the cell in the finite volume method. What confuses me is that the result is the same when I just insert the forcing term directly into the equation without fvm::Sp(f1):

            fvScalarMatrix psixEqn
            (
                fvm::laplacian(DT, psix) == f1
            );

Anyway it seems like using fvOptions is the preferred way to do things now, so maybe I should do that instead. I think some of my initial confusion was due to poor mesh quality anyway - thanks for the help!

[shouvik ghorui]

Did you get any preferred solution to this problem? I am also facing the same issue. It would be a great help if you could help me

[KTG]

fvOptions is definitely the way to go:


https://www.openfoam.com/documentati...fvoptions.html


Check out this thread for instructions on how to find all the tutorials that use fvOptions for some examples:


Looking for a simple example with fvOption in OpenFOAM

[shouvik ghorui]

Thank you for your response,

I have to implement a gradient of the C term (as given in the picture in red) in the source term of scalar transport equation.
Can you please give me some idea how to mange with this. It would be very helpful for me, I got stuck for many days.


Note : I have done

fvScalarMatrix CEqs
(
fvm::ddt(C)
+fvm::div(phiup,C)
-fvm::laplacian(gammaPsi,C)
+fvm::Sp(Un,C)
+?????

);