CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > OpenFOAM Bugs

Possible bug in meshPhiC please take a look

Register Blogs Members List Search Today's Posts Mark Forums Read

Reply
 
LinkBack Thread Tools Display Modes
Old   December 20, 2007, 09:30
Default Dear Developers Wishing yo
  #1
Senior Member
 
Join Date: Mar 2009
Posts: 248
Rep Power: 9
jaswi is on a distinguished road
Dear Developers

Wishing you all a very nice day.

I request you to take a look at this post and recheck whether the meshPhi.C has an error. I am using interFoam with mesh motion and thus need to use meshPhi(rho,U) for making fluxes absolute and relative when required. I followed the trail and below is what I found:


The file /fvc/fvcMeshPhi.H has an overloaded function meshPhi() alongwith other functions.
-----------
meshPhi(U):
-----------

tmp<surfacescalarfield> meshPhi
(
const volVectorField& U
);

---------------
meshPhi(rho,U): when rho - dimensionedScalar
---------------


tmp<surfacescalarfield> meshPhi
(
const dimensionedScalar& rho,
const volVectorField& U
);

---------------
meshPhi(rho,U): when rho - volScalarField
---------------

tmp<surfacescalarfield> meshPhi
(
const volScalarField& rho,
const volVectorField& U
);

Now the corresponding implementation of these functions in /fvc/meshPhi.C are


-----------
meshPhi(U):
-----------
tmp<surfacescalarfield> meshPhi
(
const volVectorField& vf
)
{
Info << "Calling 1st constructor" << endl;

return fv::ddtScheme<vector>::New
(
vf.mesh(),
vf.mesh().ddtScheme("ddt(" + vf.name() + ')')
)().meshPhi(vf);
}


---------------
meshPhi(rho,U): when rho - dimensionedScalar
---------------
tmp<surfacescalarfield> meshPhi
(
const dimensionedScalar& rho,
const volVectorField& vf
)
{
Info << "Calling 2nd constructor" << endl;
return fv::ddtScheme<vector>::New
(
vf.mesh(),
vf.mesh().ddtScheme("ddt(" + rho.name() + ',' + vf.name() + ')')
)().meshPhi(vf);
}




---------------
meshPhi(rho,U): when rho - volScalarField
---------------

tmp<surfacescalarfield> meshPhi
(
const volScalarField& rho,
const volVectorField& vf
)
{
Info << "Calling 3rd constructor" << endl;
return fv::ddtScheme<vector>::New
(
vf.mesh(),
vf.mesh().ddtScheme("ddt(" + rho.name() + ',' + vf.name() + ')')
)().meshPhi(vf);
}


The point is that there is no difference whether I call meshPhi(U) or meshPhi(rho,U).

This issue has been raised before. please take a look at this post:

http://www.cfd-online.com/OpenFOAM_D...es/1/3548.html

In this post by Fredrick, he at least gets an error :
-------------------------------------------------
phi += fvc::meshPhi(U).
For a compressible flow, the fluxes should be:
phi += fvc::meshPhi(rho,U)

This does not work, I got the following error:
--> FOAM FATAL ERROR : Different dimensions for +=
dimensions : [1 0 -1 0 0 0 0] = [0 3 -1 0 0 0 0]
-------------------------------------------------

and Prof. Hrvoje Jasak and Eric Lillberg had a discussion over which is the best method to correct it.


I took the clue from the post and multiplied the return by rho as suggested in that post:
*********************************************
template<class>
tmp<surfacescalarfield> meshPhi
(
const dimensionedScalar& rho,
const GeometricField<type,>& vf
)
{
//HJ, Missing rho. 13/Dec/2006
return rho*fv::ddtScheme<type>::New
(
vf.mesh(),
vf.mesh().ddtScheme("ddt(" + rho.name() + ',' + vf.name() + ')')
)().meshPhi(vf);

}
**********************************************

but attempt to recompile the library gives the following error

