Assigning anisotropic thermal conductivities in the assemble loop

[Morris Jowas]

Hello everyone,

 

I am solving the heat equation with the Laplacian expressed in terms of spatially varying thermal conductivities (k_x, k_y and k_z) as shown below:

Can I please request your suggestions on how I may go about assigning the three thermal conductivities to every cell in the assemble loop? The examples I have seen in dealii provide sufficient guidance for a homogeneous problem where k is the same in all three directions per cell.

For example, in the below loop I solve a version of the problem where there are two materials in the domain each with material properties k1 and k2 which I assign based on the cell number contained in the vector geomDef assuming that each cell is homogenous. However, if the cell is anisotropic, how can I set up the condition so that dealii knows the direction in which I am pointing so I can assign the correct thermal conductivity accordingly?

Here is the code

for (const auto& cell : dof_handler.active_cell_iterators())

{

    cell_matrix = 0;

    cell_rhs = 0;

    fe_values.reinit(cell);

    for (const unsigned int q_index : fe_values.quadrature_point_indices())

        for (const unsigned int i : fe_values.dof_indices())

        {

            for (const unsigned int j : fe_values.dof_indices())

 

                if (geomDef[curr_cell] == 1)

                {

                    cell_matrix(i, j) +=

                        (k1 * fe_values.shape_grad(i, q_index) *  // assigning k1, thermal conductivity for material 1. Help needed here for allocating the thermal conductivity correctly

                            fe_values.shape_grad(j, q_index) *

                            fe_values.JxW(q_index));

                }

                else

                {

                    cell_matrix(i, j) +=

                        (k2 * fe_values.shape_grad(i, q_index) * // assigning k2, thermal conductivity for material 2

                            fe_values.shape_grad(j, q_index) *

                            fe_values.JxW(q_index));

                }

        }

}

 

Thank you in advance for your help with this.

 

MJ

[Bruno Turcksin]

MJ,

fe_values.shape_grad() returns a Tensor<1,dim> so you can do something like this

Tensor<1,dim> k_grad;

for (int d=0; d<dim; ++d)

  k_grad[d] = k[d]*fe_values.shape_grad(i, q_index)[d];

Then you case use k_grad instead of k*fe_values.shape_grad(i, q_index)

Best,

Bruno

[Morris Jowas]

P.S. Here is the equation. Sorry it did not attach last time around: 

[Morris Jowas]

Thank you very much Bruno. I am implementing your recommendation now.

Best regards,

MJ

[Wolfgang Bangerth]

On 9/28/21 11:59 AM, Morris Jowas wrote:
> Thank you very much Bruno. I am implementing your recommendation now.

step-6 already deals with a case where the coefficient is spatially variable,
but not tensor-valued. When you've got it figured out how to use tensor-valued
coefficients, might you be interested in writing a part for the "Possibilities
for extensions" section of step-6 that explains how to use tensored
coefficients? I think that would be a really nice contribution to the library!

Best
Wolfgang

--
------------------------------------------------------------------------
Wolfgang Bangerth email: bang...@colostate.edu
www: http://www.math.colostate.edu/~bangerth/

[Morris Jowas]

That would be my pleasure Prof. I will endeavour to make a contribution once I test and validate my implementation. 

Best regards,

MJ

[Wolfgang Bangerth]

On 9/28/21 4:24 PM, Morris Jowas wrote:
> That would be my pleasure Prof. I will endeavour to make a contribution once I
> test and validate my implementation.
>

Outstanding! Feel free to ask here how to concretely submit patches for
inclusion, but in short, this is how this is done:
https://github.com/dealii/dealii/wiki/Contributing

Best
W.