getfem-users
[Top][All Lists]
Advanced

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: [Getfem-users] number of Newton-Cotes points


From: Yves Renard
Subject: Re: [Getfem-users] number of Newton-Cotes points
Date: Fri, 18 Jan 2019 09:39:31 +0100
User-agent: Mozilla/5.0 (X11; Linux x86_64; rv:60.0) Gecko/20100101 Thunderbird/60.2.1

Dear Edouard,

Yes, in Newton-Cotes integration methods, some  quadature points of the regular grid correspond to vanishing weights and these points are indeed eliminated. It would not be difficult to have an integration method with the complete grid (this is an option of the add_point fonction of getfem_integration.cc, in fact)  but this is not available for the moment.

Best regards,

Yves






Le 16/01/2019 à 20:19, Edouard Oudet a écrit :
Dear all,

I studied a simple situation : a triangular mesh and defined a MeshIm class based on Newton-Cotes methods of degree "degi". I expected for every triangle to have (degi + 1) * (degi + 2) / 2 quadrature points.

Surprisingly for me, the following illustration below shows that for a degree 4 for instance, I obtain 12 quadrature points which is different from 15.. It seems that some vertices of the triangle are not taken into account.

Is it because the weights of these points is 0.? Or perhaps I didn't call the class in the right way?

Any comments are welcome,

best,

Edouard.

/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////

import getfem
import numpy as np

def nb_imdata(degi):
    # creation mesh / fem / im
    NX = 6
    tmesh = getfem.Mesh('regular_simplices', np.arange(0, 1 + 1. / NX, 1. / NX),                                              np.arange(0, 1 + 1. / NX, 1. / NX))
    tmeshfem = getfem.MeshFem(tmesh, 1)
    tmeshfem.set_fem(getfem.Fem('FEM_PK(2, 2)'))
    tmeshim = getfem.MeshIm(tmesh, getfem.Integ('IM_NC(2,' + str(degi) + ')'))
    #print tmeshim.im_nodes()
    print 'degree ', degi, '|', \
          'nb of quad points ', tmeshim.im_nodes(2)[0].shape[0], '|', \
          'expected          ', ((degi + 1) * (degi + 2) / 2)


for  k in range(3, 10):
    nb_imdata(k)


--

  Yves Renard (address@hidden)       tel : (33) 04.72.43.87.08
  INSA-Lyon
  20, rue Albert Einstein
  69621 Villeurbanne Cedex, FRANCE
  http://math.univ-lyon1.fr/~renard

---------




reply via email to

[Prev in Thread] Current Thread [Next in Thread]