Creating highly incompressible viscou... Log Out | Topics | Search
Moderators | Register | Edit Profile

FlexPDE User's Forum » User Postings » Creating highly incompressible viscous flow « Previous Next »

Author Message
Top of pagePrevious messageNext messageBottom of page Link to this message

Jared Barber (jared_barber)
Member
Username: jared_barber

Post Number: 30
Registered: 01-2007
Posted on Wednesday, April 16, 2008 - 07:02 pm:   

Hello,

I have moving particles hitting walls. In real incompressible flow, lubrication theory tells us (for these particular particles and walls) that the particles should never hit the walls. To prove that this should be the case, I have taken the "penalty" used in "viscous.pde" (grad^2(p) = penalty*div(u)) and ramped it up to 10000 (we set the penalty to be a constant in our equations since the length scale is O(1) always). In this case, the particles don't hit the wall suggesting the reason the particles are hitting the wall in penalty = 100 case is because the incompressibility condition in these gaps (where the particles hit a wall) is not being enforced as well as it could be.

The ultimate goal, then, is to enforce more incompressibility in regions where particles are in danger of hitting walls. One way to do this is to assign penalty = 10000 everywhere for all of our runs. The problem is that this increases computational time a lot. Our simulations go from 1 min per run to 10 min per run. That is what happens if you take our code and replace penalty = 100 with penalty = 10000.

Another way to do find out what happens when penalty = 10000, however, is to take and make this into a staged problem by setting penalty = staged(100,10000). One would hope this would speed up the computation since a really good guess for the 10000 case would be obtained from the 100 case. Unfortunately, this does not seem to work as well as it could, 20 min down to 8 min (e.g.).

So, the other idea I came up with is to take the penalty and make it a function of time. That is penalty = (e.g.) 100*10^(2*t). Then you solve for time from t = 0 to 1.

Surpisingly enough this works well for our problem (which involves 42 scalar equations, a lot) taking approx 15-20 min of computation down to 4 min. I have attached a similar smaller problem (with only 1 scalar equation and a much simpler geometry...because of this just setting penalty = 10000 works faster than trying the attached time approach...this is not the case for our other problem). By looking at this problem, you can see how it works.

One problem is I don't know the accuracy of the solution returned. For instance, if I insert the correct initial values for u, v, and p (obtained from running Lube_Flow with k = .01) the values from each method (setting penalty = 10000 or using the time approach) actually agree relatively well whereas if I don't insert the correct initial values, the values can differ a bit (compare u0 for instance). I'm not sure why this is the case since there is no explicit time dependence in this problem (give a penalty, u, v, and p depend only on the value of that penalty and not on any values of any parameters from previous or subsequent time steps).

So the question is, how accurate will this time-dependent penalty approach be in the end?

As a side note, is there anyway to obtain the distance from (x,y) in the domain to the nearest boundary?

Thanks,

Jared
application/octet-streamLube_Flow
Lube_Flow.pde (2.3 k)

application/octet-streamLube_Flow_Dynamic
Lube_Flow_Dynamic.pde (2.3 k)
Top of pagePrevious messageNext messageBottom of page Link to this message

Robert G. Nelson (rgnelson)
Moderator
Username: rgnelson

Post Number: 1101
Registered: 06-2003
Posted on Friday, April 18, 2008 - 02:53 pm:   

Rather than address the time dependent problem, I went back to the beginning and tried to speed up the staged version.

If you look at the U plots, you will see there are some twists at the exit planes. This indicates a conflict between the boundary condition and the needs of the PDE. Such a conflict can cause the solution to struggle, because there is really no solution that satisfies both. I replaced your outflow condition with the "tautological" ones,
Natural(U)=sigma12 Natural(V)=sigma22 at the top
Natural(U)=-sigma12 Natural(V)=-sigma22 at the bottom.
This improves the smoothness of the solution at the exit, but there is still a whip in the contours, indicating that P=0 is not really a consistent BC either.

I made a four stage run, starting with high viscosity and low penalty, then lowering the viscosity and raising K, so that K*penalty goes from 100 to 100000.

I added histories to track U0 and dx(U) as the penalty increases. These plots show the consistency of the result at the higher penalties.

This script completes in 6 seconds on my 1.5GHz Athlon.

application/octet-streamLube_Flow2
Lube_Flow2.pde (2.5 k)

U0 vs stage
Top of pagePrevious messageNext messageBottom of page Link to this message

Jared Barber (jared_barber)
Member
Username: jared_barber

Post Number: 31
Registered: 01-2007
Posted on Sunday, May 11, 2008 - 12:00 pm:   

Hello,

I have an even more challenging problem in terms of incompressibility. I found out that using a space-dependent penalty actually seems to the job ok, for arcs approaching lines.

By the way, once again, is there any quick way to find out the d(x,y) where d(x,y) is the minimum distance from the point (x,y) to all segments of the boundary?

