Question: mesh.update() in FSI

pi06jl6 · September 30, 2014, 05:15

Quote:

Originally Posted by wuwenbo

Johan,thanks for your work for this topic about developing a useful cod to simulate the fsi in strong coupling.
I can understand most of your conclusion,but forth. Can you discrible why you say "the restart of meshPhi is not loaded correctly" detailed.
thank you!

Regarding MeshPhi: This is a resolved issue, only relevant for those that are using older version of OF, like me, since the HPCC centers lag with their updates, the issue ID : 0001106.

pi06jl6 · October 7, 2014, 13:00

Quote:

Originally Posted by pi06jl6

Case solved.

This is only a minor issue, but still relevant, is the moving mesh routine in OF limited to equidistant time step? The reason I ask is that i could see (minor) pressure oscillation whenever I switch to adaptive time stepping. From the articles I see on the subject of GLC, they assume equidistant, but from what i can see, its only a question of modify the scheme accordingly. The purpose in using adaptive is that I am currently testing pimple in order to take larger time steps for stability reasons.

pi06jl6 · October 30, 2014, 05:50

Quote:

Originally Posted by pi06jl6

Regarding MeshPhi: This is a resolved issue, only relevant for those that are using older version of OF, like me, since the HPCC centers lag with their updates, the issue ID : 0001106.

FYI: In OF2.3.0 this matter is not resolved and currently i dont know if this is a bug/feature or input mishap, but since it is not a critical error I have postponed this issue but i added to the bug report to OpenFOAM, to see if someone could shed some light upon this restart issue with DyM solvers.

pi06jl6 · November 4, 2014, 14:15

Quote:

Originally Posted by pi06jl6

Regarding GLC i have studied the code in more detail and in the sublayer and I can conclude that a recent article where it is stated that it works only for Euler is plainly wrong, just remember to activate dtPhiCorr.

Clarification, the GLC correction is done via meshPhi, however, my result goes plainly wrong when I dont activate dtPhiCorr, but due to smaller pressure oscillation due to Uf implementation, I get amplifying oscillation (ie increasing due to coupling with solid solver) which makes my solver unstable unless I add smaller damping to the solid or use coarser mesh in fluid (dissipative).

wuwenbo · November 16, 2014, 08:08

Quote:

Originally Posted by pi06jl6

Regarding MeshPhi: This is a resolved issue, only relevant for those that are using older version of OF, like me, since the HPCC centers lag with their updates, the issue ID : 0001106.

