getfem-users
[Top][All Lists]

Re: [Getfem-users] High level generic assembly procedures

 From: Marco Pischedda Subject: Re: [Getfem-users] High level generic assembly procedures Date: Fri, 11 Apr 2014 09:28:22 +0200

```Dear Yves,

thank you for your answer.  I understand how to impose the boundary condition.
I have another question: in my formulation I'm assembling a term
A(rho)*f.Test_f, where A(rho) is a symmetric matrix defined as
follows:

" [-sqr(rho(2))-sqr(rho(3)), rho(1)*rho(2), rho(1)*rho(3);"
" rho(1)*rho(2), -sqr(rho(1))-sqr(rho(3)), rho(2)*rho(3);"
" rho(1)*rho(3), rho(2)*rho(3) , -sqr(rho(1))-sqr(rho(2))]";

where rho and f are two FemVariables. When I'm assembling this term I
receive this warning message:  "WARNING: detected wrong equivalent
nodes".
What does it mean?

Then I have an expression "sin(phi)/phi", how can I insert an "if
condition" (i.e if phi<1e-8) in the expression when "phi" tends to
zero?

Marco

2014-04-08 15:29 GMT+02:00 Yves Renard <address@hidden>:
>
>
> Dear Marco,
>
> valid. The expression of an assembly string is to be evaluated at each
> Gauss point of the concerned region. An expression such as "p(0)" will
> be understand to be the component 0 of p which is not valid because the
> first component is 1.
> If you want to add some expression on a boundary (an extremity of the
> intervall [0, L] here) then it is possible but not in the same assembly
> string. In your example you have to add the assembly string "p.u" at the
> boundary x=0 and "-p.u" at the boundary x=L.
>
> However, if you need to prescribe a Dirichlet boundary condition, the
> best would be to use a corresponding model brick (if you use the model
> system).
>
>
> Yves.
>
>
>
>
>
> Le 07/04/2014 15:43, Marco Pischedda a écrit :
>> Hi,
>>
>> how can I impose the boundary condition with high level generic
>> assembly? For example I want to assemble "Grad_u.Test_p +
>> p.Grad_Test_u+p.Test_u" and then imposing the value of p on the
>> boundaries of the 1d domain, i.e x=0 and x=L. It is possible
>> i.e can I insert the boundaries conditions in the assembling string?
>>
>>
>> Marco
>>
>>
>> 2014-04-03 9:52 GMT+02:00 Marco Pischedda <address@hidden>:
>>> Dear Yves,
>>>
>>> thank you for bug correction, now it works.
>>>
>>>> Concerning your question on Grad_u for a vector field, it is usually an
>>>> order two  tensor (a matrix) except in 1D.
>>>   Ok that's good.
>>>
>>> I will let you know if there are other problems.
>>>
>>> Thank you
>>>
>>> Marco
>>>
>>> 2014-04-03 9:25 GMT+02:00 Yves Renard <address@hidden>:
>>>>
>>>> Dear Marco,
>>>>
>>>> Concerning your question on Grad_u for a vector field, it is usually an
>>>> order two  tensor (a matrix) except in 1D.
>>>> I tried to to the best to make the most operations "dimension
>>>> indepedent" so there is some permitivity in the langage (components of
>>>> size 1 are sometimes ignored).
>>>> Remember also that you can have the expression of any term with the
>>>> commant Print. For instance "Print( Grad_u).Test_p" will print the
>>>> gradient of u on each Gauss point. It does a lot of print, but at least,
>>>> you can see the format of the term.
>>>>
>>>> Yves.
>>>>
>>>>
>>>> Le 02/04/2014 17:21, Marco Pischedda a écrit :
>>>>> Hi,
>>>>>
>>>>> I have other questions:
>>>>>
>>>>> - I tried to assemble separately the following terms:
>>>>> correctly.
>>>>>   Then I want to assemble the sum of this terms, i.e: "Grad_u.Test_p +
>>>>>   but I receive the following error:
>>>>>
>>>>>    Addition or substraction of incompatible expressions or of different
>>>>> sizes
>>>>>    terminate called after throwing an instance of 'gmm::gmm_error'
>>>>>    what():  Error in getfem_generic_assembly.cc, line 3591 :
>>>>>    Error in assembly string
>>>>>    Aborted
>>>>>
>>>>> - I'm working in 1d problem with vectors as unknowns. Grad_u is
>>>>> therefore a vector or is a        tensor? When I do Grad_u.Test_p the
>>>>> result is a scalar or a vector? Test_p is a vector or a scalar?
>>>>>
>>>>>
>>>>> Marco
>>>>>
>>>>>
>>>>>
>>>>> 2014-04-01 17:51 GMT+02:00 Marco Pischedda <address@hidden>:
>>>>>> Ok thank you,
>>>>>>
>>>>>> I tried it and it works. I will let you know if there are other problems.
>>>>>>
>>>>>> Thanks
>>>>>>
>>>>>> Marco
>>>>>>
>>>>>> 2014-04-01 17:16 GMT+02:00 Yves Renard <address@hidden>:
>>>>>>> Dear Marco,
>>>>>>>
>>>>>>> All seems to me correct in your implementation. This is probably just
>>>>>>> the test line 2074 of getfem_generic_assembly.cc
>>>>>>> which is not correct for one-dimensionnal problems.
>>>>>>>
>>>>>>> I think the test should be
>>>>>>>
>>>>>>>  GA_DEBUG_ASSERT((qdim == 1 && t.sizes()[0] == N) ||
>>>>>>>                       (t.sizes()[1] == N && t.sizes()[0] == qdim) ||
>>>>>>>                       (N == 1 && t.sizes()[0] == qdim),
>>>>>>>                       "dimensions mismatch");
>>>>>>>
>>>>>>> May be you can try this. I will validate it if it works.
>>>>>>>
>>>>>>> Regards,
>>>>>>>
>>>>>>> Yves.
>>>>>>>
>>>>>>>
>>>>>>>
>>>>>>> Le 01/04/2014 16:56, Marco Pischedda a écrit :
>>>>>>>> Dear Yves,
>>>>>>>>
>>>>>>>>
>>>>>>>> I have another question:
>>>>>>>>
>>>>>>>> - I have a monodimensional problem but the unknowns are vectorial. I
>>>>>>>> set the vectorial
>>>>>>>> dimension of the unknows with mf_u.set_qdim(3) and with
>>>>>>>> mf_p.set_qdim(3).
>>>>>>>> When I define the expression " Grad_u.Test_p "  i receive the
>>>>>>>> following error:
>>>>>>>>
>>>>>>>> terminate called after throwing an instance of 'gmm::gmm_error'
>>>>>>>>   what():  Error in
>>>>>>>> ../../getfem-svn/getfem/src/getfem_generic_assembly.cc, line 2074 :
>>>>>>>> dimensions mismatch
>>>>>>>>
>>>>>>>> Grad_u should be a 2nd order tensor, while Test_p should be a vector.
>>>>>>>> I suppose the
>>>>>>>> problem is that Grad_u is interpreted as a vector while Test_p is
>>>>>>>> interpreted as a scalar.
>>>>>>>> How can I use the high level generic assembly procedures on vectorial
>>>>>>>> problems
>>>>>>>> defined on a monodimensional computational domain?
>>>>>>>>
>>>>>>>>
>>>>>>>> Marco
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>>
>>>>>>>> 2014-04-01 8:42 GMT+02:00 Yves Renard <address@hidden>:
>>>>>>>>> Dear Marco,
>>>>>>>>>
>>>>>>>>> Yes, your code is correct and it should work correctly.
>>>>>>>>> Note that you can also use the high level generic assembly with the
>>>>>>>>> model object using a generic assembly brick.
>>>>>>>>> Note also that high level assembly is only available in the svn
>>>>>>>>> repository version of Getfem
>>>>>>>>> (and a very recent version, at least r4570 because a bug have been
>>>>>>>>> corrected on coupled problems recently).
>>>>>>>>>
>>>>>>>>> Yves.
>>>>>>>>>
>>>>>>>>>
>>>>>>>>> Le 31/03/2014 14:00, Marco Pischedda a écrit :
>>>>>>>>>> Dear all,
>>>>>>>>>>
>>>>>>>>>> I want to use the high level generic assembly procedures in order to
>>>>>>>>>> use the automatic differentation for calculating the Jacobian of a
>>>>>>>>>> non-linear problem.
>>>>>>>>>>
>>>>>>>>>> For example I want to calculate the Jacobian of two vectorial
>>>>>>>>>> equations in the unknowns "u" and "p". The first vectorial equation
>>>>>>>>>> is
>>>>>>>>>> (u+p).Test_u while the second is (u+p).Test_p.
>>>>>>>>>>
>>>>>>>>>>  It is correct the following code?
>>>>>>>>>>
>>>>>>>>>>  gmm::sub_interval Iu(0, ndofs_u);
>>>>>>>>>>  gmm::sub_interval Ip(ndofs_u,ndofs_p);
>>>>>>>>>>  getfem::model_real_sparse_matrix Jac(ndofs_u+ndofs_p,
>>>>>>>>>> ndofs_u+ndofs_p);
>>>>>>>>>>  workspace.set_assembled_matrix(Jac);
>>>>>>>>>>  workspace.assembly(2);
>>>>>>>>>>
>>>>>>>>>> In the matrix Jac I have the Jacobian of the system?
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>>
>>>>>>>>>> Marco
>>>>>>>>>>
>>>>>>>>>> _______________________________________________
>>>>>>>>>> Getfem-users mailing list
>>>>>>>>>> https://mail.gna.org/listinfo/getfem-users
>>>>>>>>> --
>>>>>>>>>
>>>>>>>>>   Yves Renard (address@hidden)       tel : (33) 04.72.43.87.08
>>>>>>>>>   Pole de Mathematiques, INSA-Lyon             fax : (33)
>>>>>>>>> 04.72.43.85.29
>>>>>>>>>   20, rue Albert Einstein
>>>>>>>>>   69621 Villeurbanne Cedex, FRANCE
>>>>>>>>>   http://math.univ-lyon1.fr/~renard
>>>>>>>>>
>>>>>>>>> ---------
>>>>>>>>>
>>>>>>> --
>>>>>>>
>>>>>>>   Yves Renard (address@hidden)       tel : (33) 04.72.43.87.08
>>>>>>>   Pole de Mathematiques, INSA-Lyon             fax : (33) 04.72.43.85.29
>>>>>>>   20, rue Albert Einstein
>>>>>>>   69621 Villeurbanne Cedex, FRANCE
>>>>>>>   http://math.univ-lyon1.fr/~renard
>>>>>>>
>>>>>>> ---------
>>>>>>>
>>>>
>>>> --
>>>>
>>>>   Yves Renard (address@hidden)       tel : (33) 04.72.43.87.08
>>>>   Pole de Mathematiques, INSA-Lyon             fax : (33) 04.72.43.85.29
>>>>   20, rue Albert Einstein
>>>>   69621 Villeurbanne Cedex, FRANCE
>>>>   http://math.univ-lyon1.fr/~renard
>>>>
>>>> ---------
>>>>
>
>
> --
>
>   Yves Renard (address@hidden)       tel : (33) 04.72.43.87.08
>   Pole de Mathematiques, INSA-Lyon             fax : (33) 04.72.43.85.29
>   20, rue Albert Einstein
>   69621 Villeurbanne Cedex, FRANCE
>   http://math.univ-lyon1.fr/~renard
>
> ---------
>

```