/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Website:  https://openfoam.org                  |
|   \\  /    A nd           | Version:  9                                     |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    format      ascii;
    class       dictionary;
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application     rheoFoam;

startFrom       startTime;

startTime       0;

stopAt          endTime;

endTime         1;

deltaT          1e-7;

writeControl    adjustableRunTime;

writeInterval   0.1;

purgeWrite      0;

writeFormat     ascii;

writePrecision  10;

writeCompression compressed;

timeFormat      general;

timePrecision   10;

graphFormat     raw;

runTimeModifiable yes;

adjustTimeStep yes;

maxCo          0.1;

maxDeltaT      1e-5;

functions
{

// Monitor for convergence

    probes
    {
        // Where to load it from
        functionObjectLibs ( "libsampling.so" );

        type            probes;

        // Name of the directory for probe data
        name            probes;

        // Write at same frequency as fields
        writeControl   timeStep;
        writeInterval  50;

        // Fields to be probed
        fields
        (
            U
	    p
	    tau
        );

        probeLocations
        (
	     (-0.0001 0.00014375 -0.00001)
             (0 0.00014375 -0.00001)
	     (0.0001 0.00014375 -0.00001)
        );
    }

    QField
    {
        type coded;
        libs ("libutilityFunctionObjects.so");
        name QField;

        codeExecute
        #{
            const volVectorField& U = mesh().lookupObject<volVectorField>("U");

            // Velocity gradient
            volTensorField gradU(fvc::grad(U));

            // Symmetric (strain) and antisymmetric (rotation) parts
            volSymmTensorField D(symm(gradU));
            volTensorField W(0.5*(gradU - gradU.T()));

            // Compute invariants
            volScalarField DD(magSqr(D));
            volScalarField WW(magSqr(W));

            // Define Q with full IOobject so it can be written
            volScalarField Q
            (
                IOobject
                (
                    "Q",
                    mesh().time().timeName(),
                    mesh(),
                    IOobject::NO_READ,
                    IOobject::NO_WRITE,   // ensures it will be written
                    false                   // do not register to objectRegistry
                ),
                (DD - WW)/(DD + WW + dimensionedScalar("tiny", DD.dimensions(), SMALL))
            );
	    
	    if (mesh().time().timeOutputValue() >= mesh().time().endTime().value() - SMALL)
            {
            	Q.write();
            }
        #};
    }

    #includeFunc residuals	
}
// ************************************************************************* //
