[petsc-users] Elasticity tensor in ex52
Matthew Knepley
knepley at gmail.com
Mon Apr 21 10:24:32 CDT 2014
On Mon, Apr 21, 2014 at 10:23 AM, Miguel Angel Salazar de Troya <
salazardetroya at gmail.com> wrote:
> 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.
>
Yep, that is it.
Thanks,
Matt
> 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
>
>
--
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/512355bf/attachment-0001.html>
More information about the petsc-users
mailing list