Re: [Getfem-users] Build a transport equation with the generic procedure

[Renard Yves]

Dear Sebastien,

The finite element method is instable for the resolution of convection  
problems. At least without treatment. There is several possible  
treatments (Petrov-Galerkin, SUPG, additional diffusion term ...). Or  
you can use another method. A caracteristic Galerkin method is  
implemented in Getfem, see the command
E = gf_compute(mesh_fem MF, vec U, 'convect', mesh_fem mf_v, vec V,  
scalar dt, int nt[, string option])

But this is a very dissipative method.


I think your term should be

 D = gf.asm_volumic('a=data(mdim(#1),#2); M(#1,#1)+=
 comp(Grad(#1).Base(#1).Base(#2))(:,i,:,k).a(i,k)',mim,mf_u,mf_d,DATA)

with mf_d a scalar fem.


Yves.


sébastien janas <address@hidden> a écrit :

> Hi,
>
> I try to solve a transport equation with getfem.
>
> I have problem to deal with the term
>
> div (  RHO * Y * U )
>
> where
>
> - div is the divergence
> - RHO is the density
> - Y is a moisture content
> - U is the velocity
>
> I try to build the matrix
>
> (PSI_x^i * RHO * U_x * PSI^j ) +
> (PSI_y^i * RHO * U_y * PSI^j ) +
> (PSI_z^i * RHO * U_z * PSI^j )
>
> where - PSI_a^b are the derivative of the base functions at nodes b in
> respect to ditection a
>
>       - PSI^b are the base functions at nodes b
>
>       - U_a are the a components of velocity
>
> In my program, I have the product RHO*U in a (nb_dof x 3) array, and I
> construct my matrix like this
>
> D = gf.asm_volumic('a=data(#2); M(#1,#1)+=
> sym(comp(Grad(#1).Base(#1).Base(#1))(:,i,:,k).a(i,k))',mim,mf_u,mf_d,DATA)
>
> where mf_u is
>
> mf_u  = gf.MeshFem(m,1)
> mf_u.set_fem(gf.Fem("FEM_PK(3,1)"))
>
> and mf_d
>
> mf_d  = gf.MeshFem(m,3)
> md_d.set_fem(gf.Fem("FEM_PK(3,1)"
>
> My questions are these
>
> -> Is it the correct way to solve a transport equation ?
> -> Is anybody have an exemple of such equation ?
> -> With the exemple above, I get the following error that I don't
> understand, is there somebody who can explain it to me ?
>
>  "wrong number of indexes for the 2th argument of the reduction
> comp(Grad(#1).Base(#1).Base(#1))(:,i,:,k).a(i,k) (ranges=[0..12])"
>
> Thanks a lot for your help,
>
> Best regards,
>
> Sébastien Janas
>
>
> _______________________________________________
> Getfem-users mailing list
> address@hidden
> https://mail.gna.org/listinfo/getfem-users