[philippose]

Hello and a good day to everyone :-)!

Its amazing how quickly the year 2012 has gone by, and going by the fact that I am actually typing out a post on the forum, it looks like 21 Dec 2012 is in-fact not the "end of the world" ;-)!

I have a quick query regarding the Mesquite implementation in OpenFOAM by Sandeep and team....

Would it be possible to use Mesquite in order to perform a simulation residual based automatic refinement / coarsening of the mesh during a simulation, so that the mesh is automatically adapted as the simulation progresses?

The documentation of the Mesquite Library from Sandia Labs is not very clear about this, but it looks like the library does support residual based adaptation of meshes (http://acts.nersc.gov/events/Worksho...s/Mesquite.pdf - pages 3-5)

Have there been any thoughts or ideas in this direction? As far as I have understood, the current mesh refinement implementation in OpenFOAM (dynamicrefineFvMesh) only supports hexahedral meshes.... or am I mistaken?

Have a wonderful weekend ahead, and of course, a Merry Christmas, and Happy Holidays to the whole OpenFOAM community :-)!

Philippose

[deepsterblue]

Mesquite currently does not support topological changes. The dynamicTopoFvMesh library does. It currently adapts tetrahedral meshes (using edge refinement / collapse) based on a length-scale field, but that can easily be adapted to represent anything - including residuals.

All you would need to do is specify a scalar value at every cell, and a min/max threshold ratio based on which edges are refined.

[philippose]

Hi Sandeep,

Thanks a lot for your response. I was looking into the ITAPS project a little more in-depth, and found that there was a separate dedicated module called "MeshAdapt" for solution based mesh refinement / coarsening.

Coming back to the OpenFOAM implementation.... Wow :-)! I had not looked into the dynamicTopoFvMesh at all for some reason. Thanks a lot for the heads-up. Shall look into this option and get back to you with either results, or more questions :-)!

Is the solution-mapping in dynamicTopoFvMesh also working? In that.... are the variables mapped from the old mesh to the new mesh automatically as part of the dynamicTopoFvMesh library, or does this need to be done separately?

As of now, the one time I actually tried solution based mesh refinement, was done by:
1. Running the simulation in OpenFOAM
2. Using "simpleFoamResidual" to get the residual field.
3. Using this residual field to generate a "Mesh Size File" for Netgen.
4. Meshing the Geometry in Netgen by using the Mesh Size File for the local refinement / coarsening
5. Importing the mesh back into OpenFOAM
6. Mapping the solution from the old mesh to the new one
7. Rerunning the OpenFOAM simulation

I am quite sure this can be done completely within OpenFOAM.... hence the purpose of the post.



Now again switching back to ITAPS.... Have you worked with components of this project other than Mesquite? I was particularly interested in the EBMesh + MeshKit mix, which seems to be capable of creating Cartesian Meshes with local refinement around geometry specifically for use in "immersed boundary method" (or "Embedded Boundary" as ITAPS coins it) simulations.

Would this be relevant for OpenFOAM given the fact that the last Workshop also had a presentation and implementation of the Immersed Boundary Method in OpenFOAM?

Have a nice day!

Philippose

[deepsterblue]

All the sources are currently hosted here:

https://github.com/smenon/dynamicTopoFvMesh

I've implemented conservative solution mapping in the library, so all fields should be mapped automatically without problems. This should work with OpenFOAM-2.1.x as well - take a look at the Port-2.1.x branch.

It's possible that there's an ITAPS implementation for adaptation - I've never looked into it. Ideally, it's best to keep dependencies to a minimum.

[philippose]

Hello Sandeep,

Ah ok. I thought you were referring to the stock version of dynamicTopoFvMesh already present in the OpenFOAM-1.6-ext git repo.

I have just cloned your github repository, and finished compiling it. Everything went through fine, based on the "Install.txt" file in the repo.

Is there anything I need to know before trying to use the library? Would you happen to have a tutorial somewhere which might be useful as a "getting started" setup?

(I know that there is a circCylinder3d tutorial in the OpenFOAM-ext code)

