{"action":"created","issue":{"url":"https://api.github.com/repos/firedrakeproject/firedrake/issues/1266","repository_url":"https://api.github.com/repos/firedrakeproject/firedrake","labels_url":"https://api.github.com/repos/firedrakeproject/firedrake/issues/1266/labels{/name}","comments_url":"https://api.github.com/repos/firedrakeproject/firedrake/issues/1266/comments","events_url":"https://api.github.com/repos/firedrakeproject/firedrake/issues/1266/events","html_url":"https://github.com/firedrakeproject/firedrake/issues/1266","id":350886612,"node_id":"MDU6SXNzdWUzNTA4ODY2MTI=","number":1266,"title":"Strange bug when solving stokes","user":{"login":"sebastkm","id":13945983,"node_id":"MDQ6VXNlcjEzOTQ1OTgz","avatar_url":"https://avatars2.githubusercontent.com/u/13945983?v=4","gravatar_id":"","url":"https://api.github.com/users/sebastkm","html_url":"https://github.com/sebastkm","followers_url":"https://api.github.com/users/sebastkm/followers","following_url":"https://api.github.com/users/sebastkm/following{/other_user}","gists_url":"https://api.github.com/users/sebastkm/gists{/gist_id}","starred_url":"https://api.github.com/users/sebastkm/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/sebastkm/subscriptions","organizations_url":"https://api.github.com/users/sebastkm/orgs","repos_url":"https://api.github.com/users/sebastkm/repos","events_url":"https://api.github.com/users/sebastkm/events{/privacy}","received_events_url":"https://api.github.com/users/sebastkm/received_events","type":"User","site_admin":false},"labels":[],"state":"open","locked":false,"assignee":null,"assignees":[],"milestone":null,"comments":11,"created_at":"2018-08-15T16:59:53Z","updated_at":"2018-08-17T15:50:00Z","closed_at":null,"author_association":"NONE","body":"So I was doing a taylor remainders test when I stumbled upon this strange result.\r\n\r\nRunning this:\r\n```\r\nfrom firedrake import *\r\n\r\nmu = Constant(1.0)\r\nalphaunderbar = 2.5 * mu / (100 ** 2)\r\nalphabar = 2.5 * mu / (0.01 ** 2)\r\nq = Constant(0.01)\r\n\r\n\r\ndef alpha(rho):\r\n    \"\"\"Inverse permeability as a function of rho, equation (40)\"\"\"\r\n    return alphabar + (alphaunderbar - alphabar) * rho * (1 + q) / (rho + q)\r\n\r\nN = 10\r\nmesh = RectangleMesh(N, N, 1.5, 1)\r\nA = FunctionSpace(mesh, \"CG\", 1)\r\n\r\nU_h = VectorElement(\"CG\", mesh.ufl_cell(), 2)\r\nP_h = FiniteElement(\"CG\", mesh.ufl_cell(), 1)\r\nW = FunctionSpace(mesh, U_h*P_h)\r\n\r\n\r\ndef forward(rho):\r\n    \"\"\"Solve the forward problem for a given fluid distribution rho(x).\"\"\"\r\n    w = Function(W)\r\n    (u, p) = split(w)\r\n    (v, q) = TestFunctions(W)\r\n\r\n    F = (alpha(rho) * inner(u, v) * dx + inner(grad(u), grad(v)) * dx +\r\n         inner(grad(p), v) * dx  + inner(div(u), q) * dx)\r\n    bcs = [DirichletBC(W.sub(0).sub(1), 0, \"on_boundary\"),\r\n           DirichletBC(W.sub(0).sub(0), 1, \"on_boundary\")]\r\n\r\n    sp = {\"snes_type\": \"newtonls\",\r\n          \"snes_monitor\": None,\r\n          \"snes_linesearch_type\": \"basic\",\r\n          \"snes_linesearch_monitor\": None,\r\n          \"ksp_type\": \"preonly\",\r\n          \"pc_type\": \"lu\",\r\n          \"pc_factor_mat_solver_type\": \"mumps\",\r\n          \"mat_type\": \"aij\"}\r\n    solve(F == 0, w, bcs=bcs, solver_parameters=sp)\r\n    return w\r\n\r\n\r\nif __name__ == \"__main__\":\r\n    rho_base = interpolate(Constant(1.0/3), A)\r\n    rho_bugged = Function(A)\r\n    rho_working = Function(A)\r\n\r\n    # Bugged for {0, 1e-17}\r\n    rho_bugged.vector()[:] = rho_base.vector()[:] + 0.01/2**3 + 1e-17\r\n    # Working for {1e-16, 1}\r\n    rho_working.vector()[:] = rho_base.vector()[:] + 0.01/2**3 + 1e-16\r\n    w   = forward(rho_bugged)\r\n    (u, p) = split(w)\r\n    J = assemble(inner(u, u) * dx)\r\n    print(\"Bugged value: \", J)\r\n    w = forward(rho_working)\r\n    (u, p) = split(w)\r\n    J = assemble(inner(u, u) * dx)\r\n    print(\"Working value:\", J)\r\n```\r\nThe output is\r\n```\r\n  0 SNES Function norm 1.205469803964e+01\r\n  1 SNES Function norm 4.467216974121e+03\r\n  2 SNES Function norm 3.179445678848e+03\r\n  3 SNES Function norm 9.371457594904e+03\r\n  4 SNES Function norm 9.857177172648e+03\r\n  5 SNES Function norm 2.296564448810e+04\r\n  6 SNES Function norm 2.567023797634e+04\r\n  7 SNES Function norm 4.503989424028e+04\r\n  8 SNES Function norm 6.718281497186e+04\r\n  9 SNES Function norm 7.012727637067e+04\r\n 10 SNES Function norm 7.665161787969e+04\r\n 11 SNES Function norm 9.129369446643e+04\r\n 12 SNES Function norm 1.378040276622e+05\r\nBugged value:  2456322.5129464837\r\n  0 SNES Function norm 1.205469803964e+01\r\n  1 SNES Function norm 4.122870303636e-12\r\nWorking value: 1.499999999999998\r\n```\r\n"},"comment":{"url":"https://api.github.com/repos/firedrakeproject/firedrake/issues/comments/413908192","html_url":"https://github.com/firedrakeproject/firedrake/issues/1266#issuecomment-413908192","issue_url":"https://api.github.com/repos/firedrakeproject/firedrake/issues/1266","id":413908192,"node_id":"MDEyOklzc3VlQ29tbWVudDQxMzkwODE5Mg==","user":{"login":"wence-","id":1126981,"node_id":"MDQ6VXNlcjExMjY5ODE=","avatar_url":"https://avatars2.githubusercontent.com/u/1126981?v=4","gravatar_id":"","url":"https://api.github.com/users/wence-","html_url":"https://github.com/wence-","followers_url":"https://api.github.com/users/wence-/followers","following_url":"https://api.github.com/users/wence-/following{/other_user}","gists_url":"https://api.github.com/users/wence-/gists{/gist_id}","starred_url":"https://api.github.com/users/wence-/starred{/owner}{/repo}","subscriptions_url":"https://api.github.com/users/wence-/subscriptions","organizations_url":"https://api.github.com/users/wence-/orgs","repos_url":"https://api.github.com/users/wence-/repos","events_url":"https://api.github.com/users/wence-/events{/privacy}","received_events_url":"https://api.github.com/users/wence-/received_events","type":"User","site_admin":false},"created_at":"2018-08-17T15:50:00Z","updated_at":"2018-08-17T15:50:00Z","author_association":"CONTRIBUTOR","body":"OK, so something is funky in the C compiler as regards optimisation."}}
{"actor":{"display_login":"wence-"}}