Oceananigans.jl: Adding Immersed Boundary Capabilities
@christophernhill @glwagner @weymouth
I have been looking into what is the best way to implement an immersed boundary method here, and I hope to be able to model complex topography and possibly moving boundaries with it, for at least Dirichlet and Neumann boundary conditions, within the next year.
The immersed boundary method (IBM) is a way to deal with complex topography without a complicated mesh or coordinate transformations. Instead of conforming the mesh to the fluid domain, a cartesian grid is generated over the whole area. Then, an added forcing term incorporates the boundary conditions into the equations. Mittal and Iaccarino (2003), outline several ways we could do this, which all fall into two large categories:
Continuous Forcing: We add a forcing term to the continuous equations, discretize the new system, and solve as normal. If the IB is set to coincide with the mesh nodes, then you get a “stairstep” like boundary, whereas if the IB intersects the mesh arbitrarily, you must distribute the forcing to several nearby nodes. Either way, it may not give you a sharp enough boundary, and you get stability constraints on rigid boundaries.
Discrete Forcing: We discretize, then determine what the forcing needs to be at each time step to satisfy the boundary conditions. This method usually takes a predictor-corrector type approach as described by Balaras (2004). You would use a predictor step to interpolate the correct values at the boundary with nearby nodes. Then, you can solve for the necessary forcing in the discretized time stepper. Finally, you recompute the true solution for the next time step, with the correct forcing term. This method does not have the stability constraints and can create a sharper boundary, but implementation will involve more changes to create this routine.
I'm leaning towards the discrete forcing, since it seems like continuous just won't work well for rigid boundaries. I haven't fully worked out how the predictor-corrector idea will work with the pressure solver, etc.
WRF uses Discrete forcing in their LES IBM code, while PALM uses a masking method that involves stair step representation for complex topography.
whitleyv
All 29 comments
One extra consideration is that we also allow users to specify diffusive fluxes across boundaries. This is especially important for geophysical problems at LES scales and larger, where its appropriate to employ a "wall model" to predict momentum and tracer fluxes at boundaries where there's an unresolved or partially-resolved turbulent boundary layer (rather than prescribing a particular value or gradient). That said, I think if we are able to specify the gradients of a field across a boundary it will likely be straightforward to extend that implementation to specifying fluxes.
The paper
"Moving from momentum transfer to heat transfer – A comparative study of an advanced Graetz-Nusselt problem using immersed boundary methods" by Lu et al. (2019)
may also be relevant. Their conclusion is a bit confusing. They state
In all simulations, excellent agreement are reached between CFM and DFM results, with the deviation being below 10%.
which suggests that _accuracy_ may not be an important factor in deciding which method to use. But the next sentence is
Considering the nature of capturing the discontinuity at the fluid-solid interface, DFM might offer a more accurate result, which however requires more follow-up simulations to give a solid investigation.
which is difficult to interpret. I suppose all they can say is that their results are similar to one another, but they cannot say which one is more accurate (and perhaps it doesn't matter which method is more accurate in their case, if both methods return similar results).
That said, I think time-step considerations are really important, and seem like a good reason to choose DFM over CFM.
Balaras (2004)'s pressure equation is
where, crucially, Ω is the _computational_ domain (irrespective of the immersed boundary). So following Balaras (2004) would mean _not_ modifying the pressure solver? That's good news!
Does it make sense to first implement an algorithm that assumes the boundary coincides with the grid?
I am not super worried about the algorithmic changes require to implement a DFM immersed boundary. I think we can isolate the immersed boundary implementation from the rest of the code and interface with the time-stepping routines with a single function call that calculates the discrete forcing and applies a correction to the velocity field. Using multiple dispatch, this change to the algorithm will have no cost to simulations that don't use immersed boundaries. The main challenge I think is figuring out abstractions that make it easy to create immersed boundaries and assign boundary conditions (@ali-ramadhan and my job to figure this out) --- and we'd have this challenge for any immersed boundary implementation, whether CFM or DFM.
glwagner
on 13 Oct 2020
Greg, Doesn't Ali have a version of immersed boundary layers going?
On Mon, Oct 12, 2020, 7:50 PM Gregory L. Wagner notifications@github.com
wrote:
One extra consideration is that we also allow users to specify diffusive
fluxes across boundaries. This is especially important for geophysical
problems at LES scales and larger, where its appropriate to employ a "wall
model" to predict momentum and tracer fluxes at boundaries where there's an
unresolved or partially-resolved turbulent boundary layer (rather than
prescribing a particular value or gradient). That said, I think if we are
able to specify the gradients of a field across a boundary it will likely
be straightforward to extend that implementation to specifying fluxes.The paper
"Moving from momentum transfer to heat transfer – A comparative study of
an advanced Graetz-Nusselt problem using immersed boundary methods"
https://www.sciencedirect.com/science/article/pii/S0009250918306250 by
Lu et al. (2019)may also be relevant. Their conclusion is a bit confusing. They state
In all simulations, excellent agreement are reached between CFM and DFM
results, with the deviation being below 10%.which suggests that accuracy may not be an important factor in deciding
which method to use. But the next sentence isConsidering the nature of capturing the discontinuity at the fluid-solid
interface, DFM might offer a more accurate result, which however requires
more follow-up simulations to give a solid investigation.which is difficult to interpret. I suppose all they can say is that their
results are similar to one another, but they cannot say which one is more
accurate (and perhaps it doesn't matter which method is more accurate in
their case, if both methods return similar results).That said, I think time-step considerations are really important, and seem
like a good reason to choose DFM over CFM.Balaras (2004)'s pressure equation is
[image: image]
https://user-images.githubusercontent.com/15271942/95799302-36f2c280-0cc2-11eb-9342-a2b47a1cfdfa.pngwhere, crucially, Ω is the computational domain (irrespective of the
immersed boundary). So following Balaras (2004) would mean not
modifying the pressure solver? That's good news!Does it make sense to first implement an algorithm that assumes the
boundary coincides with the grid?I am not super worried about the algorithmic changes require to implement
a DFM immersed boundary. I think we can isolate the immersed boundary
implementation from the rest of the code and interface with the
time-stepping routines with a single function call that calculates the
discrete forcing and applies a correction to the velocity field. Using
multiple dispatch, this change to the algorithm will have no cost to
simulations that don't use immersed boundaries. The main challenge I think
is figuring out abstractions that make it easy to create immersed
boundaries and assign boundary conditions (@ali-ramadhan
https://github.com/ali-ramadhan and my job to figure this out) --- and
we'd have this challenge for any immersed boundary implementation, whether
CFM or DFM.—
You are receiving this because you are subscribed to this thread.
Reply to this email directly, view it on GitHub
https://github.com/CliMA/Oceananigans.jl/issues/1036#issuecomment-707401128,
or unsubscribe
https://github.com/notifications/unsubscribe-auth/AKXUEQXTZSGE6RBZ5E4B733SKOI4PANCNFSM4SNJ4CSA
.
johncmarshall54
on 13 Oct 2020
Greg, Doesn't Ali have a version of immersed boundary layers going?
@johncmarshall54, Ali experimented with a simple immersed boundary implemented via Oceananigans's user-defined forcing functions. The code is these 9 lines:
This implementation damps the velocity field to zero on a very fast time-scale (specified by the parameter K) within the immersed boundary. This is certainly a nice, simple immersed boundary implementation for Dirichlet / Value boundary conditions and could be a good starting point. However, I also think there's some good reasons to pursue an alternate immersed boundary implementation to what @ali-ramadhan has done:
- @ali-ramadhan's implementation doesn't obviously extend to other boundary conditions, like prescribed gradients or fluxes (the latter being crucial for the geophysical problems we're interested in). So, even if we use a continuous forcing method similar to @ali-ramadhan's implementation, we need to figure out how to enforce boundary conditions other than Dirichlet boundary conditions.
- @ali-ramadhan's implementation conforms exactly to the grid; however we would like to be able to model smoothly-varying boundaries.
- As noted by @whitleyv, @ali-ramadhan's "continuous forcing method" implementation introduces a time-step restriction due to the need to explicitly resolve the damping time-scale in the forcing function. It seems that a discrete forcing method overcomes this restriction and could prove crucial for geophysical problems that involve otherwise long time-steps.
glwagner
on 13 Oct 2020
thanks for the explanation.
On Mon, Oct 12, 2020 at 8:25 PM Gregory L. Wagner notifications@github.com
wrote:
Greg, Doesn't Ali have a version of immersed boundary layers going?
@johncmarshall54 https://github.com/johncmarshall54, Ali experimented
with a simple immersed boundary implemented via Oceananigans's user-defined
forcing functions. The code is these 9 lines:This implementation damps the velocity field to zero on a very fast
time-scale (specified by the parameter K) within the immersed boundary.
This is certainly a nice, simple immersed boundary implementation for
Dirichlet / Value boundary conditions and could a good starting point.
However, I also think there's some good reasons to pursue an alternate
immersed boundary implementation to what @ali-ramadhan
https://github.com/ali-ramadhan has done for the following reasons:
- @ali-ramadhan https://github.com/ali-ramadhan's implementation
doesn't obviously extend to other boundary conditions, like prescribed
gradients or fluxes (the latter being crucial for the geophysical problems
we're interested in). So, even if we use a continuous forcing method
similar to @ali-ramadhan https://github.com/ali-ramadhan's
implementation, we need to figure out how to enforce boundary conditions
other than Dirichlet boundary conditions.- @ali-ramadhan https://github.com/ali-ramadhan's implementation
conforms exactly to the grid; however we would like to be able to model
smoothly-varying boundaries.- As noted by @whitleyv https://github.com/whitleyv, @ali-ramadhan
https://github.com/ali-ramadhan's "continuous forcing method"
implementation introduces a time-step restriction due to the need to
explicitly resolve the damping time-scale in the forcing function. It seems
that a discrete forcing method overcomes this restriction and could prove
crucial for geophysical problems that involve otherwise long time-steps.—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
https://github.com/CliMA/Oceananigans.jl/issues/1036#issuecomment-707412244,
or unsubscribe
https://github.com/notifications/unsubscribe-auth/AKXUEQXYHJJRK3FTQF5JZTTSKOM6ZANCNFSM4SNJ4CSA
.
--
John Marshall
Earth, Atmospheric and Planetary Sciences, MIT
http://oceans.mit.edu/JohnMarshall/
johncmarshall54
on 13 Oct 2020
Sounds great! I've had a discussion with Greg about this before, and which method you need to implement really depends on the application.
- The "penalty function" boundary was the first type proposed historically, but it is completely inappropriate for unsteady flows and will really mess up a turbulent boundary layer evolution.
- Direct forcing methods and other similar approaches go in and set the value of velocity components inside the immersed solid at each time step to satisfy the velocity boundary conditions before or after the projection step. However, you have to choose between enforcing a divergence free velocity and the BC - you can't have both without modifying the pressure Poisson equation. As long as the body isn't moving, the results tend to be fine. The pressure forces will not be time accurate and so dynamic forces will suffer.
- For fluid structure interaction problems you really need to get those forces right and therefore must modify the pressure equation.
I've attach a write-up one of my students is working on which lays this out with a few simple examples. Section 3 reviews Immersed Boundary methods and 3.1 has a simple 1D FSI example.
Mr M Lauber_070e4204-4db5-451e-8b61-494f2ae9eaa2_Progression_Report__9_month_pdf_7177_0.pdf
weymouth
on 13 Oct 2020
Fun pair programming sesh with @ali-ramadhan @whitleyv lead to this!
Implementation is here:
and the script that produced the above animation:
https://github.com/CliMA/Oceananigans.jl/blob/immersed-boundary/examples/flow_around_cylinder.jl
We just did something very simple as a starting point --- hopefully more to come.
@weymouth thanks for your insights --- could make sense to schedule a meeting sometime soon to discuss next steps.
glwagner
on 30 Oct 2020
Awesome! I'm happy to sit down and chat. I've had two meetings today with
MIT groups (ocean engineering and self-assembly lab), so maybe next week...
;-)
Gabriel D Weymouth
"Computers are useless. They can only give you answers." Pablo Picasso
On Fri, Oct 30, 2020 at 3:52 PM Gregory L. Wagner notifications@github.com
wrote:
Fun pair programming sesh with @ali-ramadhan
https://github.com/ali-ramadhan @whitleyv https://github.com/whitleyv
lead to this![image: flow_around_cylinder]
https://user-images.githubusercontent.com/15271942/97726829-083c6080-1aa6-11eb-8c62-c38771eac0bf.gifImplementation is here:
and the script that produced the above animation:
https://github.com/CliMA/Oceananigans.jl/blob/immersed-boundary/examples/flow_around_cylinder.jl
We just did something very simple as a starting point --- hopefully more
to come.@weymouth https://github.com/weymouth thanks for your insights ---
could make sense to schedule a meeting sometime soon to discuss next steps.—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
https://github.com/CliMA/Oceananigans.jl/issues/1036#issuecomment-719635414,
or unsubscribe
https://github.com/notifications/unsubscribe-auth/AADSKJYHOUF5GATVQ66VU7TSNLOMHANCNFSM4SNJ4CSA
.
weymouth
on 30 Oct 2020
Algorithms.pdf gives a general outline of the algorithm for the RK3 stepping currently implemented (Algorithm 1) as well as one for the IBM implemented above (Algorithm 3), if anyone is interested!
whitleyv
on 30 Oct 2020
@glwagner @ali-ramadhan and @whitleyv that is really nice.
christophernhill
on 30 Oct 2020
I've been waiting for something like this for long now :)
Nice work @whitleyv et al!
navidcy
on 31 Oct 2020
Great stuff guys. Can we now put a ridge down our eddying channel? Is an
island possible? John
On Fri, Oct 30, 2020, 11:52 AM Gregory L. Wagner notifications@github.com
wrote:
Fun pair programming sesh with @ali-ramadhan
https://github.com/ali-ramadhan @whitleyv https://github.com/whitleyv
lead to this![image: flow_around_cylinder]
https://user-images.githubusercontent.com/15271942/97726829-083c6080-1aa6-11eb-8c62-c38771eac0bf.gifImplementation is here:
and the script that produced the above animation:
https://github.com/CliMA/Oceananigans.jl/blob/immersed-boundary/examples/flow_around_cylinder.jl
We just did something very simple as a starting point --- hopefully more
to come.@weymouth https://github.com/weymouth thanks for your insights ---
could make sense to schedule a meeting sometime soon to discuss next steps.—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
https://github.com/CliMA/Oceananigans.jl/issues/1036#issuecomment-719635414,
or unsubscribe
https://github.com/notifications/unsubscribe-auth/AKXUEQRK6RHZWBUDAPODUQLSNLOMHANCNFSM4SNJ4CSA
.
johncmarshall54
on 1 Nov 2020
Great stuff guys. Can we now put a ridge down our eddying channel? Is an island possible? John
Not quite, but we're definitely making progress toward ridges and islands in an eddying channel!
@whitleyv might be able to provide more detail but the algorithms we've been discussing seem valid for arbitrary geometries --- so it islands would be on the table, as well as three-dimensional geometries like underwater "arches".
I think @whitleyv's algorithm may easily generalize to Dirichlet boundary conditions on tracer fields, which might suffice for proof-of-concept simulations in the eddying channel context. But I think we'll want to be able to set the gradients of fields (so we can enforce insulating boundary conditions on temperature, for example) before implementing eddying channel setups with bathymetry.
As @weymouth implied we may also need to invest some effort in modifying the pressure solver for large-scale oceanography-type setups where getting form drag / pressure drag across bathymetry is important (having a hydrostatic solver might make that task easier). We'll have to look into this more when the time is ripe.
glwagner
on 1 Nov 2020
Just adding to @glwagner's comment, @whitleyv and I have discussed adding:
- Boundary interpolation, so the topography does not have to match the grid.
- Neumann boundary conditions.
- A simple log-layer wall-model.
I think we have a good sense of how to approach the above.
I'm also interested in discussing pressure solver modifications with @weymouth, we will talk on our end and then be in touch to schedule a time to chat.
Thanks!
Jacob
wenegrat
on 1 Nov 2020
Some notes regarding the velocity correction / pressure solver method:
As @weymouth mentioned, we can't get both the boundary condition and the pressure distribution correct without modifying the pressure solver
One possibility is just to "not worry too much", since:
- LES often makes approximations at the boundary anyways, and the pressure gradient error may be small compared to the error associated with the wall model or LES closure
- the error may scale with the time-step and magnitude of wall-normal pressure gradients, which may be small
One method for getting correct boundary conditions + pressure gradients is to iterate the prediction-correction algorithm that both corrects for the immersed boundary, and then projects out the divergence:
- this could work for non-moving boundaries
- if boundary motion drives the flow, this might fail
- possibly preconditioning with a coarse but accurate boundary-aware pressure solver can help (eg using multigrid solvers developed by @weymouth)
glwagner
on 13 Nov 2020
I have been doing some analysis on the immersed boundary currently implemented for the case of steady state flow over a cylinder. I'm getting the following pressure contours and just wanted to check and make sure that the pressure field I'm getting within the cylinder (ie. where there shouldn't actually be anything going on) is not an issue, at least for the case of a stationary boundary. The contours on the outside of the cylinder are correct for this scenario, but the literature usually just masks the boundary with a circle hiding what is going on within the cylinder. So I am not sure what other people are getting in their IBM models.
I am also getting a small velocity within the cylinder. I can plot the velocity within the cylinder at a couple of points compared to a rough estimate of the RK3 last stage pressure correction, and they match pretty well. I'm fairly certain that the velocity errors within the IB scale with the NH pressure correction, so as long as the time stepping is kept reasonably small, then the velocity errors within the cylinder should be negligible. Is this the correct interpretation of the velocity error due to pressure correction? @weymouth @glwagner
whitleyv
on 3 Dec 2020
That solution looks pretty good!
As far as I can tell I think it makes sense for the "full" pressure field to continue smoothly into the IB, since the pressure field is somehow a solution to Poisson's equation --- even if the RHS of the Poisson equation varies rapidly close to the IB (?)
Does the magnitude of the velocity error scale with the time-step, or resolution? Perhaps plotting the dependence of the error on some of those parameters can give us confidence that the method is working as expected.
Should we try iterating the IB correction + pressure solve to see if it reduces the velocity error, as we hypothesized it might?
glwagner
on 3 Dec 2020
Victoria, is the velocity normal to the cylinder zero? and perhaps the
tangential component too, if you are using no-slip boundary conditions. I
have a nice application of all this if you are going in 3-d. John
On Thu, Dec 3, 2020 at 2:54 PM Gregory L. Wagner notifications@github.com
wrote:
That solution looks pretty good!
As far as I can tell I think it makes sense for the "full" pressure field
to continue smoothly into the IB, since the pressure field is somehow a
solution to Poisson's equation --- even if the RHS of the Poisson equation
varies rapidly close to the IB (?)Does the magnitude of the velocity error scale with the time-step, or
resolution? Perhaps plotting the dependence of the error on some of those
parameters can give us confidence that the method is working as expected.Should we try iterating the IB correction + pressure solve to see if it
reduces the velocity error, as we hypothesized it might?—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
https://github.com/CliMA/Oceananigans.jl/issues/1036#issuecomment-738270609,
or unsubscribe
https://github.com/notifications/unsubscribe-auth/AKXUEQRQB7MCPDWK3EEILQDSS7UJBANCNFSM4SNJ4CSA
.
--
John Marshall
Earth, Atmospheric and Planetary Sciences, MIT
http://oceans.mit.edu/JohnMarshall/
johncmarshall54
on 3 Dec 2020
I don't have any experience using immersed boundary methods but I would say that if you are getting the correct behavour outside of the domain, where the fluid is physically present, then it shouldn't matter very much (if at all) what is happening in the interior.
Great to see and thanks for sharing!
francispoulin
on 3 Dec 2020
@whitleyv that plot looks very nice and pressure makes sense. Plots of pressure gradient in x and y might be good too.
It would be interesting to estimate the normal flow to the cylinder walls somehow - not entirely sure how though.....
christophernhill
on 3 Dec 2020
Does the magnitude of the velocity error scale with the time-step, or resolution? Perhaps plotting the dependence of the error on some of those parameters can give us confidence that the method is working as expected.
So far I've run two different grid sizes, but I've been changing the time step with the grid size for stability, so I'm not sure. I mean the grid spacing will make it more accurate, so it probably does a better job, but I think the pressure correction is so close to the velocities we're seeing that the time step within the pressure correction has to play a big role.
Victoria, is the velocity normal to the cylinder zero? and perhaps the tangential component too, if you are using no-slip boundary conditions. I have a nice application of all this if you are going in 3-d. John
…
The velocity roughly normal to the cylinder at the top (0-180 degrees) is close to zero. I'm not sure if my normal velocity calculation is off, but the bottom half of the cylinder is showing a lot more variation.
whitleyv
on 3 Dec 2020
@whitleyv I am not sure how small the instantaneous normal velocities should be, they look kind of big to me, but the solution looks nice - I don't have a good sense what they should be in an IBM setting. We could do a zoom with Timour Radko about how they evaluated the IBM in https://doi.org/10.1063/1.5100969 ? It has a similar FV algorithm, but a slower solver.
christophernhill
on 4 Dec 2020
From your contour plot I see that the maximum speed is about 0.5. The plot of the normal velocities has the normal velocity at the surface to go as large as 0.008 or so. It is close to 100 times smaller. Do you think this gets better with increased resolution?
It's intresting that the largest error seems to be between -90 degrees and 10 degrees. I am not sure what causes the asymmetry but it is interesting.
francispoulin
on 4 Dec 2020
Does the magnitude of the velocity error scale with the time-step, or resolution?
I think the expected behavior is that the velocity error in the object scales with the non-hydrostatic pressure gradient times the timestep of the Runge-Kutta substep. Presumably the non-hydrostatic pressure gradient in the object scales as the non-hydrostatic pressure over the length scale of the object. If the strength of the non-hydrostatic pressure field on the object of the boundary is a function of the outer-flow and object configuration, then I don't think grid refinement will help much. We should do a more rigorous validation of this at some point, but at least the dependence on the timestep gives a way to control the velocity error.
I also suspect the pattern of the normal velocities in the plot shown may be spurious, as we know the velocity gets very small in all 'solid' nodes fully inside the object. The apparent error on the boundary as a function of angle may have to do with the how the current first-pass implementation sets the solid boundary location on the c-grid.
@christophernhill thanks very much for the offer! It will be good to take you up on that down the line, but let @whitleyv and I do a bit of digging first.
wenegrat
on 4 Dec 2020
@wenegrat and @whitleyv sounds good - is the angle dependence maybe due to time dependence of flow, does it move around as the flow changes?
christophernhill
on 4 Dec 2020
Looks pretty good! A 1% violation in the BCs due to the projection seems
reasonable (as in - I would think something had gone wrong if it was 20%).
This is expected due to the issues we discussed on Zoom. While it will
certainly impact force predictions and the near wall solution slightly, the
impact should be fairly minor and might be controlled by a second BC
enforcement and projection step.
The characteristics of the pressure solution within the body are
nonphysical, and so down entirely to the implementation of the pressure
solver. My guess is that the asymmetry is due to an index sweep in the
solver smearing this error ahead of itself. Is the solution outside the
body perfectly symmetric?
Gabriel D Weymouth
"Computers are useless. They can only give you answers." Pablo Picasso
On Fri, Dec 4, 2020 at 12:08 AM Chris Hill notifications@github.com wrote:
@wenegrat https://github.com/wenegrat and @whitleyv
https://github.com/whitleyv sounds good - is the angle dependence maybe
due to time dependence of flow, does it move around as the flow changes?—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
https://github.com/CliMA/Oceananigans.jl/issues/1036#issuecomment-738463407,
or unsubscribe
https://github.com/notifications/unsubscribe-auth/AADSKJZWCFPECYN4I4ZBG2DSTAR6HANCNFSM4SNJ4CSA
.
weymouth
on 4 Dec 2020
I think the expected behavior is that the velocity error in the object scales with the non-hydrostatic pressure gradient times the timestep of the Runge-Kutta substep. Presumably the non-hydrostatic pressure gradient in the object scales as the non-hydrostatic pressure over the length scale of the object. If the strength of the non-hydrostatic pressure field on the object of the boundary is a function of the outer-flow and object configuration, then I don't think grid refinement will help much.
Thanks @wenegrat. Sounds like there might be two decent tests to confirm that the algorithm is implemented correctly: 1) run the same problem for 3-5 different time-steps, and confirm that the error reduces with 1/dt and 2) refine the grid and confirm that the error remains constant.
My guess is that the asymmetry is due to an index sweep in the
solver smearing this error ahead of itself.
Hmm... the solver uses FFTs. Perhaps @ali-ramadhan can answer, but this might be a question about the internal implementation of the FFT -- not sure.
glwagner
on 4 Dec 2020
Is the solution outside the body perfectly symmetric?
@wenegrat and @whitleyv sounds good - is the angle dependence maybe due to time dependence of flow, does it move around as the flow changes?
@christophernhill This should be a _roughly_ steady state solution for Re = 40 by the end of the simulation. Changes in velocity at that point are O(10^-5). Here is the velocity contours for the above case @weymouth. Angle-wise, 0 degrees was taken due east but the flow is not in that direction, so it should have been symmetric from [-90,90] and [90,270]. As @wenegrat mentioned, the normal velocity may be uneven due to the interpolation and derivative calculations so near the boundary.
whitleyv
on 4 Dec 2020
That solution looks pretty symmetric, but you can just subtract left from
right to get a quantitative measure. If there is a persistent asymmetry
near the body, then it would be worth examining those derivative and
interpolation routines.
Gabriel D Weymouth
"Computers are useless. They can only give you answers." Pablo Picasso
On Fri, Dec 4, 2020 at 3:34 PM Victoria Whitley notifications@github.com
wrote:
Is the solution outside the body perfectly symmetric?
@wenegrat https://github.com/wenegrat and @whitleyv
https://github.com/whitleyv sounds good - is the angle dependence maybe
due to time dependence of flow, does it move around as the flow changes?@christophernhill https://github.com/christophernhill This should be a
roughly steady state solution for Re = 40 by the end of the simulation.
Changes in velocity at that point are O(10^-5). Here is the velocity
contours for the above case @weymouth https://github.com/weymouth.
Angle-wise, 0 degrees was taken due east but the flow is not in that
direction, so it should have been symmetric from [-90,90] and [90,270]. As
@wenegrat https://github.com/wenegrat mentioned, the normal velocity
may be uneven due to the interpolation and derivative calculations so near
the boundary.[image: velocity_Re40_dx04_small]
https://user-images.githubusercontent.com/67593861/101181102-7ed6fb80-361a-11eb-92f1-041f0daa398e.gif—
You are receiving this because you were mentioned.
Reply to this email directly, view it on GitHub
https://github.com/CliMA/Oceananigans.jl/issues/1036#issuecomment-738849541,
or unsubscribe
https://github.com/notifications/unsubscribe-auth/AADSKJ4FSQI3H6QIOY6YVKLSTD6PJANCNFSM4SNJ4CSA
.
weymouth
on 4 Dec 2020
OK - different from the Oct 30 movie! It does look pretty steady, but it could still be asymptoting to fully steady. For equilibration problems sometimes animating a quantity like dv/dt with a min max scale that varies with the field is useful. It can take surprisingly long for the numerical solution to reach fully steady and that can be a way to see it. That can also visually highlight unexpected asymmetry, if there is any lurking.
Is there any sort of [ -1/(alpha.dt) u ] term within the cylinder - in flow in porous media type problems that sort of thing is sometimes added.
I think the FFT solver itself is fairly isotropic and grid unaware. In essence it just projects the RHS onto a bunch of sin and cos bases, scales the wave space solution by appropriate eigenvalues and then evaluates the result back in physical space. The FFT is done in successive directions, but the end result is in machine precision usually and I don't think I would expect to see any preferential index direction in the solution. It is possible to permute the X and Y FFT computation order to see if that impacts.
The RHS will have a lot of grid and boundary approximation in it though - but it is all local I think, so not sure it knows about different direction sweeps etc...
christophernhill
on 4 Dec 2020
Related issues
ali-ramadhan
·
5Comments
truedichotomy
·
8Comments
ali-ramadhan
·
7Comments
glwagner
·
4Comments
ali-ramadhan
·
5Comments
Most helpful comment
Sounds great! I've had a discussion with Greg about this before, and which method you need to implement really depends on the application.
I've attach a write-up one of my students is working on which lays this out with a few simple examples. Section 3 reviews Immersed Boundary methods and 3.1 has a simple 1D FSI example.
Mr M Lauber_070e4204-4db5-451e-8b61-494f2ae9eaa2_Progression_Report__9_month_pdf_7177_0.pdf