CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Software User Forums > OpenFOAM > OpenFOAM Programming & Development

if (magSqr(fvc::grad(T)) != 0)

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

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   May 22, 2014, 05:34
Default if (magSqr(fvc::grad(T)) != 0)
  #1
New Member
 
Gennaro
Join Date: May 2014
Posts: 23
Rep Power: 10
Gennaro is on a distinguished road
Dear all,

I need to solve the following equation in OpenFOAM as a part of an algorithm:

Keff = (Foam:ow(fInv/magSqr(fvc::grad(T)),scalar(1.)/scalar(3.)))/(rho*Cp);

The problem is that when grad(T) = (0 0 0), e.g on the first iteration, the ratio becomes infinite and the solver understandably crashes.

Therefore I tried to to something like:

if (magSqr(fvc::grad(T)) != 0)
{
Keff = (Foam:ow(fInv/magSqr(fvc::grad(T)),scalar(1.)/scalar(3.)))/(rho*Cp);
}
else
{
Keff = (Foam:ow(fInv/(dimensionedScalar("small", dimensionSet(0,-2,0,2,0), 1e-4)),scalar(1.)/scalar(3.)))/(rho*Cp);
}

Unfortunately this won't compile returning the error

error: no match for ‘operator==’ in ‘Foam::magSqr(const Foam::tmp<Foam::GeometricField<Type, PatchField, GeoMesh> >&) [with Type = Foam::Vector<double>, PatchField = Foam::fvPatchField, GeoMesh = Foam::volMesh]() == 0’

Please find the whole solver in attachment and a case to test it. Please note that I'm not a programmer and I'm just starting to learn C++. Forgive me if the answer is too obvious.

Thanks and best regards

Gennaro
Attached Files
File Type: gz laplacianFoamCustom.tar.gz (85.2 KB, 11 views)
Gennaro is offline   Reply With Quote

Old   May 22, 2014, 07:25
Default
  #2
Senior Member
 
Tomislav Maric
Join Date: Mar 2009
Location: Darmstadt, Germany
Posts: 284
Blog Entries: 5
Rep Power: 19
tomislav_maric is on a distinguished road
Quote:
Code:
error: no match for ‘operator==’ in ‘Foam::magSqr(const  Foam::tmp<Foam::GeometricField<Type, PatchField, GeoMesh>  >&) [with Type = Foam::Vector<double>, PatchField =  Foam::fvPatchField, GeoMesh = Foam::volMesh]() == 0’
Even though you're not a programmer, you'll have to be come one to some extent, if you intend to extend OpenFOAM(-extend). Pun intended. The first step is reading compiler errors.

That template error says that there is no 'operator==' in GeometricField that takes '0' scalar value as the second argument.

If you want to stabilize your gradient, don't use if-then-else, you can stabilize it using SMALL (1e-15 for a double scalar). For example, check out

Code:
$WM_PROJECT_DIR/src/finiteVolume/cfdTools/general/include/setDeltaT.H
and take a look at the line:

Code:
    scalar maxDeltaTFact = maxCo/(CoNum + SMALL);
To prevent division by zero that results with a floating point exception, the denominator is summed with SMALL. Try the same for the magnitude of your gradient and loose the conditional statements, this should help.
__________________
When asking a question, prepare a SSCCE.
tomislav_maric is offline   Reply With Quote

Old   May 22, 2014, 08:53
Default
  #3
New Member
 
Gennaro
Join Date: May 2014
Posts: 23
Rep Power: 10
Gennaro is on a distinguished road
Thanks Tomi,

Actually both != and == were attempts returning the same error message: no match for ‘operator==’or no match for ‘operator!='. So, the operator itself was not the problem.

Thanks a lot for your tip, it's very interesting. However, while SMALL is non-dimensional, I need a SMALL scalar with the right dimensions. Indeed, while with SMALL the solver would still compile, then when I run the case I get:

--> FOAM FATAL ERROR:
LHS and RHS of + have different dimensions
dimensions : [0 -2 0 2 0 0 0] + [0 0 0 0 0 0 0]

Any tips on how to apply dimensions to SMALL?

Thanks

Gennaro
Gennaro is offline   Reply With Quote

