simulateWithTwoWayCoupling

[Rookie]

Dear OpenLB developers,

Can you tell me the difference between simulateWithTwoWayCoupling_Mathias and simulateWithTwoWayCoupling_Davide, forwardCoupling is the force exerted by the fluid on the particle, backCoupling is the force exerted by the particle on the fluid, subSteps is the total number of simulated time steps, resetExternalField is the reset of the external force field, should I set it to false if I want to retain the original force exerted by the fluid?Please let me know if I’m wrong.

template<typename T, template<typename U> class PARTICLETYPE>
void SuperParticleSystem3D<T, PARTICLETYPE>::simulateWithTwoWayCoupling_Mathias ( T dT,
ForwardCouplingModel<T,PARTICLETYPE>& forwardCoupling,
BackCouplingModel<T,PARTICLETYPE>& backCoupling,
int material, int subSteps, bool resetExternalField, bool scale )
{
// reset external field
if (resetExternalField) {
backCoupling.resetExternalField(material);
}

for (auto pS : _pSystems) {
time_t delta = clock();
pS->_contactDetection->sort();
_stopSorting += clock() – delta;
pS->simulateWithTwoWayCoupling_Mathias(dT, forwardCoupling, backCoupling, material, subSteps, scale);
pS->computeBoundary();

}
updateParticleDistribution();
}

template<typename T, template<typename U> class PARTICLETYPE>
void SuperParticleSystem3D<T, PARTICLETYPE>::simulateWithTwoWayCoupling_Davide ( T dT,
ForwardCouplingModel<T,PARTICLETYPE>& forwardCoupling,
BackCouplingModel<T,PARTICLETYPE>& backCoupling,
int material, int subSteps, bool resetExternalField, bool scale )
{
for (auto pS : _pSystems) {
time_t delta = clock();
pS->_contactDetection->sort();
_stopSorting += clock() – delta;
pS->simulateWithTwoWayCoupling_Davide(dT, forwardCoupling, backCoupling, material, subSteps, scale);
pS->computeBoundary();
pS->eraseInactiveParticles();
}
updateParticleDistribution();
backCoupling.communicate();
} pS->computeBoundary();

}
updateParticleDistribution();
}

Regards,
Rookie

[jan]

Dear Rookie,

Since this is old legacy code, my knowledge of it is limited, but I’ll try to help.

The posted source code shows that they’re mostly the same, but the _Davide code adds a few more functions at the end. The code deletes inactive particles, communicates the back coupling, and also updates the particle distributions twice. Why this is, I do not know. It could be that eraseInactiveParticles is also called here to avoid a second loop over all particles. Also, the _Davide code somehow doesn’t include an option to reset the back coupling forces.

For now, I’d suggest you try them both and see if you see any differences in the results. Do they both do the same thing in your application? Then it’s probably more efficient to use the first one. If the results are different, then you may need the second, as this could indicate missing communications.

The parameters and your explanations seem correct. If you don’t want the reset, then yes, you should set resetExternalField to false.

Best regards,
Jan

[Rookie]

Dear Jan,

I am very glad to continue to communicate with you about this aspect of learning. I’m learning this part of the code, and the lack of explanation really bothers me. If I implement this I will be happy to share it here. But now I want to add this two coupling and how to set the parameters that it needs.

template<typename T, template<typename U> class PARTICLETYPE>
void SuperParticleSystem3D<T, PARTICLETYPE>::simulateWithTwoWayCoupling_Davide ( T dT,
ForwardCouplingModel<T,PARTICLETYPE>& forwardCoupling,
BackCouplingModel<T,PARTICLETYPE>& backCoupling,
int material, int subSteps, bool resetExternalField, bool scale )

I thought it should be set up like the addforce, so I added this line of code like this, I’m not sure if I added the right position (after supParticleSystem.addForce or after supParticleSystem.simulate).I am also not sure on what basis should I judge the choice between Local and SmoothingFunctional.

If I set up the following code, it will report an error saying that converter cannot be passed because of a constant, and I don’t know how to continue to modify it.

