I want to use the AdaptiveNonlinearVariationalSolver to adapt my mesh automatically, but my problem is nonlinear and in a mixed space. I can solve my problem generally with solve(), but now that I want to put it in the correct format for AdaptiveNonlinearVariationalSolver, it doesn’t work. Currently the error is: ValueError: too many values to unpack (expected 2). I gathered that it is because I define it as a Nonlinear problem (instead of MixedNonlinear), and that’s not correct. However, I can’t figure out the right way to fix it.
This is my MWE:
from dolfin import * from mshr import * import matplotlib.pyplot as plt import numpy as np mesh = RectangleMesh(Point(0, 0), Point(0.2, 1), 8, 8) n = FacetNormal(mesh) # Define Taylor--Hood function space W V = VectorElement("CG", triangle, 2) Q = FiniteElement("CG", triangle, 1) S = FiniteElement("CG", triangle, 1) W = FunctionSpace(mesh, MixedElement([V, Q, S])) # Define Function and TestFunction(s) w = TrialFunction(W) (u, p, T) = split(w) (v, q, s1) = split(TestFunction(W)) # Define the parameters and bcs mu = Constant(1.0) Qc2 = Constant(1.0) # this will be a function of x Ta = Constant(1.0) # this will be a function of x Pe = Constant(27.0) Bi = Constant(11.6) u_in = Constant(-2.0) u_c = Constant(-1.0) inflow = 'near(x, 1.0) && x<=0.1' wall = 'near(x, 0.2)' centre = 'near(x, 0.0)' outflow = 'near(x, 0.0)' bcu_inflow = DirichletBC(W.sub(0), (0.0, u_in), inflow) bcu_wall = DirichletBC(W.sub(0), (0.0, u_c), wall) bcu_outflow = DirichletBC(W.sub(0), (0.0, u_c), outflow) bcu_symmetry = DirichletBC(W.sub(0).sub(0), Constant(0.0), centre) bcT_inflow = DirichletBC(W.sub(2), 0.0, inflow) bcs = [bcu_inflow, bcu_wall, bcu_outflow, bcu_symmetry, bcT_inflow] x = SpatialCoordinate(mesh) def epsilon(v): return sym(as_tensor([[v.dx(0), v.dx(1), 0], [v.dx(0), v.dx(1), 0], [0, 0, 0]])) # stress tensor def sigma(v, p): return 2*mu*epsilon(v)-Id(p) def Id(p): return as_tensor([[p, 0, 0], [0, p, 0], [0, 0, p]]) def cond(v): return sym(as_tensor([[v.dx(0), v.dx(1)], [v.dx(0), v.dx(1)]])) def sigmabc(v, p): return 2*mu*cond(v) - p*Identity(2) facet_f = MeshFunction("size_t", mesh, mesh.topology().dim() - 1) # FACET function CompiledSubDomain('near(x, 1.0) && x<=0.1').mark(facet_f, 0) CompiledSubDomain('near(x, 1.0) && x>=0.1').mark(facet_f, 1) CompiledSubDomain('near(x, 0.2)').mark(facet_f, 2) # wall CompiledSubDomain('near(x, 0.0)').mark(facet_f, 3) # outflow CompiledSubDomain('near(x, 0.0)').mark(facet_f, 4) CompiledSubDomain('near(x, 1.0) && x>0.1 && x<0.2').mark(facet_f, 5) # Create the measure ds = Measure("ds", domain=mesh, subdomain_data=facet_f) # Define the variational form f = Constant((0, -1)) a = (inner(sigma(u, p), epsilon(v)) - div(u) * q + (dot(u, grad(T)) * s1 + ( 1 / Pe) * inner(grad(s1), grad(T)))) * dx() - (1 / Pe) * (-Bi * s1 * T * ds(2)) L = (- dot(f, v) - (1/Pe) * Qc2*s1) * dx() - (1/Pe) * (s1 * Bi * Ta * ds(2)) F = a + L # For nonlinear problem, write F = 0 w = Function(W) # Define goal functional (quantity of interest) M = inner(w, w)*dx() # Define error tolerance tol = 1.e-6 # Solve equation a = L with respect to u and the given boundary # conditions, such that the estimated error (measured in M) is less # than tol problem = NonlinearVariationalProblem(F, w, bcs) # The problem is here #problem = MixedNonlinearVariationalProblem(F, w, bcs) # This doesn't work solver = AdaptiveNonlinearVariationalSolver(problem, M) solver.parameters["error_control"]["dual_variational_solver"]["linear_solver"] = "umfpack" #parameters["refinement_algorithm"] = "plaza_with_parent_facets" solver.solve(tol) solver.summary() # Plot solutions (u0, p0, T0) = w.root_node().split() (u1, p1, T1) = w.leaf_node().split() R = w.leaf_node().function_space() plot(R.mesh())
I think the error is because of the term (dot(u, grad(T)) * s1 , but I don’t know how else to write the variational formulation. Any help would be greatly appreciated.
NOTE: I am aware that my problem for u and p is decoupled from my problem for T, so I could just solve for u and p, and then use that solution to solve for T separately. However, I am working my way up to the fully coupled problem, which is why I have written in that way.