Velocity distribution function sampling issue in DSMC

flotus1 · December 12, 2014, 05:31

This is rather a mathematical question:
When sampling velocity distribution functions in in a Direct Simulation Monte Carlo (DSMC) we first have to choose a small velocity interval $\Delta c$ .
Now during the sampling process, all molecules with velocities between
$c - \frac{\Delta c}{2}$ and $c + \frac{\Delta c}{2}$
are counted with a velocity of

.

So after the simulation instead of the values of the distribution function for a certain velocity

we actually get the mean value over the velocity interval $\Delta c$ :
$f_\text{DSMC}(c) = \frac{1}{\Delta c} \int_{c - \frac{\Delta c}{2}}^{c + \frac{\Delta c}{2}} f(c) dc.$

Now here is the question since this sampling method is only second order accurate in $\Delta c$ .
Is there a way to get a better approximation for the actual value of the velocity distribution function of the velocity

?
I mean a other than reducing the velocity interval $\Delta c$ which would result in a higher amount of statistical scatter.
Maybe by using the neighbor values $f_\text{DSMC}(c - \Delta c)$ and $f_\text{DSMC}(c + \Delta c)$ ?

Edit:
Nevermind, I found a solution that results in a fourth order accurate interpolation of the velocity distribution function.
As soon as I get gnuplot to plot my results I will post the procedure here.

FMDenaro · December 12, 2014, 12:23

Quote:

Originally Posted by flotus1

This is rather a mathematical question:
When sampling velocity distribution functions in in a Direct Simulation Monte Carlo (DSMC) we first have to choose a small velocity interval $\Delta c$ .
Now during the sampling process, all molecules with velocities between
$c - \frac{\Delta c}{2}$ and $c + \frac{\Delta c}{2}$
are counted with a velocity of

.

So after the simulation instead of the values of the distribution function for a certain velocity

we actually get the mean value over the velocity interval $\Delta c$ :
$f_\text{DSMC}(c) = \frac{1}{\Delta c} \int_{c - \frac{\Delta c}{2}}^{c + \frac{\Delta c}{2}} f(c) dc.$

Now here is the question since this sampling method is only first order accurate in $\Delta c$ .
Is there a way to get a better approximation for the actual value of the velocity distribution function of the velocity

?
I mean a other than reducing the velocity interval $\Delta c$ which would result in a higher amount of statistical scatter.
Maybe by using the neighbor values $f_\text{DSMC}(c - \Delta c)$ and $f_\text{DSMC}(c + \Delta c)$ ?

Edit:
Nevermind, I found a solution that results in a fourth order accurate interpolation of the velocity distribution function.
As soon as I get gnuplot to plot my results I will post the procedure here.

Hello,
I am not sure to get correctly your question but the average value when c is centered between +dc and -dc is a second order approximation of the pointwise value.
It is similar to what happens in FVM.
Higher order approximations require suitable reconstruction/deconvolution techniques

flotus1 · December 12, 2014, 12:39

You are right, I tend to confuse orders of accuracy.
The thing is that for the purpose of the simulation involved, the second order accurate approximation was not accurate enough.

If the remedy is so obvious I will rather not post my solution because I feel like I reinvented the wheel now.

FMDenaro · December 12, 2014, 12:53

Quote:

Originally Posted by flotus1

You are right, I tend to confuse orders of accuracy.
The thing is that for the purpose of the simulation involved, the second order accurate approximation was not accurate enough.

If the remedy is so obvious I will rather not post my solution because I feel like I reinvented the wheel now.

you can find some techniques in the field of FVM since the flux reconstruction is based on pointwise values, that is F=F(U), whereas the FV update is expressed in terms of the average value U_bar. Therefore, at each time step, some functional relation U=U(U_bar) has to be prescribed. U=U_bar is order 2. Accuracy higher than 2 can be obtained in some way.
A way can be developing in Taylor series the pointwise function and integrate each terms.

flotus1 · December 13, 2014, 06:05

For the sake of completeness, here is what I did.
The sampled DSMC value equals the following integral:

$f_\text{DSMC}(c) = \frac{1}{\Delta c} \int_{c - \frac{\Delta c}{2}}^{c + \frac{\Delta c}{2}} f(c) dc = \frac{1}{\Delta c} \left[ F\left(c+\frac{\Delta c}{2} \right) + F \left( c-\frac{\Delta c}{2} \right) \right].$

