Question #105900 on DOLFIN changed: https://answers.launchpad.net/dolfin/+question/105900
Anders Logg proposed the following answer: On Tue, Mar 30, 2010 at 07:21:25AM -0000, Garth Wells wrote: > Question #105900 on DOLFIN changed: > https://answers.launchpad.net/dolfin/+question/105900 > > Garth Wells proposed the following answer: > > On 30/03/10 02:32, Marie Rognes wrote: > > Chris Richardson wrote: > >> New question #105900 on DOLFIN: > >> https://answers.launchpad.net/dolfin/+question/105900 > >> > >> I want to use RT or BDM elements in a mixed-poisson type approach to > >> solve Darcy flow. > >> > >> Because the "pressure" is discretised on P0 elements, I cannot define > >> it on the boundary, so I have included a surface source term in the > >> RHS. Is this the correct approach? > > > > Yes, this looks correct. > > > >> Unfortunately, although it works on the left hand end of my domain, > >> the right hand end gives a slight discontinuity in slope. > >> > > > > > > I suspect that this is simply an artifact associated with plotting DG_0 > > functions. > > Could well be. plot(..) performs a projection onto continuous P1 > elements in the background. I think it does interpolation. Projection should be fine. -- Anders > The VTK file output format supports > cell-wise data (i.e. DG0), so try > > File("pressure.pvd") << p > > and view the result with ParaView. > > Garth > > > A similar "dip" appears with this simple script: > > > > from dolfin import * > > mesh = UnitSquare(20, 20) > > V = FunctionSpace(mesh,"DG",0) > > f = Expression("- 2.0*x[0] + 1.0", degree=1) > > g = project(f, V) > > plot(g, interactive=True) > > > > > > > > > > > >> Any help much appreciated, > >> Chris > >> > >> > >> mesh = UnitSquare(20,20) > >> n= FacetNormal(mesh) > >> > >> R = FunctionSpace(mesh,"BDM",1) > >> W = FunctionSpace(mesh,"DG",0) > >> V = R+W > >> (u,p) = TrialFunctions(V) > >> (v,q) = TestFunctions(V) > >> > >> # boundary conditions > >> def boundary0(x): > >> return x[1]<DOLFIN_EPS def boundary1(x): > >> return x[1]>1.0-DOLFIN_EPS > >> > >> u0 = Constant((0.0,0.0)) > >> bc0 = DirichletBC(V.sub(0), u0, boundary0) > >> bc1 = DirichletBC(V.sub(0), u0, boundary1) > >> bcs=[bc0,bc1] > >> > >> f=Expression("-1.0*(x[0]==0.0)+1.0*(x[0]==1.0)") > >> > >> # 'Permeability' k - gaussian in y > >> k=Expression("1.0+exp(-40*pow(x[1]-0.5,2))") > >> > >> a = (k*dot(v,u) - div(v)*p + q*div(u))*dx > >> L = dot(v,f*n)*ds > >> > >> # Compute solution > >> problem = VariationalProblem(a, L,bcs) > >> (u, p) = problem.solve().split() > >> > >> # Plot solution > >> plot(p, interactive=True) > >> plot(u, interactive=True) > >> > >> > > > > > > _______________________________________________ > > Mailing list: https://launchpad.net/~dolfin > > Post to : [email protected] > > Unsubscribe : https://launchpad.net/~dolfin > > More help : https://help.launchpad.net/ListHelp > -- You received this question notification because you are a member of DOLFIN Team, which is an answer contact for DOLFIN. _______________________________________________ Mailing list: https://launchpad.net/~dolfin Post to : [email protected] Unsubscribe : https://launchpad.net/~dolfin More help : https://help.launchpad.net/ListHelp

