
[Sponsors] 
July 19, 2013, 15:48 
Setting normal component of a vector to zero

#1 
Member
Chris
Join Date: Aug 2012
Location: Calgary, Alberta, Canada
Posts: 73
Rep Power: 5 
In my efforts to implement a thermal diffusion model I've run into another snag. Essentially what I think is happening is I have species diffusing through the wall. I suspect this is because I was severely underpredicting the wall concentration when I was giving no special treatment to the diffusion coefficient at the wall.
Initially I tried just setting the diffusion coefficient at the wall boundary to zero, but that seems to lead to an over prediction of the wall concentration. I suspect it isn't working because that doesn't allow for species to diffuse away from the wall. So what I need to do is check the direction of the species transport and if it is into the wall set the wall normal component to zero, otherwise let it happen... here's what I have so far: Code:
volVectorField JTherm = DT*fvc::grad(T); forAll(JTherm.boundaryField(),patchi) { if (isA<wallFvPatch>(mesh.boundary()[patchi])) { const vectorField& surfaceNormal = mesh.boundary()[patchi].nf(); forAll(JTherm.boundaryField()[patchi],celli) { scalar dotprod = surfaceNormal[celli] & JTherm[celli]; if(dotprod < 0) { //set normal component to zero Edit: well that was a bit silly, two cross products should do exactly this. 

July 22, 2013, 03:23 

#2 
Member
Artem Shaklein
Join Date: Feb 2010
Location: Russia, Izhevsk
Posts: 43
Rep Power: 8 
Why don't you consider setting up diffusive flux of your variable T to zero by using zeroGradient (dT/dn = 0) boundary condition?


July 22, 2013, 03:29 

#3 
Member
Chris
Join Date: Aug 2012
Location: Calgary, Alberta, Canada
Posts: 73
Rep Power: 5 
You're right that would work for cases where T (temperature) has a zero gradient at the wall, however this is not true for the case I'm trying to run (fixed wall temperature based on experiments).


July 22, 2013, 05:36 

#4 
Member
Artem Shaklein
Join Date: Feb 2010
Location: Russia, Izhevsk
Posts: 43
Rep Power: 8 
Hello, ChrisA.
Ok, I get your point) When you computes JTherm = DT*fvc::grad(T); nonzero diffusive flux is already in cell adjacent to boundary. And JTherm boundary will not be used in further computations (strictly speaking, it will be used, but will make changes only for variable fields with calculated bc types). Why you are trying to take gradient of vector? Gradient means changing of your variable across each of three dimensions. So, you want to get 3x3 tensor matrix using gradient operation (∂T1/∂x, ∂T1/dy, ∂T1/∂z, ∂T2/∂x....)? Anyway, if you want to compute grad(T) without wall diffusive flux affection, just dig dipper: fvc::grad(volField) with Gauss scheme is computed in gaussGrad.C (you can find file by "locate gaussGrad.C"). To be more precise, function to search for is gradf. This part of function operates with boundaries Code:
forAll(mesh.boundary(), patchi) { const labelUList& pFaceCells = mesh.boundary()[patchi].faceCells(); const vectorField& pSf = mesh.Sf().boundaryField()[patchi]; const fvsPatchField<Type>& pssf = ssf.boundaryField()[patchi]; forAll(mesh.boundary()[patchi], facei) { igGrad[pFaceCells[facei]] += pSf[facei]*pssf[facei]; } } Code:
forAll(mesh.boundary()[patchi], facei) { igGrad[pFaceCells[facei]] += pSf[facei]*pssf[facei]; } 

July 22, 2013, 09:32 
one additional question...

#5 
Member
Andreas Ruopp
Join Date: Aug 2009
Location: Stuttgart / Germany
Posts: 31
Rep Power: 8 
@ARTem:
Many thanks for the explanation. One question. After: Code:
... forAll(mesh.boundary()[patchi], facei) { igGrad[pFaceCells[facei]] += pSf[facei]*pssf[facei]; } } the next function call is: Code:
gGrad.correctBoundaryConditions(); Code:
forAll(vsf.boundaryField(), patchi) { if (!vsf.boundaryField()[patchi].coupled()) { vectorField n = vsf.mesh().Sf().boundaryField()[patchi] /vsf.mesh().magSf().boundaryField()[patchi]; gGrad.boundaryField()[patchi] += n * ( vsf.boundaryField()[patchi].snGrad()  (n & gGrad.boundaryField()[patchi]) ); } } For what does this part stand for? Best regards, Andy 

July 22, 2013, 12:21 

#6 
Member
Chris
Join Date: Aug 2012
Location: Calgary, Alberta, Canada
Posts: 73
Rep Power: 5 
Before I start I'd just like to note that I don't 100% understand how OF handles boundaries so please forgive any ignorance included below.
This is good discussion, probably useful for those who are having issues setting species diffusion to zero at boundaries (mass flow inlets for example) among other things but I think you are missing what I'm trying to do, which is mostly my fault. I should have been more clear with what I'm trying to implement. I'm trying to use a species thermal diffusion model, which is a species tendency to travel along temperature gradients. This is why you see a grad(T) in my code. T is temperature and its "transport", for lack of a better word, of course gets handled in the energy equation. With the model I end up with a mass flux, JTherm, based on the temperature gradient. The next step in the code is to compute the transport of species Yi by taking the divergence of JTherm. While manipulating the way the gradient is taken at the wall will fix the species transport it will impact the energy transport in a nonphysical way (I think?). Now if I understand the first part of what you said right (the bit about diffusive flux being in the cell adjacent to the boundary). I don't think this is the case because when I change DT at the boundary (example set it to zero or the internal field value) the solution changes. If what you're saying is true wouldn't the solution be independent of DT's boundary value? I'm also not sure where I'm taking the gradient of a vector? Perhaps inherent in the confusion as to what T is in this case? If that is the case I apologize for not being clearer from the start. Last edited by ChrisA; July 22, 2013 at 15:15. 

July 23, 2013, 03:03 

#7  
Member
Artem Shaklein
Join Date: Feb 2010
Location: Russia, Izhevsk
Posts: 43
Rep Power: 8 
Hello, andyru.
Don't forget about functions with same name but taking different arguments. Code:
gGrad.correctBoundaryConditions(); and Code:
template<class Type> void Foam::fv::gaussGrad<Type>::correctBoundaryConditions ( const GeometricField<Type, fvPatchField, volMesh>& vsf, GeometricField < typename outerProduct<vector, Type>::type, fvPatchField, volMesh >& gGrad ) 2ChrisA. Sorry for grad of a vector stuff, my mess) So, you have an equation ∂(ρ Yi)/∂t + div(ρ Yi U)  div( (µ/Sc + µ_t/Sc_t) grad(Yi) ) = div ( DT grad(T) ) ? Strange, but DT = 0 at walls should do the job. Quote:
You can try to do next (with DT=0 at walls). Consider just 1d case with 2 walls. Firstly, set, e.g., T_left = 400, T_right = 100. When steadystate is achieved, switch values between T_left and T_right and continue a simulation. Compare intermediate and final solutions. Is there a flux from cell near to wall in domain with DT = 0 at walls? About gradT and energy equation. If you create new set of functions and will call fvc::gradZeroBC(T), with original functions untouched, all will be fine. 

July 23, 2013, 12:36 

#8 
Member
Chris
Join Date: Aug 2012
Location: Calgary, Alberta, Canada
Posts: 73
Rep Power: 5 
ARTem,
That is the general idea behind my species transport equation yes. It's known as Soret diffusion and is often neglected due to its relatively high computation cost but it plays a non negligible role in the diffusion of notably different sized molecules when a temperature gradient is present (that's an oversimplification of when its nonnegligible but that's the case for what I'm working on). I see what you mean now with the grad operator. DT=0 at the walls ended up severely over predicting the wall concentration. When the DT value is interpolated to the internal faces does it use the wall value at all in that interpolation? I ask because I think the over prediction is because at some point things want to start diffusing away from the wall. The equation is in such a form that the sign will flip if the concentration of the diffusing species gets high enough, so what was once diffusing to the wall will diffuse away till steady state and whatnot. I should also note that I've been suspecting the boundary as the source of my issue as the concentrations away from the wall are lining up quite well with simulations in previous work on the same problem. 

July 24, 2013, 03:01 

#9  
Member
Artem Shaklein
Join Date: Feb 2010
Location: Russia, Izhevsk
Posts: 43
Rep Power: 8 
Quote:
create 1d mesh with 2 nodes (or 3d with 8 nodes  2 in each directions  to be sure by 100%). Create volField DT with, e.g. value "1", and in one boundary (or in all boundaries) set DT to be "0". Then call for Code:
surfaceScalarField DTf("DTf", fvc::interpolate(DT)); DTf.write(); 

Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Setting Opening Boundary to Move Normal to Itself  NCle  CFX  10  April 20, 2012 07:54 
Cells with t below lower limit  Purushothama  CDadapco  2  May 31, 2010 21:58 
my stl surface is seen as just a line  rcastilla  OpenFOAM Meshing & Mesh Conversion  2  January 6, 2010 02:30 
NACA0012 geometry/design software needed  Franny  Main CFD Forum  13  July 7, 2007 15:57 
Warning 097  AB  CDadapco  6  November 15, 2004 05:41 