The sum of this value and its two neighbors is

$f_\text{DSMC}(c-\Delta c) + f_\text{DSMC}(c) + f_\text{DSMC}(c+\Delta c) = \frac{1}{\Delta c} \int_{c - \frac{3 \Delta c}{2}}^{c + \frac{3 \Delta c}{2}} f(c) dc$ $= \frac{1}{\Delta c} \left[ F \left( c+\frac{3 \Delta c}{2} \right) + F \left( c-\frac{3 \Delta c}{2} \right) \right].$

Using a 4th order central difference formula for the first derivative we can reconstruct the value of

:

$f(c) = \frac{d}{dc} F(c)$ $= \frac{1}{\Delta c} \left[ \frac{1}{24}F \left( c-\frac{3 \Delta c}{2} \right) - \frac{9}{8}F \left( c-\frac{\Delta c}{2} \right) + \frac{9}{8}F \left( c+\frac{\Delta c}{2} \right) - \frac{1}{24}F \left( c+\frac{3 \Delta c}{2} \right) \right]$ $+ \mathcal{O} \left( \Delta c^4 \right).$

Here we can substitute the DSMC values:

$f(c) = \frac{d}{dc} F(c) = -\frac{1}{24} \left[ f_\text{DSMC}(c-\Delta c) + f_\text{DSMC}(c) + f_\text{DSMC}(c+\Delta c) \right] + \frac{9}{8} f_\text{DSMC}(c)$ $+ \mathcal{O} \left( \Delta c^4 \right)$

which results in

$f(c) = \frac{d}{dc} F(c) = -\frac{1}{24} \left[ f_\text{DSMC}(c-\Delta c) + f_\text{DSMC}(c+\Delta c) \right] + \frac{13}{12} f_\text{DSMC}(c)$ $+ \mathcal{O} \left( \Delta c^4 \right).$

This does exactly what I want and effectively removes the sampling bias I observed in the simulations.

FMDenaro · December 13, 2014, 06:16

If I understand, you are using reconstruction via primitive function, something like the reconstruction used in the ENO/WENO scheme.

You can also proceed this way:

f(c) = f(0) + f'(0)*c + f''(0)*c^/2 + f'''(0)*c^3/6 + ....

Then, integrate analytically each terms between +Deltac/2 and -Deltac/2.
You see that terms like c, c^3 produce zero value of the integral.
This way you have a relation between the average value (LHS) and pointwise values (RHS). You can truncate at suitable order that you need

FMDenaro · December 13, 2014, 06:18

PS: of course, I suppose that f is a continuous function of c, having all regular derivatives .....

December 12, 2014, 05:31	Velocity distribution function sampling issue in DSMC	#1
flotus1 Super Moderator Alex Join Date: Jun 2012 Location: Germany Posts: 3,399 Rep Power: 46	This is rather a mathematical question: When sampling velocity distribution functions in in a Direct Simulation Monte Carlo (DSMC) we first have to choose a small velocity interval $\Delta c$ . Now during the sampling process, all molecules with velocities between $c - \frac{\Delta c}{2}$ and $c + \frac{\Delta c}{2}$ are counted with a velocity of . So after the simulation instead of the values of the distribution function for a certain velocity we actually get the mean value over the velocity interval $\Delta c$ : $f_\text{DSMC}(c) = \frac{1}{\Delta c} \int_{c - \frac{\Delta c}{2}}^{c + \frac{\Delta c}{2}} f(c) dc.$ Now here is the question since this sampling method is only second order accurate in $\Delta c$ . Is there a way to get a better approximation for the actual value of the velocity distribution function of the velocity ? I mean a other than reducing the velocity interval $\Delta c$ which would result in a higher amount of statistical scatter. Maybe by using the neighbor values $f_\text{DSMC}(c - \Delta c)$ and $f_\text{DSMC}(c + \Delta c)$ ? Edit: Nevermind, I found a solution that results in a fourth order accurate interpolation of the velocity distribution function. As soon as I get gnuplot to plot my results I will post the procedure here. Last edited by flotus1; December 12, 2014 at 12:42.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
[snappyHexMesh] How to define to right point for locationInMesh	Mirage12	OpenFOAM Meshing & Mesh Conversion	7	March 13, 2016 14:07
get velocity gradient with a boundary function	Chris Lucas	OpenFOAM Programming & Development	4	August 27, 2012 03:11
Moving Mesh Velocity Issue: Mesh velocity does not equal displacement	Doginal	CFX	2	September 8, 2011 12:02
channelFoam for a 3D pipe	AlmostSurelyRob	OpenFOAM	3	June 24, 2011 13:06
[blockMesh] BlockMesh FOAM warning	gaottino	OpenFOAM Meshing & Mesh Conversion	7	July 19, 2010 14:11

