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/)
-   -   Modifying the laplacian operator (http://www.cfd-online.com/Forums/openfoam-solving/65292-modifying-laplacian-operator.html)

mlawson June 10, 2009 15:45

Modifying the laplacian operator
 
Hi everyone,

I would like to modify the laplacian so the derivative is only calculated in one direction - that is I would like to compute only d^2T/dx^2 and not d^2T/dx^2 + d^2T/dy^2 + d^2T/dz^2

Then I want to solve a system of equations using this modified laplacian. I tried to do this by dotting the result of the laplacian operator with a unit normal vector in the x-direction. A simple example would be:

solve
(
fvm::ddt(T)
- fvm::laplacian(DT, T) & xDir // where xDir is a dimentionedVector with the units [0 1 -1 0 0 0 0] (1 0 0)
);


Unfortunately, I get the following error when trying to compile my solver.

---------------------------------------
1dScalarTransportFoam.C: In function ‘int main(int, char**)’:
1dScalarTransportFoam.C:67: error: no match for ‘operator&’ in ‘Foam::operator-(const Foam::tmp<Foam::fvMatrix<Type> >&, const Foam::tmp<Foam::fvMatrix<Type> >&) [with Type = double](((const Foam::tmp<Foam::fvMatrix<double> >&)((const Foam::tmp<Foam::fvMatrix<double> >*)(& Foam::fvm::laplacian(const Foam::dimensioned<Type2>&, Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>&) [with Type = double, GType = double](((Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>&)(& T))))))) & xDir’
---------------------------------------

I assume I get this error because the & operator is not defined for a fvScalarMatrix.

Does anyone know a method modifying the laplacian operator so that the derivative is only calculated in one direction?

-Mike

ziemowitzima August 3, 2011 17:27

Hi,
I am trying to do the same thing.
Did you manage to do this ?
Method you described is wrong, because laplacian(DT,T) is not a vector , but scalar:
Txx + Tyy + Tzz, so you cant make dot product operation for scalar and vector...
I did it in explicit way, but this not the best way...

If you have some updates how to solve this problem, pleas let me know.

best
-ZM

mlawson August 3, 2011 17:54

You are correct the result of the laplacian operator is a scalar in this case. What should have asked is, is there a way to dot the laplacian operator itself with a unit vector before it is applied to the variable T (i.e (d/dx^2 + d/dy^2 + d/dz^2) dot (1 0 0) * T).

Either way, I never did solve this problem in a robust way. I ended up tricking the solver into only computing diffusion in one direction by making my two dimensional domain very long in one of the non-empty directions, the y-dir in my case, so the d/dy^2 term was negligible with respect to the d/dx^2 term. I know this is not a good way to solve the problem, but it ended up working for the specific case I was working on.

Good luck and please let me know if you find a solution to the problem.

-Mike

ziemowitzima August 3, 2011 18:23

in explicit way u can do this like that:
volVectorField gradT=fvc::grad(T);
volScalarField gradTx = gradT & vector(1,0,0);
volScalarField Tx = gradT.component(0); // T_x
volVectorField gradT2 = fvc::grad(Tx); // (T_xx, T_xy)
volScalarField Txx = gradT2.component(0); // T_xx

akidess August 4, 2011 04:32

Michael, it sounds like you want anisotropic diffusion - why not just pass a diffusion coefficient tensor to laplacian instead of a scalar? I don't think you need any modifications to the code to do that (but I haven't tried it myself).

mlawson August 4, 2011 11:32

Hi Anton,

Mathematically that makes sense, but I'm unsure how to convert the scalar field T into a tensor in OpenFOAM. Do you have an idean about how this would be done?

-Mike

akidess August 4, 2011 11:49

T stays a scalar field, you only have to change the diffusion coefficient into a tensor (you wouldn't do anything else on pen and paper). Redefine DT from dimensionedScalar to dimensionedTensor (see e.g. http://www.foamcfd.org/Nabla/guides/...sGuidese5.html) and recompile.

ziemowitzima August 4, 2011 13:13

Hi Akidess,
Thanks for your suggestion.
It works just as it should.
Human being is learning all its life... :)

To Mlawson:
x - direction diffusion:
DTe DTe [ 0 2 -1 0 0 0 0 ] (3e-02 0 0 0 0 0 0 0 0 );

x,y and - directions diffusion:
DTe DTe [ 0 2 -1 0 0 0 0 ] (3e-02 0 0 0 3e-02 0 0 0 3e-02 );

Thanks again!

mvoss September 7, 2011 12:51

hey there.

could anybody please give me a hint on where to find some sort of documentation for a part of the FSI-solver (icoFSIFoam??) like mentioned above. i am espacially interessted in how the pde´s have been derived for the solid/stress part.

e.g. fvm::laplacian (2*mu, + lambda, Usolid, "laplacian (DU,U)")

what is that third part good for? I wasn´t able to get any doc. for a third part in laplacian(..).
Is this some sort of comment??

Thanks in advance,
neewbie

ziemowitzima September 7, 2011 13:05

Third part "laplacian(DU,U)" is just needed for fvScheme file, where laplacian(DU,U) is defined. Here it just says that variable DU = 2*mu + lambda.

ZMM

mvoss September 7, 2011 13:06

okay.
found it. i guess.
it´s a label for that part of the eq. to hook it up with the schemes in fvSchemes.

ziemowitzima September 7, 2011 13:30

I think that all needed information you can find in the solver description:l
inear-elastic, small-strain deformation of a solid body, with optional thermal
diffusion and thermal stresses.
For example in solver solidDisplacementFoamfact there are only two PDE solved (I dont know where FSIsolver is):

fvm::ddt(T) == fvm::laplacian(DT, T)

and

fvm::d2dt2(D)
==
fvm::laplacian(2*mu + lambda, D, "laplacian(DD,D)")
+ divSigmaExp

first equation is just standart heat equation,
second equation is just linear-elastic equation with optional thermat stress (+ divSigmaExp) influence (calculated by solving first equation).
This are standart equations so they derivation can be easly found in basic PDE books, or online:
http://en.wikipedia.org/wiki/Linear_elasticity

ZMM



if (thermalStress)


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