I find that "mesh.movePoints(prev)" can not make the mesh to oldTime, if you use a solid whose motion has be gived (the motion of the solid body isn't affected by fluid), and subiteration more than two times, you will find that every subiteration will return different result, and this is incorrect. So, if any body can tell me how can I get the oldTime mesh, I think I can resove the problem that the presure oscillation.

pi06jl6 · November 19, 2014, 11:55

Quote:

Originally Posted by wuwenbo

I find that "mesh.movePoints(prev)" can not make the mesh to oldTime, if you use a solid whose motion has be gived (the motion of the solid body isn't affected by fluid), and subiteration more than two times, you will find that every subiteration will return different result, and this is incorrect. So, if any body can tell me how can I get the oldTime mesh, I think I can resove the problem that the pressure oscillation.

I though I was clear on this, same conclusions, well, If you could figure that out, i am open ears.

However, I resolved the pressure oscillation: Dont try use movePoints twice in a subiteration, just move the points and assume that the movement is sufficient smaller than the cell size, adjust the time step so that this assumption is valid.

The only issue I have left is that I am using non-conservative transfer load, unfortunately this can lead to pressure noise proportional to the velocity speed..

wuwenbo · November 20, 2014, 06:38

you mean that movePoints can not be using twice as mesh.update(). From these days investigation, I have same view as you. So, I think the openfoam is not suit to simulate stong coupled fsi, I want to chang other software.

pi06jl6 · November 21, 2014, 03:55

Quote:

Originally Posted by wuwenbo

you mean that movePoints can not be using twice as mesh.update(). From these days investigation, I have same view as you. So, I think the openfoam is not suit to simulate stong coupled fsi, I want to chang other software.

Ahem, actually, there is strongly coupled FSI solvers using OF with several published articles. My software works excellent. Just take sufficient small time step, move boundary points only during subiterations, and then apply movePoints at the end of the cycle. I known, this is an engineering solution, but it works! Good luck with your project!

pi06jl6 · January 23, 2015, 15:59

Important update: From the code migration from OF2.1 to OF 2.3, I discovered a minor bug in my code (order of function calls), of minor character that resolved the smaller oscillations mentioned above and it even allows my to use mesh.update() in each subcycle (!), as a matter of fact, now the code were validated with backward with respect to still tank experiment. It is not an achievement, it is a very easy and fast case, but useful for debug purpose, since you remove dependency in the flow and limits the study to the flow generated around the body. This case revealed the bug, which i otherwise wouldnt be able not to see in my experiments with turbulent flow. Outline of the still tank experiment: put zero-gradient velocity on exits to a channel, then apply force (harmonic/constant) on the bluff body. This presented correct added mass ( for multiple cases with mass ratio 1:1 to 1:20 , E=1e6-1e11, for 1:5,1:10,1:20 aspect ratio cantilever/flap with constant front area, do not go below 1:3 ). Comparing results with above trick with moving boundary points, one can see that it introduces negative mass that sometimes for high boundary velocity can affect the solution, however of more importance, the required time step often becomes very small and thus limits its use.

Thank you all for valuable feedback, in post as well by email.

pi06jl6 · January 24, 2015, 15:15

Clarification on this thread.

1. The update() is changed significantly between OF.2.1.x and OF.2.2.x. If I had started with OF.2.2 and used still tank case, I might even not been posting this thread. So while reading the thread, note to which version of the OF the discussion concerns.

2. Moving the boundary points during subcycle and update only at the final iteration dont trigger the necessary update of all components in the fields, but it works for sufficient small time steps and for cases with small velocity gradients.

3. The double movePoint trick proposed by another research didnt work me, but as I figured out later, this i caused by that I did that in a later version of OF than used for that tricks. Hence, could be worthwhile doing this prior OF.2.2.

4. To save one iteration, place update prior fluid step, the only reason i make the update after is that i needed for the AMI having a field for given mesh in "correct" format for the first iteration of the simulation, but i guess it is relative obvious from the view of other existing FSI solvers (using node-node stitch). I have removed this awkward dependency in my current solver.

5. In OF.2.3, a major change in solvers for DyMFoam, when you cut a loop be careful whether you should make the phi absolute/relative for next step.

6. Restart bug mentioned is for backward only, works fine for Euler. I have reported but unfortunately it was required to supply a case setting, but i couldnt provide my solver, it is a mixture between different opensource components.

7. Use the benefit with Total-Lagrangian and create the map function between solid and fluid only once (ie dont re-evaluate nearest neighbor) if you not are using node to node interface between solid/fluid interface, since otherwise the AMI can introduce first order error which triggers strongly with backward. A solution to this is by introduce damping.

8. Turn off the skin friction unless you have a proper wall model, current investigation (my current verification step, not finished) show improper drag effect for high boundary velocities.

wuwenbo · May 5, 2015, 00:51

Hello, pi06jl6! I have do something else for my phD. But I am still very interested in this issue.
In your clarification, you mentioned that "phi absolute/relative" should modified in DyMFoam in OF.2.3. Can you tell me how can i modify it, or sent the part of your code to compute the fluid domain to me generously.
thank you very much!

pi06jl6 · May 5, 2015, 10:35

Quote:

Originally Posted by wuwenbo

Hello, pi06jl6! I have do something else for my phD. But I am still very interested in this issue.
In your clarification, you mentioned that "phi absolute/relative" should modified in DyMFoam in OF.2.3. Can you tell me how can i modify it, or sent the part of your code to compute the fluid domain to me generously.
thank you very much!

I am currentlly writing articles/thesis, so i dont share code until I am finished, but of course I can give you details in how to modify with a snippet code.

wuwenbo · May 5, 2015, 22:00

surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU));

volVectorField HbyA("HbyA", U);
HbyA = rAU*UEqn().H();

if (pimple.nCorrPISO() <= 1)
{
UEqn.clear();
}

surfaceScalarField phiHbyA
(
"phiHbyA",
(fvc::interpolate(HbyA) & mesh.Sf())
+ rAUf*fvc::ddtCorr(U, Uf)
);

if (p.needReference())
{
fvc::makeRelative(phiHbyA, U);
adjustPhi(phiHbyA, U, p);
fvc::makeAbsolute(phiHbyA, U);
}

