CFD Online Discussion Forums

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)
-   OpenFOAM (https://www.cfd-online.com/Forums/openfoam/)
-   -   access to velocity gradient for Lagrangian particles (https://www.cfd-online.com/Forums/openfoam/74452-access-velocity-gradient-lagrangian-particles.html)

jiejie March 31, 2010 19:28

access to velocity gradient for Lagrangian particles
 
Dear all

I am kind of get stuck at finding the velocity gradient at a particular position for the Lagrangian particle dynamics. My question is how I should calculate the velocity gradient if I have the velocity vector "U_". Is there some utility like "curl" that I can use to calculate it?

Thank you very much.

su_junwei April 1, 2010 14:30

Hi jie
It seems not. And it seems that the gradient is only reasonable for a Euler field. A possible way is calculating the gradient from an Eulerian field converted from particle lagrangian velocity field
Junwei

jiejie April 1, 2010 21:33

Quote:

Originally Posted by su_junwei (Post 252787)
Hi jie
It seems not. And it seems that the gradient is only reasonable for a Euler field. A possible way is calculating the gradient from an Eulerian field converted from particle lagrangian velocity field
Junwei

Yes I am trying to get the gradient for Euler Field and use it to calculate the particle movement. Do you have any idea how to calculate the Euler field gradient?

su_junwei April 1, 2010 23:05

Oh sorry.

I am confused last night. Just use fvc::grad(U_) and make an interpolation on the particle location like this

volTensorField gradU=fvc::grad(U_);
autoPtr<interpolation<tensor> >gradUInterpolator_ =
interpolation<tensor>::New(interpolationSchemes, gradU);
//interpolationSchemes is a word for assigning interpolation schemes("cell", "cellPoint" or "cellPointFace")

tensor gradUatPoint = gradUInterpolator->interpolate(ball.position(), ball.cell());

If the gradient vary little across the space, interpolation can be neglected, and use the gradient at the cell center to represent the gradient at the particle location. In this situation, just use the following code
tensor gradUatPoint =gradU[ball.cell()];
ball is the particle defined in OpenFOAM.

Hope this helps

regards, Junwei Su

sankarv April 6, 2010 08:54

Sampling tool for Lagrangian Particles
 
Dear su_junwei

Is there a sampling tool for Lagrangian particles ? i.e, I would
like to collect particle statistics such as mean particle size, mean
particle velocity etc., at an arbitrary surface inside the domain
DURING the simulation ?

I saw from another post that there is an utility called
functionObject for doing this for Eulerian fields. I am looking for the
lagrangian particles.

Please let me know

Thanks
Vaidya

su_junwei April 6, 2010 09:31

Hi Vaidya

I don't think there is a tool in OpenFOAM which can sample particles on a arbitrary surface in the solution domain. I didn't find one. Actually, at a certain time point, there may be no particles hitting the surface concerned exactly. If you want to get a distribution(of particle size, velocity, etc) on a certain plane, I think you'd better make a conversion from Lagrangian field to Eulerian field and then use the sampling tool in OpenFOAM.

The sample utility in OpenFOAM is located at
../OpenFOAM/OpenFOAM-1.6/src/sampling

Regards, Junwei

sankarv April 6, 2010 10:01

Thanks Su-junwei. I will keep looking. Hopefully the developers Hrvoje Jasak or Henry Weller might have some idea. But unfortunately no body else seems to know

Thanks
Vaidya

jiejie April 6, 2010 22:00

Quote:

Originally Posted by su_junwei (Post 252836)
Oh sorry.

volTensorField gradU=fvc::grad(U_);
autoPtr<interpolation<tensor> >gradUInterpolator_ =
interpolation<tensor>::New(interpolationSchemes, gradU);
//interpolationSchemes is a word for assigning interpolation schemes("cell", "cellPoint" or "cellPointFace")

tensor gradUatPoint = gradUInterpolator->interpolate(ball.position(), ball.cell());

Hi JunWei

I am sorry about the confusion caused, the "U_" is the parcel velocity and it is a vector but not a vector field. Hence, I got some error complaining about "volTensorField gradU=fvc::grad(U_);"

I assume that maybe I could calculate the velocity gradient with the velocity filed U and interpolate the velocity gradient at the parcel position.

Thanks a lot

su_junwei April 6, 2010 22:48

Hi Jie
What force do you want to calculate for a particle, I haven't encountered a force for granular flow. Are you trying to implement a meshless method?
I usually make a conversion from the lagrangian field to euler field when calculating the Euler properties of granular fields.
If you don't want a conversion, you can search all the particles around the particle concerned, calculate the gradient and do a average. This is common in meshless method, the formulation can be found in text book about meshless method.

Regards, Junwei

jiejie April 6, 2010 23:55

Hi JunWei

I am trying to calculate the drag and lift force acting on the particle so that I need the velocity gradients from the Eulerian field at the position of the particle.

I will try to create a volVectorField GradU that contains all the velocity gradients and interpolate them on to the position of the particle.

One quick question about the createFields.H:

Info << "Reading field U\n" << endl;
volVectorField U
(
IOobject
(
"U",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);

It creates and initialises the volVectorField U and U is calculated from other mechanisms, right?

So I can just do the followings to create GradU?

Info << "Reading field GradU\n" << endl;
volVectorField GradU
(
IOobject
(
"GardU",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);

However, there are 9 components of the velocity gradients. Should I create and initialise each of them as GradU1, GradU2 and etc?

Or should I use the "voltensorFields"?

volTensorField gradU(fvc::grad(U));

I am pretty new to OpenFOAM so I am still learning myself.

Thanks

su_junwei April 7, 2010 00:29

Quote:

Originally Posted by jiejie (Post 253462)
Hi JunWei

volVectorField GradU
(
IOobject
(
"GardU",
runTime.timeName(),
mesh,
IOobject::MUST_READ,
IOobject::AUTO_WRITE
),
mesh
);

Thanks


Hi Jie

you can do it like the following code (you don't have to initialize it using the file)

volTensorField GradU
(
IOobject
(
"GardU",
runTime.timeName(),
mesh,
IOobject::NO_READ,
IOobject::AUTO_WRITE
),
fvc::grad(U)
);

PS: Don't forget recalculate gradU before using it.

Regards,Junwei

jiejie April 7, 2010 00:40

That means at each point of the tensor field of GradU, it should have 9 velocity gradient components as GradU.xx(), GradU.xy(), GradU.xz() ...?

su_junwei April 7, 2010 01:10

exactly

Junwei

jiejie April 8, 2010 22:44

Hi JunWei

I managed to get gradU calculated and updated during the simulation. However, I have problem to access gradU. All my particle dynamics is defined in a file called "mySolidParticle.C". It gives some error of "mySolidParticle.C:83: error: ‘gradU’ was not declared in this scope" when I tried to compile my solver after I added in the following code:

autoPtr<interpolation<tensor> >gradUInterpolator_ =
interpolation<tensor>::New(cell, gradU);

How should I access the tensor field gradU other than the main "SOLVER.C" (SOLVER.C contains the main mechanism for solver)?

Thank you

su_junwei April 9, 2010 23:49

Hi Jie
You don't have to declare the gradU in createField.H. You just get the value where you use it. In your case, declare it before the interpolator for instance like this

volTensorField gradU=fvc::grad(U_);
autoPtr<interpolation<tensor> >gradUInterpolator_ =
interpolation<tensor>::New(cell, gradU);

If you U_ is not accessible. you can lookup it through ObjectRegistry using the the function
const volVectorField &U=[db or simular, you can get ObjectRegistry].lookupObject<volVectorField>("U");
you can also lookup gradU, if it is declared elsewhere createFields.H for instance.
const volTensorField & gradU=[db...].lookupObject<volTensorField>("gradU");

Regards, Junwei

jiejie April 10, 2010 00:34

Thanks JunWei

I will give the lookup it through ObjectRegistry a try.

jiejie April 20, 2010 01:49

Hi JunWei

Now I am trying to write out the gradU for each particle. So I add "gradU.write();" in the "solidparticleIO.C". This should create a file at /xxx/lagrangian/defaultCloud/gradU where xxx is the time step at which OpenFOAM output the data.

At /xxx/lagrangian/defaultCloud/, there are files of d, positions, U, gradU where d is the diameter of the particles, position is the locations of the particles, U is the velocity at the location of the particles, gradU should be the velocity gradients at the location of the particles.

In my d and U fields, they contain the corresponding information of all the particles in the flow field. However, it only contains the velocity gradients information for one particle at the location where it releases the particle in the flow field in gradU.

Would you have any idea about this?

sankarv April 20, 2010 08:58

how to set correct mass flow rate for lagrangian particles ?
 
Dear Junwei

I am running a simulation where a simple nozzle is injecting compressible gas + solid particles jet into an open atmosphere. I am using OF 1.6.x

I know the mass flow rate of the particles, density of the particles and the particle diameter distribution from the experiments.

I am using the patchInjection option to inject particles from the inlet patch. The code works without a problem

But I am not quite sure how to set the mass flow rate of the particles correctly in the kinematicCloud1Properties file

My simple question is this:

Can you please tell me what is the relationship between
parcelsPerSecond, massTotal, volumeflowRate, SOI and duration ?

Please help me

Thanks
Vaidya

su_junwei April 20, 2010 11:33

Quote:

Originally Posted by jiejie (Post 255369)
Hi JunWei

In my d and U fields, they contain the corresponding information of all the particles in the flow field. However, it only contains the velocity gradients information for one particle at the location where it releases the particle in the flow field in gradU.

Would you have any idea about this?

Perhaps. you forgot to update gradU every time step. Just recalculate gradU every time step. add the following line to time loop
gradU=fvc::grad(U) ;

regards,Junwei

su_junwei April 20, 2010 12:00

Quote:

Originally Posted by sankarv (Post 255447)
Dear Junwei
Can you please tell me what is the relationship between
parcelsPerSecond, massTotal, volumeflowRate, SOI and duration ?

Please help me

Thanks
Vaidya

Hi Vaidya
The meaning of these keywords may be

parcelsPerSecond: The total number of particles should be injected a second after start starttime of Injection(SOI)

massTotal: the total mass of particles should be injected between SOI and SOI+duration.

volumeFlowRate: the total volume of the particles should be injected a second after SOI

SOI: start of injection

duration: The injection time segment. Injecting particle from SOI to SOI+duration

regards,Junwei

jiejie April 20, 2010 19:50

Quote:

Originally Posted by su_junwei (Post 255478)
Perhaps. you forgot to update gradU every time step. Just recalculate gradU every time step. add the following line to time loop
gradU=fvc::grad(U) ;

regards,Junwei

Hi JunWei

I updated gradU every time step as I have gradU in the data set for each time step. However, the there is no gradU stored in the lagrangian folder. As I have the position data and gradU field data, I want to use the particle position information to check the gradU stored in the lagrangain particle is correct.

su_junwei April 20, 2010 23:05

Hi Jie
I cannot quite following you. Do you want to store the gradU at particle centers into the Lagrangian folder or gradU has been stored but the data was not correct.
If they was not stored, just add a new variable for your particle class and add this IO utility for it like other properties in particle class. The IO utility is implemented in the following function

1) constructor function like this

Foam::solidparticle::solidparticle
(
const Cloud<solidparticle>& cloud,
Istream& is,
bool readFields
)
2) readFields
3) writeFields
4) operator<<

If gradU has been written into lagrangian folder but not correct. I think they are not updated. Note that the volTensorField gradU you declare in your createFields.H is not the field of gradU in particle center. The first is the gradU in the patchcenters and cellCenters, whereas the second is the gradu in several discrete points.

regards, Junwei

jiejie April 20, 2010 23:49

Hi JunWei

Sorry about the confusion. I tried to store the gradU at particle centers into the Lagrangian folder, but it only contains the fist particle.

I added the following lines in the constructor:

######################################
Foam::solidParticle::solidParticle
(
const Cloud<solidParticle>& cloud,
Istream& is,
bool readFields
)
:
Particle<solidParticle>(cloud, is, readFields)
{
if (readFields)
{
if (is.format() == IOstream::ASCII)
{
d_ = readScalar(is);
is >> U_;
//- ** velocity gradients
is >> gradU_;
}
else
{
is.read
(
reinterpret_cast<char*>(&d_),
//sizeof(d_) + sizeof(U_)

//- **velocity gradients
sizeof(d_) + sizeof(U_) + sizeof(gradU_)
);
}
}

// Check state of Istream
is.check("solidParticle::solidParticle(Istream&)") ;
}
#####################################

and I asked it to write out the gradU

####################################
void Foam::solidParticle::writeFields(const Cloud<solidParticle>& c)
{
Particle<solidParticle>::writeFields(c);

label np = c.size();

IOField<scalar> d(c.fieldIOobject("d", IOobject::NO_READ), np);
IOField<vector> U(c.fieldIOobject("U", IOobject::NO_READ), np);

//- **velocity gradients
IOField<tensor> gradU(c.fieldIOobject("gradU", IOobject::NO_READ), np);

label i = 0;
forAllConstIter(Cloud<solidParticle>, c, iter)
{
const solidParticle& p = iter();

d[i] = p.d_;
U[i] = p.U_;

//- **velocity gradients
gradU[i] = p.gradU_;

i++;
}

d.write();
U.write();

//- **velocity gradients
gradU.write();
}
####################################

in the lagrangain folder, gradU only contains:

FoamFile
{
version 2.0;
format ascii;
class tensorField;
location "101.3";
object gradU;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

1301{(0 0 0 0 0 0 0 0 0)}

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

Instead of having all zeros in 1301{(0 0 0 0 0 0 0 0 0 0)}, it should have the gradU for all 1301 particles. I have no problem with U.

jiejie April 21, 2010 01:19

Hi JunWei

I have attached my solidparticleIO.C here.

It looks like it counts the number of particles in the gradU in lagrangian folder, but all the gradU for each particle is not there.

su_junwei April 21, 2010 05:02

Hi Jie
1301{(0 0 0 0 0 0 0 0 0 0)} means there are 1301 particles with the same value. tensor::zero. The possible reason is that the gradU of the particle location was not updated. The attachment you provided is right I think.
Regards, Junwei

jiejie April 22, 2010 20:27

Hi JunWei

I am still looking at how to update the gradU for each of the particle and I came across couple question as follows:

1. if gradU has 1301{(0 0 0 0 0 0 0 0 0 0)} at each time step, it means all the gradU stored in the particles are zero?

2. In the particle.C:

template<class ParticleType>
00411 void Foam::Particle<ParticleType>::transformPosition(const tensor& T)
00412 {
00413 position_ = transform(T, position_);
00414 }

It says this transforms the physical properties of the particle.
But in parcel.C

void Foam:: parcel::transformProperties(const tensor& T)
00689 {
00690 U_ = transform(T, U_);
00691 }

It says this transforms the position and physical properties of the particle.
How come it uses the velocity U_ to transform the position of the particle and it use position_ to transform the position of the particle?

Thank you very much.

su_junwei April 24, 2010 04:19

1) if gradU has 1301{(0 0 0 0 0 0 0 0 0 0)} at each time step, it means all the gradU stored in the particles are zero?
Yes
2) transform can be used in rotation. T is rotational tensor.

