|
[Sponsors] |
April 11, 2012, 11:16 |
Track vortices using Swak4FOAM
|
#1 |
Senior Member
Vincent RIVOLA
Join Date: Mar 2009
Location: France
Posts: 283
Rep Power: 18 |
Hi all,
I am simulating the propagation in space and time of some vortices. I would like to be able to trace the center of these vortices into my field and then be able to compute the circulation on a closed surface around this center. So, first step I need to trace the center of my vortex. It looks like defining some functions in the controlDict and using the libsimpleSwakFunctionObjects could be an option but I am not sure how to do it. To start easy I was planning to find the center of my vortex by looking for the maximum vorticity in my field. I wrote those lines in my dictionary: Code:
functions ( defineFollowVortex { type addGlobalVariable; outputControl timeStep; outputInterval 1; globalScope followVortex; globalVariables { yThres { valueType scalar; value 0.2; } } }; createVortexCenter { type createSampledSet; outputControl timeStep; outputInterval 1; setName vortexCenter; set { type cloud; axis y; points( (0.5 2 0.0) ); } } vorticityMax { type swakExpression; verbose true; valueType set; setName vorticityMax; expression "(pos.y()>yThres) ? magVorticity : 0"; accumulations ( max ); set { type swakRegistryProxy; axis y; setName vortexCenter; } interpolate true; interpolationType cellPoint; } vortexHeight { type swakExpression; verbose true; valueType set; setName fillingHeight; globalScopes ( "fillDam" ); expression "(magVorticity >= 0.99*vorticityMax) ? pos().y : 0"; accumulations ( max ); set { type swakRegistryProxy; axis y; setName vortexCenter; } interpolate true; interpolationType cellPoint; } ); Code:
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // Create time [0] [0] [0] --> FOAM FATAL IO ERROR: [0] "ill defined primitiveEntry starting at keyword 'functions' on line 60 and ending at line 124" [0] [0] file: /CFDBB/vincent/OpenFOAM/vincent-2.1.0/run/tutorials/compressible/rhoPimpleDyMFoamWithSources/LES/AfterMeeting120330/testVortexHeight/processor0/../system/controlDict at line 124. [0] [0] From function primitiveEntry::readEntry(const dictionary&, Istream&) [0] in file lnInclude/IOerror.C at line 132. [0] FOAM parallel run exiting |
|
April 11, 2012, 15:53 |
|
#2 | |
Assistant Moderator
Bernhard Gschaider
Join Date: Mar 2009
Posts: 4,225
Rep Power: 51 |
Quote:
Anyway. What might help you is minPosition/maxPosition (this is currently only available in the development-version of swak. You'll have to get that via Mercurial): they take an expression and return the position where that expression has its minimum/maximum (you'll have to guess from the names which does which) |
||
April 12, 2012, 03:08 |
|
#3 |
Senior Member
Vincent RIVOLA
Join Date: Mar 2009
Location: France
Posts: 283
Rep Power: 18 |
Dear Bernhard,
Thank you for your reply. Indeed it looks like minPosition and maxPosition could be a good starting point. I went to the wiki in order to see how to download the development version. I got it and apparently managed to install in on my openfoam version 2.1. I looked at the examples and did not find one including the min or max Position function. Could you tell me if something like this should work: Code:
vorticityMax { type swakExpression; verbose true; valueType pointField; setName vorticityMax; expression "maxPosition(magVorticity)"; } Or do you have an example? |
|
April 12, 2012, 09:08 |
|
#4 |
Senior Member
Vincent RIVOLA
Join Date: Mar 2009
Location: France
Posts: 283
Rep Power: 18 |
I got it working. At least I can track one vortex center with this code:
Code:
functions ( maxVortPosY { type swakExpression; valueType internalField; expression "maxPosition(mag(vorticity)).y"; accumulations ( max ); verbose true; } maxVortPosX { type swakExpression; valueType internalField; expression "maxPosition(mag(vorticity)).x"; accumulations ( max ); verbose true; } ); |
|
April 12, 2012, 15:52 |
|
#5 |
Assistant Moderator
Bernhard Gschaider
Join Date: Mar 2009
Posts: 4,225
Rep Power: 51 |
Great. Sorry for a lack of an example. I tested that feature with another case (that I can't publish). But if you have a simple example I'll be happy to include it into swak
More will be difficult with the min/max-approach. Unless you "mask out" the region of the first vortex Sorry. PyFoamPlotWatcher can currently only do plots over time. Thanks for the reminder on that missing feature (but of course you can easily plot the data written by it easily with any program you like) |
|
April 13, 2012, 01:56 |
|
#6 | |
Senior Member
Vincent RIVOLA
Join Date: Mar 2009
Location: France
Posts: 283
Rep Power: 18 |
Quote:
I am not sure how to handle that. I would need to say something like "track the maximum vorticity only where vorticity.z is positive/negative" but I'm not sure if there is a way to do it. Nothing to be sorry of, you already do a lot! When you say that I can plot the data written by PyFoamPlotWatcher, are you talking about the file in the different directories swakExpression_*? or is there another way to use the data analysed by PyFoamPlotWatcher in order to create the graph I want? Last edited by vinz; April 13, 2012 at 02:32. |
||
April 16, 2012, 02:38 |
|
#7 |
Senior Member
Vincent RIVOLA
Join Date: Mar 2009
Location: France
Posts: 283
Rep Power: 18 |
Dear Bernhard,
do you think it is possible to create a subset of my mesh based on the vorticity value or on a coordinate for instance and then look for a maximum only in that subset? |
|
April 16, 2012, 04:50 |
|
#8 | |
Assistant Moderator
Bernhard Gschaider
Join Date: Mar 2009
Posts: 4,225
Rep Power: 51 |
Quote:
vortPos1=maxPosition(vorticity) maskField=(mag(pos()-vortPos1)<rCut) ? 0 : 1 vortPos2=maxPosition( maskField ? vorticity : -1000) maskField=(mag(pos()-vortPos2)<rCut) ? 0 : maskField vortPos3= .... Of course the criteria for masking could be something more advanced (the problem here is that rCut being to small might give you "vortex centers" that are in the other vortex and being too large might "hide" another vortex center) Alternative would be "dynamic" cellSets based on an expression, but that feature is not yet implemented (partly because I'm not sure whether the parts of the OF-code where this could be interesting can "cope" with cellSets that change their size) |
||
April 16, 2012, 10:32 |
|
#9 |
Senior Member
Vincent RIVOLA
Join Date: Mar 2009
Location: France
Posts: 283
Rep Power: 18 |
Ok I see, not straight forward but brilliant!
I manage to get something which looks interesting with the following lines: Code:
mask { type swakExpression; valueType internalField; variables ( "maskField=(vorticity.z>0 || pos().y<0.3) ? 0 : 1;" ); storedVariables ( { name maskField; initialValue "0"; } ); expression "maskField"; accumulations ( max ); verbose true; } maxVortPos { type swakExpression; valueType internalField; variables ( "vortPos2= maxPosition((maskField>0) ? magVorticity : -1000);" ); storedVariables ( { name vortPos2; initialValue "vector(-0.4,1,0)"; } ); expression "vortPos2"; accumulations ( max ); verbose true; } maxVortPosY { type swakExpression; valueType internalField; expression "vortPos2.y"; accumulations ( max ); verbose true; } maxVortPosX { type swakExpression; valueType internalField; expression "vortPos2.x"; accumulations ( max ); verbose true; } I tried the following to create such a sphere : Code:
distField { type expressionField; outputControl timeStep; outputInterval 1; fieldName distanceFromVortex; expression "mag(pos()-vortPos2)"; autowrite true; } surfaceIsoDist { type swakExpression; valueType surface; surfaceName isoVortexDist; surface { type isoSurface; isoField distanceFromVortex; isoValue 0.15; interpolate true; } verbose true; expression "mag(U)"; accumulations ( min max average ); } Then I am wondering what type is expected at this place? Any idea? |
|
April 16, 2012, 14:04 |
|
#10 | |
Assistant Moderator
Bernhard Gschaider
Join Date: Mar 2009
Posts: 4,225
Rep Power: 51 |
I think it is the other way round: straight forward (brute force) but not brilliant/elegant
Quote:
a) putting vortPos2 in an expression field b) using a global variable (see slide 35ff in http://openfoamwiki.net/images/d/da/...Leoben2011.pdf). Something like (untested) Code:
calculateVortex { type calculateGlobalVariables; outputControl timeStep; outputInterval 1; valueType internalField; toGlobalNamespace vortexLoc; toGlobalVariables ( vortexPos2 ); variables ( "vortPos2= maxPosition((maskField>0) ? magVorticity : -1000);" ); } Code:
globalScopes ( vortexLoc ); |
||
April 19, 2012, 02:47 |
|
#11 |
Senior Member
Vincent RIVOLA
Join Date: Mar 2009
Location: France
Posts: 283
Rep Power: 18 |
Thanks again for your help Bernhard!
I managed to get it working, and to track my two vortices. Here is an example of what I was trying to do: |
|
April 19, 2012, 07:23 |
|
#12 | |
Assistant Moderator
Bernhard Gschaider
Join Date: Mar 2009
Posts: 4,225
Rep Power: 51 |
Quote:
BTW: looking at the picture it occurred to me that by searching for the minimum and the maximum of "vorticity*pos().x" the "masking" could have been avoided |
||
April 19, 2012, 10:52 |
|
#13 |
Senior Member
Vincent RIVOLA
Join Date: Mar 2009
Location: France
Posts: 283
Rep Power: 18 |
Honestly I have never user the wikki so I have no idea how to include this
However, it is not urgent, it's an ongoing work so I cannot share everything right now. The contract ends on july and I should have other interesting (hopefully) thing to show afterwards and place it on the wikki. Regarding the other possibility to track the vortices, I think it is less general. And although it might work in this specific case, I have other cases with wind comming from the left or where one vortex might turn around the other and I think your last proposal would not work anymore. Thanks again for your help. |
|
April 23, 2012, 19:18 |
|
#14 | ||
Assistant Moderator
Bernhard Gschaider
Join Date: Mar 2009
Posts: 4,225
Rep Power: 51 |
Quote:
It's quite easy. Really Quote:
But I guess there is no general solution that will always work |
|||
|
|
Similar Threads | ||||
Thread | Thread Starter | Forum | Replies | Last Post |
Why my track file is 0 K after the simulation is saved | calmdown911 | STAR-CCM+ | 0 | December 19, 2017 10:22 |
Hairpin vortices in turbulent channel flow using DNS | Nitin Minocha | Main CFD Forum | 2 | November 21, 2016 17:52 |
[swak4Foam] Rotating Wall Velocity and Swak4Foam | mpanch13114 | OpenFOAM Community Contributions | 0 | September 14, 2016 04:55 |
[swak4Foam] and twoPhaseEulerFoam | mnikku | OpenFOAM Community Contributions | 1 | February 19, 2016 05:21 |
[swak4Foam] fails in parallel with -otherTime? | Phicau | OpenFOAM Community Contributions | 3 | June 26, 2013 13:00 |