Dear Yves, Thank you very much for your reply. I checked with what you proposed, I got the same result. I still do not know how J(u)-J(uh) is different from J(u-uh) at certain points whilst J is linear.
Best regards, Phuoc On Wed, Apr 18, 2018 at 5:28 PM, Yves Renard <[email protected]> wrote: > > Dear Phuoc, > > If you just want to compute J(u) = \int_{\Omega} div(u) dx, then I would > say that the more straigthforward computation is > > gf.asm_generic(mim,0,"Div_u",OMEGA > ,"u",False,mfu,md.variable('u') > ,"t",True,mfef, np.zeros(mfef.nbdof())) > > However, if mfef is a Lagrange finite element, what you wrote will also do > the job, may be except that it sums all the components of course. > > Best regards, > > Yves. > > ----- Original Message ----- > From: "Huu Phuoc BUI" <[email protected]> > To: [email protected] > Sent: Wednesday, April 18, 2018 1:03:59 PM > Subject: [Getfem-users] Generic assembly of GetFEM++ > > Dear GetFEM++ users, > > I am working on adaptive refinement of a linear elasticity problem using a > posteriori error estimate. > > At each adaptive refinement step, I need to compute some quantity of > interest defined on a subdomain $\Omega$. Let's call this quantity of > interest J(u) = \int_{\Omega} div(u) dx. > > I use 'asm_generic' of GetFEM++ to compute this quantity. The python code > looks like: > > QoI = "(Div_u)*Test_t" > QoI_asm = gf.asm_generic(mim,1,QoI,OMEGA > ,"u",False,mfu,md.variable('u') > ,"t",True,mfef, np.zeros(mfef.nbdof())) > QoI_asm_elem = QoI_asm [ QoI_asm.size - mfef.nbdof() : QoI_asm.size ] > qoi = abs(np.sum(QoI_asm_elem)) > > with > mfef = gf.MeshFem(m,1) > mfef.set_fem(gf.Fem('FEM_PK_DISCONTINUOUS(2,{d})'.format(d=0))) > > I compute then the relative error of the quantity of interest > (J(u)-J(u_h))/J(u), with u is the solution computed from very fine mesh, > u_h is the solution of the adaptive mesh. > > What I got is that (J(u)-J(u_h))/J(u) does not converge well under mesh > refinement. Secondly, J(u)-J(u_h) differs from J(u-u_h) for adaptive > refinement case, which is not acceptable since J is linear. For the uniform > refinement case, they are however identical. > > > > > I checked that the region of interest OMEGA is refined, and is updated > correctly at each refinement step. > > I do not know where the problem comes from. Do you think the generic > assembly code I wrote is correct? Any hint would be very helpful and > appreciated. > > Best, > Phuoc >