Junwei

jiejie April 27, 2010 01:39

figuring out how to update gradU.

Thanks for your time JunWei

jiejie May 10, 2010 03:01

One quick question: if I need to do an inner product of a tensor "T" and a vector "a", should it be "T & a" instead of "a & T"?

ehsan December 13, 2010 14:34

How to find gradU from tauMC
 
Dear All,

tauMC is calculated via:

volTensorField tauMC("tauMC",mu*dev2(fvc::grad(U)().T()));

Could someonhe help me to find gard(U) from tauMC (revers of above relation)

the notation gard(U)=tauMC/mu gives error.

Thanks

francescomarra January 26, 2011 10:57

Further variables in lagrangian folder
 
Dear All,

After running dieselFoam, I would be like to observe further properties of the lagrangian field, as the Weber number of each particle, similarly to the way the diameter d or the density rho of each particle is reported in the lagrangian folder at each time step. However I am a newbe with C++ and I do not understand how I can adopt the function We available in the parcelFunctions.C to compute and store this quantity for each particle. Is this a simple task or does it require a long coding effort ?
Any help will be very gratefully appreciated.

Thank you in advance,

Franco

ycui December 2, 2016 04:56

Dear Su Jun Wei,

I was trying to use your way to calculate the velocity gradient tensor.

It works fine with one core. But error comes for parallel computing when try to access fvc::grad().

Do you have any idea?

Best,

ycui

Quote:

Originally Posted by su_junwei (Post 252836)
Oh sorry.

I am confused last night. Just use fvc::grad(U_) and make an interpolation on the particle location like this

volTensorField gradU=fvc::grad(U_);
autoPtr<interpolation<tensor> >gradUInterpolator_ =
interpolation<tensor>::New(interpolationSchemes, gradU);
//interpolationSchemes is a word for assigning interpolation schemes("cell", "cellPoint" or "cellPointFace")

tensor gradUatPoint = gradUInterpolator->interpolate(ball.position(), ball.cell());

If the gradient vary little across the space, interpolation can be neglected, and use the gradient at the cell center to represent the gradient at the particle location. In this situation, just use the following code
tensor gradUatPoint =gradU[ball.cell()];
ball is the particle defined in OpenFOAM.

Hope this helps

regards, Junwei Su



All times are GMT -4. The time now is 10:05.