How to add another foil in IBLevelSet/ex1

[Han Han]

Hello all,

Is it possible to add another foil in this example?

Best,

Han

[Amneet Bhalla]

Yes, it should be possible to add more structures in IBLevelSet.

@Boyce: We are trying to validate a self-propelling airfoil with IBLevelSet. Our code gives results different than what is reported in the literature. We had given Amin an FE mesh to simulate the same case with IIM, but he is facing stability issues with that method. I’m wondering if you or somebody in your group can take a look into the IIM setup. That way we can validate if IBLevelSet and IIM are giving the same results or not. There are a lot of users wanting to do airfoil simulations in IBAMR. 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/06308afd-6fa7-45d9-9fe6-e06f691807d5n%40googlegroups.com.

--

--Amneet 

[Boyce Griffith]

On Dec 13, 2021, at 12:11 PM, Amneet Bhalla <mail2...@gmail.com> wrote:

Yes, it should be possible to add more structures in IBLevelSet.

@Boyce: We are trying to validate a self-propelling airfoil with IBLevelSet. Our code gives results different than what is reported in the literature. We had given Amin an FE mesh to simulate the same case with IIM, but he is facing stability issues with that method. I’m wondering if you or somebody in your group can take a look into the IIM setup. That way we can validate if IBLevelSet and IIM are giving the same results or not. There are a lot of users wanting to do airfoil simulations in IBAMR. 

One thing to confirm is that the boundary conditions and domain sizes are consistent --- I know that they can make important differences for very simple problems, like flow past a cylinder or sphere, and I would anticipate the same would hold for airfoils.