while (pimple.correctNonOrthogonal())
{
fvScalarMatrix pEqn
(
fvm::laplacian(rAUf, p) == fvc::div(phiHbyA)
);

pEqn.setReference(pRefCell, pRefValue);

pEqn.solve(mesh.solver(p.select(pimple.finalInnerI ter())));

if (pimple.finalNonOrthogonalIter())
{
phi = phiHbyA - pEqn.flux();
}
}

#include "continuityErrs.H"

// Explicitly relax pressure for momentum corrector
p.relax();

U = HbyA - rAU*fvc::grad(p);
U.correctBoundaryConditions();
fvOptions.correct(U);

{
Uf = fvc::interpolate(U);
surfaceVectorField n(mesh.Sf()/mesh.magSf());
Uf += n*(phi/mesh.magSf() - (n & Uf));
}

// Make the fluxes relative to the mesh motion
fvc::makeRelative(phi, U);

this is the "pEqn" in the OF 2.3.0. if the bold line ( "fvc::makeRelative(phi, U)" ) needed here. In the "icoFsiElasticNonLinULSolidFoam" of OF-ext, it is discarded.

wuwenbo · May 5, 2015, 22:10

// Calculate absolute flux from the mapped surface velocity

phi = mesh.Sf() & Uf;

if (mesh.moving() && correctPhi)
{
#include "correctPhi.H"
}

// Make the flux relative to the mesh motion
fvc::makeRelative(phi, U);

if (mesh.moving() && checkMeshCourantNo)
{
#include "meshCourantNo.H"
}

// --- Pressure-velocity PIMPLE corrector loop
while (pimple.loop())
{
#include "UEqn.H"

// --- Pressure corrector loop
while (pimple.correct())
{
#include "pEqn.H"
}

if (pimple.turbCorr())
{
turbulence->correct();
}
}

This is the code after "mesh.update()" in the "pimpleDyMFoam" of OF 2.3.0. Is the bold code located at right place.

It should be supplemented that the code in previous tip is also in the "pimpleDyMFoam" of OF 2.3.0.

Thank you for your kind! your help will be the biggest boost for me to resolve this problem.

pi06jl6 · May 6, 2015, 05:59

Everything seems correct, as I told you, there is a significant change between distributions, my error was using code snippets from different distributions, using pimpleDyMFoam as it is in current distribution OF2.3 works fine for me now. The point in above post is with respect that you have to keep in mind that at phi is made absolute at start, and it is relative at the end of the loop, so if you using phi after fluid step for any reason, keep in mind to use the right state of phi.

Three major error caused significant problems in my calculations:

1. Mixed codes from different distributions
2. Interpolation at every subiteration
3. A bug in my own code in updating variables between time steps.

Keep in mind the largest contributions to the errors were the feature in older version: swapping of the pointers in update. This is not the case with OF2.3.

Correction:

At the end of the pimpleDyMFoam i added phi = mesh.Sf() & Uf; Then at the beginning of the pimpleDymFoam I saved the starting variables, so that the only variable that changes during subiteration is the displacement of the solid. This gives a slightly more stable coupling.

wuwenbo · May 6, 2015, 07:29

Thank you!

I think i understand your point. You don't change your mesh in subiteration except final step, and start the your every subiteration with the same fluid variation. Is it true. So, I think you don't using the "movingWallVelocity" in you moving wall boundary. I had try this method too, and the code still can not work.

I will try again with your suggestion, if i work out, i will tell you.

Thank you vary much for your help.

wuwenbo · May 6, 2015, 10:23

can you tell me how can i "saved the starting variables".

thank you!

wuwenbo · May 6, 2015, 10:36

I use this code:

Usave = U;
psave = p;

subiteration
{
p==psave;
U==Usave;
}

is this right?

pi06jl6 · May 7, 2015, 05:40

Quote:

Originally Posted by wuwenbo

is this right?

Correctly!

pi06jl6 · May 7, 2015, 05:48

Quote:

Originally Posted by wuwenbo

... "movingWallVelocity" ...

That works for me. Out of curiosity, what solver for moving the mesh are you using?

