from fenics import *
%matplotlib inline
# Create mesh and define function space
mesh = UnitSquareMesh (8 , 8 )
V = FunctionSpace ( mesh , "Lagrange", 1 )
# Define boundary condition
u0 = Function ( V )
bc = DirichletBC (V , u0 , "x[0] < DOLFIN_EPS || x[0] > 1.0 - DOLFIN_EPS")
u = TrialFunction ( V )
v = TestFunction ( V )
f = Expression ("10*exp ( -( pow(x[0] - 0.5, 2) + pow (x[1] - 0.5, 2)) / 0.02)", degree=1 )
g = Expression ("sin(5*x[0])", degree =1 )
a = inner ( grad ( u ) , grad ( v ) )*dx
L = f*v*dx + g*v*ds
plot(project(f, V))
u = Function(V)
M = u*dx
tol = 1.0e-5
solver_parameters = {"error_control":
{"dual_variational_solver":
{"linear_solver": "cg"}}}
solve(a == L, u, bc, tol=tol, M=M,solver_parameters = solver_parameters)
plot(u.root_node() , title =" Solution on initial mesh ")
interactive
plot(u.leaf_node() , title =" Solution on final mesh ")
interactive()
file = File("u.pvd")
file << u.leaf_node()