symbolic differentiation with respect to dependent function

[Vinayak Vijay]

Hi,

I am trying to use symbolic differentiation to obtain the Kirchhoff stress tensor using the material model given in step 44. I have attached the minimal example herewith, which reproduces the following error:

In line 63, the differentiation of the energy function (psi_sd) with respect to left cauchy green tensor (b_sd) results in  the Kirchhoff stress tensor tau_sd = zero tensor. I believe this is due to the fact that tau_sd is stored as a function of F_sd (symbolic deformation gradient) and is therefore unable to recognise the dependence with b_sd (symbolic left cauchy green tensor).

How can i resolve this issue? Thanks in advance.

Regards

Vinayak

[Wolfgang Bangerth]

On 7/17/24 09:15, Vinayak Vijay wrote:
>
> In line 63, the differentiation of the energy function (psi_sd) with respect
> to left cauchy green tensor (b_sd) results in  the Kirchhoff stress tensor
> tau_sd = zero tensor. I believe this is due to the fact that tau_sd is stored
> as a function of F_sd (symbolic deformation gradient) and is therefore unable
> to recognise the dependence with b_sd (symbolic left cauchy green tensor).
>
> How can i resolve this issue? Thanks in advance.

If I understand your question right, you have a function
f(x)
and you are computing
d/dy f(x)
which is indeed zero because you have not said anything about how x and y are
related. If x can be expressed as a function of y, say x=g(y), then you need
to tell the autodifferentiation machinery that your function f is really
f(g(y))
i.e., everywhere where you previously used x, you need to replace it by g(y)
to make it clear that f is a function of y.

Does this help?

Best
W.

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

[Vinayak Vijay]

Dear prof. Bangerth,

Thanks for your response. Two things: 1) I am using symbolic differentiation. 2) Here's a simpler implementation of my problem in line with your example:

const Differentiation::SD::Expression x("x");

const Differentiation::SD::Expression y = std::pow(x,2);

std::cout << "y: " << y << std::endl;

const Differentiation::SD::Expression g = std::pow(y,2);

std::cout << "g: " << g << std::endl;

const Differentiation::SD::Expression dg_dy = Differentiation::SD::differentiate(g,y);

std::cout << "dg_dy: " << dg_dy << std::endl;

const Differentiation::SD::Expression dg_dx = Differentiation::SD::differentiate(g,x);

std::cout << "dg_dx: " << dg_dx << std::endl;

Output:

y: x**2
g: x**4
dg_dy: 0
dg_dx: 4*x**3

Here, I have a function g=g(y(x)). I wish to know its derivative with respect to x and y. Manually, we would obtain both using the chain rule. 

But here, it looks like the information about the composition of the functions is not being stored somehow, maybe because it reduces the expressions to the independent variables ("x" in this case). What change should I make to obtain dg_dy? The independent variable remains "x" of course.

Regards

Vinayak

[Wolfgang Bangerth]

On 7/18/24 00:22, Vinayak Vijay wrote:
>
> But here, it looks like the information about the composition of the functions
> is not being stored somehow, maybe because it reduces the expressions to the
> independent variables ("x" in this case). What change should I make to obtain
> dg_dy? The independent variable remains "x" of course.

I do not actually know enough about SD to help with this. The way I read the
code is that you declare x to be an independent variable, whereas y is not.
Intuitively, I would have expected the underlying library to raise an
exception when you try to compute the derivative with regard to anything other
than an independent variable.

I think you'll have to find out what underlying SD library your example is
using, and then read through the documentation of that library to see how to
approach the issue. Alternatively, you can of course define a function g_of_x
and another g_of_y, and depending on whether you want to compute the
derivative with regard to x or y, you choose one or the other.

[Vinayak Vijay]

Dear prof. Bangerth,

Thanks. I have decided to implement my problem using the alternate approach. However, I am now facing another issue, somewhat different from the above but related to the Differentiation::SD class.

I am trying to compute the fourth-order elasticity tensor in the spatial description using the formula Jc = 4*b*d2W_db_db*b (as in step-44). I am using symbolic differentiation to compute all of the quantities in the formula. Finally, I would like to have Jc of the type SymmetricTensor<4, dim, NumberType>. This is where I run into a problem, which I explain below.

