real a = 1; real c = 6*a; real rs = (c/2 + a)*4; border b1(t = 0, rs) {x = t; y = 0;}; border b2(t = 0, pi/2) {x = rs*cos(t); y = rs*sin(t);}; border b3(t = rs, c/2 + a) {x = 0; y = t;}; border b4(t = pi/2, -pi/2) {x = a*cos(t); y = c/2 + a*sin(t);}; border b5(t = c/2 - a, 0) {x = 0; y = t;}; mesh Th = buildmesh(b1(10) + b2(10) + b3(10) + b4(10) + b5(5)); fespace Vh(Th, P2); Vh u, v; problem problem1(u, v) = int2d(Th) ((dx(u)*dx(v) + dy(u)*dy(v))*x) + on(b4, u=1) + on(b2, u=0); for (int i = 0; i < 3; i++) { problem1; Th = adaptmesh(Th, u); }; problem1; plot(u, fill = true, value = true); real capacitance = int2d(Th) ((dx(u)*dx(u) + dy(u)*dy(u))*2*pi*x); capacitance = 1/(1/(2*capacitance) + 1/(4*pi*rs)); cout << "capacitance/(4*pi*a) = " << capacitance/(4*pi*a) << endl;