Philippose

[deepsterblue]

The recent OpenFOAM workshop had a tutorial session with a document / presentation by Kyle Mooney. I guess you could look at those.

[philippose]

Hello again Sandeep,

A Good day to you!

So... I have been slowing putting together a modified version of simpleFoam in order to write out a field with the required mesh size derived from the residuals.

I have a couple of questions regarding the dynamicMeshDict settings for the field based refinement. The documentation from the 7th OpenFOAM Workshop by Kyle Mooney does not really touch on this topic.

I looked into the refinement settings used by interDyMFoam, and though the settings in your library are similar, there are some differences and hence, some questions... :-)!

1. What is implied by the setting: fieldLengthScale?
---- Is this something like a "zero" level for the fieldLengthScale at which no mesh changes would happen (neither refinement nor coarsening)?

2. What does the following setting do: meanScale?
---- This could also be interpreted as a "zero" level, or is this used for something else?

I am assuming that the "lower" and "upper" refinement levels are just the min and max for the lengthScale, which would act like saturations on the lower and upper ends for the mesh edge lengths.... right?

Now.... there are also the settings: minLengthScale and maxLengthScale. How do these play a role in the field based refinement?

Also... are the length scales specified in the lengthScale field used directly as the corresponding edge lengths for the remeshing, or is there some kind of scaling / mapping between the field values and the actual edge lengths?

Wishing you a great Saturday ahead!

Philippose

[deepsterblue]

Philippose,

It's been a while, so I had to look at the sources again to see what was going on.

1. For fieldLengthScale, the intention is to specify a mesh length-scale corresponding to the zero-level. This would imply a mesh edge-length at the "alpha" interface in interFoam, for example.

2. For meanScale, the intention is to specify a mesh length-scale that exists at regions sufficiently far away from the zero-level. This would imply a mesh edge-length at the alpha=0 or alpha=1 regions in interFoam.

I suppose these options don't really sound too intuitive, now that you mention it, but I couldn't really think of anything better at the time.

The lower and upper refinement levels are supposed to be equivalent to the one in interDyMFoam, where an integer value is used to define each level of refinement. In this case, there's no strict "unrefinement" operation that recovers the original mesh, but an edge collapse method, which roughly achieves the same effect.

The min/maxLengthScale values correspond directly to mesh edge lengths. It would probably help to take a look at the lengthScaleEstimator class under the dynamicTopoFvMesh directory to see how things are done.

[philippose]

Hello Sandeep,

A Happy New Year to you :-)! Hope you had a great Christmas break.

I apologise for the delay in replying to your post. I must have missed it. Its only today when I tried to find this post again to ask you something more that I saw that you had already replied to my earlier post.

Thanks a lot for the information. I looked into the lengthEstimator code before seeing your post, and I have some thoughts / questions regarding the implementation of the field based refinement.

1. From the code I understand that the field which drives the length scales is a "volScalarField", and hence, the length scale is defined on a per cell basis rather than a per point basis. This is great because this is what I had in mind earlier.

2. In the part where the field is actually applied to the length scales, I think there is a difference in the way you have implemented it compared to what I had in mind...... Basically, you are checking to see if the average field value of the current cell and the neighbour cell is between the specified "lower" and "upper" refinement levels, and if so, the length scale for that cell is set to the "lengthScale" value specified in the dyamicMeshDict.

To see if I have understood it right.... suppose:
... lengthScale = 0.1
... lowerRefinementLevel = 0.01
... upperRefinementLevel = 0.35
... Average field value = 0.2

In this case, the length scale for that cell and its neighbour would be set to 0.1 right?

That would imply, that the actual length scale of the mesh is not proportional to the field... rather... it is either the already existent mesh length scale (i.e. unmodified), or 0.1.

What I had in mind (and I could be wrong... so please correct me if I am) is, that I could specify a field, calculated based on say... the error residual during the simulation, and the lengthScaleEstimator basically proportionally maps the values in the field to the min and max length scales, giving a resulting length scale which mirrors the values in the field.

