CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM (http://www.cfd-online.com/Forums/openfoam/)
-   -   help for a modified probe file (wave elevation vs time ) (http://www.cfd-online.com/Forums/openfoam/64888-help-modified-probe-file-wave-elevation-vs-time.html)

zhajingjing May 27, 2009 11:42

help for a modified probe file (wave elevation vs time )
 
Hello everyone,
I want to monitor the free surface elevation vs time ,and I've try hard to search the relevant posts. My current idea is to average the gamma values on those cells at x in order to get the wave elevation vs time at x. So I just modify the "patchAverage"file and try to make it as much simple as possible. After "wmake", a exe file has been produced,but when I use it ,error occur ,and I don't know where's wrong.Please help me !


This's the modified file:

Application
freeSurfaceElevation

Description
Calculates the average of the specified field over the specified x point.

\*---------------------------------------------------------------------------*/

#include "fvCFD.H"

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

// Main program:

int main(int argc, char *argv[])
{
timeSelector::addOptions();
argList::validArgs.append("fieldName");

//argList::validArgs.append("x"); //to avoid error I get rid of second arg

# include "setRootCase.H"
# include "createTime.H"
instantList timeDirs = timeSelector::select0(runTime, args);
# include "createMesh.H"

word fieldName(args.additionalArgs()[0]);

//word x(args.additionalArgs()[1]);
//cells to be probed,obtainde from the value of x




forAll(timeDirs, timeI)
{
runTime.setTime(timeDirs[timeI], timeI);
Info<< "Time = " << runTime.timeName() << endl;

IOobject fieldHeader
(
fieldName,
runTime.timeName(),
mesh,
IOobject::MUST_READ
);

// Check field exists
if (fieldHeader.headerOk())
{
mesh.readUpdate();

/*label patchi = mesh.boundaryMesh().findPatchID(patchName);How to access any part of datas you want??? that is how to access a cell and get its ID by the coordinate ,how to access a cell by its ID?

if (patchi < 0)
{
FatalError
<< "Unable to find patch " << patchName << nl
<< exit(FatalError);
}*
/




if (fieldHeader.headerClassName() == "volScalarField")
{
Info<< " Reading volScalarField " << fieldName << endl;
volScalarField field(fieldHeader, mesh);

scalar sumField = 0;
labelList cellLabels;
pointField cellProbe;

for(int i=0;i<30;i++) // I just choose 30 points each of which stays in a cell at x=5 according to my case. There are 30 cells at x with uniform grid
{
cellProbe[i]=point(5,-0.6+i*0.8/30+0.8/60,0);
cellLabels[i]=mesh.findCell(cellProbe[i]);
sumField+=field[cellLabels[i]];
}



Info<< fieldName << " over probe "
<< " = "
<< sumField/30 << endl;
}
else
{
FatalError
<< "Only possible to average volScalarFields "
<< nl << exit(FatalError);
}
}
else
{
Info<< " No field " << fieldName << endl;
}

Info<< endl;
}

Info<< "End\n" << endl;

return 0;
}

// ************************************************** *********************** //



This's the error message:

[jingjing@jingjing stokesFirst]$ freeSurfaceElevation gamma
/*---------------------------------------------------------------------------*\
| ========= | |
| \\ / F ield | OpenFOAM: The Open Source CFD Toolbox |
| \\ / O peration | Version: 1.5 |
| \\ / A nd | Web: http://www.OpenFOAM.org |
| \\/ M anipulation | |
\*---------------------------------------------------------------------------*/
Exec : freeSurfaceElevation gamma
Date : May 28 2009
Time : 06:44:37
Host : jingjing
PID : 9971
Case : /home/jingjing/OpenFOAM/jingjing-1.5/run/othercases/surfacewaveBC/waveBC2/stokesFirst
nProcs : 1

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 0

Time = 0
Reading volScalarField gamma
#0 Foam::error:printStack(Foam::Ostream&) in "/home/jingjing/OpenFOAM/OpenFOAM-1.5/lib/linux64GccDPOpt/libOpenFOAM.so"
#1 Foam::sigSegv::sigSegvHandler(int) in "/home/jingjing/OpenFOAM/OpenFOAM-1.5/lib/linux64GccDPOpt/libOpenFOAM.so"
#2 __restore_rt at sigaction.c:0
#3 main in "/home/jingjing/OpenFOAM/jingjing-1.5/application/bin/linux64GccDPOpt/freeSurfaceElevation"
#4 __libc_start_main in "/lib64/libc.so.6"
#5 Foam::regIOobject::readIfModified() in "/home/jingjing/OpenFOAM/jingjing-1.5/applications/bin/linux64GccDPOpt/freeSurfaceElevation"
segment fault

Please give me some hint about what's wrong.Writing a appropriate probe file is the most direct way to get the wave elevation vs time ,but it seems that I'm on a long road .Hope any suggestions!
Thank you for reading this post!

jingjing

zhajingjing May 27, 2009 19:04

Quote:

Originally Posted by zhajingjing (Post 217336)
scalar sumField = 0;
labelList cellLabels;
pointField cellProbe;

for(int i=0;i<30;i++) // I just choose 30 points each of which stays in a cell at x=5 according to my case. There are 30 cells at x with uniform grid
{
cellProbe[i]=point(5,-0.6+i*0.8/30+0.8/60,0);
cellLabels[i]=mesh.findCell(cellProbe[i]);
sumField+=field[cellLabels[i]];
}

I just modified above into :

scalar sumField = 0;
label cellLabel;
scalar x=5;
scalar z=0.05;

for(int i=0;i<30;i++)
{
scalar y=-0.6+i*0.8/30;

point cellProbe(x,y,z);

cellLabel=mesh.findCell(cellProbe);

if(cellLabel==-1)
{
FatalError<<"select x again please"
<< nl << exit(FatalError);

}

/* Info<< "cellProbe point are"
<< cellProbe << endl;
Info<< "cellLabel are "
<< cellLabel << endl;*/


sumField+=field[cellLabel];
}


I got the data I want. Because I'm not quite clear about the data types ,mistakes happened.I have modified sevral times until it runs well.I have to get familiar with those basic and important class types further.

Thanks.

zhajingjing November 24, 2009 04:32

obstacle for updating the postprocessing tool
 
2 Attachment(s)
I want to improve my postprocessing tool for my wave tank to make it suitable for a wave tank with non-uniform mesh.In order to get the average value of a field with non-uniform mesh,I have to calculate the size of each cell (the heigt,length,width of a rectangular cell) as a weight factor,but I don't know how to get this kind of geomertric information of a cell.
One method is to get the coordinate value of the vertices of a cell,then make one subtract another.But I
can't go on at the detail.
Code:


  forAll(mesh.cellCentres,cellid)
    {
 
        if(abs(mesh.cellCentres[cellid].x()-x)<=deltaX/2)  //deltaX:the length of a cell, choose out the cells at the x point along the 2D wave tank
               
      {
        ncells++;
        sumField+=field[cellid]*deltaY;//deltaY: the height of a cell
        weight+=deltaY;
 
        Info<<ncells<<"cell No"<< cellid <<"field value"<<field[cellid]<<endl;
      }                       
 
    } 
 
  sumField=sumField/weight;  //the average value of the gamma field at the x point along the wave tank

The question is , how to calculate the deltaX,deltaY according to the cellid that has been chosen out.

in the following short script:
Code:

const labelList patchCells = this->patch().faceCells();
const fvPatchScalarField& gamma=patch().lookupPatchField<volScalarField,scalar >("gamma");
const pointField& points = mesh.points();
forAll(patchCells,i)
{
 labelList cellpoints=mesh.cellPoints(patchCells[i]);
 label size=cellpoints.size();
 for(int j=0,j<size,j++)    //
 {   
  points[cellpoints[j]].y()
 
 ...
 gamma[patchCells[i]]*...
    waveAtWavemaker_=...

I know something about how to access the coordinate value of a vertex of a cell,but I don't know how to access two vertices at the same time, if I could, of course I will get the height of a cell(2D)

Is there anyone can help me?Thanks


two attachments: the old version of the wave tank postprocessing tool,
and another is a new version waiting to be improved to postprocess the non-uniform mesh wave tank.

zhajingjing March 5, 2010 12:38

a new version!
 
2 Attachment(s)
In order to get the wave elevation(VS time or X)," just integrate the gamma value vertically at one point along the wave tank.", in my opinion, is the most direct and convenient method.
The following snippet is for the above:

after compiling , just use a command like this:
waveElevationVStime gamma n 5 > x5
waveElevationVSx gamma n 20 > t20

the number " 5 " means to extract the time history of the wave elevation at x=5 of the wave tank
"20" means to extract the wave profile along the wave tank at the moment of 20s
the parameter "n"/"y" is a switch ,just for more details.If you want more details,choose "y".

a simple conversion is needed when draw the profile:
elevation = averageGamma*tankHigh-waterDepth

Right now it seems enough and convenient for me. My wave tank is a rectangular with mesh vertically non-uniform and horizontally uniform. the "select X ..."statement sometimes may pops out, but it seems no bad effect. It seems also useful for wave tank with simple mesh motion.
I'm glad if it's helpful for someone.
If you find a problem, please notify me. Let's try to make it more perfect.

nuovodna September 8, 2010 10:32

it works
 
Thanks zhajingjing!! Your app works well!

owenjiang August 5, 2011 08:25

Thanks for your work. Recently I've read your master thesis and I'd like to improve your attach codes. Then I will write the FSI codes including multibody dynamics based on OF and your work.
My Email: owen9020@126.com.


All times are GMT -4. The time now is 09:26.