It should be very easy to make comparisons with IBFE, using either filled or unfilled domains. (See this paper for what we think are best practices: https://arxiv.org/abs/2105.14536.) With the IIM implementation, there could possibly be issues with the sharp corner at the trailing edge. It might also make sense to try to use only velocity penalization with IIM for a stationary boundary. We surely can investigate, although I am not sure about the precise timeframe.

On Mon, Dec 13, 2021 at 8:58 AM Han Han <stacea...@gmail.com> wrote:

Hello all,

Is it possible to add another foil in this example?

Best,

Han

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/06308afd-6fa7-45d9-9fe6-e06f691807d5n%40googlegroups.com.

--

--Amneet 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ3bhoLxuJMqf3M3xV0uRcH2ecy0545vnRRp-JqFJiEcAw%40mail.gmail.com.

[Amneet Bhalla]

On Mon, Dec 13, 2021 at 9:35 AM Boyce Griffith <boy...@gmail.com> wrote:

On Dec 13, 2021, at 12:11 PM, Amneet Bhalla <mail2...@gmail.com> wrote:

Yes, it should be possible to add more structures in IBLevelSet.

@Boyce: We are trying to validate a self-propelling airfoil with IBLevelSet. Our code gives results different than what is reported in the literature. We had given Amin an FE mesh to simulate the same case with IIM, but he is facing stability issues with that method. I’m wondering if you or somebody in your group can take a look into the IIM setup. That way we can validate if IBLevelSet and IIM are giving the same results or not. There are a lot of users wanting to do airfoil simulations in IBAMR. 

One thing to confirm is that the boundary conditions and domain sizes are consistent --- I know that they can make important differences for very simple problems, like flow past a cylinder or sphere, and I would anticipate the same would hold for airfoils.

We had done some additional cases to test the sensitivity. The variations and trends in our results and those reported in https://doi.org/10.1017/jfm.2020.955 are quite far off to be attributed to just the domain size and boundary conditions. 

It should be very easy to make comparisons with IBFE, using either filled or unfilled domains. (See this paper for what we think are best practices: https://arxiv.org/abs/2105.14536.) With the IIM implementation, there could possibly be issues with the sharp corner at the trailing edge. It might also make sense to try to use only velocity penalization with IIM for a stationary boundary. We surely can investigate, although I am not sure about the precise timeframe.

Attached are the ME mesh and input file for the domain size, boundary conditions, and other physical parameters of the problem considered in Fig. 4. We had also tried cases of Fig. 7, but those also differ from ours. 

 

On Mon, Dec 13, 2021 at 8:58 AM Han Han <stacea...@gmail.com> wrote:

Hello all,

Is it possible to add another foil in this example?

Best,

Han

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/06308afd-6fa7-45d9-9fe6-e06f691807d5n%40googlegroups.com.

--

--Amneet 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ3bhoLxuJMqf3M3xV0uRcH2ecy0545vnRRp-JqFJiEcAw%40mail.gmail.com.

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/E8B74010-7598-40D6-8C06-B1E7A14311D2%40gmail.com.

--

--Amneet 

[Amneet Bhalla]

Hi Han,

Please find the attached code for simulating two airfoils in the domain using IBLevelSet. Link this application to IBAMR master using the Makefile included here. This code will make its way to IBAMR GitHub once Boyce et al. validate the rigid body dynamics of a single foil using a different IB method. However, I'm reasonably confident that the dynamics I obtain are correct.  

The airfoil geometry and simulation parameters correspond to Fig 4 of this paper:

https://doi.org/10.1017/jfm.2020.955

The code outputs the rigid body dynamics and hydrodynamic force and torque on each airfoil in rbd.foilA, rbd.foilB, hydro_force_torque.foilA, hydro_force_torque.foilB, respectively.

Also, feel free to explore the code and simulate another test case from a different paper (to validate what I've implemented is correct).  

Thanks, 

--Amneet 

 

--

--Amneet 

[Han Han]

Thank you so much for the information! 

[Amneet Bhalla]

FYI: the airfoils can be visualized in VisIt using their zero contours. See attached images. In this case the lower airfoil is flapping with a pitch amplitude of 10 degrees, whereas the upper one is pitching with an amplitude of 20 degrees (and hence moving faster) 

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/dd3cd50b-c2d3-44cb-a785-06cd4e039faan%40googlegroups.com.

--

--Amneet 

[Amneet Bhalla]

One cool thing about this setup is that the two airfoils are far away from each other and their dynamics are pretty much independent at this Reynolds number (= 200). This could be useful in carrying out large parametric study, where one can do two or more airfoils with different kinematics within a single simulation :-)  

--

--Amneet 

[Han Han]

Hi Amneet,

Sorry for the late reply. Thank you again for your information, I am very interested in this.

I recently did some verification and found that there are numerical oscillations in the force-time curve as shown in the attached figure.

Have you encountered this phenomenon? Any suggestions to eliminate it?

Best,

Han

[Amneet Bhalla]

Hi Han,

Is the figure showing the total hydrodynamic force (viscous + pressure) or a part of it? Also, is this a different case than the one I did in the code I sent you?

The oscillations in the force is due to the direct evaluation of the hydrodynamic stress tensor on the surface of the moving body. Although there are other ways of computing the forces on the body, as described in here (https://doi.org/10.1016/j.jcp.2017.06.047), most of the oscillations can be eliminated by doing a simple moving point averaging. For example, force at point t can be replaced by average of force at {t-2 dt, t - dt, t, t+ dt,  t+2dt} as a post-processing step. The actual rigid body dynamics is however not affected by this operation, which uses the original force data in the second law of motion.

Thanks,

—Amneet

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/46cd6619-a607-4e27-8824-eeefe3743b2en%40googlegroups.com.

--

--Amneet 

[Amneet Bhalla]

One thing I forgot to mention is that you may want to run this code at higher solver tolerance than the default one of 1e-6. This is due to the implicit FSI solver. This is the command line I typically use for IBLevelSet codes:

mpiexec -np 4 ./main2d input2d -stokes_ksp_rtol 1e-8 -adv_diff_rtol 1e-10 -adv_diff_pc_type none

--

--Amneet 

[Han Han]

Hi Amneet,

Thank you for the quick reply and your suggestions! The figure shows the "hydro_force_pressure[1]" and it is from the same code as you sent to me.

I tried to run the command you mentioned but it seems to be an error. I am not sure if I did it right.

Best,

Han

[Amneet Bhalla]

On Wed, Dec 29, 2021 at 9:33 AM Han Han <stacea...@gmail.com> wrote:

Hi Amneet,

Thank you for the quick reply and your suggestions! The figure shows the "hydro_force_pressure[1]" and it is from the same code as you sent to me.

Ok, so this is the pressure force in the y-direction. The total hydrodynamic force in the y-direction will be  hydro_force_pressure[1] + hydro_force_viscous[1]. The viscous part is quite small for aerofoil cases at moderate to high Reynolds numbers, but not sure about their ratio at Re = 200. 

I tried to run the command you mentioned but it seems to be an error. I am not sure if I did it right.

I don't think there is anything wrong with what you did. Looks like an MPI problem to me. What command did you use earlier: mpirun or mpiexec? Does that command work without the additional solver flags? I would make sure the mpirun or mpiexec comes from the same MPI you used to build IBAMR. Also, the serial C++ (to be super sure serial C and fortran compilers as well) compiler and the MPI C++ (and MPI C/Fortran) compiler are the same. For example, you can check this by typing

$ g++ --version 

$ mpic++ --version

Both should return the same underlying serial C++ compiler. 

 

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/913315dd-f574-4b11-a602-0e4a3e851d4dn%40googlegroups.com.

--

--Amneet 

[Amin Kolahdouz]

Hi Amneet, 

Sorry it took me so long to finally get to this. I spent some time yesterday and properly set up the problem (as I understood it's a combination of unconstrained translational motion + prescribed rotation around a pivot point that is not COM). Taking the thinnest airfoil case with f=1 and theta_0=20 degree (buoyant), I actually seem to be getting the same results as in the JFM paper (figure 4). Please see attached for some figures and the code, and  let me know of any question.

Btw, to run the IIM code at the time being, one needs to compile it against either of these 2 branches:

git clone https://github.com/AminKolah/IBAMR.git

git checkout iim-tests-2

Or

git clone https://github.com/drwells/IBAMR.git

git checkout iim-1

PS. I'd recommend to run with a tight solver tolerance -stokes_ksp_rtol 1.0e-10 -ksp_rtol 1.0e-10

[Amin Kolahdouz]

and here is the airfoil geometry (should be placed in the same directory as the example code).

Best,

Amin

[Amneet Bhalla]

Hi Amin,

Thank you for your efforts in this matter. This will help to validate this problem further. I quickly checked the JFM results, our IBLevelSet results, and your IIM results attached in the previous email. Actually your results for Ux are in between JFM and those from IBLevelSet. Also, it appears that Ux in your case is still increasing and has not yet settled fully (?). However, for Uy and Y the two results are very close. The comparison of the results from the  three studies are attached herewith.

I will be on vacation/travel but will run your code with different airfoil parameters and mass ratios to investigate it further, when I get a chance.    

Again appreciate your efforts in undertaking this task!

--Amneet

 

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/f681e7e7-760e-45ba-a0aa-c190abd17119n%40googlegroups.com.

--

--Amneet 

[Amneet Bhalla]

(The negative sign in Y range for Ux comparison got hidden under the left image. Check this file instead). 

--

--Amneet 

[Boyce Griffith]

Are these grid converged results?

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ3jVhQcuri9UipQ7E%3D%3Da%2BQzgzw0D6FgM3T5yOJ6xc2UTw%40mail.gmail.com.
<IIM_IBLevelSet_JFM_comparison.pdf>

[Amin Kolahdouz]

@Amneet: Sure thing.  These results look nice, and both seem pretty consistent with the JFM results too (at least for U_y and y data)... I remember though you mentioned some sort of inconsistency that was worth investigating with our methods?  I admit the agreement for U_x is less nice; mine actually continues to stay in the range of 0.55 +/- 0.1 (it does not increase). The mean value in your result appears to be about 0.8, and the paper has reported a mean value of about 0.4?

@Boyce: I haven't done a rigorous convergence but I ran with a couple of mesh spacings at least and two half time-steps to make sure the results stay roughly consistent throughout. I also added a flag further tightening the tolerance for the deviations between the two ILE representations. With the attached input file I get under 0.1*h discrepancies which is great for such thin and sharp geometry, though it uses a smaller time step... It may be worth investigating using direct forcing too. One last thing about this particular example that makes it nice and unique is that it only seems to work with discontinuous L2_Lagrange family and not regular Lagrange, no matter how small the time-step.

Best,

Amin 

[Amin Kolahdouz]

Attached is the input file with the tighter tolerance of ILE representations i.e. more accurate result.

[Amneet Bhalla]

@Boyce: I think the results are converged. In IBLevelSet we took atleast 10 cells per diameter of the circular head. The number of cells in the chord is ~100 or more. This should be enough for FSI at Re = 200, based on my prior experience. I think Amin also took the same mesh resolution.

@Amin: The mean in the plot is between 0.3 - 0.35. The thick black dots make it look like around 0.4. However, we were not able to match any of other cases with the JFM paper also, e.g., at different mass ratio and frequency of pitching.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/03D60339-15D4-4E68-A374-70CD26290180%40gmail.com.

--

--Amneet 

[Amneet Bhalla]

One thing I will do is plot the level set geometry and the FE mesh and see if they are the same or not. We used constructive solid geometry operators on primitive shapes to create the airfoil analytically. 

--

--Amneet 

[Amin Kolahdouz]

Yeah, I think that'd be a good thing to check, because I believe I tried to set everything else based on your input file (note that I'm rescaling the .e by 0.1 in the code and rotate it 180 degrees) 

As for Fig 4(d), from what I see (eye-balling the curve) the u_x mean for the case with f=2 and theta_0=10 degree is around 0.3, while f=1 and theta_0=20 degree (our case) is closer to 0.4, I might be wrong though...I agree it'd be nice to check a few other parameter sets as you mentioned. 

Best,

Amin

[Amneet Bhalla]

The geometries match pretty well. The zero-contour of the level set (black) and the (scaled, translated, and rotated) FE mesh (in red) are plotted. VisIt fudges the trailing edge of the zero-contour. Matlab, however, does a better job in showing the zero-contour of the signed distance function, which we had verified separately.

Also, attached are the IBLevelSet results for the two cases considered in Fig. 7. The mass ratios are 1 and 20, while f = 2, theta = 10, and b/c = 0.1 for both of them. The Ux velocities in the paper are significantly lower to what I obtain. Although, I plan on running your code in the next year [ ;-)  ], you could also try these more contrasting cases before me.  

  

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/aec648d1-95cb-422a-8562-53cb02324daen%40googlegroups.com.

--

--Amneet 

[Yadong Zeng]

Hi Prof. Amneet Bhalla,

Just a quick question. This is a prescribed motion, which means the mass ratio does not affect the results. Why the multiple curves in the last figure are different?

Maybe I missed sth......

Best.

To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users+unsubscribe@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/aec648d1-95cb-422a-8562-53cb02324daen%40googlegroups.com.

--

--Amneet 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.

To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users+unsubscribe@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ3sGQUD9jRL7rLUgbwsVf9HuVYY05gdK7RSnYRDfN5Jpg%40mail.gmail.com.

--

Yadong (Jordan) Zeng

Ph.D. Candidate

Department of Mechanical Engineering  

University of Minnesota 

[Amneet Bhalla]

Hi Yadong,

The prescribed motion is only for the rotational/pitching degree of freedom, whereas the two translational (x and y) are free. The mass of the body affects the free DoFs (through second law of motion).

Thanks

—Amneet

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/aec648d1-95cb-422a-8562-53cb02324daen%40googlegroups.com.

--

--Amneet 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.

To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ3sGQUD9jRL7rLUgbwsVf9HuVYY05gdK7RSnYRDfN5Jpg%40mail.gmail.com.

--

Yadong (Jordan) Zeng

Ph.D. Candidate

Department of Mechanical Engineering  

University of Minnesota 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAFZuioBnH2vtqzmv0uN7znd4Ab0U5r_9YaCP428XRqpm3vszsQ%40mail.gmail.com.

--

--Amneet 

[Yadong Zeng]

Ha, I got it.  Thanks so much.

Now I am still doing my testing. 

To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users+unsubscribe@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/aec648d1-95cb-422a-8562-53cb02324daen%40googlegroups.com.

--

--Amneet 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.

To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users+unsubscribe@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ3sGQUD9jRL7rLUgbwsVf9HuVYY05gdK7RSnYRDfN5Jpg%40mail.gmail.com.

--

Yadong (Jordan) Zeng

Ph.D. Candidate

Department of Mechanical Engineering  

University of Minnesota 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.

To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users+unsubscribe@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAFZuioBnH2vtqzmv0uN7znd4Ab0U5r_9YaCP428XRqpm3vszsQ%40mail.gmail.com.

--

--Amneet 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.

To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users+unsubscribe@googlegroups.com.
To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ3YYHsBJ1zvD6%2BK5eBep%3Dpu34L016iA_t9TtMyvYx_a6A%40mail.gmail.com.

[Han Han]

Hi Amneet,

There are still numerical oscillations in the figure (Figure 1) of the total hydrodynamic force. After doing a moving point averaging, the situation has improved, but it is still not smooth (Figure 2).

I successfully ran the codes you sent to me, and the result does not seem to be a significant improvement. I'm still thinking of ways to make it smooth.

@Amin: I tried to run your codes, but I failed at this step, could you please give me some suggestions (Figure 3).

Best,

Han

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/aec648d1-95cb-422a-8562-53cb02324daen%40googlegroups.com.

--

--Amneet 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.

To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ3sGQUD9jRL7rLUgbwsVf9HuVYY05gdK7RSnYRDfN5Jpg%40mail.gmail.com.

--

Yadong (Jordan) Zeng

Ph.D. Candidate

Department of Mechanical Engineering  

University of Minnesota 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.

To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAFZuioBnH2vtqzmv0uN7znd4Ab0U5r_9YaCP428XRqpm3vszsQ%40mail.gmail.com.

--

--Amneet 

--
You received this message because you are subscribed to the Google Groups "IBAMR Users" group.

To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ3YYHsBJ1zvD6%2BK5eBep%3Dpu34L016iA_t9TtMyvYx_a6A%40mail.gmail.com.

[Amneet Bhalla]

Hi Han,

I think these oscillations are at an acceptable level. Note that it is difficult to get rid of all the oscillations in an IB method from the direct evaluation of the stress tensor on the surface of the body. 

You may want to experiment with the width of the moving point average filter. Another option would be to perform an FFT (fast Fourier transform) on the force data to remove the high frequency noise. You can use FFT on the moving point average filtered data as well. 

Thanks,

--Amneet 

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/c10e6ecb-c2dc-4404-a8ee-b10ebdcbe5e6n%40googlegroups.com.

--

--Amneet 

[Boyce Griffith]

On Jan 4, 2022, at 11:06 PM, Amneet Bhalla <mail2...@gmail.com> wrote:

Hi Han,

I think these oscillations are at an acceptable level. Note that it is difficult to get rid of all the oscillations in an IB method from the direct evaluation of the stress tensor on the surface of the body. 

I don't know what this would correspond to with IBLevelSet, but as I mentioned in the other thread, we have tried to systematically investigate computing these kinds of forces using IBFEMethod. We find that there can be large impacts of the structural mesh spacing and choice of kernel function; see preprints here: https://arxiv.org/abs/2105.14536 and here: https://arxiv.org/abs/2111.09958. In general, we find that we get lower quality integrated forces if the Lagrangian point spacing is too fine with respect to the background Cartesian grid, but we can avoid these kinds of oscillations by choosing particular kernel functions (BSPLINE_3 seems to be the most robust) and making the structural mesh relatively coarse with respect to the background Cartesian grid.

You may want to experiment with the width of the moving point average filter. Another option would be to perform an FFT (fast Fourier transform) on the force data to remove the high frequency noise. You can use FFT on the moving point average filtered data as well. 

Thanks,

--Amneet 

On Sat, Jan 1, 2022 at 10:55 PM Han Han <stacea...@gmail.com> wrote:

Hi Amneet,

There are still numerical oscillations in the figure (Figure 1) of the total hydrodynamic force. After doing a moving point averaging, the situation has improved, but it is still not smooth (Figure 2).

I successfully ran the codes you sent to me, and the result does not seem to be a significant improvement. I'm still thinking of ways to make it smooth.

@Amin: I tried to run your codes, but I failed at this step, could you please give me some suggestions (Figure 3).

Best,

Han

<han_1.jpg><han_2.jpg><han_3.jpg>

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ1MwtJvYju0rf-KjSSFxeLec3QuY9vd_L1rwjpXvQ-50Q%40mail.gmail.com.

[Amneet Bhalla]

In IBLevelSet there is no Lagrangian mesh and no IB kernel. It is free from Eulerian-Lagrangian interactions. 

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/925718AF-FADC-4065-9F62-029834C39446%40gmail.com.

--

--Amneet 

[Han Han]

Thanks for your advice!

[Boyce Griffith]

On Jan 5, 2022, at 10:08 AM, Amneet Bhalla <mail2...@gmail.com> wrote:

In IBLevelSet there is no Lagrangian mesh and no IB kernel. It is free from Eulerian-Lagrangian interactions.

There is a representation of the structure (e.g. its geometry + rigid body DOFs) that is coupled to the Navier-Stokes solver somehow. The details of how that coupling works may be important.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/CAMETWJ0QoXsJTDeO%2BSzrQLM3GCzza-vLwD_JViWw%2BhX-eAgwpw%40mail.gmail.com.

[Boyce Griffith]

On Jan 5, 2022, at 11:16 AM, Boyce Griffith <boy...@gmail.com> wrote:

On Jan 5, 2022, at 10:08 AM, Amneet Bhalla <mail2...@gmail.com> wrote:

In IBLevelSet there is no Lagrangian mesh and no IB kernel. It is free from Eulerian-Lagrangian interactions.

There is a representation of the structure (e.g. its geometry + rigid body DOFs) that is coupled to the Navier-Stokes solver somehow. The details of how that coupling works may be important.

BTW, my other point is that noisy forces are not an inherent part of IB formulations. 😁

[Amneet Bhalla]

Hi Amin, 

Your code links against the branch you mentioned earlier. Thanks for suggesting that. However, I'm not able to figure out how to change the frequency and angle of rotation. It seems that the kinematics are baked into the code and not adjustable via the input file. Can you instruct me what changes to makes for :

1) Mass (I think this can be changed via the RHOS variable in the input file?)

2) Angle of rotation 

