CFD Online Discussion Forums

CFD Online Discussion Forums (
-   OpenFOAM Running, Solving & CFD (
-   -   Modifying motionU in the development version (

philippose February 1, 2007 07:48

Hello, I was creating a sol
I was creating a solver in which I need to modify the contents of the motionU file (change the velocity on the patches based on calculations).

Using the development version I downloaded from Frank Bos' website, I get an error while trying to compile the code at the point where I do the following:

motionU.boundaryField()[patchID] == vector(0.0,0.0,calcVelocity)

The error states that a vector<double> assignment is not possible with a tetFemMatrix.

I tried compiling the same code with the standard OpenFOAM 1.3 version, and it worked without any problems. But in the standard version, there seems to be a bug with the way moving meshes are implemented.... always resulting in fluid velocities much lower than expected. This bug was not present in the development version.

In case a direct dump of the compiler output is required, I can send one as soon as I get back home this evening.

Any idea what might be wrong?

Have a nice day!


hjasak February 1, 2007 08:00

Highly dubious: how did the te
Highly dubious: how did the tetFemMatrix get into this story? motionU is a tetFemField.

In any case, I do this in the fluid-structure interaction solver all the time so it should/does/will work. In case of trouble, you can refCast the patch before assignment:

fixedValueTetPolyPatchVectorField& motionUFuidPatch =

motionUFuidPatch ==


philippose February 1, 2007 13:56

Hello Hrv, Good Evening!
Hello Hrv,
Good Evening!

Finally got back home... so.. for a little more detail... I created the motionU field using the following code which is called outside the runtime loop (I included tetPointFields.H at the beginning of the solver code):

tetPointVectorField& motionU = const_cast<tetpointvectorfield&>
mesh.objectRegistry::lookupObject<tetpointvectorfi eld>

After this definition, I do not do any other modifications to this field till the following lines:

label patchIndex = mesh.boundaryMesh().findPatchID(forcePatches[fPI]);
motionU.boundaryField()[patchIndex] == vector(0.0,0.0,calcVelocity);

In the above code snip, the variable "calcVelocity " is of type "scalar".

My $WM_DECOMP_INC is -DFACE_DECOMP and $WM_DECOMP_LIBS is -lfaceDecompFiniteElement.

As a reminder... I dont use the original OpenFOAM 1.3 version, but the precompiled development version 11.09.2006. When I try compiling the code with "wmake", I get the following output:

Making dependency list for source file turbForceFoam.C

SOURCE=turbForceFoam.C ; g++ -m32 -Dlinux -DDP -Wall -W -Wno-unused-parameter -Wold-style-cast -O3 -DNoRepository -ftemplate-depth-30 -I/home/philippose/OpenFOAM/OpenFOAM-1.3/src/dynamicFvMesh/lnInclude -I/home/philippose/OpenFOAM/OpenFOAM-1.3/src/dynamicMesh/lnInclude -I/home/philippose/OpenFOAM/OpenFOAM-1.3/src/finiteVolume/lnInclude -I/home/philippose/OpenFOAM/OpenFOAM-1.3/src/meshTools/lnInclude -I/home/philippose/OpenFOAM/OpenFOAM-1.3/src/turbulenceModels -I/home/philippose/OpenFOAM/OpenFOAM-1.3/src/transportModels -I/home/philippose/OpenFOAM/OpenFOAM-1.3/src/tetDecompositionFiniteElement/lnInc lude -DFACE_DECOMP -I/home/philippose/OpenFOAM/OpenFOAM-1.3/src/OpenFOAM/lnInclude -IlnInclude -I. -fPIC -pthread -c $SOURCE -o Make/linuxGcc4DPOpt/turbForceFoam.o
/home/philippose/OpenFOAM/OpenFOAM-1.3/src/OpenFOAM/lnInclude/PointPatchField.H: In instantiation of 'Foam::PointPatchField<foam::tetpolypatchfield,> >':
/home/philippose/OpenFOAM/OpenFOAM-1.3/src/tetDecompositionFiniteElement/lnInclu de/tetPolyPatchField.H:62: instantiated from 'Foam::tetPolyPatchField<foam::vector<double> >'
solveForceBalance.H:56: instantiated from here
/home/philippose/OpenFOAM/OpenFOAM-1.3/src/OpenFOAM/lnInclude/PointPatchField.H: 368: error: invalid use of undefined type 'struct Foam::tetFemMatrix<foam::vector<double> >'
/home/philippose/OpenFOAM/OpenFOAM-1.3/src/tetDecompositionFiniteElement/lnInclu de/tetPolyPatchField.H:52: error: declaration of 'struct Foam::tetFemMatrix<foam::vector<double> >'
/home/philippose/OpenFOAM/OpenFOAM-1.3/src/finiteVolume/lnInclude/readTimeContro ls.H: In function 'int main(int, char**)':
/home/philippose/OpenFOAM/OpenFOAM-1.3/src/finiteVolume/lnInclude/readTimeContro ls.H:40: warning: unused variable 'maxCo'
make: *** [Make/linuxGcc4DPOpt/turbForceFoam.o] Error 1

The same code copied and pasted into an original OpenFOAM 1.3 system running in 64bit mode compiled without any errors.

Hope we can sort out the problem... I am personally quite excited about this solver... :-)!

I thank you for your help, and have a nice day!


philippose February 1, 2007 14:51

Hello again, I just tried us
Hello again,
I just tried using the reference casting method you suggested.... I get the same errors during compilation....

"invalid use of undefined type struct Foam::tetFemMatrix<foam::vector<double> >"


hjasak February 1, 2007 14:52

#include "tetFemMatrix.H" i
#include "tetFemMatrix.H"

in the offending file.


philippose February 1, 2007 15:05

Hello Hrv, Thank you very m
Hello Hrv,

Thank you very much for the really quick response (as usual :-)!)...

I included "tetFemMatrix", and the solver compiled without any more errors.

Since I keep fiddling around with new solvers, it would be interesting to know why the compilation worked in the original version of OpenFOAM 1.3 while it required this additional include directive in the development version... could it be that you are calling "tetFemMatrix" in one of the other include files which was later removed?

Now I can go ahead and test the solver and continue with the next phase :-)! I am basically trying to solve the rigid body motion equations along with the turbFoam code in order to move the spool of a valve based on flow forces.... And.. I am trying to keep it as general as possible with dictionary based parameter selection, etc...

