Dear Alessandra,
The model object in Getfem do not support Petrov-Galerkin methods (it is assumed that u and Test_u are described on the same FEM). So that you will be able to obtain the stiffness matrix of your problem with the assembly functions but will have to prescribe the boundary conditions and solve on your own (the system in a priori not necessarily a square one in your case ...) To obtain a finite element cut by an interface, the simplest is indeed to use a level-set and the object mesh_level_set and mesh_fem_level_set objects as in crack.cc, yes. Then, you can obtain obtain the Petrov-Galerkin stiffness matrix with the high-level assembly with a assembly string as "Grad_u.Grad_Test_v" where u is a variable on the cut fem and v a variable on the uncut one. You can prescribe the jump to be zero on your interface with a multiplier but you have to ensure a LBB condition. Nitsche's method can be a better choice in that situation. Of course, the method as a interest only if the interface do not fit with the element faces. In the conformal case, you end-up with the same space for your primal variable and the test functions. Best regards, Yves ----- Mail original ----- De: "Alessandra Arrigoni" <[email protected]> À: "getfem-users" <[email protected]> Envoyé: Samedi 28 Décembre 2019 17:45:48 Objet: Petrov-Galerkin method for interface problems Dear GetFem++ users, I am trying to solve a simple elliptic interface problem on the unit square: the domain features an interface at x=0.5 that defines two rectangular areas (Omega0 and Omega1) with different diffusion coefficients (a0 and a1). The problem data contain the usual source term, the boundary conditions and two additional interface conditions q0 (on the jump of the solution u0-u1) and q1 (on the jump of its conormal derivative). According to the weak formulation I am given, the functional spaces for the solution u (with restrictions u0 and u1 on the two subdomains) and the test functions are different: for the test functions I am using the classic P1 basis functions on the whole unit square, while for the solution I am using the functions that belong to the P1 space on each subdomain, and discontinuous on the interface. To construct this space I am thinking of simply doubling the degrees of freedom on the interface line and associating one of them to u0 and the other to u1. Both the Dirichlet boundary conditions and the interface condition q0 should be imposed in a direct way, that is, following the idea of the method add_Dirichlet_condition_with_multipliers(md, mim, varname, degree, region, dataname = std::string()) linked to the Dirichlet condition brick; regarding the interface condition, the rows of the matrix should be modified by setting the interface dofs associated to u0 and u1 to 1 and -1 (respectively) and by filling them with zeros. My question is the following: is there a way to construct such a finite element method with the high level procedures contained in GetFem++ (latest version)? At the moment I am trying to manually modify the global number of dofs and the stiffness matrix. In the documentation, I read about the mesh::region and the possibility of defining new bricks (which I might use to enforce the interface conditions as in the aforementioned method for the Dirichlet BC). Alternatively I am thinking of using the level_set to define the interface, but then (even looking at the examples such as tests/crack.cc), I can't understand how to specify different discrete spaces for the solution and the test functions when the bricks are added to the model. Can I exploit one of these ideas (and how) or are there better ways to achieve what I am looking for? As a last remark, this method is expected to be used in the future for the linear elasticity problem, so I would be grateful if you can help me with an implementation that works even in the case of a vectorial problem. Thank you for your time and your support, Kind regards, Alessandra Arrigoni