openfoam@taifun:~/OpenFOAM/OpenFOAM-1.4.1/src/finiteVolume> wmake libso
SOURCE=finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtSchemes.C ; g++ -m64 -Dlinux64 -DDP -Wall -Wno-strict-aliasing -Wextra -Wno-unused-parameter -Wold-style-cast -march=opteron -O3 -DNoRepository -ftemplate-depth-40 -I/home/openfoam/OpenFOAM/OpenFOAM-1.4.1/src/triSurface/lnInclude -I/home/openfoam/OpenFOAM/OpenFOAM-1.4.1/src/meshTools/lnInclude -IlnInclude -I. -I/home/openfoam/OpenFOAM/OpenFOAM-1.4.1/src/OpenFOAM/lnInclude -fPIC -c $SOURCE -o Make/linux64GccDPOpt/EulerDdtSchemes.o
SOURCE=finiteVolume/fvc/fvcMeshPhi.C ; g++ -m64 -Dlinux64 -DDP -Wall -Wno-strict-aliasing -Wextra -Wno-unused-parameter -Wold-style-cast -march=opteron -O3 -DNoRepository -ftemplate-depth-40 -I/home/openfoam/OpenFOAM/OpenFOAM-1.4.1/src/triSurface/lnInclude -I/home/openfoam/OpenFOAM/OpenFOAM-1.4.1/src/meshTools/lnInclude -IlnInclude -I. -I/home/openfoam/OpenFOAM/OpenFOAM-1.4.1/src/OpenFOAM/lnInclude -fPIC -c $SOURCE -o Make/linux64GccDPOpt/fvcMeshPhi.o
finiteVolume/fvc/fvcMeshPhi.C: In function 'Foam::tmp<foam::geometricfield<double,> > Foam::fvc::meshPhi(const Foam::volScalarField&, const Foam::volVectorField&)':
finiteVolume/fvc/fvcMeshPhi.C:84: error: conversion from 'Foam::tmp<foam::field<double> >' to non-scalar type 'Foam::tmp<foam::geometricfield<double,> >' requested
make: *** [Make/linux64GccDPOpt/fvcMeshPhi.o] Error 1


Please take a look and comment whether the bug is still there or it has been corrected and I am making a wrong assumption.

Thanks alot for your attention.
with best regards

jaswinder
jaswi is offline   Reply With Quote

Old   December 20, 2007, 09:49
Default Yes it has an error. I have c
  #2
Senior Member
 
Join Date: Mar 2009
Posts: 854
Rep Power: 13
henry is on a distinguished road
Yes it has an error. I have corrected it in our development version and here is the fix for 1.4.1: fvcMeshPhi.C which I haven't tested but it should be OK, let me know it is works.

Henry
henry is offline   Reply With Quote

Old   December 20, 2007, 10:12
Default In all honesty, this is not qu
  #3
Senior Member
 
Hrvoje Jasak
Join Date: Mar 2009
Location: London, England
Posts: 1,758
Rep Power: 21
hjasak will become famous soon enough
In all honesty, this is not quite as clear-cut as it seems. Firstly, I need the rho as the name of the argument in order to capture the correct ddt scheme. However, what do you expect meshPhi() to return: as simple volume swept by the face in motion consistent with the ddt term, allowing you to manipulate it afterwards.

The conclusion from the previous round was that meshPhi() should return the mesh flux - if you don't like it, just multiply by phi interpolated in the way you want.

So much from my side,

Hrv
__________________
Hrvoje Jasak
Providing commercial FOAM/OpenFOAM and CFD Consulting: http://wikki.co.uk
hjasak is offline   Reply With Quote

Old   December 20, 2007, 12:03
Default Thanks a lot "Professors" for
  #4
Senior Member
 
Join Date: Mar 2009
Posts: 248
Rep Power: 9
jaswi is on a distinguished road
Thanks a lot "Professors" for a quick response. :-). I made the changes and compiled the library. it now gives this error.


--> FOAM FATAL ERROR : Different dimensions for -=
dimensions : [0 3 -1 0 0 0 0] = [1 0 -1 0 0 0 0]


From function dimensionSet::operator-=(const dimensionSet& ds) const
in file dimensionSet/dimensionSet.C at line 183.

when i do the following

return fvc::interpolate(rho)*fv::ddtScheme<vector>::New
(
vf.mesh(),
vf.mesh().ddtScheme("ddt(" + rho.name() + ',' + vf.name() + ')')
)().meshPhi(vf);

then flux[m^3/s] is equated to [(kg/m^3)*(m^3/sec)]
and this dimensionally not correct. The interesting point is that the post i referred to mentions this same error before the correction was made.