For example, if:

... minLengthScale = 0.1e-03
... maxLengthScale = 0.25e-03

... lowerRefinementLevel = 0.2
... upperRefinementLevel = 10.0

... Average field value = 2.0

Then, the resulting length scale for that cell (and the neighbour if applicable) should be:
lengthScale = [(2.0 - 0.2)/(10.0 - 0.2)*(0.25e-03 - 0.1e-03)] + 0.1e-03

And field values above and below the upper and lower refinement levels should be saturated to the upper and lower limits.

What do you think? I hope that makes sense?

Do you think something like this would be a good addition to the current implementation of the remeshing library? And if so, how complicated do you think the implementation would be?

Wishing you a great day ahead!

Philippose

[deepsterblue]

Philippose,

What you're saying makes sense. Again, you're free to use any approach that you find convenient. If you look at the dynamicTopoFvMesh class, you have a DynamicList<scalar> member called lengthScale_, which carries a value specified for each cell in the mesh, and this field that is checked during mesh manipulation. In this case, I fill values in them using the lengthEstimator class (check the dynamicTopoFvMesh::calculateLengthScale function for details), but you could easily fill it using any other volScalarField (will need some minor modifications to the function I mentioned).

What I propose is this - if you can find a scaling function that works for you, create a volScalarField at the top-level in your solver, and add that to the objectRegistry using an IOobject. Once that is done, you can add a dictionary entry to dynamicTopoFvMesh which takes the name of your registered field, and fills values for lengthScale_ from its internal field, rather than calling lengthScaleEstimator for it. Of course, the responsibility of maintaining that field at every time-step is yours. Since it is a registered volScalarField, the mesh will automatically resize it when the mesh changes, so you can be guaranteed that it is properly sized.

Savvy?

[philippose]

Hello again Sandeep,

Cool :-)! Thanks a lot for those tips. I was under the impression that the lengthScaleEstimator was a mandatory part of your system, which cannot be replaced as such.

What I have done so far, is to create a modified version of simpleFoam which writes an additional "volScalarField" containing what in future would be an indication of the length scale for the mesh based on the simulation.

I could change that to make it not just an indicator, but the actual required length scale, which could then be passed on as a field and read in by your library.

I need to have a more detailed peek into the dynamicTopoFvMesh class to determine what the flow of information is like.

Will a change like the one you suggested break any other parts of the library? Or is the lengthScale completely encapsulated and as of now modified only via the lenghScaleEstimator?

In case I have to modify the library, what do you think is the best way to do so? Create a local git branch?


In a way I am actually surprised that no one has wanted a similar solution based re-meshing system for tetrahedral meshes within OpenFOAM so far..... Hmmmm....


Am pretty sure I will be back in a short time with more questions :-)!

Have a nice day ahead!

Philippose

[deepsterblue]

The change I proposed should work. Since the lengthScale_ field is dynamic, values are typically copied or averaged when new cells are added, so I'm assuming that's not a problem for you.

If you're interested in modifying the library, create a fork on GitHub, and once you've tested changes and have it working reliably, initiate a Pull Request, and I can merge it. You have several possible branches to work with - either 1.6-ext or 2.1.x. Pick whichever suits you best.

[philippose]

Hello again Sandeep,

Hope you are having a great weekend!

As you might have already seen, I forked your dynamicTopoFvMesh project in Github and have been working on modifying the code.

I basically finished everything, and tried running a test case using an actual flow solver (variant of simpleFoam) rather than just moveDynamicMesh.

Looks like there might be some issue. In both, the original code, as well as my modified version, the simulation crashes with a segmentation fault when a remeshing is performed (I set the remeshing interval to something like 100, so though the mesh quality and length scales are calculated in each time-step, the actual remeshing only happens every 100 iterations).

After poking around a bit with status outputs at various points in the code, I have found, that the code crashes at the following location:

dynamicTopoFvMesh/fieldMapping/conservativeMapFields.H

between lines 69 and 88.

