Dear Yves Renard,
Thank you so much for your answers. The alternative to the Ghost penalty is
a good idea to solve with P_k FEMs. For the remark about the slice, I don't
use it to calculate the error, but just to export the solution in the
physical domain. To calculate the error I use exactly the integration
method with MeshLevelSet as follow:
ls = gf.LevelSet(m,2,'(x-{X0})*(x-{X0}) + (y-{Y0})*(y-{Y0}) -
{R}*{R}'.format(X0=x0, Y0=y0, R=r))
mls = gf.MeshLevelSet(m)
mls.add(ls)
mls.adapt()
mim = gf.MeshIm('levelset',mls,'INSIDE', gf.Integ('IM_TRIANGLE(6)'))
mim.set_integ(6)
#After creating the model and solve it, I calculate the error as follows:#
L2error = gf.compute(mfu, Uap-Uex, 'L2 norm', mim)
H1error = gf.compute(mfu, Uap-Uex, 'H1 norm', mim)
where mfu is given by:
mfu = gf.MeshFem(m, 1)
mfu.set_fem(gf.Fem('FEM_PK(2,1)'))
But I still not able to get the optimal rate of convergence (I get
approximately the same error for NX=16, 32, 64,etc. But if I take NX=400
for example, I get an approximate solution very similar to the exact
solution with a very small error where the maximum value=2.1e-06).
Thank you for the help.
Le jeu. 27 déc. 2018 à 09:54, Yves Renard <[email protected]> a
écrit :
>
> Dear Yassine,
>
> Concerning the Ghost penalty, I am not sure that penalization of normal
> derivative is not sufficient even for P2 or P3 method. I think it should be
> sufficient.
>
>
> Note also that an alternative to Ghost penalty is to directly interpolate
> the normal derivative in a neighbour element for elements with a small
> intersection with the domain as we proposed in the paper
>
> J. Haslinger, Y. Renard. A new fictitious domain approach inspired by the
> extended finite element method. Siam J. on Numer. Anal., 47(2):1474--1499,
> 2009.
> (http://math.univ-lyon1.fr/~renard/papers/2009_fictif.pdf)
>
> Best regards,
>
> Yves
>
>
> ----- Original Message -----
> From: "yves renard" <[email protected]>
> To: "Yassine ZAIM" <[email protected]>
> Cc: "getfem-users" <[email protected]>
> Sent: Wednesday, December 26, 2018 10:02:11 PM
> Subject: Re: [Getfem-users] Calculate the error in the physical domain
>
> Dear Yassine,
>
> Once again, the slice operations are only for graphical post-processing.
> You cannot compute an error with them. You have to use cut integration
> methods instead.
>
> Best regards,
>
> Yves
>
> ----- Original Message -----
> From: "Yassine ZAIM" <[email protected]>
> To: "yves renard" <[email protected]>
> Cc: "getfem-users" <[email protected]>
> Sent: Wednesday, December 19, 2018 5:00:01 PM
> Subject: Re: [Getfem-users] Calculate the error in the physical domain
>
> Dear Yves Renard and getfem++ users,
> Thank you so much for your answer. I am inspired by the
> "demo_fictitious_domains.py"
> example to achieve my test and the proposed method is what I did at the
> beginning. But I get the same value of the error for different values of NX
> and I can't get the optimal order. For this, I am saying maybe the outside
> domain of Omega which creates the problem, for this I asked how to restrict
> the computation of error to Omega. I enclosed my code if you can see it
> quickly, maybe you will find out the mistake.
> I have another question about the Ghost penalty. For the P_2 element, I
> guess that I have to work with the hessian "hess_u" and "hess_Test_u" But
> how we can implement it for P_3, P_4,... spaces.
> Thank you so much for your help.
>
> Le mer. 19 déc. 2018 à 14:39, Yves Renard <[email protected]> a
> écrit :
>
> >
> > Dear Yassine,
> >
> > The slice operation is only for graphical post-processing. You cannot use
> > a sliced solution into computations because the obtained vector represent
> > the value of the solution on a cut mesh representing the computed
> > intersection. The obtained vector is no longer a dof vector on the
> meshfem
> > mfu. The standard way to compute this error is to produce a cut
> integration
> > method, and just compute the L2 or H1 norms in a standard way with it.
> You
> > can produce cut integration method with the MeshLevelSet and
> > MeshIM('levelset') objects. An example of use in the python test program
> > "demo_fictitious_domains.py"
> >
> > Best regards,
> >
> > Yves
> >
> >
> > Le 18/12/2018 à 18:54, Yassine ZAIM a écrit :
> >
> > Dear getfem++ users,
> > I am trying to learn how to implement the fictitious domain method for
> the
> > simple problem of Poisson. My domain of interest Omega is a circle of
> > center (*x0=0,y0=0*) and radius *r=0.35*. I made the resolution in the
> > fictitious domain (the square [-0.5, 0.5]^2) and after that, I used the
> > slice to restrict my solution to the physical domain.
> > sl = gf.Slice(('comp',('ball', +1, [x0, y0], r)), m, 5)
> > sl.export_to_vtk('App_Solution.vtk', 'ascii', mfu, Uap)
> > sl.export_to_vtk('Exc_Solution.vtk', 'ascii', mfu, Uex)
> > sl.export_to_vtk('Error-Ex-App.vtk', 'ascii', mfu, Uap-Uex)
> > When I see the approximate and the exact solution in addition to the
> error
> > I can say that I get a good approximate solution. But I am trying to get
> > the optimal order of convergence (h and h^2) like in the papers "J.
> > Haslinger and Y. Renard" or "E. Burman and P. Hansbo". So I interpolate
> the
> > exact and approximate solution in the slice (domain of interest).
> > Ue = gf.compute(mfu, Uex, 'interpolate on', sl)
> > U = gf.compute(mfu, Uap, 'interpolate on', sl)
> > And after that I tried to calculate the error for different values of
> > NX=[16,32,..]
> > L2error = gf.compute(mfu, U-Ue, 'L2 norm', mim)
> > H1error = gf.compute(mfu, U-Ue, 'H1 norm', mim)
> > such that :
> > mim = gf.MeshIm('levelset',mls,'inside', gf.Integ('IM_TRIANGLE(6)'))
> > and mfu.set_fem(gf.Fem('FEM_PK(2,1)'))
> > But I get the following error in the "compute" function (of the errors) :
> > return getfem('compute', mf, U, what, *args)
> > RuntimeError: (Getfem::InterfaceError) -- The trailing dimension of
> > argument 2 (an array of size 4670) has 4670 elements, 289 were expected.
> > My question is how to calculate the error just in the interesting domain
> > (physical domain, slice sl in my case).
> > Thank you in advance.
> >
> >
> > --
> > *ZAIM Yassine *
> > *PhD in Applied Mathematics*
> >
> >
> > --
> >
> > Yves Renard ([email protected]) tel : (33) 04.72.43.87.08
> > INSA-Lyon
> > 20, rue Albert Einstein
> > 69621 Villeurbanne Cedex, FRANCE
> > http://math.univ-lyon1.fr/~renard
> >
> > ---------
> >
> >
>
> --
> *ZAIM Yassine *
> *PhD in Applied Mathematics*
>
--
*ZAIM Yassine *
*PhD in Applied Mathematics*
