dolfin team mailing list archive

[Martin Sandve Alnæs <martinal@xxxxxxxxx>]

Kyle,
you will get faster answers if you post
the error messages you get in addition
to the program. :)

Martin

On 3 May 2011 18:27, Kyle <question155708@xxxxxxxxxxxxxxxxxxxxx> wrote:
> Question #155708 on DOLFIN changed:
> https://answers.launchpad.net/dolfin/+question/155708
>
> Kyle requested more information:
> After reading Johan's post I was able to get past the nonlinear problem
> but now my program crashes. Could anyone help me with this?
>
> Here is my program:
>
> #!/usr/bin/python
>
> from dolfin import *
> import numpy, sys
>
> # MESHING
> mesh = Rectangle(0.0, -0.5, 1.0, 0.5, 10, 10, "right/left")
>
> # FUNCTION SPACE
> scalar = FunctionSpace(mesh, "CG", 1) #Pressure
> vector = VectorFunctionSpace(mesh, "CG", 2) #Velocity
> system = vector * scalar #Mixed Function Space
>
> # BOUNDARIES
> right = compile_subdomains('x[0] == 1.0')
> left = compile_subdomains('x[0] == 0.0')
> top = compile_subdomains('x[1] == 0.5')
> bottom = compile_subdomains('x[1] == -0.5')
>
> # INFLOW VELOCITY BC FOR TOP
> bc1 = DirichletBC(system.sub(0), Constant((1.0, 0.0)), top)
>
> # NO-SLIP BC FOR TOP-BOTTOM
> noslip = Constant((0.0, 0.0))
> bc0 = DirichletBC(system.sub(0), noslip, left)
> bc3 = DirichletBC(system.sub(0), noslip, bottom)
> bc2 = DirichletBC(system.sub(0), noslip, right)
>
> # STORAGE OF BCs
> bcs = [bc1, bc2, bc0, bc3]
>
> # VARIATIONAL PROBLEM
> ts = TestFunction(system)
> ds = TrialFunction(system)
> s = Function(system)
>
> v, q = split(ts)
> dv, dq = split(ds)
> u , p = split(s)
>
> f = Constant((0, 0))
> nu = 0.001
>
> #F = (nu*inner(grad(v), grad(u)) - div(v)*p + q*div(u) - inner(v,f))*dx # Stokes
> F = (nu*inner(grad(v), grad(u)) - div(v)*p + q*div(u) + inner(dot(u,grad(u)),v) - inner(v,f))*dx # Navier-Stokes
> J = derivative(F,s,ds)
>
> # SET UP PDE
> problem = VariationalProblem(F, J, bcs)
> #problem.is_nonlinear=True
>
> # SOLVE PDE
> (U, P) = problem.solve().split()
>
> # FILES FOR PARAVIEW
> ufile_pvd = File("velocity.pvd")
> ufile_pvd << U
> pfile_pvd = File("pressure.pvd")
> pfile_pvd << P
>
> # PLOTS
> plot(U)
> plot(P)
> interactive()
>
> --
> 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     : dolfin@xxxxxxxxxxxxxxxxxxx
> Unsubscribe : https://launchpad.net/~dolfin
> More help   : https://help.launchpad.net/ListHelp
>