auto dragModel = std::make_shared<SchillerNaumannDragModel<T, DESCRIPTOR, PARTICLE>>(converter);
LaddForwardCouplingModel<T, DESCRIPTOR, PARTICLE> forwardCoupling(converter, sLattice, superGeometry, dragModel);
int overlap = 2. * converter.getConversionFactorLength();
LocalBackCouplingModel<T, DESCRIPTOR, PARTICLE> backCoupling(converter, sLattice, superGeometry, overlap);
int subSteps = converter.getLatticeTime(particleMaxPhysT);
supParticleSystem.simulateWithTwoWayCoupling_Mathias ( dT, forwardCoupling, backCoupling, 1, subSteps, false, true);

Best regards,
Rookie

[jan]

Dear Rookie,

what does the error say exactly?

I assume that the UnitConverter is declared constant in your app. However, the drag model expects a reference to a non-constant UnitConverter: SchillerNaumannDragModel(UnitConverter<T, Lattice>& converter);. The constructor of LaddForwardCouplingModel has the same problem. You could try to replace UnitConverter<T, Lattice>& converter with UnitConverter<T, Lattice> const& converter in the corresponding constructors or remove the const keyword from the line where you create the UnitConverter object.

Best regards,
Jan

[Rookie]

Dear Jan,

I also wanted to remove const, but I don’t know if it will affect the rest of the converter. Like you said, I changed this place to：UnitConverterFromResolutionAndLatticeVelocity<T, DESCRIPTOR> converter

And then I changed it based on the error report：int overlap = 2;The original code had such a definition： supParticleSystem.setOverlap(2. * converter.getConversionFactorLength());

I understand it as overlapping area, can you explain the meaning of this in detail, when I analyze the data in the fluid part, the two layers of data in the middle of the Z-axis will be repeated, I initially thought it was caused by parallel Settings, I don’t know what the correlation is, and I set it equal to 2 correctly?

Best regards,
Rookie

[jan]

Dear Rookie,

removing the const keyword shouldn’t have any impact as long as you don’t try to manipulate the values that are stored in it.

Could you please clarify what you replaced with int overlap = 2; and why?

The overlap defines the size of the area in which the communication of data happens. For communication, the blocks must so that storage for the necessary neighborhood data is provided and accessible. Usually an overlap of 2 is sufficient.

Best regards,
Jan

[Rookie]

Dear Jan,

I understand your explanation, I also wanted to use 2 lattices at the beginning, so I set this at the beginning;int overlap = 2. * converter.getConversionFactorLength(); and saw the error report that it must be greater than 1, so I changed it to 2, here it is 2 lattices.

My result is not out yet, I will tell you the first time after it is out, you can also check if there is any problem with my above parameter Settings, I have chosen a local way for coupling, I am not sure whether the local is more stable, and then I will try other results to share here.

Best regards,
Rookie

[Rookie]

Dear Jan,

I’m sorry to bother you again, can you point me to the literature on these two couplings, I’m trying to understand these based on the literature and the code.

simulateWithTwoWayCoupling_Mathias
simulateWithTwoWayCoupling_Davide

Best regards,
Rookie

[jan]

Dear Rookie,

Henn has used this code extensively, so I suggest you check out his dissertation (DOI: 10.5445/IR/1000061601).

For more in-depth discussions, you should consider attending the next occurrence of our annual [Spring School](https://www.openlb.net/spring-school-2024/).

Best regards,
Jan

[Rookie]

Dear Jan,

Thank you for the suggestion, I have downloaded this paper before, I will read this paper and other papers published by Henn, the problem I am currently experiencing is that the fluid and particles cannot be simulated at the same time, the original code simulates the fluid first and then the particles. After modifying it to reach the time step where I add the particles it reports an error. I would love to attend spring school to solve my problem, but please forgive me for not being able to attend.

terminate called after throwing an instance of ‘std::domain_error’
what(): read only access to data which is not available locally
[banana:883562] *** Process received signal ***
[banana:883562] Signal: Aborted (6)
[banana:883562] Signal code: (-6)

————————————————————————–
Primary job terminated normally, but 1 process returned
a non-zero exit code. Per user-direction, the job has been aborted.
————————————————————————–
————————————————————————–
mpirun noticed that process rank 1 with PID 0 on node banana exited on signal 6 (Aborted).
————————————————————————–

Best regards,
Rookie

[jan]

Dear Rookie,

if the Spring School doesn’t work for you, then maybe becoming a member of the consortium (https://www.openlb.net/consortium/) is an alternative option, because unfortunately I can’t provide such detailed support here.

Best regards,
Jan

[Author]

Posts