One more idea is that as I am using mesh motion with interfoam . In interfoam's createField we declare and define a surface scalar field "rhoPhi" which has the same dimensionality what the modified meshPhi(rho,U) returns. But when we correct the fluxes we use statements like
phi+=meshPhi(rho,U) and phi-=meshPhi(rho,U). May be there is something wrong with that approach. I have very less experience with interfoam solver thus cannot say it with confidence.

I have been loosing my sleep over it . the loss of liquid phase volume fraction due to the mesh motion has really boggled me.

Please comment what shall be done to correct this problem.

With Best Regards
Jaswinder
jaswi is offline   Reply With Quote

Old   December 20, 2007, 12:16
Default Did you recompile interFoam?
  #5
Senior Member
 
Join Date: Mar 2009
Posts: 854
Rep Power: 13
henry is on a distinguished road
Did you recompile interFoam?
henry is offline   Reply With Quote

Old   December 20, 2007, 12:21
Default Thanks again for the discussio
  #6
Senior Member
 
Join Date: Mar 2009
Posts: 248
Rep Power: 9
jaswi is on a distinguished road
Thanks again for the discussion.


Yes I did that. To make sure that it does I

wclean it and the wmake it. But the error is still there.
jaswi is offline   Reply With Quote

Old   December 20, 2007, 12:25
Default Scrap that, there is no need t
  #7
Senior Member
 
Join Date: Mar 2009
Posts: 854
Rep Power: 13
henry is on a distinguished road
Scrap that, there is no need to recompile interFoam.

I removed mesh motion from the standard interFoam in the 1.4.1 release because it wasn't entirely consistent and had various modes of failure. I have written a completely new interDyMFoam for the 1.5 release which includes all kinds of mesh change both geometric and topological, an early version of which created the animation on the front page of our web-site. We hope to release 1.5 in the 1st quarter 2008.

Henry
henry is offline   Reply With Quote

Old   December 20, 2007, 12:26
Default Good Evening Henry I looked
  #8
Senior Member
 
Join Date: Mar 2009
Posts: 248
Rep Power: 9
jaswi is on a distinguished road
Good Evening Henry

I looked through the interFoam again. Maybe we should correct both phi and rhoPhi for the moving mesh.

The correction block can have something like:


//Make fluxes absolute
phi+=meshphi(U)
rhoPhi+=meshPhi(rho,U)


//Make fluxes relative
phi-=meshphi(U)
rhoPhi-=meshPhi(rho,U)

where and when required in the mesh motion procedure. It is completely a arrow in the dark and I might ne completely wrong.

Please Comment

With Best Regards
Jaswinder
jaswi is offline   Reply With Quote

Old   December 20, 2007, 12:30
Default Hi Henry Can I request for
  #9
Senior Member
 
Join Date: Mar 2009
Posts: 248
Rep Power: 9
jaswi is on a distinguished road
Hi Henry

Can I request for a Christmas Gift :-)

Please If i can have that solver now because I have been loosing sleep since a week over this issue and you know how the advisors some time react when there are no results and one can't explain what is wrong.

my email is singhj@in.tum.de

Please let me know if I can have it.

With Best Regards
Jaswinder
jaswi is offline   Reply With Quote

Old   December 20, 2007, 12:43
Default The new solver relies on a man
  #10
Senior Member
 
Join Date: Mar 2009
Posts: 854
Rep Power: 13
henry is on a distinguished road
The new solver relies on a many new developments and will not compile with 1.4.1, sorry. We will get 1.5 out as soon as we can.

Henry
henry is offline   Reply With Quote

Old   December 20, 2007, 13:03
Default Thanks alot for considering it
  #11
Senior Member
 
Join Date: Mar 2009
Posts: 248
Rep Power: 9
jaswi is on a distinguished road
Thanks alot for considering it.

It seems I have to explore other options then, as I am running out of time.

By the way, is there some mass transfer modeling coming up as well. I just can't wait to have the new version. The more i understand this stuff the more I am surprised what all it can do.

Thank you for the great work.

Wish you a Merry Christmas and a Happy New Year.
Happy Holidays
Jaswinder
jaswi is offline   Reply With Quote

Reply

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are On
Pingbacks are On
Refbacks are On



All times are GMT -4. The time now is 03:52.