Now, the new problem is that I want to keep corners from hitting lines. I have attached one time step from my full-fledged code in a zipped file so that what I am talking about can be seen.

In the simulation, a "corner" of one object is very close to a "wall" of a different object. In my code it is realistic to assume that the "corner" will oppose approach to the "wall" (experiencing forces on the order of 1/(distance to wall)^3 pushing it away from the wall). As such there are a select few equations with very large coefficients (that are trying to enforce slow approach) while all of the rest have O(1) coefficients. I want to solve this system but it seems FlexPDE is having difficulties doing this for me.

Way #1: Try to solve the system straightforward. (To do this, just open "Norm_Flow.pde" and run it.) Problem: this takes in excess of 20 minutes (I have yet to have had the patience to let it finish.) In addition, if I coax it by using some nice transfer data (including initial mesh and good initial guess for u, v, and p), I once got it to solve in 16 minutes. Even then it seems that while all the scalar equations should equal 0, they actually equal O(1000). In addition, the corner is still approaching the wall much faster than it should be.

Way #2: Turn on the corner forces using a staged variable ("lub_on" in my code, just change lub_on in c_vs_defs.txt file to be staged(0,1) or staged(0,0.5,1) or anything like that, same results). Problem, the nodal equations corresponding to the corner node and the colliding wall nodes are not 0 but rather O(1000000).

So the question is, how do I force the equations governing the movement of the corner and the wall near the corner to be 0? It seems if I use staging, those equations are allowed to be O(1000000). If I don't use staging, it seems it tries not to let it be O(1000000) but it takes a long long time for it to get an answer which is still not that good in the end. Is there some sort of request I can make to enforce those special equations more than the rest? (And why does staging allow the residual on such scalar equations to be so big when unstaged runs don't?)

Please let me know.

Also, I am fairly certain that I have no small errors lurking beneath the mounds of include files I have included but I can always try to make a simpler example to look at if anyone would think that would be more helpful.

Hoping for a couple suggestions,
Jared
application/x-tar
For_Posting_FlexPDE.tar (3358.7 k)
Top of pagePrevious messageNext messageBottom of page Link to this message

Robert G. Nelson (rgnelson)
Moderator
Username: rgnelson

Post Number: 1119
Registered: 06-2003
Posted on Sunday, May 11, 2008 - 04:52 pm:   

As I believe I mentioned before, I think the process of scaling equations is not working well when there are many GLOBAL equations. Frankly, we never expected anyone to have so Many....

We have discovered that the problem behaves much better if you turn matrix scaling off and use a much smaller overshoot:
SELECT SCALING=OFF OVERSHOOT=1e-6

The conjugate-gradient linear system solver iterates until the solution updates are resolved to ERRLIM*OVERSHOOT. If some components of the solution are much larger or smaller than the others, requiring a tighter convergence helps resolve the smaller components.

With these changes, the problem you sent completes in 18 minutes on our 3GHz Pentium, with
global Ux_x and Vx_x values like 1e-9.

Top of pagePrevious messageNext messageBottom of page Link to this message

Jared Barber (jared_barber)
Member
Username: jared_barber

Post Number: 32
Registered: 01-2007
Posted on Sunday, May 11, 2008 - 09:52 pm:   

Hello,

I am very sorry but I completely forgot to turn the "problem" on. I'm not sure what possessed me. Anyway, I have attached a new version with the particular "lub_on" parameter set to 1 rather than 0. I think if you run this you may discover different results. (Hopefully along the lines of what I spoke of before.)

Jared
application/x-tar
For_Posting_FlexPDE.tar (3389.4 k)
Top of pagePrevious messageNext messageBottom of page Link to this message

Marek Nelson (mgnelson)
Moderator
Username: mgnelson

Post Number: 34
Registered: 07-2007
Posted on Monday, May 12, 2008 - 09:00 pm:   

I can't explain why, but it appears that by introducing the 'lub_on' components the matrix is ill-conditioned enough that FlexPDE does not converge. By gradually introducing the 'lub_on' factor in stages, as one might do in a non-linear problem, FlexPDE is able to converge.

I staged the problem with lub_on = STAGED(0, 0.01, 0.1, 0.25, 0.5, 1) and the problem ran in about 30 min.
(I also used the previously mentioned SELECT SCALING=off)

I didn't do this, but you can use HISTORY plots to display the changes to your global variables as the stages are applied. This would help you to pick more meaningful staged values of 'lub_on'. I just picked a few arbitrarily.

Also, now that the problem converges, you may want to decrease your error limit. You are using errlim=0.01 which is fairly loose.

Add Your Message Here
Post:
Username: Posting Information:
This is a private posting area. Only registered users and moderators may post messages here.
Password:
Options: Enable HTML code in message
Automatically activate URLs in message
Action:

Topics | Last Day | Last Week | Tree View | Search | Help/Instructions | Program Credits Administration