Surprisingly, it seems to reaching this location whether I have mapping enabled or not (via the "skipMapping" option). I assume that the actual case of whether a mapping is happening or not takes place one level deeper, using either a dummy mapper, or an actual mapping routine.

Just to make it a little clearer.... both the code versions are working fine when the simulation is run with "moveDynamicMesh". It only crashes when executed with an actual fluid solver.

And... I have not made any changes to this part of the code. However, I guess that is more or less clear because the same problem occurs in the original code.


Would it be possible for you to take a quick look into the original code? Or am I missing something?

(P.S.... The actual modification of the length scales seem to be working quite well.... I am still calling "calculateLengthScale". Once things are working you can see whether the approach I have used is acceptable to you or not)

Have a great Sunday!

Philippose

[deepsterblue]

I can take a look. Post a test case at a location that I can access.

[philippose]

Hello Sandeep,

I have copied the test case and the modified version of simpleFoam into my google drive. You can access them using the following link:

https://docs.google.com/folder/d/0B0syP1Z1w2s1YzZiY2UzMTQtNDBhOC00MmRmLWIzZWYtYjJjM TVlYjExZTk3/edit

In the folder, the two relevant files are:

1. mySimpleFoam.tgz - Modified version of simpleFoam to handle dynamic meshes, as well as writing out a field which indicates the length scale

2. 001_Benchmark_Blende_v03_Remeshing_001.tgz - This is a simple test case of an orifice with an inlet and outlet.

In order to run the case, you need to compile "mySimpleFoam", and use this solver.

You can use the original version of the dynamicTopoFvMesh library from your Guthub account without any modifications.

As it is set up now, the remeshing happens in the 10th iteration, and thats when the simulation quits with a segmentation fault, irrespective of whether mapping is enabled or not.


Hope this helps. In case you have trouble accessing the shared files, please let me know.

Philippose

[deepsterblue]

Quick question - did you try this with the Port-2.1.x branch? If not, could you do that and tell me if you get the same error? Unfortunately, I updated my OpenSUSE machine to 12.2, and 1.6-ext no longer compiles with the stock gcc-4.7.1 version, but 2.1.x already has all the source changes in place.

[philippose]

Hello Sandeep,

Nope, I have never used the official OpenFOAM releases yet. Always stuck to the OpenFOAM-ext versions.

I dont know if you have the time / inclination..... I had made a branch on the OpenFOAM-ext Git repository called "hotfix/Gcc47" which I recently updated (merged from the Master) to allow OpenFOAM-1.6-ext to be compiled with GCC-4.7.x (works on my Fedora 17 64-bit VM, which I am currently using).

However, in the meantime I can try to test the library and the case on a version of OpenFOAM-2.1.x once I have the base OpenFOAM-2.1.x version up and running on my system.

Did you try the case on OF-2.1.x with the "mySimpleFoam" solver from me? I dont know if it directly compiles on a 2.1.x version without any changes to the source.

Technically speaking, you can also try the case using "pimpleDyMFoam" in 2.1.x.... you will need to make the relevant changes to "fvSolution" and "fvSchemes" though.

Philippose

[deepsterblue]

I've actually merged the hotfix/gcc47 into master and am compiling right now... I've also tried compiling with 2.1.x, but the errorEstimation directory doesn't exist, and don't know where its been relocated to.

It would be interesting to know if this has any issues in OpenFOAM-2.1.x, just to eliminate the possible compiler-related issues.

[deepsterblue]

I actually ran it with with master version just now, and everything worked okay... Here's a few screenshots at 0, 10 and 20 (attached), and relevant bits of the log.

Granted, the mesh quality is not that great, but you can play around with that.

[philippose]

Wow... thanks a lot for the information!

Oh ohhhh :-O! Now thats a little worrying.... So basically you ran "mySimpleFoam" on the case I provided you, without any modifications to the case? And it just worked?

Hmmm.... I put in a lot of "Info" statements, and converged on the location of the crash to be the function call: conservativeMapVolFields

Any ideas what the problem on my end could be?