[Daniele111]

Thanks you Martin, it works fine. I'll go to Germany and I'll offer you a good beer! Do you get correct syntax also for turbolent case this time in the code? For single core I wrote code coping it from OpenFoam utility.
Daniele

[MartinB]

Beer is fine

If you need more data to be written, just post your single core code, I'll have a look at it. I don't use turbulence for my simulations here, so I'm not sure if it will work with turbulent flow. But if you can provide a simple test case for simpleFoam we can check it...

See you in Germany

Martin

[Daniele111]

I changed fot single core your old code in this way:


#include "RASModel.H"
// write matlab file
{
Info << "Writing wall shear stress for matlab" << endl;

word wallPatchName = "bottom";
Info << "Searching for patch " << wallPatchName << endl;

// Check, if the patch name exists:
label patchID = mesh.boundaryMesh().findPatchID(wallPatchName);
if (patchID < 0)
{
Info << "No patch found with name \"" << wallPatchName << "\"!"
<< endl;
}
else
{
/* for laminar case:
Info<< "Reading viscosity nu\n" << endl;
dimensionedScalar nu_
(
laminarTransport.lookup("nu")
);
*/
// calculating wall shear rate


#include "createFields2.H"
mesh.readUpdate();
volSymmTensorField Reff(RASModel->devReff());
volVectorField tau
(
IOobject
(
"tau",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
mesh,
dimensionedVector
(
"tau",
Reff.dimensions(),
vector::zero
)
);

forAll(tau.boundaryField(), patchi)
{
tau.boundaryField()[patchi] =
(
-mesh.Sf().boundaryField()[patchi]
/mesh.magSf().boundaryField()[patchi]
) & Reff.boundaryField()[patchi];
}


// Pointer to the patch we want to have shear stress for:
const polyPatch& cPatch = mesh.boundaryMesh()[patchID];
label nFaces = cPatch.size();

pointField faceCenterCoords(nFaces, vector(0,0,0));
pointField faceAreas(nFaces, vector(0,0,0));
pointField wallShearStress(nFaces, vector(0,0,0));

label runFaces = 0;

// collecting the interesting data
for (int i=0; i<nFaces; i++)
{
faceCenterCoords[runFaces] = mesh.Cf().boundaryField()[patchID][i];
faceAreas[runFaces] = mesh.Sf().boundaryField()[patchID][i];
wallShearStress[runFaces] = tau.boundaryField()[patchID][i];
runFaces++;
}

// make a cute filename
fileName casePath = cwd();//casePath.name();
fileName caseName = casePath.name();
fileName myFile = "tau.dat";
Info << "Data is written to file " << myFile << endl;

// create a new file
OFstream *myStream = NULL;
myStream = new OFstream(myFile);

for (int i=0; i<nFaces; i++)
{
*myStream << i << " "
// coordinates on face centers:
<< mesh.Cf().boundaryField()[patchID][i].component(0) << " "
<< mesh.Cf().boundaryField()[patchID][i].component(1) << " "
<< mesh.Cf().boundaryField()[patchID][i].component(2) << " "

// wall shear stress as a vector:
<< tau.boundaryField()[patchID][i].component(0) << " "
<< tau.boundaryField()[patchID][i].component(1) << " "
<< tau.boundaryField()[patchID][i].component(2) << " "

// magnitude of wall shear stress:
<< mag(tau.boundaryField()[patchID][i]) << endl;
}

}
}

[MartinB]

I merged your version and my parallel one, so can you try this:

#include "RASModel.H"
// write matlab file
{
    Info << "Writing wall shear stress for matlab!" << endl;

    word wallPatchName = "bottom";
    Info << "Searching for patch " << wallPatchName << endl;

    // Check, if the patch name exists:
    label patchID = mesh.boundaryMesh().findPatchID(wallPatchName);
    if (patchID < 0)
    {
        Info << "No patch found with name \"" << wallPatchName << "\"!"
             << endl;
    }
    else
    {
        // calculating tau

#       include "createFields2.H"
        mesh.readUpdate();
        volSymmTensorField Reff(RASModel->devReff());
        volVectorField tau
        (
            IOobject
            (
                "tau",
                runTime.timeName(),
                mesh,
                IOobject::NO_READ,
                IOobject::AUTO_WRITE
            ),
            mesh,
            dimensionedVector
            (
                "tau",
                Reff.dimensions(),
                vector::zero
            )
        );

        forAll(tau.boundaryField(), patchi)
        {
            tau.boundaryField()[patchi] =
            (
                -mesh.Sf().boundaryField()[patchi]
                /mesh.magSf().boundaryField()[patchi]
            ) & Reff.boundaryField()[patchi];
        }



        // Pointer to the patch we want to have shear stress for:
        const polyPatch& cPatch = mesh.boundaryMesh()[patchID];
        label nFaces = cPatch.size();

        pointField faceCenterCoords(nFaces, vector(0,0,0));
        pointField faceAreas(nFaces, vector(0,0,0));
        pointField wallShearStress(nFaces, vector(0,0,0));

        label runFaces = 0;

        // collecting the interesting data
        for (int i=0; i<nFaces; i++)
        {
            faceCenterCoords[runFaces] = mesh.Cf().boundaryField()[patchID][i];
            faceAreas[runFaces] = mesh.Sf().boundaryField()[patchID][i];
            wallShearStress[runFaces] = tau.boundaryField()[patchID][i];
            runFaces++;
        }

        // for parallel
        reduce(nFaces, sumOp<label>());
        Info << "Number of global faces: " << nFaces << endl;


        //  make a cute filename
        fileName casePath = cwd();//casePath.name();
        fileName caseName = casePath.name();
        fileName myFile = "tau.dat";
        Info << "Data is written to file " << myFile << endl;

        // create a new file
        OFstream *myStream = NULL;
        if (Pstream::master())
        {
            myStream = new OFstream(myFile);
        }

        // slaves send data to master, master writes data to file
        if (Pstream::parRun())
        {
            if (Pstream::myProcNo() != Pstream::masterNo())
            {
                OPstream toMaster(Pstream::blocking, Pstream::masterNo());
                toMaster << faceCenterCoords;
                toMaster << faceAreas;
                toMaster << wallShearStress;
            }
            else // master part
            {
                // first write own data
                for (int i=0; i<faceCenterCoords.size(); i++)
                {
                    *myStream << i << ";"
                        << faceCenterCoords[i].component(0) << "; "
                        << faceCenterCoords[i].component(1) << "; "
                        << faceCenterCoords[i].component(2) << "; "

                        << wallShearStress[i].component(0) << "; "
                        << wallShearStress[i].component(1) << "; "
                        << wallShearStress[i].component(2) << "; "

                        << mag(wallShearStress[i]) << endl;
                }

                label run = faceCenterCoords.size();

                // then receive slave data and write
                pointField bufferFaceCenterCoords;
                pointField bufferFaceAreas;
                pointField bufferWallShearStress;


                for (int slave=Pstream::firstSlave(); slave <= Pstream::lastSlave(); slave++)
                {
                    IPstream fromSlave(Pstream::blocking, slave);
                    fromSlave >> bufferFaceCenterCoords;
                    fromSlave >> bufferFaceAreas;
                    fromSlave >> bufferWallShearStress;

                    for (int i=0; i<bufferFaceCenterCoords.size(); i++)
                    {
                        *myStream << run << ";"
                            << bufferFaceCenterCoords[i].component(0) << "; "
                            << bufferFaceCenterCoords[i].component(1) << "; "
                            << bufferFaceCenterCoords[i].component(2) << "; "

                            << bufferWallShearStress[i].component(0) << "; "
                            << bufferWallShearStress[i].component(1) << "; "
                            << bufferWallShearStress[i].component(2) << "; "

                            << mag(bufferWallShearStress[i]) << endl;
                        run++;
                    }
                }
            }
        }
        else // single processor only
        {
            for (int i=0; i<nFaces; i++)
            {
                *myStream << i << "; "
                        // coordinates on face centers:
                        << mesh.Cf().boundaryField()[patchID][i].component(0) << "; "
                        << mesh.Cf().boundaryField()[patchID][i].component(1) << "; "
                        << mesh.Cf().boundaryField()[patchID][i].component(2) << "; "

                        // wall shear stress as a vector:
                        << tau.boundaryField()[patchID][i].component(0) << "; "
                        << tau.boundaryField()[patchID][i].component(1) << "; "
                        << tau.boundaryField()[patchID][i].component(2) << "; "

                        // magnitude of wall shear stress:
                        << mag(tau.boundaryField()[patchID][i]) << endl;
            }
        }

    }
}

Martin