3) Frequency of oscillation.

Thanks, 

--Amneet 

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/f681e7e7-760e-45ba-a0aa-c190abd17119n%40googlegroups.com.

--

--Amneet 

[Amneet Bhalla]

There is an error while reading the mesh. I'm using libMesh v1.6.2

Error opening ExodusII mesh file: airfoil_bc_0.1_TRI3_0.1.e

Stack frames: 9

0: libMesh::print_trace(std::ostream&)

1: libMesh::MacroFunctions::report_error(char const*, int, char const*, char const*)

2: libMesh::ExodusII_IO_Helper::open(char const*, bool)

3: libMesh::ExodusII_IO::read(std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&)

4: libMesh::NameBasedIO::read(std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&)

5: libMesh::UnstructuredMesh::read(std::__cxx11::basic_string<char, std::char_traits<char>, std::allocator<char> > const&, void*, bool, bool)

6: ./main2d_exe() [0x48ccff]

7: __libc_start_main

8: ./main2d_exe() [0x4a071f]

[0] ../libMesh/src/mesh/exodusII_io_helper.C, line 523, compiled nodate at notime

Abort(1) on node 0 (rank 0 in comm 0): application called MPI_Abort(MPI_COMM_WORLD, 1) - process 0

--

--Amneet 

[Han Han]

