Contrib/groovyBC

A library that introduces a boundary-condition groovyBC. This boundary condition is basically a mixed-BC where value, gradient and valueFraction are specified as expressions instead as fields.

It can be used to set non-uniform boundary-conditions without programming.

Valid versions:

Note: Active development of groovyBC now takes place in swak4Foam so please have a look there for more up-to-date versions. Support for this version of groovyBC is discontinued with 2.0.

1 About

This library is provided as-is and is a permanent Beta-Version (but it works fine for me, if it doesn't for you: tell me about it and I'll see what the problem is)

1.1 Words of warning

This library can save you the work to program your own boundary-conditions, but

• you should be familiar with the C++ expression syntax
• it makes it easier to 'shoot yourself in the foot' (do stupid things)
• especially for large cases a custom-made library might be more efficient

(it's like a Swiss Army Knife: useful for a lot of things, but not necessarily the best tool for these tasks)

1.2 Pre-requisites

To compile this utility Bison has to be installed. The compilation is known to work with version 2.4.1. It may probably work with older versions. Check with

bison -V

on the command line before trying to compile it and report the version number when reporting a problem.

2 Usage

2.1 Compilation and installation

First, download all of the source files from the bottom of this page (see section six)

Then, just do

wmake libso

in the directory of the sources. The library will be compiled and installed in a place (%FOAM_USER_LIBBIN) where it is usable

2.2 Older versions of bison

For systems with older versions of bison, the compilation sometimes fails. A makeshift bypass could be implemented as follows:

1. first compile the library (as described above) on a machine with the correct version of bison (call this machine A)
2. copy the entire directory of the generated files to the machine with the old version of bison (call this machine B)
3. On machine B, edit the Make/files by changing

   PatchValueExpressionParser.yy to PatchValueExpressionParser.C

   PatchValueExpressionLexer.ll to PatchValueExpressionLexer.C

4. go back to the parent directory on machine B and execute the command:

wmake libso

2.3 Additions to controlDict

A function object is activated by adding an entry like this to the system/controlDict of a case:

libs ( "libgroovyBC.so" ) ;

If you are experiencing problems with paraFoam try this combination.

libs ( "libOpenFOAM.so" "libgroovyBC.so" ) ;

2.4 Usage of the boundary condition

Just set the type of a patch to groovyBC

2.5 Parameters in the patch

valueExpression 

String with the value to be used if a Dirichlet-condition is needed. Defaults to zero

value 

is used if no "valueExpression" is given. value is also used for the first timestep/iteration if "valueExpression" is specified. If "valueExpression" is specified without setting "value", 0 is taken for the first timestep/iteration. (might cause a Floating Point Exception)

gradientExpression 

String with the gradient to be used if a Neumann condition is needed. Defaults to zero

fractionExpression 

Determines whether the face is Dirichlet (1) or Neumann (0). Defaults to 1

variables 

List with temporary variables separated by a semicolon. May make the writing of expressions shorter. Defaults to empty. Names defined here "shadow" fields of the same name

timelines 

List with subdictionaries that specify interpolation tables over time. See the original timeVaryingUniform-condition. Currently only scalars are allowed. The parameter name specifies the name under which this may be accessed. The name "shadows" fields of the same name

3 Expression syntax

The most complete documentation of the expression syntax is the source file for the Bison-grammar (*.yy and *.ll). Sorry.

3.1 These C++ operators are implemented:

+,-,*,/ 

Arithmetic operators. Can be used for vectors and scalars (only if useful. For instance: vectors can't be added to scalars)

 % 

A modulo-operator. Somehow differently defined from the standard C++ %-operator: the value of $(x-\epsilon)$ % $x$ is $-\epsilon$ (not $x-\epsilon$ - if you don't understand what I mean, please test for yourself - for most applications of this operator this is in my opinion the more practical implementation)

&,^ 

The vector operators as defined by OpenFOAM

<,>,<=,>=,!=,== 

Comparison operators (only defined for scalars)

&&,|| 

Logical Operators

Conditional operator

The conditional operator ( test ? val1 : val2) is defined for scalars and vectors

Operator precedence should be the same as for C++.

All the fields that are defined on the patch can be used. If the field is also the field that this is a BC for the old value is used.

3.2 These pseudo-variables are defined:

pi 

Guess ;)

These functions are defined:

pow,log,exp,sqr,sqrt,sin,cos,tan 

Only defined for scalars

  example: 
  x^y = pow(x,y)

mag 

defined for scalars and vectors

min,max 

minimum and maximum of a scalar field

sum 

sum of the values for a scalar field

3.3 These pseudo-functions are defined:

toPoint 

Takes a scalar (or vector) field defined on the faces and interpolates it to the points

toFace 

Takes a scalar (or vector) field defined on the points and interpolates it to the faces

pos 

Vector field with the face-centers

pts 

Vector field with the vertices

normal 

Vector field with the face normals

Sf 

field with the surface vectors (mag will give you the area)

rand 

Scalars-field with random numbers from [0,1]

randNormal 

Random-number scalar field that is Gauss-distributed

deltaT 

a field that returns the time-step

time 

a field that returns the current time

3.4 Advanced functions you're probably not going to need are:

Cn 

neighbour cell centers of the patch

delta 

Return cell-centre to face-centre vector

weights 

Return patch weighting factors

3.5 Functions that need another field are:

snGrad 

Gradient of that field

internalField 

internal values of that field

neighbourField 

neighbour values for a coupled patch

3.6 Accessing fields from other patches

groovyBC provides the ability to access fields from remote patches (other patches present within the current simulation). When evaluating the expression in run-time, the state of the fields on the remote patches from the previous time-step or iteration are utilised. The syntax for this feature is:

  <variable name>@<patch name> = <expression>

Note: Any field names used within the <expression> refer to fields local to the remote patch on which the variable is being defined.

Note: New versions of groovyBC included with swak4Foam replace the syntax <variable name>@<patch name> = <expression> with <variable name>{<patch name>} = <expression>

A variable defined as specified above can then be utilised like an ordinary variable within other groovyBC expressions defined for the local patch.

  Example:
   outlet1
   {
       type groovyBC;
       valueExpression "0.5*(pInlet+pOutlet2)";
       variables "pInlet@inlet=sum(p*mag(Sf()))/sum(mag(Sf()));pOutlet2@outlet2=p;";
       value uniform 100010;
   } 

In the above example, the variable pInlet has been defined on the remote patch inlet, and accesses the pressure and surface normal fields of the inlet patch. Likewise, the variable pOutlet2 has been defined on the remote patch outlet2, and accesses the pressure field of the outlet2 patch.

Both these remote variables are then used within the valueExpression of the local patch outlet1 to determine the pressure field to be imposed on this patch.

It is very important to note, that no form of patch-to-patch interpolation of the fields is performed. If the number of faces on the local and remote patches are not identical, entire fields cannot be mirrored from remote patches to local patches. In the example above, the number of faces on the patch inlet is not identical to the number of faces on the local patch outlet1, and hence an area weighted average of the pressure field is used. On the other hand, the number of faces on the remote patch outlet2 matches that of the local patch. In this case, the remote pressure field can be directly mapped onto the local patch.

4 Usage Examples

4.1 Demo-Cases

There are four demo-cases included in the repository in a folder Demos. All of these cases only require a blockMesh.

Don't try to find experimental results for these cases. They are what the name of the folder says: Demos

The three cases are

4.1.1 pulsedPitzDaily

To be run with oodles. Runs the pitzDaily-case with a parabolic inlet condition that pulsates

   inlet           
   {
     type            groovyBC;
     variables "yp=pts().y;minY=min(yp);maxY=max(yp);para=-(maxY-pos().y)*(pos().y-minY)/(0.25*pow(maxY-minY,2))*normal();";
     valueExpression "10*(1+0.5*sin(500*time()))*para";
     value           uniform (10 0 0);
   }

on the outlet the regular inletOutlet-BC is emulated (just as a demonstration)

   outlet          
   {
     type            groovyBC;
     valueExpression "vector(0,0,0)";
     gradientExpression "vector(0,0,0)";
     fractionExpression "(phi > 0) ? 0 : 1";
     value           uniform (0 0 0);
   }

4.1.2 wobbler

A solidDisplacementFoam case where a bar that is fixed on one side is pushed according to a predefined timeline on a fraction of the boundary (the rest is zero gradient). On the relevant patch the boundary conditon for D is

 forced
 {
   type groovyBC;
   value uniform (0 0 0);
   //    valueExpression "vector(0,0.1*sin(time()),0)";
   timelines (
       {
           name impulse;
           outOfBounds clamp;
           fileName "$FOAM_CASE/impulse.data";
       }
   );
   valueExpression "-impulse*normal()";
   gradientExpression "vector(0,0,0)";
   fractionExpression "(time()<5) ? ((pos().x>0.45 && pos().x<0.55) ? 1 : 0) : 0";
 }

4.1.3 circulatingSplash

A interDyMFoam-case that crashes halfway through the calculation (but looks quite nice until then). Whoever wants to fix that case is welcome to do so

gamma is set on top on a small circular part that moves around during time:

   atmosphere      
   {
       type            groovyBC;
       variables "r1=0.25*(max(pts().x)-min(pts().x));r2=r1*0.2*(1-0.5*cos(54*time()));";
       valueExpression "(sqrt(pow(pos().x-r1*sin(10*time()),2)+pow(pos().z-r1*cos(15*time()),2))<r2) ? 1 : 0";
       value           uniform 0;
   }

The Velocity is set depending on gamma:

   atmosphere
   {
       type            groovyBC;
       valueExpression "-(gamma+internalField(gamma))*0.5*normal()";
       value           uniform (0 0 0);
   }

4.1.4 average-t-junction

This is a turbFoam case, which showcases the ability of groovyBC to evaluate and access fields from remote patches (other patches present within the simulation case).

The case consists of a T-junction, with one inlet and two outlets. The total-pressure at the patch inlet is varied as a function of time. The patch outlet2 has a fixed pressure, whereas the pressure at the patch outlet1 is set using groovyBC to always be the average of the pressures at the inlet at outlet2 patches.

The corresponding boundary conditions for p are:

   boundaryField
   {
       inlet
       {
           type timeVaryingTotalPressure;
           p0 100040; // only used for restarts
           outOfBounds clamp;
           fileName "$FOAM_CASE/constant/p0vsTime";
    	
           U U;
           phi phi;
           rho none;
           psi none;
           gamma 1;
           value uniform 100040;
       }
   
       outlet1
       {
           type groovyBC;
           valueExpression "0.5*(pInlet+pOutlet2)";
           variables "pInlet@inlet=sum(p*mag(Sf()))/sum(mag(Sf()));pOutlet2@outlet2=p;";
           value uniform 100010;
       }
   
       outlet2
       {
           type fixedValue;
           value uniform 100000;
       } 
   
       defaultFaces
       {
           type zeroGradient;
       }
   }

4.2 User-Cases

4.2.1 groovyWaveTank

An example of using groovyBC to generate 2nd-order Stokes waves. It is designed to run with interFoam. For a wave with amplitude = 0.1, and wavelength = 5, gamma is set on the inlet:

   inlet
   {
       type            groovyBC;
       valueExpression "(pos().z<=A*cos(-w*time())+0.5*k*A*A*cos(2*(-w*time()))) ? 1 : 0";
       variables       "l=5;A=0.1;g=vector(0,0,-9.81);k=2*pi/l;w=sqrt(k*mag(g));";
       timelines       ();
   }

The corresponding velocity field U is set using:

   inlet
   {
       type            groovyBC;
       valueExpression "(pos().z<=A*cos(-w*time())+0.5*k*A*A*cos(2*(-w*time()))) ? vector( A*w*exp(k*pos().z)*cos(-w*time()), 0, A*w*exp(k*pos().z)*sin(-w*time())) : vector(0,0,0)";
       variables       "l=5;A=0.1;g=vector(0,0,-9.81);k=2*pi/l;w=sqrt(k*mag(g));";
       timelines       ();
   }

A test case can be downloaded from: groovyWaveTank.tgz. As usual, run blockMesh, setFields, and interFoam.

4.2.2 Swirl on inlet

Here's an example of how you can specify a swirl on the inlet. You only need to specify 2 parameters, Un and rpm, where

Un is the normal velocity and rpm is the rotational speed.

   inlet
   {
       type            groovyBC;
       variables       "rpm=8000.0;Un=55.8;c=sum(pos()*mag(Sf()))/sum(mag(Sf()));n=sum(normal())/mag(sum(normal()));p=pos()-c;r=mag(p)+1.0e-10;R=max(r);xt=vector(n.y,-n.x,0);xT=xt/mag(xt);yt=vector(-n.x*n.z,-n.y*n.z,n.x*n.x+n.y*n.y);yT=yt/mag(yt);";
       valueExpression "-Un*normal() + (rpm*pi/30)*((p & yT)*xT - (p & xT)*yT)";
   }

c is the centre of the patch, n is the averaged patch normal. You can easily extend it to have a variation on the patch normal component, for instance a parabolic inlet would only require the addition of *(1-pow(r/R,2)) to Un*normal(). The only time this will not work is when n is equal to (0,0,1), then you need another transformation, for instance this one (which will not work when swirl is around (1,0,0))

   inlet
   {
       type            groovyBC;
       variables       "rpm=8000.0;Un=55.8;c=sum(pos()*mag(Sf()))/sum(mag(Sf()));n=sum(normal())/mag(sum(normal()));p=pos()-c;r=mag(p)+1.0e-10;R=max(r);xt=vector(0,n.z,-n.y);xT=xt/mag(xt);yt=vector(n.y*n.y+n.z*n.z,-n.x*n.y,-n.x*n.z);yT=yt/mag(yt);";
       valueExpression "-Un*normal() + (rpm*pi/30)*((p & yT)*xT - (p & xT)*yT)";
   }

4.2.3 Simple 2-Way coupling of patches

This is an example of using groovyBC for creating a simple 2-Way coupling between two patches, such that it acts as more or less, one continuously connected mesh during the simulation. As the title specifies, it is a very limited and simplistic approach to the complex concept of full 2-Way patch to patch coupling, but is sufficient in a lot of situations.

Some applications where this approach would make sense (feel free to add or detract from this list):

• Cases where only the mean values of the variables at the respective patches are of importance, and not the actual spatial distribution (unless both patches have exactly the same number of faces of course).
• Cases where physically connecting up the two patches would add on an unnecessarily large number of cells resulting in a bloated up mesh.
• Situations where field values need to be modified by a simple mathematical model of some device connected between the two patches which does not necessarily need to be simulated using CFD.
• A way of showing off some of the cool features of groovyBC ;-)
• and so on...

The example uses the simpleFoam solver, and consists of three completely unconnected mesh regions say, A, B and C. Each mesh region contains two patches as listed below:

   Region A: inlet ; interface11
   Region B: interface12 ; interface21
   Region C: interface22 ; outlet

The pressure (p) and velocity (U) fields of the following sets of patches are coupled to each other using groovyBC, with the aim of emulating a continuous mesh without any breaks:

   interface11 :: interface12
   interface21 :: interface22

The system uses pressure driven flow, and the pressure at the patch inlet is ramped up from 10 bar to 35 bar over 50 iterations (the ramp is used for a softer start), and held constant thereafter. The pressure at the patch outlet is held constant at 10 bar. The case uses the k-epsilon turbulence model.

Here is an excerpt of the pressure boundary conditions:

   interface11
   {
       type            groovyBC;
       variables       "p_int12@interface12=sum(p*mag(Sf()))/sum(mag(Sf()));p_int11=sum(p*mag(Sf()))/sum(mag(Sf()));p_relax=0.3;";
       valueExpression "(p_int11 + p_relax*(p_int12 - p_int11))";
       value           $internalField;
   }
   interface12
   {
       type            zeroGradient;
   }
   
   interface21
   {
       type            groovyBC;
       variables       "p_int22@interface22=sum(p*mag(Sf()))/sum(mag(Sf()));p_int21=sum(p*mag(Sf()))/sum(mag(Sf()));p_relax=0.3;";
       valueExpression "(p_int21 + p_relax*(p_int22 - p_int21))";
       value           $internalField;
   }
   interface22
   {
       type            zeroGradient;
   }

At patch interface11, the area weighted average of the pressure field on patch interface12 is used. However, instead of directly imposing the value of the pressure, it is introduced through a relaxation factor which was required to ensure stability of the simulation. A similar set up can be seen also for the patch interface21.

The corresponding velocity boundary conditions are:

   interface11
   {
       type            zeroGradient;
   }
   interface12
   {
       type            groovyBC;
       variables       "Q_int11@interface11=-1*sum(phi);Q_int12=sum(phi);U_relax=0.3;";
       valueExpression "(Q_int12 + U_relax*(Q_int11 - Q_int12))/sum(mag(Sf()))*normal()";
       fractionExpression "(phi > 0) ? 0 : 1"
       value           $internalField;
   }
   
   interface21
   {
       type            zeroGradient;
   }
   interface22
   {
       type            groovyBC;
       variables       "Q_int21@interface21=-1*sum(phi);Q_int22=sum(phi);U_relax=0.3;";
       valueExpression "(Q_int22 + U_relax*(Q_int21 - Q_int22))/sum(mag(Sf()))*normal()";
       fractionExpression "(phi > 0) ? 0 : 1"
       value           $internalField;
   }

Here, the coupling is done in the opposite direction, with the total flow at the patch interface11 being imposed on the patch interface12 again, through a relaxation factor in order to ensure stability as the simulation evolves. A similar set up can be observed also for the patch interface22. In addition, it should be noted, that in this particular case, an inletOutlet boundary condition has been emulated.

The complete simulation case can be downloaded from: groovyBC_2Way_coupling_001.tgz (let me know if the link does not work) This case does not use blockMesh, and the mesh has already been imported and set up. Either use the Allrun script, or simply run simpleFoam from within the case root.

As always, feedback (good, bad and ugly) is always welcome :-)

5 Technical

• the BC only works for fvPatchFields but can "read" point and face-fields too
• it also works for pointPatchFields. Please note that constants must be enclosed in toPoint, otherwise they will be defined as face-values

5.1 Other

The parser-part can be used in other boundary conditions (that are not covered by the mixed-BC).

x, X, y, Y, z, Z can NOT be used as variables, since they are symbols for the parser.

5.2 Known Bugs

This utility is still under development. Currently there are no known known bugs (but propably some unknown ones)

5.3 Bug reporting

Bugs can be reported with the this bug-tracker

6 Download

Valid versions: - For the very latest version, see swak4Foam.

Valid versions: - The latest independent sources can always be downloaded via Subversion:

svn checkout svn://svn.code.sf.net/p/openfoam-extend/svn/trunk/Breeder_1.6/libraries/groovyBC

Valid versions: - Older sources are also available via Subversion:

svn checkout svn://svn.code.sf.net/p/openfoam-extend/svn/trunk/Breeder_1.5/libraries/groovyBC

7 Plans

• allow expressions for tensors
• allow "permanent" variables
• implement the BC for face-fields

8 History

• 2009-02-19: Initial upload

--Bgschaid 23:50, 19 February 2009 (CET)