CFD Online Discussion Forums

CFD Online Discussion Forums (http://www.cfd-online.com/Forums/)
-   OpenFOAM Running, Solving & CFD (http://www.cfd-online.com/Forums/openfoam-solving/)
-   -   Updating Boundary Conditions Each Iteration (http://www.cfd-online.com/Forums/openfoam-solving/88762-updating-boundary-conditions-each-iteration.html)

thomas99 May 25, 2011 08:24

Updating Boundary Conditions Each Iteration
 
Hello All,

I have a situation where I would like to calculate some boundary conditions from the internal field. The problem is that, for each time step, the boundary conditions should be iteratively updated as part of the solution. I have searched the forum and noticed that many people have recommended the groovyBC or swak4Foam add ons for recalculating a boundary condition from the flow field but these tools seem to do this at the end of a time step and not as part of the iterative process. Can I use one of these tools to somehow calculate the BCs iteratively?

Alternatively, I noticed there is a BC built into openFoam called calculated that might work but I can't really find too much information about this in the forum and the tutorials don't really make it clear to me exactly how it works.

Does anyone have any idea how I can iteratively calculate boundary conditions?

Regards,
Tom

marupio May 25, 2011 09:02

I don't think groovyBC does iterative solutions. You may need to write your own boundary condition, and put the iterative solver into it. See

http://openfoamwiki.net/index.php/Ho...dary_condition

for an overview of writing new boundary conditions.

thomas99 May 25, 2011 09:31

Hi Marupio,

Thanks for the reply and the link.

I was hoping that I wouldn't have to do that but I'll give it a try.

Best,
Tom

bigphil May 25, 2011 09:35

Hi thomas99,

As marupio said you could write your own boundary condition,

alternatively, you could just include a header file which updates the boundary contact inside your solver solution loop,
for example I sometimes use the "fixedGradient" boundary condition in my solid stress solver to apply a traction, this means I must update this fixedGradient every solution iteration.
The header file might look something like the following (where U is the volVectorField I solve for):
Code:

label patchID = mesh.boundaryMesh().findPatchID("patch_of_interest_name");
if
(
  U.boundaryField()[patchID].type()
  == fixedGradientFvPatchVectorField::typeName
)
{
          fixedGradientFvPatchVectorField& Upatch =
            refCast<fixedGradientFvPatchVectorField>(U.boundaryField()[patchID]);
         
          Upatch.gradient() = (.......calculate something here.....);
}


Hope it helps,
Philip

gschaider May 25, 2011 10:03

Quote:

Originally Posted by thomas99 (Post 309182)
Hello All,

I have a situation where I would like to calculate some boundary conditions from the internal field. The problem is that, for each time step, the boundary conditions should be iteratively updated as part of the solution. I have searched the forum and noticed that many people have recommended the groovyBC or swak4Foam add ons for recalculating a boundary condition from the flow field but these tools seem to do this at the end of a time step and not as part of the iterative process. Can I use one of these tools to somehow calculate the BCs iteratively?

Depends on what you mean by iterative. groovyBC is recalculated every time other BCs are recalculated too (sol also in a loop inside the time-loop)

Quote:

Originally Posted by thomas99 (Post 309182)
Alternatively, I noticed there is a BC built into openFoam called calculated that might work but I can't really find too much information about this in the forum and the tutorials don't really make it clear to me exactly how it works.

Forget this one. It means "I was calculated once, I won't change again. And I don't know if I'm Dirichlet or Neuman so I will behave like a Dirichlet"

thomas99 May 26, 2011 04:20

Hi Philip and Bernhard,

Thanks for the information. I think the inclusion of a header file is a good idea and I'll try that first.

Regards,
Tom

mathslw August 29, 2012 12:21

Quote:

Originally Posted by thomas99 (Post 309300)
Hi Philip and Bernhard,

Thanks for the information. I think the inclusion of a header file is a good idea and I'll try that first.

Regards,
Tom

Hi Tom,

Do you succeed in updating the boundary in the loop?
I want to update the boundary condition after each iteration for steady state calculation.
Thanks!

Wei

marupio August 29, 2012 12:26

I'm not infront of my FOAM machine right now, but I seem to recall an issue where the boundary condition would update only once per timestep. I was doing multiple iterations per timestep, and I had to call some function to force it to update. I forget the details... but if you encounter the situation where it only seems to update once, then let me know, I'll dig for the answer.

-Dave

mathslw August 29, 2012 12:30

Quote:

Originally Posted by marupio (Post 379349)
I'm not infront of my FOAM machine right now, but I seem to recall an issue where the boundary condition would update only once per timestep. I was doing multiple iterations per timestep, and I had to call some function to force it to update. I forget the details... but if you encounter the situation where it only seems to update once, then let me know, I'll dig for the answer.

-Dave

Hi Dave,

I am using the steady state solver, so I want to update the valure after each iteration.
Do you know about that?

Thanks!

Wei

mathslw August 29, 2012 13:08

Quote:

Originally Posted by bigphil (Post 309195)
Hi thomas99,

As marupio said you could write your own boundary condition,

alternatively, you could just include a header file which updates the boundary contact inside your solver solution loop,
for example I sometimes use the "fixedGradient" boundary condition in my solid stress solver to apply a traction, this means I must update this fixedGradient every solution iteration.
The header file might look something like the following (where U is the volVectorField I solve for):
Code:

label patchID = mesh.boundaryMesh().findPatchID("patch_of_interest_name");
if
(
  U.boundaryField()[patchID].type()
  == fixedGradientFvPatchVectorField::typeName
)
{
          fixedGradientFvPatchVectorField& Upatch =
            refCast<fixedGradientFvPatchVectorField>(U.boundaryField()[patchID]);
         
          Upatch.gradient() = (.......calculate something here.....);
}


Hope it helps,
Philip

Hi Philip,

What should I do if I used the fixedValue and want to update it?
Thanks!

Wei

bigphil August 29, 2012 13:29

Quote:

Originally Posted by mathslw (Post 379357)
What should I do if I used the fixedValue and want to update it?
Thanks!

Hi Wei,

For a fixedValue boundary condition, you can do the following (U is the volVectorField which is solved):
Code:

// find ID of patch
label patchID = mesh.boundaryMesh().findPatchID("patch_of_interest");
// check patch has been found
if(patchID == -1)
{
  FatalError << "patch not found!" << exit(FatalError);
}

// set value on patch to what ever you want
U.boundaryField()[patchID] == vector(1,2,3);

Note the double equal signs, this is how you set the value on a patch. You don't need to cast the patch if it is fixedValue.

Hope it helps,
Philip

mathslw August 29, 2012 14:01

Quote:

Originally Posted by bigphil (Post 379363)
Hi Wei,

For a fixedValue boundary condition, you can do the following (U is the volVectorField which is solved):
Code:

// find ID of patch
label patchID = mesh.boundaryMesh().findPatchID("patch_of_interest");
// check patch has been found
if(patchID == -1)
{
  FatalError << "patch not found!" << exit(FatalError);
}
 
// set value on patch to what ever you want
U.boundaryField()[patchID] == vector(1,2,3);

Note the double equal signs, this is how you set the value on a patch. You don't need to cast the patch if it is fixedValue.

Hope it helps,
Philip

Hi Philip,

Your suggestion is very helpful!

If I use U.boundaryField()[patchID] == vector(1,2,3) , the boundary would be uniform.
Can I use a loop to set non-uniform boundary conditions?
Besides, to set the new boundary, I have to make use of the velocity on the boundary face, velocity on the boundary cell and velocity gradient on the boundary cell. Do you know how to access these values?

Many thanks!

Wei

bigphil August 30, 2012 08:50

Hi Wei,

You can set non-uniform boundary values based on anything you want:
Code:

// find ID of patch
label patchID = mesh.boundaryMesh().findPatchID("patch_of_interest");
// check patch has been found
if(patchID == -1)
{
  FatalError << "patch not found!" << exit(FatalError);
}
 
// set value on patch
forAll(U.boundaryField()[patchID], facei)
{
  vector faceOldVel = U.boundaryField()[patchID][facei];
  vector faceCellOldVel =
        U.internalField()[mesh.boundaryMesh()[patchID].faceCells()[facei]];
  vector faceSnGrad = U.boundaryField()[patchID].snGrad()[facei];
 
  // set whatever you want here
  U.boundaryField()[patchID][facei] == faceOldVel + faceCellOldVel + faceSnGrad;
}

Best regards,
Philip

mathslw August 30, 2012 12:01

Quote:

Originally Posted by bigphil (Post 379492)
Hi Wei,

You can set non-uniform boundary values based on anything you want:
Code:

// find ID of patch
label patchID = mesh.boundaryMesh().findPatchID("patch_of_interest");
// check patch has been found
if(patchID == -1)
{
  FatalError << "patch not found!" << exit(FatalError);
}
 
// set value on patch
forAll(U.boundaryField()[patchID], facei)
{
  vector faceOldVel = U.boundaryField()[patchID][facei];
  vector faceCellOldVel =
        U.internalField()[mesh.boundaryMesh()[patchID].faceCells()[facei]];
  vector faceSnGrad = U.boundaryField()[patchID].snGrad()[facei];
 
  // set whatever you want here
  U.boundaryField()[patchID][facei] == faceOldVel + faceCellOldVel + faceSnGrad;
}

Best regards,
Philip

Hi Philip,

Thanks for your reply!
I will try it in my program to see how it is works!

Best Regards!

Wei

mathslw August 30, 2012 12:13

Quote:

Originally Posted by bigphil (Post 379492)
Hi Wei,

You can set non-uniform boundary values based on anything you want:
Code:

// find ID of patch
label patchID = mesh.boundaryMesh().findPatchID("patch_of_interest");
// check patch has been found
if(patchID == -1)
{
  FatalError << "patch not found!" << exit(FatalError);
}
 
// set value on patch
forAll(U.boundaryField()[patchID], facei)
{
  vector faceOldVel = U.boundaryField()[patchID][facei];
  vector faceCellOldVel =
        U.internalField()[mesh.boundaryMesh()[patchID].faceCells()[facei]];
  vector faceSnGrad = U.boundaryField()[patchID].snGrad()[facei];
 
  // set whatever you want here
  U.boundaryField()[patchID][facei] == faceOldVel + faceCellOldVel + faceSnGrad;
}

Best regards,
Philip

Hi Philip,

As far as I know, the U.boundaryField()[patchID].snGrad()[facei] returns the patch-normal gradient. How to obtain the velocity gradient on the patch cell, which is a tensor. Can I define a
tensorField gradU = fvc::grad(U),
then use the
tensor faceCellGradientU =
gradU.internalField()[mesh.boundaryMesh()[patchID].faceCells()[facei]];
to obtain the velocity gradient on the patch cell?

Many thanks!

Wei

mathslw September 9, 2012 17:19

Quote:

Originally Posted by bigphil (Post 379492)
Hi Wei,

You can set non-uniform boundary values based on anything you want:
Code:

// find ID of patch
label patchID = mesh.boundaryMesh().findPatchID("patch_of_interest");
// check patch has been found
if(patchID == -1)
{
  FatalError << "patch not found!" << exit(FatalError);
}
 
// set value on patch
forAll(U.boundaryField()[patchID], facei)
{
  vector faceOldVel = U.boundaryField()[patchID][facei];
  vector faceCellOldVel =
        U.internalField()[mesh.boundaryMesh()[patchID].faceCells()[facei]];
  vector faceSnGrad = U.boundaryField()[patchID].snGrad()[facei];
 
  // set whatever you want here
  U.boundaryField()[patchID][facei] == faceOldVel + faceCellOldVel + faceSnGrad;
}

Best regards,
Philip

Hi Philip,

Sorry to bother u again.
If I use

U.boundaryField()[patchID] == vector(0.8, 0, 0);

I could update the boundary condition;
But if I use

forAll(U.boundaryField()[patchID],i)
{
U.boundaryField()[patchID][i] == vector(0.8, 0, 0);
}

The boundary would not be updated.

Do you know the reason?

Thanks!

Wei

bigphil September 10, 2012 10:03

Hi Wei,

To answer your first question, to calculate the full gradient (tensor) at the boundary face you use the fvc::grad operator as you have shown:
tensor faceCellGradientU =
gradU.internalField()[mesh.boundaryMesh()[patchID].faceCells()[facei]];

However, be careful because most fvc::grad schemes do not calculate the total gradient fully at boundary faces, they essentially calculate the normal gradient like snGrad and then they set the tangential gradient to be the same as the tangential gradient at the centre of the boundary cell.

As regards the forAll loop not updating the boundary I am not sure but here is a work-around:
Code:

vectorField newPatchValues(U.boundaryField()[patchID].size(), vector::zero);

forAll(U.boundaryField()[patchID],i)
{
    newPatchValues[i] = vector(0.8, 0, 0);
 }

U.boundaryField()[patchID][i] == newPatchValues;

Also, did you try chaining the "==" to "=" inside the forAll loop?

mathslw September 10, 2012 13:07

Quote:

Originally Posted by bigphil (Post 381068)
Hi Wei,

To answer your first question, to calculate the full gradient (tensor) at the boundary face you use the fvc::grad operator as you have shown:
tensor faceCellGradientU =
gradU.internalField()[mesh.boundaryMesh()[patchID].faceCells()[facei]];

However, be careful because most fvc::grad schemes do not calculate the total gradient fully at boundary faces, they essentially calculate the normal gradient like snGrad and then they set the tangential gradient to be the same as the tangential gradient at the centre of the boundary cell.

As regards the forAll loop not updating the boundary I am not sure but here is a work-around:
Code:

vectorField newPatchValues(U.boundaryField()[patchID].size(), vector::zero);
 
forAll(U.boundaryField()[patchID],i)
{
    newPatchValues[i] = vector(0.8, 0, 0);
 }
 
U.boundaryField()[patchID][i] == newPatchValues;

Also, did you try chaining the "==" to "=" inside the forAll loop?

Hi Philip,

You are right, I tried to change the "==" to "=" inside the forAll loop and it works.
In your code, you write "U.boundaryField()[patchID][i] == newPatchValues", is it U.boundaryField()[patchID]== newPatchValues? without the [i]?

Many thanks!

Wei

bigphil September 10, 2012 14:17

Quote:

Originally Posted by mathslw (Post 381094)
In your code, you write "U.boundaryField()[patchID][i] == newPatchValues", is it U.boundaryField()[patchID]== newPatchValues? without the [i]?

Yes, sorry that was a typo.

Best regards,
Philip

ancolli June 26, 2015 17:49

Quote:

Originally Posted by bigphil (Post 379492)
Hi Wei,

You can set non-uniform boundary values based on anything you want:
Code:

// find ID of patch
label patchID = mesh.boundaryMesh().findPatchID("patch_of_interest");
// check patch has been found
if(patchID == -1)
{
  FatalError << "patch not found!" << exit(FatalError);
}
 
// set value on patch
forAll(U.boundaryField()[patchID], facei)
{
  vector faceOldVel = U.boundaryField()[patchID][facei];
  vector faceCellOldVel =
        U.internalField()[mesh.boundaryMesh()[patchID].faceCells()[facei]];
  vector faceSnGrad = U.boundaryField()[patchID].snGrad()[facei];
 
  // set whatever you want here
  U.boundaryField()[patchID][facei] == faceOldVel + faceCellOldVel + faceSnGrad;
}

Best regards,
Philip


I have the following error when i try to compile the code for a scalar T:

DCP.C: In function ‘int main(int, char**)’:
DCP.C:75:66: error: no match for ‘operator[]’ (operand types are ‘Foam::tmp<Foam::Field<double> >’ and ‘Foam::label {aka int}’)
newPatchValues[i] = -T.boundaryField()[patchID].snGrad()[i];
^
make: *** [Make/linux64GccDPOpt/DCP.o] Error 1


any idea?


All times are GMT -4. The time now is 12:57.