Thanks once again... and have a nice day!


hjasak February 1, 2007 15:06

Not sure if you've seen my mes
Not sure if you've seen my message: edit turbForceFoam.C and before main() add

#include "tetFemMatrix.H"

You may need to add some include paths in Make/options, but I doubt it - it should just work.


hjasak February 1, 2007 15:36

Do you really want to know? :-
Do you really want to know? :-)

OK, let's see. On the kind of meshes we are working with, you can specify the boundary conditions:

- on the surface, meaning "per face"
- on exposed (boundary) vertices, meaning "per vertex"

The first is typical for the FVM and surface forces in FEM (which unfortunately get split up into nodal values, while the second is a FEM0 kind.

So, I first wrote the FEM solvers and implemented the boundary conditions. One of the job of a boundary condition is to enforce the constraints (FEM) or diagonal/source contribution on the discretisation matrix. Therefore, a boundary condition needs to know about the matrix right?

I then took on the point interpolation. This formally also has boundary conditions (think: symmetry plane) and needs to make sure boundary values are "interpolated" from face values only rather than extrapolated from the cells. Clearly, the functionality and behaviour of b.c.-s is exactly identical to my FEM solver and we should re use the code. Thus, what used to be a FEM boundary condition is nicely templated and generalised and used for point interpolation.

Clear so far?

However, the FEM solvers in OpenFOAM have got a matrix, but point interpolation is just that: an interpolation. Therefore, patch constraints need to be enforced but I don't care about the matrix functionality. Because the code is templated, the stuff never got instantiated and we lived with it happily for a number of years.

Then came gcc-4.1.0 (or was it 4.0.2?) and we got found out: how can you have matrix dependent behaviour in the boundary condition when in point interpolation there's no matrix?

The simplest solution is to give up, copy the whole lot and have two kinds of b.c.-s. However:
- this is cheating: OpenFOAM is the ONLY code I know that actually shares discretisation, matrices etc between the FVM and FEM solver (as it should)
- this is below my I.Q. grade.

Therefore, the FVM and FEM matrix gets a typedef for a kind of boundary condition interaction it needs, point interpolation gets a dummy matrix with typedef only (do nothing) and the generalised templated code is still there, waiting for me to implement full FEM. Good, huh?

Finally, the issue of the include: the boundary condition in FEM is now a straight typedef, but a matrix depends on the b.c. and vice versa. This makes for elliptic dependency and I cannot include one header into another. Remember, this is ALL in templated code.

Therefore, in order for the compiler to instantiate the template (which is what the compiler will do from the top level), you need to let it know which b.c. you wish to use (you have inherited this one) AND what is the associated matrix: hence the additional include.

I could sort all this out by explicitly instantiating all the necessary templates with some tricks I did recently, but let's save this for my next trick. :-) In any case, only true C++ experts would be able to follow such a dodge...

The point of OpenFOAM is to encourage you to do some brain-work, and there's always an expert to help.



P.S. I strongly suspect you may need to read this a few times to pick it all up
P.P.S. Copyright, Hrvoje Jasak

philippose February 1, 2007 16:25

Oops... One learns the hard wa
Oops... One learns the hard way what questions to ask without thinking twice... what questions to ask in a cautious manner and what questions to bury as soon as it pops into ones head :-)!!

I guess I just happened to stumble on a question of type 3 :-)!

The thing is, though I do C and some C++ programming, the very foundations of C++... the concepts of inheritance, polymorphism, virtualisation and the other ghosts which lurk around, have always managed to escape being captured by me due to a lack of the need to understand them (possibly.. till now :-)!)

I normally stick to C, and some simple C++. Since my field of work hasnt required me to get into hardcore software development, I havent really sat down with C++ in a structured manner yet :-)! (I can see you frown :-)!) Who knows... maybe OpenFOAM will provide the driving force :-)!

As you rightly mentioned.... I have already read that two times... and am kind of still scratching my head :-)! Shall give it a couple more shots though... not used to giving up easily :-)! Was up quite late last night trying to figure out what was going on with that code which refused to compile!!

And about OpenFOAM encouraging me to do some brain work... thats an understatement :-)! I am sure you will bloat with pride if I tell you that OpenFOAM has helped me understand hydraulics to a level not easily achievable in such a short time frame :-)! The immense possibilities one has when the source code can be modified is unbelievable... and the amount of insight one gains about the physical system while setting up... running and interpreting the simulations is just amazing :-)!

I thank you, the rest of the team and the participants of the forum for all the help and support :-)!



All times are GMT -4. The time now is 13:59.