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

twoPhaseEulerFoam BC calculated from wall shear and alpha

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

Reply
 
LinkBack Thread Tools Search this Thread Display Modes
Old   October 6, 2014, 08:08
Default twoPhaseEulerFoam BC calculated from wall shear and alpha
  #1
Member
 
Jan Potac
Join Date: Aug 2009
Posts: 30
Rep Power: 16
potac is on a distinguished road
I would like to define point-wise mesh deformation boundary condition for my dynamic mesh twoPhaseEulerFoam (OF 2.3) where the mesh displacement is calculated from values of concentration at cells adjacent to the wall and wall shear stress. For simplification, let's say the displacements = alpha*tau.
I have found angularOscillatingDisplacement source code to be probably suitable template to introduce my variables in.
Firstly, I have focused on calculation of wall shear stress. I found discussion on correcting the calculation of wall shear stress
Since my C++ knowledge is limited, I have been struggling to implement the code to obtain the shear stress. So far, I have started with code below but not sure it is the right start.
Code:
void surfaceDisplacementPointPatchVectorField::updateCoeffs()
{
    if (this->updated())
    {
        return;
    }
// * * * * * * * * * * * * * * * Wall shear  * * * * * * * * * * * * * //
    
    const fvMesh& mesh,
    const Time& runTime,
    const volVectorField& U,
    volVectorField& shear2,
    volVectorField& wallShearStress,
    volVectorField& nn
    
    #include "createPhi.H"

    singlePhaseTransportModel laminarTransport(U, phi);

    autoPtr<incompressible::RASModel> model
    (
        incompressible::RASModel::New(U, phi, laminarTransport)
    );

    const volSymmTensorField Reff(model->devReff());

    forAll(wallShearStress.boundaryField(), patchI)
    {
        wallShearStress.boundaryField()[patchI] =
        (
           -mesh.Sf().boundaryField()[patchI]
           /mesh.magSf().boundaryField()[patchI]
        ) & Reff.boundaryField()[patchI];
    }

    forAll(nn.boundaryField(), patchI)
    {
        nn.boundaryField()[patchI] =
           -mesh.Sf().boundaryField()[patchI]
           /mesh.magSf().boundaryField()[patchI];
    }

    forAll(shear2.boundaryField(), patchI)
    {
        shear2.boundaryField()[patchI]=wallShearStress.boundaryField()[patchI]-(wallShearStress.boundaryField()[patchI] & 
        nn.boundaryField()[patchI])*nn.boundaryField()[patchI];
    }
    
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // 
    const polyMesh& mesh = this->dimensionedInternalField().mesh()();
    const Time& t = mesh.time();

    scalar angle = angle0_ + amplitude_*sin(omega_*t.value());
    vector axisHat = axis_/mag(axis_);
    vectorField p0Rel(p0_ - origin_);

    vectorField::operator=
    (
        p0Rel*(cos(angle) - 1)
      + (axisHat ^ p0Rel*sin(angle))
      + (axisHat & p0Rel)*(1 - cos(angle))*axisHat
    );

    fixedValuePointPatchField<vector>::updateCoeffs();
}
In next step, I would like to call boundary field alpha. But cannot find an example easy to follow and getting lost in code lines.
Any suggestions and hints would be helpful.
J.
potac is offline   Reply With Quote

Reply

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 21:45.