[wenjie]

Hi all,

I think I just found a possible error in the DsmcCloud::collisions() in DsmcCloud.C .

There, in the number of collision pairs, it should be nC * nC instead of nC* (nC-1) according to Bird's book (p219).

Does anybody agree or disagree with me?

-Wenjie

[wyldckat]

Greetings Wenjie,

I had a look at the source code and took me a while to spot the piece of code you were referring to:

            scalar selectedPairs =
                collisionSelectionRemainder_[cellI]
              + 0.5*nC*(nC - 1)*nParticle_*sigmaTcRMax*deltaT
               /mesh_.cellVolumes()[cellI];

            label nCandidates(selectedPairs);

            collisionSelectionRemainder_[cellI] = selectedPairs - nCandidates;

To me, the use of "nC-1" is because the current particle isn't counted in possible collisions. Nonetheless, I have no idea which exact book you're referring to
It's possible that the implementation done in OpenFOAM is either based in a more recent paper or version of said book, or it could be the other way around, where they could have based themselves only on the very first edition of the book.

In addition, according to this git blame: https://github.com/wyldckat/OpenFOAM...ud/DsmcCloud.C - that specific part of the code you're referring to exists since the first commit for OpenFOAM 1.6.x, so there is the possibility that this is still pending some additional revision.

Best regards,
Bruno

[wenjie]

Originally Posted by wyldckat

Greetings Wenjie,

I had a look at the source code and took me a while to spot the piece of code you were referring to:

            scalar selectedPairs =
                collisionSelectionRemainder_[cellI]
              + 0.5*nC*(nC - 1)*nParticle_*sigmaTcRMax*deltaT
               /mesh_.cellVolumes()[cellI];

            label nCandidates(selectedPairs);

            collisionSelectionRemainder_[cellI] = selectedPairs - nCandidates;

To me, the use of "nC-1" is because the current particle isn't counted in possible collisions. Nonetheless, I have no idea which exact book you're referring to
It's possible that the implementation done in OpenFOAM is either based in a more recent paper or version of said book, or it could be the other way around, where they could have based themselves only on the very first edition of the book.

In addition, according to this git blame: https://github.com/wyldckat/OpenFOAM...ud/DsmcCloud.C - that specific part of the code you're referring to exists since the first commit for OpenFOAM 1.6.x, so there is the possibility that this is still pending some additional revision.

Best regards,
Bruno

Hi Bruno,

Yes, that's the code that I was talking about. Indeed by using nC*(nC-1) instead of nC*nC, you are not counting the current particle, but it is explained in the Bird's book "Molecular Gas Dynamics and the DIrect Simulation of Gas FLows" on page 219 that the number of choices should be N^2 /2 instead of N(N-1) /2 for large number of gas molecules and this error becomes significant for small values of N. And the whole rest part is about how to increase the probability and decrease the fraction.

So I still think that nC*nC should be the correct one to use here. What do you think?

Best,
Wenjie

[DSMC123]

It has to be N*(N-1) for the no time counter model to be correct. See section 7 in this document from Graham Bird: http://www.gab.com.au/dsmc07notes.pdf

In addition, a quick look at Bird's 1D source code confirms that he uses the N*(N-1). This is from line 5989 of DS1.f90:

ASEL=0.5D00*ICCELL(2,N)*(ICCELL(2,N)-1)*WFC*FNUM*CCELL(4,N)*DTC/AAA+CCELL(2,N)

Where ICCELL is the number of particles in the cell.

[wenjie]

Originally Posted by DSMC123

It has to be N*(N-1) for the no time counter model to be correct. See section 7 in this document from Graham Bird: http://www.gab.com.au/dsmc07notes.pdf

In addition, a quick look at Bird's 1D source code confirms that he uses the N*(N-1). This is from line 5989 of DS1.f90:

ASEL=0.5D00*ICCELL(2,N)*(ICCELL(2,N)-1)*WFC*FNUM*CCELL(4,N)*DTC/AAA+CCELL(2,N)

Where ICCELL is the number of particles in the cell.

HI Bruno,

Thanks for the document!. Now it's clear : )