Hi Amneet,

I met the same error before and I moved airfoil_bc_0.1_TRI3_0.1.e to the output folder. Now I successfully ran the code.

BTW, if I want to modify the shape of the airfoil for more tests, what steps are needed?

Best,

Han

[Amneet Bhalla]

Hi Han, 

Thanks! I'm also able to resolve the error by using the .e file attached in this email thread. Amin's attachment in the other thread produced the  error I reported earlier. 

Anyways, I simulated the Fig. 7 case of the JFM paper with IIM method (it is still running) and compared it against the IBLevelSet method. Specifically, I took the mass ratio = 20, frequency = 2, and theta_max = 10.  Two observations:

1) Both methods are predicting Ux velocity that is much higher than what is reported in the JFM paper. IIM is predicting even higher Ux than IBLevelSet.

2) IIM method has a lot of small scale oscillations in both velocity and hydrodynamic forces. Perhaps these can be reduced by re-tuning the penalty parameters. I do not have the expertise to fix this and perhaps Amin is the best person to help you out in this regard. However, if the IIM method requires re-tuning the penalty parameters for different mass ratios and other airfoil parameters, then .... :-)

The comparison of the two methods is attached here. One thing is clear that the JFM paper results are not reliable.  IBLevelSet and IIM  methods are very different and they both are predicting higher Ux---it cannot be a coincidence. I am still waiting for Yadong to report the results with a third method to definitely say that the JFM results are wrong.

