[celestial]

Hi

I am looking at the include file solveEnergy.H of the conjugateHeatFoam solver where the temperature for a fluid mesh adjacent to a solid mesh is to be solved using the diffusion equation. The code creates one

global matrix equation Ax =b encompassing both meshes

and my BC at their interface is

k1 * grad T1 = k2 * grad T2

So obviously that BC has to be part of the overall matrix equation

But I can't find that BC listed anywhere under the FvPatchFields subfolders


Any ideas as to how I should proceed ?

Regards to all

[celestial]

I forgot to mention that this thread is only meaningful to those familiar with foam-extend.

Even with only one mesh to work with, the next level of difficulty after Dirichlet and Neumann BCs is the Robin BC. Now mixed BC is not a bona fide general Robin BC but you cast it into such a form because luckily there is no check for bounds 0<f<1 . Although, if you are interested in non-uniform coefficients, you have no choice but to use the groovyBC version of mixed BC.

Now, since my last post, I could not find anything more sophisticated than that mixed BC. There does not seem to be any readily available function to modify the matrix coefficients when considering the conductivity flow BC at the interface between 2 meshes:

k1 * grad T1 = k2 * grad T2

I don t know how hard it is to implement.

For the Poisson equation instead of the diffusion equation, minimizing its functional shows that the BC is implicitly part of the finite element method solution. You don t need to add it. ( e.g. see software Elmer)

[celestial]

To make it simpler, I would be happy to discretize all these BCs as

Xp = e * Xc + f

for one mesh where X is a scalar field and Xp is its value at the patch and Xc is its value at the center of the cell with that patch as one of the surfaces enclosing that cell.

For the solving of

AX = b

I would like to know where that matrix A and vector b are located in the code and also the coordinates (i,j) of matrix A where the coefficient e will be stored
and thus the coordinate i of vector b where coefficient f will be stored
(which gives the coordinate (i,i) in A where the coefficient 1 of Xp is stored)

Of course if there are two meshes and the patch is at their interface, the matrix A encompassing both meshes will now be larger.

Moreover my BC now becomes

Xp = e * Xc + f + g * Xn

where Xn represents the same thing as Xc but for the other mesh and I will also need to know the coordinates (i,j') of matrix A where the coefficient g will be stored.


It could very well be that only the experts who wrote the code can answer that query.

Have a nice weekend to all

[celestial]

To sum up, unless somebody proves me wrong, OpenFOAM was never designed to modify the coefficients of the matrix equation to be solved even for so standard a boundary condition as

k1*grad T1 = k2*grad T2 + S

(adjacent meshes with different heat conductivities and a source term S at their interface. I am interested in k1, k2, S each a function of space and time but the BC does not even exist for constant k1, k2 and S )

I fear that adding this feature may require considerable programming effort.

[dlahaye]

Not sure about S <> 0.

In case that S = 0 , k1 grad T1 = k2 grad T2 is part of the formulation that e.g. chtMultiRegionFoam solves. Older versions of chtMultiRegionFoam iterate between subdomains. Newer version implement a monolithic solver.

For S<>0, I would look into line sources first.

[celestial]

Thank you very much Domenico but I think the monolithic solver you mention is actually

conjugateHeatFoam

available in foam-extend.

However, it is not enough to use patch regionCoupling in a dictionary. I will need to modify the source code - probably in regionCouplingFvPatchField.C.

This will be hard as one has to first understand the source code and the fact that there are hardly any comments accompanying the code for sure does not make the task easier ( I was taught more than 3 decades ago to always add relevant comments to code).

[dlahaye]

Dear Celestial,

Thank you you for pointing me to conjugateHeatFoam.

I imagine that this coupled solver (at https://sourceforge.net/p/foam-exten.../solveEnergy.H for future reference) does away with the need of implementing the interface conditions. This would similar to the Elmer solution and what we do for finite elements.

Remains the question what you want to implement.

Is it your wish to implement a line source of the solid/gas interface? Are you happy to diffuse the support of the source over a few cells to start from? Can you employ fvOptions to this end?

Cheers. Domenico.

[celestial]

I am actually solving am electrostatics problem - so the Poisson equation over 2 meshes with the BC

e1 * grad E1 = e2 * grad E2 + Surface Charge

at the interface where E stands for electric field and e for electrical permittivity.

I was justing recasting the same equation in a static mechanical engineering context so that most CFD folks would be very likely familiar with heat flow.

It is true that in Elmer you just specify the surface charge and the user does not need to write any source code for the above BC since it is already built-in.

Looking at solveEnergy.H though, I can see it does not implement anything about BCs.

Which means I will probably have to write it as I alluded to in my last post.
I am not very familiar with fvOptions but anyhow I don't want to shift the BC to inside one the mesh and close to the interface like you were suggesting.