interesting, if you just add a
print(mim)
before calling the second assembly, then it works correctly. This is a bug,
somewhere we have forgotten to update mesh_im when mesh is changed.
On Thu, Jul 21, 2022 at 11:02 AM Lorenzo Ferro <lorenzo.i...@gmail.com>
wrote:
> Hello Konstantinos,
> I've realized for this troubleshooting a simplified example, where I
> generate a 1D mesh corresponding to the radius of a circular section.
> Then I've added the mim and I've calculated the corresponding area both
> analytically and numerically with asm_generic considering the axisymmetric
> configuration.
> I removed convexes and points from the mesh with the del_convex and
> del_point methods and then I generated other elements on the same mesh
> object.
> The numerical area calculation turns out to be wrong.
>
> ####### CODE START #######
> import getfem as gf
> import numpy as np
>
> RoundRadius=0.01
> NumberOfElements = 5
> nodes = np.arange(0, RoundRadius+RoundRadius/NumberOfElements,
> RoundRadius/NumberOfElements)
> mesh = gf.Mesh("cartesian", nodes)
>
> im_degree = 3
> mim = gf.MeshIm(mesh, im_degree)
>
> RoundSection = np.pi*R**2
> RoundSection_asm = gf.asm_generic(mim, 0, '2*pi*X' , -1)
> print ("STARTING AREA CALCULATION")
> print ("analitical section =", RoundSection, "VS numerical section = ",
> RoundSection_asm)
>
> # Change completely the mesh
> mesh.del_convex(mesh.cvid())
> mesh.del_point(mesh.pid())
>
> NumberOfElements = 10
> nodes = np.arange(0, RoundRadius+RoundRadius/NumberOfElements,
> RoundRadius/NumberOfElements)
> mesh.add_point([nodes])
> for x in range(NumberOfElements):
> mesh.add_convex(gf.GeoTrans('GT_QK(1,1)'), [[nodes[x], nodes[x+1]]])
> RoundSection_asm = gf.asm_generic(mim, 0, '2*pi*X' , -1)
> print ("AREA CALCULATION AFTER MESH CHANGE")
> print ("analytical section =", RoundSection, "VS numerical section = ",
> RoundSection_asm)
>
> # Let's try to apply another mesh method
> mesh.set_pts([nodes])
> RoundSection_asm = gf.asm_generic(mim, 0, '2*pi*X' , -1)
> print ("AREA CALCULATION AFTER MESH.SET_PTS METHOD")
> print ("analytical section =", RoundSection, "VS numerical section = ",
> RoundSection_asm)
>
> ####### CODE FINISH#######
>
> Thank you,
> Lorenzo
>
> Il giorno mer 20 lug 2022 alle ore 11:43 Konstantinos Poulios <
> logar...@googlemail.com> ha scritto:
>
>> Dear Lorenzo,
>>
>> It would be very useful if you could provide a minimal test case that
>> demonstrates this issue. Strange that md.set_time() helps updating the
>> model correctly.
>>
>> Best regards
>> Kostas
>>
>> On Mon, Jul 18, 2022 at 4:36 PM Lorenzo Ferro <lorenzo.i...@gmail.com>
>> wrote:
>>
>>> Dear All,
>>>
>>> I'm experimenting with mesh changes (a sequence like del_convex ->
>>> del_point -> add_point -> add_convex -> set_region) and relevant
>>> MeshFEM automatic adaptation with the Python interface.
>>>
>>> I found out that as soon as I update mesh and regions, if after that I
>>> try to run a test assembly method for example just to calculate the total
>>> area of the mesh elements (via gf.asm_generic) I achieve a completely wrong
>>> result,
>>> but,
>>> after just running a Model method like "md.set_time(0.)" and then
>>> running immediately "gf.asm_generic" then I achieve the correct result!
>>>
>>> Is this a bug? or is there some logic behind this behaviour?
>>>
>>> Thank you,
>>> Lorenzo
>>>
>>