January 23, 2015, 15:59		#29
pi06jl6 Member Johan Lorentzon Join Date: Mar 2009 Location: Lunds University, Sweden Posts: 78 Rep Power: 23	Important update: From the code migration from OF2.1 to OF 2.3, I discovered a minor bug in my code (order of function calls), of minor character that resolved the smaller oscillations mentioned above and it even allows my to use mesh.update() in each subcycle (!), as a matter of fact, now the code were validated with backward with respect to still tank experiment. It is not an achievement, it is a very easy and fast case, but useful for debug purpose, since you remove dependency in the flow and limits the study to the flow generated around the body. This case revealed the bug, which i otherwise wouldnt be able not to see in my experiments with turbulent flow. Outline of the still tank experiment: put zero-gradient velocity on exits to a channel, then apply force (harmonic/constant) on the bluff body. This presented correct added mass ( for multiple cases with mass ratio 1:1 to 1:20 , E=1e6-1e11, for 1:5,1:10,1:20 aspect ratio cantilever/flap with constant front area, do not go below 1:3 ). Comparing results with above trick with moving boundary points, one can see that it introduces negative mass that sometimes for high boundary velocity can affect the solution, however of more importance, the required time step often becomes very small and thus limits its use. Thank you all for valuable feedback, in post as well by email. wyldckat likes this. Last edited by pi06jl6; February 10, 2015 at 05:18.

January 24, 2015, 15:15		#30
pi06jl6 Member Johan Lorentzon Join Date: Mar 2009 Location: Lunds University, Sweden Posts: 78 Rep Power: 23	Clarification on this thread. 1. The update() is changed significantly between OF.2.1.x and OF.2.2.x. If I had started with OF.2.2 and used still tank case, I might even not been posting this thread. So while reading the thread, note to which version of the OF the discussion concerns. 2. Moving the boundary points during subcycle and update only at the final iteration dont trigger the necessary update of all components in the fields, but it works for sufficient small time steps and for cases with small velocity gradients. 3. The double movePoint trick proposed by another research didnt work me, but as I figured out later, this i caused by that I did that in a later version of OF than used for that tricks. Hence, could be worthwhile doing this prior OF.2.2. 4. To save one iteration, place update prior fluid step, the only reason i make the update after is that i needed for the AMI having a field for given mesh in "correct" format for the first iteration of the simulation, but i guess it is relative obvious from the view of other existing FSI solvers (using node-node stitch). I have removed this awkward dependency in my current solver. 5. In OF.2.3, a major change in solvers for DyMFoam, when you cut a loop be careful whether you should make the phi absolute/relative for next step. 6. Restart bug mentioned is for backward only, works fine for Euler. I have reported but unfortunately it was required to supply a case setting, but i couldnt provide my solver, it is a mixture between different opensource components. 7. Use the benefit with Total-Lagrangian and create the map function between solid and fluid only once (ie dont re-evaluate nearest neighbor) if you not are using node to node interface between solid/fluid interface, since otherwise the AMI can introduce first order error which triggers strongly with backward. A solution to this is by introduce damping. 8. Turn off the skin friction unless you have a proper wall model, current investigation (my current verification step, not finished) show improper drag effect for high boundary velocities. wyldckat likes this. Last edited by pi06jl6; January 25, 2015 at 11:05.

May 6, 2015, 05:59		#35
pi06jl6 Member Johan Lorentzon Join Date: Mar 2009 Location: Lunds University, Sweden Posts: 78 Rep Power: 23	Everything seems correct, as I told you, there is a significant change between distributions, my error was using code snippets from different distributions, using pimpleDyMFoam as it is in current distribution OF2.3 works fine for me now. The point in above post is with respect that you have to keep in mind that at phi is made absolute at start, and it is relative at the end of the loop, so if you using phi after fluid step for any reason, keep in mind to use the right state of phi. Three major error caused significant problems in my calculations: 1. Mixed codes from different distributions 2. Interpolation at every subiteration 3. A bug in my own code in updating variables between time steps. Keep in mind the largest contributions to the errors were the feature in older version: swapping of the pointers in update. This is not the case with OF2.3. Correction: At the end of the pimpleDyMFoam i added phi = mesh.Sf() & Uf; Then at the beginning of the pimpleDymFoam I saved the starting variables, so that the only variable that changes during subiteration is the displacement of the solid. This gives a slightly more stable coupling. Last edited by pi06jl6; May 7, 2015 at 05:50.

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
initial file for fsi modeling?	mortazavi	CFX	0	July 21, 2009 02:06
1-way FSI Vs. 2-way FSI	k_buz	Main CFD Forum	0	May 12, 2009 04:10
fsi	sun	Siemens	8	January 19, 2009 23:17
FSI on a Rotating Structural Model	drgolf	CFX	9	August 30, 2007 12:08
Restart of FSI simulation	V. Kumar	CFX	3	July 20, 2006 13:23