Since all the quantities in the rhs are SymmetricTensors, the operator*() recognises the multiplication as a double contraction. However, we want it to be a single contraction. So I must first convert 'b' to the type Tensor and then apply the formula. This will lead to Jc being of the type Tensor and not SymmetricTensor. Now there's a function symmetrize() to symmetrize the type Tensor<2,dim>, but not for higher-order tensor types. One could manually do this, but is there a better way?  If not, would it be a good idea to extend the symmetrize() function to Tensor<4,dim> since this type is often used to represent many quantities, such as the elasticity tensors?

Thanks & regards

Vinayak

[Wolfgang Bangerth]

On 7/24/24 10:34, Vinayak Vijay wrote:
>
> Since all the quantities in the rhs are SymmetricTensors, the
> operator*() recognises the multiplication as a double contraction.
> However, we want it to be a single contraction. So I must first convert
> 'b' to the type Tensor and then apply the formula. This will lead to Jc
> being of the type Tensor and not SymmetricTensor. Now there's a function
> symmetrize() to symmetrize the type Tensor<2,dim>, but not for
> higher-order tensor types. One could manually do this, but is there a
> better way?  If not, would it be a good idea to extend the symmetrize()
> function to Tensor<4,dim> since this type is often used to represent
> many quantities, such as the elasticity tensors?

Vinayak,
I would accept a patch for such a function, though perhaps under a
different name. The issue with the term 'symmetrize()' is that there are
multiple ways in which one could interpret it for tensors of rank 4. Do
you want the result to satisfy
A_{ijkl} = A_{ijlk} = A_{jikl} = A_{jilk}
or do you want it to be
A_{ijkl} = A_{klij}
or something else? The same discussion of course applies to the name
"symmetric tensor of rank 4", which specifically uses the first of these
definitions above -- i.e., it maps symmetric tensors of rank 2 to
symmetric tensors of rank 2.

Best
W.

[Jean-Paul Pelteret]

Hi Vinayak,

I wrote a reply a few days ago, but apparently I managed to address it to Wolfgang rather than the mailing list. Oops... I know that you've made some progress since then but since maybe these details could be of some further interest I've copied it verbatim. I think that the last sentence aligns with what you've already done.

-----------------

To answer your earlier question [from 18 July] directly, I think that the feature that you're looking for is called make_symbolic_function() for scalars, or make_[vector, tensor, symmetric_tensor]_of_symbolic_functions() for higher dimension entities. These allow you to generically express symbolic functions with other functions as arguments, and differentiate correctly.

That said, I don't think that this is what you want to do. There is an example in the tests for step-44, namely tests/symengine/step-44-sd-quadrature_level_0[1,2,3].cc that exploits the relationship between det(b) and det(F) = sqrt(det(b)) = sqrt(det(C)) to compute the kinetic variables and their derivatives. I think that you'd be able to solve your problem most easily if you follow suite. Otherwise you'd really have to differentiate with respect to F, which you state as being the true parameterisation of the energy, and then push forward the Piola-Kirchoff stress tensor and its linearisation (e.g. using the functions in the namespace Physics::Transformations if any of them seem like good candidates).

-----------------

Best,
Jean-Paul

[Vinayak Vijay]

Prof. Bangerth,

I agree that one will have to first define "symmetry" for higher-order tensors. Perhaps one solution, at least for the 4th-order tensors, could be to have a function 'symmetrize()' (or some other name) with additional arguments specifying the type of symmetry, i.e. the so-called - minor symmetry (A_{ijkl} = A_{ijlk} = A_{jikl} = A_{jilk}) or the major symmetry (A_{ijkl} = A_{klij}) or both. I'll try to write a patch.

JP,

Thanks! I decided to change my parameterization for now.

Regards

Vinayak

[Wolfgang Bangerth]

On 7/25/24 09:35, Vinayak Vijay wrote:
>
> I agree that one will have to first define "symmetry" for higher-order
> tensors. Perhaps one solution, at least for the 4th-order tensors, could be to
> have a function 'symmetrize()' (or some other name) with additional arguments
> specifying the type of symmetry, i.e. the so-called - minor symmetry (A_{ijkl}
> = A_{ijlk} = A_{jikl} = A_{jilk}) or the major symmetry (A_{ijkl} = A_{klij})
> or both. I'll try to write a patch.

Great, we look forward to it! I think the additional argument is a good approach.