[petsc-users] Elasticity tensor in ex52

Mon Apr 21 10:23:25 CDT 2014

Thank you. Now it makes sense. Just to confirm though, back to the weak
form I wrote above in index notation:

phi_{i,j} C_{ijkl} u_{k,l}

The g3 indices are equivalent to the indices of this weak form in this
manner:

ic = i
jc = k
id = j
jd = l

Is this correct? This is what I understand when you talk about the
components and the derivatives of the components.

On Mon, Apr 21, 2014 at 9:48 AM, Matthew Knepley <knepley at gmail.com> wrote:

> 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
>

-- 
*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
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.mcs.anl.gov/pipermail/petsc-users/attachments/20140421/0a96fe1c/attachment.html>