Old   May 22, 2014, 09:48
Default
  #4
Senior Member
 
Tomislav Maric
Join Date: Mar 2009
Location: Darmstadt, Germany
Posts: 284
Blog Entries: 5
Rep Power: 19
tomislav_maric is on a distinguished road
Quote:
Actually both != and == were attempts returning the same error message: no match for ‘operator==’or no match for ‘operator!='. So, the operator itself was not the problem.
Nope, this means that the GeometricField does not implement any of the two requested operators at the call sites, which results in an undefined reference error.

To give a scalar a dimension, include "dimensionedType.H" in the library code. Do not insert it into 'Make/files' for compilation - it's a class template.

Then you can instantiate a dimensionedScalar:

Code:
dimensionedScalar smallGrad ("small", gradT.dimensions(), SMALL);
Also, don't call 'fvc::grad' all over the place, it will execute the same gradient calculations over and over, unless you have set the flag for caching gradients, which is not set by default. Compute the gradient field and store it under gradT, then you can get its dimensions as shown above.
__________________
When asking a question, prepare a SSCCE.
tomislav_maric is offline   Reply With Quote

Old   May 23, 2014, 11:12
Default
  #5
New Member
 
Gennaro
Join Date: May 2014
Posts: 23
Rep Power: 10
Gennaro is on a distinguished road
Thanks Tomi,

your help is incredibly valuable. Congratulations!

I gave it a try, but it doesn't work.

With reference to the case in attachment, I added
Code:
#include "dimensionedType.H
in laplacianFoamCustom.C, trying both inside and outside int main.

Then I added
Code:
dimensionedScalar smallGrad ("small", gradT.dimensions(), SMALL);
in createFields.H.

Finally, I changed the formula in Eeqn.H to
Code:
Keff = (Foam::pow(fInv/magSqr(fvc::grad(T) + small),scalar(1.)/scalar(3.)))/(rho*Cp);
, trying both small and SMALL.

Can you please advise?

Thanks

Genn
Gennaro is offline   Reply With Quote

Old   May 26, 2014, 08:46
Default
  #6
New Member
 
Gennaro
Join Date: May 2014
Posts: 23
Rep Power: 10
Gennaro is on a distinguished road
Hi Tomi, any news?

The error is:

Code:
createFields.H:123:40: error: ‘gradT’ was not declared in this scope
meaning that
Code:
#include "dimensionedType.H"
is not enough?

I look forward to hearing from you

Thanks

Genn


Quote:
Originally Posted by tomislav_maric View Post
Nope, this means that the GeometricField does not implement any of the two requested operators at the call sites, which results in an undefined reference error.

To give a scalar a dimension, include "dimensionedType.H" in the library code. Do not insert it into 'Make/files' for compilation - it's a class template.

Then you can instantiate a dimensionedScalar:

Code:
dimensionedScalar smallGrad ("small", gradT.dimensions(), SMALL);
Also, don't call 'fvc::grad' all over the place, it will execute the same gradient calculations over and over, unless you have set the flag for caching gradients, which is not set by default. Compute the gradient field and store it under gradT, then you can get its dimensions as shown above.
Gennaro is offline   Reply With Quote

Old   February 4, 2016, 11:32
Default
  #7
Member
 
Avdeev Evgeniy
Join Date: Jan 2011
Location: Togliatty, Russia
Posts: 69
Blog Entries: 1
Rep Power: 19
j-avdeev will become famous soon enough
Send a message via Skype™ to j-avdeev
Quote:
Originally Posted by Gennaro View Post

Code:
createFields.H:123:40: error: ‘gradT’ was not declared in this scope

Hi, Genn.

I am not sure, that it optimal way to get gradT, but you can try
Code:
surfaceVectorField gradT = fvc::interpolate(fvc::grad(T));
surfaceScalarField sumGradT = gradT.component(vector::X)+gradT.component(vector::Y)+gradT.component(vector::Z);
j-avdeev is offline   Reply With Quote

Reply

Tags
algorithm, c++, grad(t), operator

Thread Tools Search this Thread
Search this Thread:

Advanced Search
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 Off
Pingbacks are On
Refbacks are On



All times are GMT -4. The time now is 09:37.