November 20, 2014, 06:38		#27
wuwenbo New Member wuwenbo Join Date: Jan 2013 Posts: 17 Rep Power: 13	you mean that movePoints can not be using twice as mesh.update(). From these days investigation, I have same view as you. So, I think the openfoam is not suit to simulate stong coupled fsi, I want to chang other software.

May 5, 2015, 00:51		#31
wuwenbo New Member wuwenbo Join Date: Jan 2013 Posts: 17 Rep Power: 13	Hello, pi06jl6! I have do something else for my phD. But I am still very interested in this issue. In your clarification, you mentioned that "phi absolute/relative" should modified in DyMFoam in OF.2.3. Can you tell me how can i modify it, or sent the part of your code to compute the fluid domain to me generously. thank you very much!

May 5, 2015, 22:00		#33
wuwenbo New Member wuwenbo Join Date: Jan 2013 Posts: 17 Rep Power: 13	surfaceScalarField rAUf("rAUf", fvc::interpolate(rAU)); volVectorField HbyA("HbyA", U); HbyA = rAUUEqn().H(); if (pimple.nCorrPISO() <= 1) { UEqn.clear(); } surfaceScalarField phiHbyA ( "phiHbyA", (fvc::interpolate(HbyA) & mesh.Sf()) + rAUffvc::ddtCorr(U, Uf) ); if (p.needReference()) { fvc::makeRelative(phiHbyA, U); adjustPhi(phiHbyA, U, p); fvc::makeAbsolute(phiHbyA, U); } while (pimple.correctNonOrthogonal()) { fvScalarMatrix pEqn ( fvm::laplacian(rAUf, p) == fvc::div(phiHbyA) ); pEqn.setReference(pRefCell, pRefValue); pEqn.solve(mesh.solver(p.select(pimple.finalInnerI ter()))); if (pimple.finalNonOrthogonalIter()) { phi = phiHbyA - pEqn.flux(); } } #include "continuityErrs.H" // Explicitly relax pressure for momentum corrector p.relax(); U = HbyA - rAUfvc::grad(p); U.correctBoundaryConditions(); fvOptions.correct(U); { Uf = fvc::interpolate(U); surfaceVectorField n(mesh.Sf()/mesh.magSf()); Uf += n(phi/mesh.magSf() - (n & Uf)); } // Make the fluxes relative to the mesh motion fvc::makeRelative(phi, U); this is the "pEqn" in the OF 2.3.0. if the bold line ( "fvc::makeRelative(phi, U)" ) needed here. In the "icoFsiElasticNonLinULSolidFoam" of OF-ext, it is discarded.

May 5, 2015, 22:10		#34
wuwenbo New Member wuwenbo Join Date: Jan 2013 Posts: 17 Rep Power: 13	// Calculate absolute flux from the mapped surface velocity phi = mesh.Sf() & Uf; if (mesh.moving() && correctPhi) { #include "correctPhi.H" } // Make the flux relative to the mesh motion fvc::makeRelative(phi, U); if (mesh.moving() && checkMeshCourantNo) { #include "meshCourantNo.H" } // --- Pressure-velocity PIMPLE corrector loop while (pimple.loop()) { #include "UEqn.H" // --- Pressure corrector loop while (pimple.correct()) { #include "pEqn.H" } if (pimple.turbCorr()) { turbulence->correct(); } } This is the code after "mesh.update()" in the "pimpleDyMFoam" of OF 2.3.0. Is the bold code located at right place. It should be supplemented that the code in previous tip is also in the "pimpleDyMFoam" of OF 2.3.0. Thank you for your kind! your help will be the biggest boost for me to resolve this problem.

May 6, 2015, 07:29		#36
wuwenbo New Member wuwenbo Join Date: Jan 2013 Posts: 17 Rep Power: 13	Thank you! I think i understand your point. You don't change your mesh in subiteration except final step, and start the your every subiteration with the same fluid variation. Is it true. So, I think you don't using the "movingWallVelocity" in you moving wall boundary. I had try this method too, and the code still can not work. I will try again with your suggestion, if i work out, i will tell you. Thank you vary much for your help.

May 6, 2015, 10:23		#37
wuwenbo New Member wuwenbo Join Date: Jan 2013 Posts: 17 Rep Power: 13	can you tell me how can i "saved the starting variables". thank you!

May 6, 2015, 10:36		#38
wuwenbo New Member wuwenbo Join Date: Jan 2013 Posts: 17 Rep Power: 13	I use this code: Usave = U; psave = p; subiteration { p==psave; U==Usave; } is this right?