Matrix manipulation
1 Attachment(s)
Hello,
I am in need of help (OpenFOAM 5.0), unfortunately I do not know much where to start and so I am here in the forum. I need to implement a function in the k-Epsilon model, called the Peterlin function (attached image). 'L' is a constant and 'C' is a matrix. How could I write the code and tell OpenFOAM that I want the elements of the main diagonal (in the division)? And how will it interpret that the process should be in steps (first the C11 ... then the C22 ...)? Understanding this is very important to me. |
If C is a square matrix, then you can loop over the elements to get the diagonal as a field like this:
Code:
scalarField diag(C.n()); Code:
scalarField diag(3); |
Quote:
I have one more question and I hope you can help me. Look this: Quote:
The point I want your attention is, this code will give me 6 values (since the matrix is symmetric), consequently fvMatrix == fvMatrix. Therefore, I can not enter a scalar value into this equation, right?! I ask this because I have another question and maybe you can help me too. When I analyze the code: Quote:
So... I'm inferring the logic I wrote above (fvMatrix == fvMatrix). Now I have a '' == scalar - fvMatrix_ ", or no?! ...and another detail is, when I call k_(), I'm calling a scalar or a matrix? What does fvScalarMatrix do? All this questioning is because I tried to make this modification in the code: Quote:
All these questions are very important to me, believe me... if you can help me, it will make a big difference. |
fvScalarMatrix constructs a matrix for the scalar equation 'k'. So in the 'k' equation, 'G' is scalar or Production Term (Sij.Sij, Sum of squares of the components of rate of strain tensor).
Coming to the second equation which you post, here you are constructing a solver matrix for Tensor variable, hence all the terms should be of Tensor type (dimensional balance) Now coming to the modification you are trying to do for the k-Epsilon, I think you will be adding some terms to k equation from Cij, and you will get Cij from other equation. For this, I will suggest you to write the 'k-e' equation without any notation form and see the indexes of Cij coming in. The thing is 'k-e' are basically scalar transport equation, so the term which your adding should be scalar. Now, to extract the scalar terms from the tensor, have a look in the programming manual of Openfoam, you will find appropriate function. - gtarang |
Let me try to answer your questions in a systematic manner to the best of my knowledge.
Quote:
Quote:
Quote:
Quote:
So, your modification, Code:
// Stress transport equation USV |
1 Attachment(s)
Hi usv001,
...thanks for your reply. I want to insert the equation below. You can use this link as a guide. The first terms of the equation are div(phi(), A_) and laplacian(nut(), A_). The rest I'm having trouble writing. As I said, the link serves as a reference. I understand that to simplify the process, I could call Ckk from tr(C) and unfortunately I do not know if this would affect the results (maybe) but, the logic you showed me, I understand how complex to my knowledge of C++ and I do not know how to implement it. If you can help me. So, Ckk is one of my biggest problems (...and this also involves the div and the laplacian) and also the second term after equality (whose are scalar values).. I was thinking of multiplying k_ * I, but, I dont know... or something like this: (scalar) * I. EDIT1: How would you write the S² equation? volSymmTensorField twoD = twoSymm(fvc::grad(U()); ? |
Hi usv001,
I've made some great progress but I still need some help because I do not know if what I'm doing is correct. The solver compiles normally and when I ask to execute it it presents the floating point error, follows a piece of the log: Quote:
Quote:
I need to learn how to make the product of the tensor strain rate (Sij), because in the equation it is squared. I believe that if I do this correctly, the values Cxy ... Cyz ... will reappear in the log. |
Quote:
If you could help me, I would be very happy. |
You could try something like the following to solve for individual diagonal components:
Code:
// labelList corresponding to diagonal components in symmTensor USV |
1 Attachment(s)
Quote:
Apparently the off-diagonal are not solved, according to the equations attached below. I understand that by using Ckk it is only referring to the elements of the main diagonal. Taking advantage of it, is it possible to insert these equations directly into the solver simpleFoam and thus, to simultaneously use the turbulence models already implemented in the OF? |
If I want to calculate a dumping function that is a function of uTau, how do I insert the uTau equation into the kEpsilon.C file??
|
Hello Guilherme,
Quote:
USV |
Quote:
I wrote this snippet of code: Code:
Dxy_ Unfortunately, however, the gradient of xy on the wall is giving zero and this causes the solver to give a floating-point error. log: LINK Could someone help me with how to solve this problem? Best. *By the way, if I want to enter the value of Dxy in a tau (symmetric) equation, how would you recommend that I do? |
Hello Guilherme,
Your code only seems to be transferring the internalField only. To transfer the boundary values, try the following: Code:
forAll(D.boundaryField(), patchi) Please post the equation for tau in code form so that people know what exactly you are trying to do because Dxy_ (a volScalarField) could be introduced many ways, as a source term, as a coefficient multiplying another term, etc. in an fvSymmTensorMatrix. USV |
1 Attachment(s)
Quote:
I've tried writing the code that sent me up this way: Code:
volSymmTensorField D = symm(gU); Code:
error: no match for ‘operator==’ (operand types are ‘const Foam::fvPatchField<double>’ and ‘Foam::tmp<Foam::Field<double> >’) LINK I'm trying to write equations 8 and 11 described in this article. These are the set of the conformation tensor (equation 18) that compose the tau equation (equation 19). The way I've written both the last two is below. Code:
volScalarField fP = (L_-3.)/(L_-tr(C_)); It seems to me that M is related to: volSymmTensorField D = symm(fvc::grad(U())); ...but I dont know. Another way in which it is possible to write equation 8 is one that was presented by a precursor author to the study presented. I attached the equation as an image. So, for not knowing how to write in a way (eq. 8), I'm limiting myself to writing the problem in a 2D way, using the attached equation, it uses the squared gradUxy (this is possible since it's a scalar). And that's why I asked you for help. This leads me to question that if I make this adaptation, the equation of the symmTensor does not necessarily need to be fvSymmTensorMatrix, it would only be fvScalarMatrix, since it will solve only the xy component and consequently it will only give me tauxy (which is what I need for a problem 2D). I want to remember that NLT is in function of M and both are not demonstrated in the code that I attached because not yet how to write them. *Incidentally, equations 22 and 24 described in the article are also dependent on M, but in this case, I BELIEVE that they are only dependent on the xy component of the U gradient. For me, the ideal would be to write the equation 8 and thus, to be able to model a 3D problem. All this argument is in the hope that you can help me. PS: Incidentally, if I can ask for more detail, in case I want to write gradU², that is, sqr(fvc::grad(U())), how would it be done? |
All times are GMT -4. The time now is 20:32. |