Based upon these tests, you can use the IBLevelSetMethod code I sent you for airfoil simulations. It is producing correct results and has an acceptable/minor level of oscillations in the hydrodynamic forces. The best thing about IBLevelSet is that you do not require to generate any mesh (as you asked in your last email). 

Thanks,

--Amneet

    

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/c918a914-6db1-413a-8124-1488d8c302cbn%40googlegroups.com.

--

--Amneet 

[Amin Kolahdouz]

Hi Amneet,

I'm a bit surprised you're getting these oscillations with the input file input2d.sdt that I shared a couple of weeks ago. There could be minor oscillations in the first sec or so but the plots at a later time should be all smooth (check back those comparison files in the pdf attachment you shared yourself a while ago). If you change the frequency to significantly larger values though, I agree re-tuning penalty parameters could be necessary, though it should not be bad at all.

As for the parameter set, I agree I hastily hard coded some parameters. Here is an updated version introducing a separate class for the parameter set that can be determined via the input file. 

Please let me know if I can be of any other help with your investigation.

Thanks,

Amin

[Amneet Bhalla]

Hi Amin,

Yes, I used the input2d.sdt file and had made a similar change to the main.cpp file for imposing the new kinematics (freq = 2 and theta_max = 10). See:

double freq = 2.0;

        double T_iter = (M_PI/18)*sin(2 * M_PI* freq* current_time);

        double T_M = -(M_PI/18)*sin(2 * M_PI*freq* (current_time - dts));

        double W_iter = (2 * M_PI * M_PI/18 * freq) * cos(2 * M_PI* freq* current_time); 

