dolfin team mailing list archive
-
dolfin team
-
Mailing list archive
-
Message #11789
Re: Nonlinear Coupled PDEs, Newton Solver
Hi Hatef,
It sounds as if you are interested in the Navier-Stokes equations, or
similar? If so, you may want to check out the Unicorn solver:
http://www.fenics.org/wiki/Unicorn
Direct any Unicorn questions to: unicorn-dev@xxxxxxxxxx
Best,
Johan
> <p>Hi</p><p>I installed the new version of dolfin (0.9.0) as you wanted
> and am trying to model a coupled set of PDEs which are nonlinear. As I
> told in my previous emails, there are currenlty not much documents on
> nonlinear problems. so based on what I have seen on the demos I came up
> with the following code which does not seem to be working. I would
> appreciate if anyone could suggest any material or hint on modeling these
> kind of problems in which there are more that one variable. in my case
> there is a vector variable (Velocity) and a scalar one (Pressure) and
> Nonlinearity is only due to the term (u.grad)u and I am trying to use
> Newton's solver for nonlinear iterations.</p><p><br /><br
> />Thanks</p><p>Hatef</p><p>//////////////////////////////////////////////////////////////////////////////////////////////////////////////////////<br
> />//Noslip function<br /><br />class VelNoslip : public Function<br
> /> {<br /> public:<br /><br /> void
> eval(real* values, const real* x) const<br /> {<br
> /> <br
> />
> values[0]=0.0 ;<br
> />
> values[1]=0.0 ;<br /> <br
> /> }<br /><br
> /> };<br /><br
> />////////////////////////////////////////////////////////////////////////////////////<br
> />//-Noslip BC subdomain<br />class Noslip_boundary_end : public
> SubDomain<br /> {<br /> bool inside(const real* x,
> bool on_boundary) const<br /> {<br
> /> return <br
> />
> (x[0]>1.0-DOLFIN_EPS) && on_boundary;<br
> /> <br /> }<br />
> };<br /><br />class Noslip_boundary_side : public SubDomain<br />
> {<br /> bool inside(const real* x, bool on_boundary)
> const<br /> {<br /> return
> <br
> />
> (x[0]>DOLFIN_EPS && x[0]<1.0-DOLFIN_EPS) &&
> on_boundary;<br /> <br />
> }<br /> };<br /><br /><br />// define subdomain for inflow<br
> />class Inflow_boundary : public SubDomain<br /> {<br
> /> bool inside(const real* x, bool on_boundary) const<br
> /> {<br /> return <br
> />
> x[0]< DOLFIN_EPS && on_boundary;<br
> /> <br /> }<br />
> };<br /><br
> />////////////////////////////////////////////////////////////////////////////////////////////<br
> />int main()<br />{<br /><br /><br />//Read Mesh and Subdomain marker and
> define function space<br />UnitSquare mesh(32,32);</p><p>//Define
> fuction space<br /> NSWEFunctionSpace V(mesh);<br /> SubSpace
> velocity(V,0);<br /> SubSpace pressure(V,1);</p><p><br /><br
> />MeshFunction <unsigned int> sub_domains(mesh,
> mesh.topology().dim() -1);<br /><br /><br />// Mark all facets as sub
> domain 3<br /> sub_domains = 3;<br />// Mark end Noslip subdomain<br
> /> Noslip_boundary_end noslip_boundary_end;<br />
> noslip_boundary_end.mark(sub_domains,0);<br /><br /><br />// mark side
> noslip subdomain<br /> Noslip_boundary_side noslip_boundary_side;<br
> /> noslip_boundary_side.mark(sub_domains,1);<br /><br /><br />//mark
> Inflow subdomain<br /> Inflow_boundary inflow_boundary;<br />
> inflow_boundary.mark(sub_domains,2);<br /> //pseudo time<br />
> real t = 0.0;<br /> <br /> //solution vector<br />
> Function v(V);<br /> Function U(velocity);<br /> Function
> P(pressure); <br /> <br /> Function Un(velocity);<br />
> Un.vector().zero();<br /> Function Pn(pressure);<br />
> Pn.vector().zero();<br /> <br /> <br /> Constant
> f(2,0.0);<br /> VelNoslip velnoslip;<br /> Constant
> inflow(0.1);<br /> <br /> <br /><br /> <br /> <br
> /> //set up Forms<br /> NSWEBilinearForm a(V,V) ;<br />
> a.U0=U;<br /> NSWELinearForm L(V);<br /> L.Un=Un; L.Pn=Pn ;
> L.F=f; L.U0=U; L.P0=P;<br /> <br /> <br /> //set up
> boundary condition<br /> <br /> DirichletBC bc0(velocity,
> velnoslip,sub_domains,0);<br /> DirichletBC bc1(velocity,
> velnoslip,sub_domains,1);<br /> DirichletBC bc2(velocity,
> velnoslip,sub_domains,2);<br /> DirichletBC bc3(pressure,
> inflow,sub_domains,2);<br /> //collect all the BCs<br />
> Array<BoundaryCondition*>
> bcs(&bc0,&bc1,&bc2,&bc3);<br /> <br />
> //parameters for time stepping<br /><br /> real T = 0.005;<br
> /> real k = 0.001;<br /> <br /> <br /> <br
> /> File file_2("Pressure.pvd");<br /> File
> file_1("Velocity.pvd");<br /><br /> <br /><br />
> //time stepping<br /> Progress prg("Time-stepping");<br
> /> while ( t < T )<br /> { <br /> <br />
> VariationalProblem nlpde(a,L,bcs,true);<br /> <br />
> nlpde.solve(v);<br /> <br /> U = v[0];<br /> P =
> v[1];</p><p><br /> // Save the solution to file<br /> file_2
> << P;<br /> file_1 << U;<br /> <br /> //
> Move to the next interval </p><p> Un.vector()=U.vector() ; <br
> /> Pn.vector()=P.vector() ; <br /> prg = t / T;<br /> t
> += k;<br /> <br /> }<br /><br />}</p>
> _______________________________________________
> DOLFIN-dev mailing list
> DOLFIN-dev@xxxxxxxxxx
> http://www.fenics.org/mailman/listinfo/dolfin-dev
>
References
