Dear Kostas, I used ' mumps ' as a linear solver. I also tried to use ' cg/ildlt ' but it did not succeed to solve the system , with some warnings about the pivot is too small . What do you think? Thank you!
Best regards, Phuoc On Thu, Apr 19, 2018, 20:20 Konstantinos Poulios <[email protected]> wrote: > Dear Phuoc > > Which linear solver are you using? > > Best regards > Kostas > > On Thu, Apr 19, 2018 at 3:30 PM, Huu Phuoc BUI <[email protected] > > wrote: > >> 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 >>> >> >> >