The parameter changes from your case to this one are: rhos = 1 to rhos = 20, freq = 1 to freq = 2, and theta_max = 20 to theta_max = 10. I did not change any penalty parameters. Apart from these large oscillations, the Navier Stokes solver stopped converging after t = 1.5 s or so.

9809 KSP unpreconditioned resid norm 3.525848442169e+00 true resid norm 3.399182161987e+00 ||r(i)||/||b|| 4.788796628937e-05

9810 KSP unpreconditioned resid norm 3.399171189179e+00 true resid norm 3.525087598833e+00 ||r(i)||/||b|| 4.966173275084e-05

9811 KSP unpreconditioned resid norm 3.399170646010e+00 true resid norm 3.524404608146e+00 ||r(i)||/||b|| 4.965211072018e-05

9812 KSP unpreconditioned resid norm 3.399161189072e+00 true resid norm 3.521161049138e+00 ||r(i)||/||b|| 4.960641518607e-05

9813 KSP unpreconditioned resid norm 3.399131756034e+00 true resid norm 3.516100517833e+00 ||r(i)||/||b|| 4.953512199229e-05

9814 KSP unpreconditioned resid norm 3.399125352772e+00 true resid norm 3.514171383796e+00 ||r(i)||/||b|| 4.950794418854e-05