December 12, 2014, 12:39		#3
flotus1 Super Moderator Alex Join Date: Jun 2012 Location: Germany Posts: 3,399 Rep Power: 46	You are right, I tend to confuse orders of accuracy. The thing is that for the purpose of the simulation involved, the second order accurate approximation was not accurate enough. If the remedy is so obvious I will rather not post my solution because I feel like I reinvented the wheel now.

December 13, 2014, 06:05		#5
flotus1 Super Moderator Alex Join Date: Jun 2012 Location: Germany Posts: 3,399 Rep Power: 46	For the sake of completeness, here is what I did. The sampled DSMC value equals the following integral: $f_\text{DSMC}(c) = \frac{1}{\Delta c} \int_{c - \frac{\Delta c}{2}}^{c + \frac{\Delta c}{2}} f(c) dc = \frac{1}{\Delta c} \left[ F\left(c+\frac{\Delta c}{2} \right) + F \left( c-\frac{\Delta c}{2} \right) \right].$ The sum of this value and its two neighbors is $f_\text{DSMC}(c-\Delta c) + f_\text{DSMC}(c) + f_\text{DSMC}(c+\Delta c) = \frac{1}{\Delta c} \int_{c - \frac{3 \Delta c}{2}}^{c + \frac{3 \Delta c}{2}} f(c) dc$ $= \frac{1}{\Delta c} \left[ F \left( c+\frac{3 \Delta c}{2} \right) + F \left( c-\frac{3 \Delta c}{2} \right) \right].$ Using a 4th order central difference formula for the first derivative we can reconstruct the value of : $f(c) = \frac{d}{dc} F(c)$ $= \frac{1}{\Delta c} \left[ \frac{1}{24}F \left( c-\frac{3 \Delta c}{2} \right) - \frac{9}{8}F \left( c-\frac{\Delta c}{2} \right) + \frac{9}{8}F \left( c+\frac{\Delta c}{2} \right) - \frac{1}{24}F \left( c+\frac{3 \Delta c}{2} \right) \right]$ $+ \mathcal{O} \left( \Delta c^4 \right).$ Here we can substitute the DSMC values: $f(c) = \frac{d}{dc} F(c) = -\frac{1}{24} \left[ f_\text{DSMC}(c-\Delta c) + f_\text{DSMC}(c) + f_\text{DSMC}(c+\Delta c) \right] + \frac{9}{8} f_\text{DSMC}(c)$ $+ \mathcal{O} \left( \Delta c^4 \right)$ which results in $f(c) = \frac{d}{dc} F(c) = -\frac{1}{24} \left[ f_\text{DSMC}(c-\Delta c) + f_\text{DSMC}(c+\Delta c) \right] + \frac{13}{12} f_\text{DSMC}(c)$ $+ \mathcal{O} \left( \Delta c^4 \right).$ This does exactly what I want and effectively removes the sampling bias I observed in the simulations.

December 13, 2014, 06:16		#6
FMDenaro Senior Member Filippo Maria Denaro Join Date: Jul 2010 Posts: 6,764 Rep Power: 71	If I understand, you are using reconstruction via primitive function, something like the reconstruction used in the ENO/WENO scheme. You can also proceed this way: f(c) = f(0) + f'(0)c + f''(0)c^/2 + f'''(0)*c^3/6 + .... Then, integrate analytically each terms between +Deltac/2 and -Deltac/2. You see that terms like c, c^3 produce zero value of the integral. This way you have a relation between the average value (LHS) and pointwise values (RHS). You can truncate at suitable order that you need

December 13, 2014, 06:18		#7
FMDenaro Senior Member Filippo Maria Denaro Join Date: Jul 2010 Posts: 6,764 Rep Power: 71	PS: of course, I suppose that f is a continuous function of c, having all regular derivatives .....