[petsc-users] Elasticity tensor in ex52

[Matthew Knepley]

On Sun, Apr 20, 2014 at 4:17 PM, Miguel Angel Salazar de Troya <
salazardetroya at gmail.com> wrote:

> I understand. So if we had a linear elastic material, the weak form would
> be something like
>
> phi_{i,j} C_{ijkl} u_{k,l}
>
> where the term C_{ijkl} u_{k,l} would correspond to the term f_1(U,
> gradient_U) in equation (1) in your paper I mentioned above, and g_3 would
> be C_{ijkl}. Therefore the indices {i,j} would be equivalent to the indices
> "ic, id" you mentioned before and "jc" and "jd" would be {k,l}?
>
> For g3[ic, id, jc, jd], transforming the four dimensional array to linear
> memory would be like this:
>
> g3[((ic*Ncomp+id)*dim+jc)*dim+jd] = 1.0;
>
> where Ncomp and dim are equal to the problem's spatial dimension.
>
> However, in the code, there are only two loops, to exploit the symmetry of
> the fourth order identity tensor:
>
>   for (compI = 0; compI < Ncomp; ++compI) {
>     for (d = 0; d < dim; ++d) {
>       g3[((compI*Ncomp+compI)*dim+d)*dim+d] = 1.0;
>     }
>   }
>
> Therefore the tensor entries that are set to one are:
>
> g3[0,0,0,0]
> g3[1,1,0,0]
> g3[0,0,1,1]
> g3[1,1,1,1]
>
> This would be equivalent to the fourth order tensor \delta_{ij}
> \delta_{kl}, but I think the one we need is \delta_{ik} \delta_{jl},
> because it is the derivative of a matrix with respect itself (or the
> derivative of a gradient with respect to itself). This is assuming the
> indices of g3 correspond to what I said.
>

I made an error explaining g3, which is indexed

  g3[ic, jc, id, jd]

I thought this might be better since, it is composed of dim x dim blocks. I
am not opposed to changing
this if there is evidence that another thing is better.

   Matt


> Thanks in advance.
>
> Miguel
>
>
> On Apr 19, 2014 6:19 PM, "Matthew Knepley" <knepley at gmail.com> wrote:
>
>> On Sat, Apr 19, 2014 at 5:25 PM, Miguel Angel Salazar de Troya <
>> salazardetroya at gmail.com> wrote:
>>
>>> Thanks for your response. Now I understand a bit better your
>>> implementation, but I am still confused. Is the g3 function equivalent to
>>> the f_{1,1} function in the matrix in equation (3) in this paper:
>>> http://arxiv.org/pdf/1309.1204v2.pdf ? If so, why would it have terms
>>> that depend on the trial fields? I think this is what is confusing me.
>>>
>>
>> Yes, it is. It has no terms that depend on the trial fields. It is just
>> indexed by the number of components in that field. It is
>> a continuum thing, which has nothing to do with the discretization. That
>> is why it acts pointwise.
>>
>>    Matt
>>
>>
>>> On Sat, Apr 19, 2014 at 11:35 AM, Matthew Knepley <knepley at gmail.com>wrote:
>>>
>>>> On Fri, Apr 18, 2014 at 1:23 PM, Miguel Angel Salazar de Troya <
>>>> salazardetroya at gmail.com> wrote:
>>>>
>>>>> Hello everybody.
>>>>>
>>>>> First, I am taking this example from the petsc-dev version, I am not
>>>>> sure if I should have posted this in another mail-list, if so, my
>>>>> apologies.
>>>>>
>>>>> In this example, for the elasticity case, function g3 is built as:
>>>>>
>>>>> void g3_elas(const PetscScalar u[], const PetscScalar gradU[], const
>>>>> PetscScalar a[], const PetscScalar gradA[], const PetscReal x[],
>>>>> PetscScalar g3[])
>>>>> {
>>>>>   const PetscInt dim   = spatialDim;
>>>>>   const PetscInt Ncomp = spatialDim;
>>>>>   PetscInt       compI, d;
>>>>>
>>>>>   for (compI = 0; compI < Ncomp; ++compI) {
>>>>>     for (d = 0; d < dim; ++d) {
>>>>>       g3[((compI*Ncomp+compI)*dim+d)*dim+d] = 1.0;
>>>>>     }
>>>>>   }
>>>>> }
>>>>>
>>>>> Therefore, a fourth-order tensor is represented as a vector. I was
>>>>> checking the indices for different iterator values, and they do not seem to
>>>>> match the vectorization that I have in mind. For a two dimensional case,
>>>>> the indices for which the value is set as 1 are:
>>>>>
>>>>> compI = 0 , d = 0      ----->     index = 0
>>>>> compI = 0 , d = 1      ----->     index = 3
>>>>> compI = 1 , d = 0      ----->     index = 12
>>>>> compI = 1 , d = 1      ----->     index = 15
>>>>>
>>>>> The values for the first and last seem correct to me, but they other
>>>>> two are confusing me. I see that this elasticity tensor (which is the
>>>>> derivative of the gradient by itself in this case) would be a four by four
>>>>> identity matrix in its matrix representation, so the indices in between
>>>>> would be 5 and 10 instead of 3 and 12, if we put one column on top of each
>>>>> other.
>>>>>
>>>>
>>>> I have read this a few times, but I cannot understand that you are
>>>> asking. The simplest thing I can
>>>> respond is that we are indexing a row-major array, using the indices:
>>>>
>>>>   g3[ic, id, jc, jd]
>>>>
>>>> where ic indexes the components of the trial field, id indexes the
>>>> derivative components,
>>>> jc indexes the basis field components, and jd its derivative components.
>>>>
>>>>    Matt
>>>>
>>>>
>>>>> I guess my question is then, how did you vectorize the fourth order
>>>>> tensor?
>>>>>
>>>>> Thanks in advance
>>>>> Miguel
>>>>>
>>>>> --
>>>>> *Miguel Angel Salazar de Troya*
>>>>> Graduate Research Assistant
>>>>> Department of Mechanical Science and Engineering
>>>>> University of Illinois at Urbana-Champaign
>>>>> (217) 550-2360
>>>>> salaza11 at illinois.edu
>>>>>
>>>>>
>>>>
>>>>
>>>> --
>>>> What most experimenters take for granted before they begin their
>>>> experiments is infinitely more interesting than any results to which their
>>>> experiments lead.
>>>> -- Norbert Wiener
>>>>
>>>
>>>
>>>
>>> --
>>> *Miguel Angel Salazar de Troya*
>>> Graduate Research Assistant
>>> Department of Mechanical Science and Engineering
>>> University of Illinois at Urbana-Champaign
>>> (217) 550-2360
>>> salaza11 at illinois.edu
>>>
>>>
>>
>>
>> --
>> What most experimenters take for granted before they begin their
>> experiments is infinitely more interesting than any results to which their
>> experiments lead.
>> -- Norbert Wiener
>>
>


-- 
What most experimenters take for granted before they begin their
experiments is infinitely more interesting than any results to which their
experiments lead.
-- Norbert Wiener
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20140421/f1df8a20/attachment-0001.html>