9815 KSP unpreconditioned resid norm 3.399124708903e+00 true resid norm 3.515255109079e+00 ||r(i)||/||b|| 4.952321180214e-05 

I think my main objective to use the IIM code was to check if the IBLevelSet is producing the correct results or not. Based upon the IIM results I obtained so far, I am confident that what I did is indeed correct. I'm confident that the IIM issues for this case can be easily resolved as well. However, I won't be able to do that myself. 

Thanks,

--Amneet   

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/51ea0aa2-0e6f-41ee-b645-f939ac33134dn%40googlegroups.com.

--

--Amneet 

[Amneet Bhalla]

Hi Li-Ming, 

Thanks for sending us the JFM reference for validating the airfoil FSI:

https://doi.org/10.1017/jfm.2020.955  

Thinking more about the results reported in the JFM paper, Fig. 7 seems paradoxical. In particular, consider mass ratios 10 and 40 (whose trends are more clearly visible in the plot). For about the first 20 seconds both bodies accelerate at the same rate. This does not feel right to me. A body that is heavier (in this case 4x heavier) has more inertia and should accelerate slower than the lighter body in the beginning. Furthermore, although it is difficult to see in the plot, it appears that all of the heavier bodies accelerate more than the lightest body (with mass ratio 1) for the first 5 seconds. 

What do you think? 

Thanks, 

--Amneet 

 

[X A (Mike Lin)]

Hi,

I would like to follow up this thread as I am attempting to validate the IBlevelSet.

I would like to enquire if there is any existing paper in terms of pitching / heaving rigid foils based on IBlevelSet, ConstraintIB, or other IBAMR modules. Also, have we compared the results against any experimental data?

With best wishes,

Mike Lin

[Amneet Bhalla]

Hi Mike, 

We published a paper on IBLevelSet for pitching airfoils recently:

https://link.springer.com/article/10.1007/s00162-023-00664-z

In that we performed a grid convergence test, but did not compare it with a prior study. The only study we know of that did similar geometry and motion is this JFM paper 

https://doi.org/10.1017/jfm.2020.955  

However, we suspect the results reported there are incorrect(?). Additionally, Amin Kolahdouz simulated the same airfoil case using his sharp interface immersed interface method (IIM), and his IIM results matched quite well with those obtained using IBLevelSet. The comparison was reported in one of the email threads on the IBAMR mailing list. Additionally, I had checked with a researcher doing FSI with immersed finite element methods, and they mentioned that it will take them some "substantial effort" to set this case in their in-house code, and this avenue was not pursued. Validating it against using a commercial code may be another possibility. 

I would be very happy to know if you can somehow validate IBLevelSetMethod for airfoil motion with any other published study and/or a different code.

Thanks, 

--Amneet   

[X A (Mike Lin)]

Hi Prof. Bhalla,

Thanks for the info. I have read through the relevant threads. 

For the convenience of validation, I would like to double-check that if the configuration of IBLevelSet/ex4 is as listed below:

