[wallace]

… where the source term depends on the current velocity field. How can I make sure that it converges at each time-step in an unsteady simulation?

I'm using pimpleFoam and I've made a new fvOption which is similar to actuationDiskSource. I understand that it works as follows:

Essentially, the class has a public function AddSup where eqn.source() is used to set the source term, where eqn would appear to reference UEqn when the IOOptionList is created by the solver, simple-/pimpleFoam. Using eqn.psi(), the fvOption class also has access to the current velocity, and mesh_ is a reference to the mesh database provided by the base class option. So it's straightforward to make a momentum source term, which depends on the current velocity field, in my case sampled at different points.

In pimpleFoam, the UEqn is defined with == fvOptions(U) on the r.h.s., and then the solution loop is passed through. I hope I'm wrong, but from pisoFoam.C, I don't get the impression that fvOptions(U) is called again, (and addSup recalculates the source term based on the intermediate velocity field solution), during each time-step.

Perhaps this is a question on how literally to interpret the top-level solver syntax (if you haven't had the time to trace back through everything). I hope I'm not getting terminology wrong here but in

    
    3 tmp<fvVectorMatrix> UEqn
    4 (
    5     fvm::ddt(U)
    6   + fvm::div(phi, U)
    7   + turbulence->divDevReff(U)
    8  ==
    9     fvOptions(U)
   10 );

I assume that the l.h.s terms are implicitly included in the matrix, while the r.h.s is explicit. Therefore, for a velocity dependent source term in an unsteady solution, there would need to be an inner loop for each time-step to converge the source term.

Assuming I'm right so far, my next step would be to modify s-/pimpleFoam to repeat the solution loop at each step. On my first attempt (for simplicity, with simpleFoam) I'm confused by the results:

Clearly it's not as simple adding a loop,

               for (int i = 0; i < 2; i++){}

ultimately a while loop to be controlled by some measure of convergence of the calculated source term, around the following:

   59         // --- Pressure-velocity SIMPLE corrector
   60         {
   61             #include "UEqn.H"
   62             #include "pEqn.H"
   63         }
   64 
   65         turbulence->correct();
   66 
   67         runTime.write();

With no source term, I might have expected the second loop to converge immediately, but it appears to do almost as much work as the first loop (same no. of iterations, very similar initial residual).

I guess my questions are more implicit than explicit here
Basically, I'd be glad if someone "in-the-know" could confirm my understanding of the way the fvOptions are adding source terms. I'm quite confident about that part, but then I'm much less sure how they are being dealt with in the solver. As my simpleFoam loop test didn't do what I expected, I'm unsure of how to make a time-step inner loop in pimpleFoam to converge the source term.

[pbachant]

Did you ever make any progress on this?

[pruthvi1991]

Hell everybody,

I'm working on a similar problem, trying to implement a solver in the non-inertial frame for a flapping airfoil study. I haven't tried using fvOptions but I explicitly added the body force term on the RHS of Ueqn. However I'm not confident that my implementation is correct. Here is a thread about it.

Can you guys check it out and give your opinion on this issue?

Thanks,
Pruthvi.

[potentialFoam]

I would like to add a user-defined momentum source term to the Navier-Stokes equations.
The source depends on the actual flow field, it should be used in the solvers simpleFoam, pimpleFoam, ... and for OF 3.0.0 (or 2.4.0).

I suppose, a reasonable approach would be the use of a 'fvOptions'-dictionary file with the 'vectorCodedSource'-type.

Do you think, this is the best solution?
If yes, where do I find information on how to write this 'fvOptions'-file?

The OF-tutorials provide no help for this source type.

[jherb]

google is your friend:
https://www.google.de/search?q=openf...term+fvOptions
http://www.sourceflux.de/blog/adding...ons-fvoptions/
http://www.sourceflux.de/blog/softwa...ign-fvoptions/

[potentialFoam]

Thanks for your fast answer! Unfortunately, Google had not helped me.

http://www.sourceflux.de/blog/adding...ons-fvoptions/
-> deals only with specific 'MRFSource'-option and gives no further information on a 'codedSource'-type option,

http://www.sourceflux.de/blog/softwa...ign-fvoptions/
-> is the subsequent post and deals only with the modification of a special solver.

My problem concerns still the 'vectorCodedSource'-type, especially the formulation of the 'code...'-parameters, or more exact the definition of the momentum source term.

Here is my actual 'fvOptions'-file:

FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvOptions;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

momentumSource
{
    type            vectorCodedSource;
    active          on;
    
    vectorCodedSourceCoeffs
    {
        fieldNames  (U); // ??

        selectionMode   cellZone; //cellSet;//cellZone; ??
        cellZone           sourceBox;

    redirectType  sourceTime; //velocitySource; //none;
    
    codeInclude
    #{

    #};
        
    codeCorrect
    #{
        Pout<< "**codeCorrect**" << endl;
    #};
    
        codeAddSup
    #{
        vectorField& USource = eqn.source();
        USource = (1 1 1); // Definition of the user-defined momentum source term
    #};
        
    codeSetValue
    #{
        Pout<< "**codeSetValue**" << endl;
    #};
        
    // Dummy entry. Make dependent on above to trigger recompilation
    code
    #{
        $codeInclude
        $codeCorrect
        $codeAddSup
        $codeSetValue
    #};
    }
    
    sourceTimeCoeffs
    {
    // Dummy entry
    }
  
}

Besides, I am wondering whether the variables
- fieldNames
- redirectType
are defined correctly.

Regards,
Peter

[pbachant]

If you're looking to add the source to the velocity field, you have the fieldNames entry correct. I'm not sure what redirectType means. The documentation is usually in the header files, but there is no explanation for that entry: https://github.com/OpenFOAM/OpenFOAM.../CodedSource.H

[potentialFoam]

Dear Foamers,

first: thanks Pete!

Unfortunately, I still cannot introduce a momentum source into OF.
And both google and the source code seem to be no help.

At first, I would be quite happy to introduce a constant momentum source like e.g. (1 1 1).

I think the

 codeAddSup{}

function inside the vectorCodedSource-type should do it.
This is my 'fvOptions'-file:
(I use OF 2.4.0 or 3.0.0.)

FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvOptions;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * 

momentumSource
{
    type            vectorCodedSource;
    active          on; // true;

    vectorCodedSourceCoeffs
    {
    fieldNames (U);

    selectionMode   cellZone;
    cellZone        sourceBox;

    redirectType    none;
//     redirectType    velocitySource;
//     redirectType    sourceTime;
//     redirectType    ramp;
    
    codeInclude
    #{

    #};
        
    codeCorrect
    #{
        Pout<< "**codeCorrect**" << endl;
    #};
    
        codeAddSup
    #{
        ....... // ### Momentum Source ###
        Pout << "***codeAddSup***" << endl;
    #};
        
    codeSetValue
    #{
        Pout<< "**codeSetValue**" << endl;
    #};
        
    // Dummy entry
    #{
//         $codeInclude
//         $codeCorrect
        $codeAddSup
//         $codeSetValue
    #};

    }
    
//     sourceTimeCoeffs
    noneCoeffs
    {
    // Dummy entry
    }
  
}

for ### Momentum Source ### I tried the following combinations:

vectorField& Usource = eqn.source();
Usource += vector(10000, 0, 0);

-> Runs, but there is no source term?!
and

volVectorField& Usource = eqn.source();
Usource = (10000 0 0);

-> error: invalid initialization of reference of type ‘Foam::volVectorField& {aka Foam::GeometricField<Foam::Vector<double>, Foam::fvPatchField, Foam::volMesh>&}’ from expression of type ‘Foam::Field<Foam::Vector<double> >’

I cannot solve this tiny nasty problem! And would be really glad for some help!

[pbachant]

Well your second version doesn't work because you are trying to set a volVectorField& to a scalar. The first version looks like it should work. You should check that your cellZone is created properly by reading the topoSet log or inspecting in ParaView. Note that the newest version of the momentum sources in OpenFOAM use cellSets rather than cellZones. You could also try selecting all cells to make sure the source is being added.

[potentialFoam]

Thanks again, Pete!
You were right, I did a mistake in 'topoSet'. So the upper code snippet works.
Thanks a lot!