freefem++-doc 3.19.1-1 / usr / share / doc / freefem++ / examples / examples++-load / IpoptVI.edp

//  Solver - Delta u = f , u < g et u = 0 on Gamma 
load "ff-Ipopt";
int nn=20;
mesh Th=square(nn,nn);
fespace Vh(Th,P2);

func f = 1.; //rhs function
real r=0.03,s=0.1; //some parameters for g
func g = r - r/2*exp(-0.5*(square(x-0.5)+square(y-0.5))/square(s));


macro Grad(u) [dx(u),dy(u)]//
varf vP(u,v) = int2d(Th)(Grad(u)'*Grad(v)) - int2d(Th)(f*v);

matrix A = vP(Vh,Vh,solver=CG);
real[int] b = vP(0,Vh);

// Warning the boundary condition are given with lb and ub on border
// equal value ... 
Vh u=0;
Vh lb=-1.e19;
varf vGamma(u,v) = on(1,2,3,4,u=1);
real[int] onGamma=vGamma(0,Vh);
Vh ub=g;
ub[] = onGamma ? 0. : ub[];
lb[] = onGamma ? 0. : lb[];

IPOPT([A,b],u[],lb=lb[],ub=ub[]);

plot(u,wait=1);