Geometry: identical to the one used in JFM paper https://doi.org/10.1017/jfm.2020.955  

Re = 200 = U_x * c * rho / mu

Chord length: c = 1

Translational Velocity at x direction: U_x = 1, 

Fluid density: rho =1

So, mu = 1 / Re = 1 / 200

solid-fluid mass ratio = 1

Degrees of Freedom: unconstrained to translate on both x, y directions, while the rotation is prescribed as theta(t) = theta_0 sin(2 * pi * f * t)

Amplitude of pitching degree: theta_0 = 20 degrees

pitching frequency: f = 1

Centre of Rotation: center of the semi-circle foil head (at the black dot)

chord thickness: b = 0.1

radius of head circle: r = b/2 = 0.05

periodic boundary conditions at all four boundaries.

comput. domain: X * Y = 60 * 20 with initial rotation centre position at  (X_0, Y_0) = (52, 10)

If it is my understanding is correct, I will go on to setup the case up in another numerical model.

Also, the mesh density used in ex4 is already the mesh-converged one, right? If not I will need to do the mesh convergence test on ex4.

Thanks,

Mike Lin

[X A (Mike Lin)]

Also, the JFM paper says that the flow is incompressible and viscous. So I guess the ex4 setup with fixed density and viscosity as well? 

I am asking because I see the input2d of ex4 says that "RHO_IS_CONST = FALSE" .

Best,

Mike

[Amneet Bhalla]

I will look into the setup later today or tomorrow, and will confirm it with you.

Regarding  RHO_IS_CONST = FALSE, this is just to indicate that we need a variable coefficient Navier Stokes solver. Specifically the penalization method requires a variable coefficient solver. The density is set programmatically equal to a uniform one by SetFluidProperties.cpp (or something like that) file in the example folder.

--

You received this message because you are subscribed to the Google Groups "IBAMR Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/f9b3853d-52bd-46cf-8639-70cc44ebf09en%40googlegroups.com.

[Amneet Bhalla]

Hi Mike, 

From the paper and IBLevelSet/ex4, I can confirm that what you wrote corresponds to one of the cases simulated in the JFM paper. 

Thanks, 

--

You received this message because you are subscribed to the Google Groups "IBAMR Users" group.
To unsubscribe from this group and stop receiving emails from it, send an email to ibamr-users...@googlegroups.com.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/3c8f293b-2cda-4193-a4c2-803250615541n%40googlegroups.com.

--

--Amneet 

[X A (Mike Lin)]

Hi Prof. Bhalla,

Thanks for the confirmation. We are now attempting the validation.

Best,

Mike

[X A (Mike Lin)]

Hi Prof Bhalla,

I have been attempting to test the ex4 against another case of prescribed heaving foil. For that purpose, the foil has to carry out up-and-down forced oscillation (heaving)

Is it only necessary to:

1. modify the imposed kinematic, to prescribe the velocity of the oscillation in:

https://github.com/IBAMR/IBAMR/blob/6d3b35333d2b1541838fd40c04ec52ba46b8257d/examples/IBLevelSet/ex4/example.cpp#L193

As:

U_com[1] = -heave_amp * M_PI * heave_freq * std::sin(2 * M_PI * heave_freq * data_time);

2. set degree of freedom to "(0,0,0)" in:

https://github.com/IBAMR/IBAMR/blob/6d3b35333d2b1541838fd40c04ec52ba46b8257d/examples/IBLevelSet/ex4/example.cpp#L430

As:

free_dofs << 0, 0, 0;

I have tried and the results look as expected.

Hopefully I am not missing any point? 

Best,

Mike Lin

[Amneet Bhalla]

Yes, these are the two required changes. What is meant by expected results? Do they match with some paper. I guess in this case one can only compare wake dynamics or hydrodynamic forces.

To view this discussion on the web visit https://groups.google.com/d/msgid/ibamr-users/2caf2312-3e3a-402b-aad9-29fe0e32c72dn%40googlegroups.com.

[X A (Mike Lin)]

Hi thanks,

So far I have only obtained "expected results", meaning the foil is heaving as I expected while trajectory, force, the pressure and wake vortex looking sensible.

I will compare it against the data (Yes, probably hydrodynamic forces) in some exisiting paper, after dx / dt / domain size convergence is achieved.

Best,

Mike

[X A (Mike Lin)]

I am doing mesh convergence, because I found that refining both dx and dt by 4 times in the original ex4 caes, can cause 20% difference in Ux.