# How to use Lagrange Multipliers in FreeFEM **URL:** \ **Category:** General Discussion\ **Created:** [July 22, 2020, 7:11pm UTC](https://community.freefem.org/t/how-to-use-lagrange-multipliers-in-freefem/514 "2020-07-22T19:11:18Z")\ **Posts on this page:** 9\ **Page:** 1
**Author:** ![kevinmo08](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/kevinmo08/32/217_2.png) [@kevinmo08](https://community.freefem.org/u/kevinmo08)\ **Post date:** [July 22, 2020, 7:11pm UTC](https://community.freefem.org/t/how-to-use-lagrange-multipliers-in-freefem/514/1 "2020-07-22T19:11:19Z")
Hi all, I am very new to the finite element field, and I want to try to solve the a variant of the fokker-plank equation with periodic boundary conditions. To do this, I need to find a way to impose an actual concentration in the system. Since I have a Robin boundary condition (no flux) at the surface of an inner wall, but periodic everywhere else, the problem does not have a unique solution. I’d like to impose some type of constraint that makes the integral over all space equal a constant or something like that to impose a unique solution. From my reading, it seems like the way to do that is through a lagrange multiplier, but I am not sure how to implement that in FreeFEM++. For Example, here is how I would write and solve my equations when I had a neumann and a dirichlet boundary conditions and so it was unique (I still was periodic on just the z direction.) ``` fespace Vh(Th3 , P1, periodic=[[5,x,y],[6,x,y]]); Vh u, v // probability and test function. problem FPeqn(prob3d,testfunc3d) = int3d(Th3)( "bilinear term") + int3d(Th3)("linear term") + int2d(Th3,innerlabel)("robin boundary term") // Robin BC + on(outerlabel, v = 1/(2*pi)) // Dirichlet BC ; FPeqn ``` But from the brief example they have on the documentation [https://doc.freefem.org/examples/finite-element.html#examplelagrangemultipliers](https://doc.freefem.org/examples/finite-element.html#examplelagrangemultipliers) , the way that they solve an equation with Lagrange Multipliers is very different than what I am used to. They use "varf"s instead of “problem” and they also make a block matrix. I guess my question is, what is this block of code that they have. And how would I go about building that for a 3d system? ``` // Parameters func f = 1 + x - y; // Mesh mesh Th = square(10, 10); // Fespace fespace Vh(Th, P1); int n = Vh.ndof; int n1 = n+1; Vh uh, vh; // Problem varf va (uh, vh) = int2d(Th)( dx(uh)*dx(vh) + dy(uh)*dy(vh) ) ; varf vL (uh, vh) = int2d(Th)(f*vh); varf vb (uh, vh) = int2d(Th)(1.*vh); matrix A = va(Vh, Vh); real[int] b = vL(0, Vh); real[int] B = vb(0, Vh); // Block matrix matrix AA = [[ A, B], [B', 0] ]; set(AA, solver=sparsesolver); real[int] bb(n+1), xx(n+1), b1(1), l(1); b1 = 0; // Builds the right hand side block bb = [b, b1]; // Solve xx = AA^-1 * bb; // Set values [uh[],l] = xx; // Display cout << " l = " << l(0) << " , b(u, 1) =" << B'*uh[] << endl; ``` ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 22, 2020, 7:15pm UTC](https://community.freefem.org/t/how-to-use-lagrange-multipliers-in-freefem/514/2 "2020-07-22T19:15:52Z")
Maybe the documentation from [FEniCS](https://fenicsproject.org/docs/dolfin/1.6.0/python/demo/documented/neumann-poisson/python/documentation.html) is cleaner? It solves exactly the same problem, so the variational formulations are the same. If you want to switch to 3D, switch `mesh` to `mesh3` and `int2d` to `int3d`. ---
**Author:** ![kevinmo08](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/kevinmo08/32/217_2.png) [@kevinmo08](https://community.freefem.org/u/kevinmo08)\ **Post date:** [July 22, 2020, 8:05pm UTC](https://community.freefem.org/t/how-to-use-lagrange-multipliers-in-freefem/514/3 "2020-07-22T20:05:49Z")
That does help a bit thank you. It seems like the section ```auto matrix A = va(Vh, Vh); real[int] b = vL(0, Vh); real[int] B = vb(0, Vh); // Block matrix matrix AA = [[ A, B], [B', 0] ]; set(AA, solver=sparsesolver); real[int] bb(n+1), xx(n+1), b1(1), l(1); b1 = 0; // Builds the right hand side block bb = [b, b1]; // Solve xx = AA^-1 * bb; // Set values [uh[],l] = xx; ``` Is unique to FreeFEM. Is this block of code something that needs to be changed or is this just how I would extract the solutions to it? ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 22, 2020, 8:11pm UTC](https://community.freefem.org/t/how-to-use-lagrange-multipliers-in-freefem/514/4 "2020-07-22T20:11:46Z")
I don’t think you need to change that for your application. The formulation in FreeFEM is really similar to what you would do in MATLAB, by defining a block matrix and concatenating vertically vectors or extracting scalars. ---
**Author:** ![kevinmo08](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/kevinmo08/32/217_2.png) [@kevinmo08](https://community.freefem.org/u/kevinmo08)\ **Post date:** [July 28, 2020, 5:03pm UTC](https://community.freefem.org/t/how-to-use-lagrange-multipliers-in-freefem/514/5 "2020-07-28T17:03:00Z")
Thank you for your help with this! After tinkering with it for a few days, I have another question if it isn’t too much of a bother. I am trying to figure out what these sections are doing. ``` varf vL (uh, vh) = int2d(Th)(f*vh); varf vb (uh, vh) = int2d(Th)(1.*vh); ``` I am under the impression here that f is the non homogeneous term in the original equation. If my equation is homogeneous, then can I set f = 0, or can I just delete this line with vL? (if so, what do I do with b = vL(0, Vh); and bb = [b, b1] later on?) Also, what is the vb term representing? Is that what the solution integrates to? For example, if my constraint is that int3d(u) = Constant, would I put that constant in the vb term above or would I put it something below like b1 = Constant? ``` matrix A = va(Vh, Vh); real[int] b = vL(0, Vh); real[int] B = vb(0, Vh); // Block matrix matrix AA = [[ A, B], [B', 0] ]; set(AA, solver=sparsesolver); real[int] bb(n+1), xx(n+1), b1(1), l(1); b1 = 0; ``` I’m thinking I won’t have to change the block matrix stuff above, except maybe the b1 thing? ``` // Builds the right hand side block bb = [b, b1]; // Solve xx = AA^-1 * bb; // Set values [uh[],l] = xx; ``` ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 28, 2020, 7:52pm UTC](https://community.freefem.org/t/how-to-use-lagrange-multipliers-in-freefem/514/6 "2020-07-28T19:52:45Z")
`f` is the source term, `vb` discretizes the constraint `int(u) = b1`. ---
**Author:** ![kevinmo08](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/kevinmo08/32/217_2.png) [@kevinmo08](https://community.freefem.org/u/kevinmo08)\ **Post date:** [July 28, 2020, 8:42pm UTC](https://community.freefem.org/t/how-to-use-lagrange-multipliers-in-freefem/514/7 "2020-07-28T20:42:43Z")
So if there is no source term, I can get away with f = 0? Thank you! ---
**Author:** ![prj](https://avatars.discourse-cdn.com/v4/letter/p/ecae2f/32.png) [@prj](https://community.freefem.org/u/prj)\ **Post date:** [July 29, 2020, 6:02am UTC](https://community.freefem.org/t/how-to-use-lagrange-multipliers-in-freefem/514/8 "2020-07-29T06:02:29Z")
Why couldn’t you? Try it and you’ll see for yourself… ---
**Author:** ![frederichecht](https://yyz2.discourse-cdn.com/flex030/user_avatar/community.freefem.org/frederichecht/32/15_2.png) [@frederichecht](https://community.freefem.org/u/frederichecht)\ **Post date:** [January 11, 2021, 9:15pm UTC](https://community.freefem.org/t/how-to-use-lagrange-multipliers-in-freefem/514/9 "2021-01-11T21:15:25Z")