Re: [firedrake] imposing boundary conditions for Hermite element MixedFunctionSpace
On Thu, Jul 17, 2025 at 11:06 AM Smail Merabet <merabetsmail@gmail.com> wrote:
The lower is a constant obstacle psi =-0.03 I have the plot using freefem++
1. You cannot leave off the Cc if you want an answer from the Firedrake people 2. Did you plot your lower bound vector? The reason I ask is that I am absolutely sure that SNESVI will force the solution to be above the lower bound you input. Therefore I think you are inputting a lower bound different from what you expect. Thanks, Matt
On Thu, Jul 17, 2025 at 4:02 PM Matthew Knepley <knepley@gmail.com> wrote:
On Thu, Jul 17, 2025 at 10:49 AM Smail Merabet <merabetsmail@gmail.com> wrote:
Dear Matthew, thanks for reply, I tried this:
psi = Constant(-0.03)
# Create the full function in mixed space
obs = Function(Vh, name="Obstacle")
upp = Function(Vh, name="Upper")
obs_sub1 = obs.sub(1)
upp_sub1 = upp.sub(1)
# Interpolate into the subfunctions
obs_sub1.interpolate(conditional(lt(psi, 0), 0, psi))
I do not understand the Firedrake side, but did you plot your lower bound? I have no idea what conditional gives back.
Thanks,
Matt
upp_sub1.interpolate(Constant(1e10))
But the result doesn't satisfy the constraint $w \geq \psi$
Many Thanks
On Thu, Jul 17, 2025 at 1:22 PM Matthew Knepley <knepley@gmail.com> wrote:
On Thu, Jul 17, 2025 at 7:52 AM Smail Merabet <merabetsmail@gmail.com> wrote:
Dear firedrake users, I try to solve a variational inequality using snes-firedrake. I have two unknown functions say (theta,w) so the functionspace is Vh = MixedFunctionSpace([Th,Wh]) where Th = FunctionSpace(mesh, "Lagrange", 1) Wh = FunctionSpace(mesh, "Lagrange", 2)
the inequality constraint concerns only $w$
I set : obstacle = Function(Vh.sub(1)).interpolate((Constant(-0.03))) upper = Function(Vh.sub(1)).interpolate(Constant(1e10))
The SNES is operating on the full space, but you have interpolated constraints on the subspace. You have to inject those to full space vectors before setting the bounds.
Thanks,
Matt
*solver.solve(bounds=(obstacle, upper)) *
but this line causes the following error:
File "/Users/maths/Desktop/Firedrake/VI/Timoshenko-obs.py", line 64, in <module> solver.solve(bounds=(obstacle, upper)) File "petsc4py/PETSc/Log.pyx", line 188, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func File "petsc4py/PETSc/Log.pyx", line 189, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func File "/Users/maths/firedrake/src/firedrake/firedrake/adjoint_utils/variational_solver.py", line 89, in wrapper out = solve(self, **kwargs) File "/Users/maths/firedrake/src/firedrake/firedrake/variational_solver.py", line 312, in solve self.snes.setVariableBounds(lb, ub) File "petsc4py/PETSc/SNES.pyx", line 1978, in petsc4py.PETSc.SNES.setVariableBounds petsc4py.PETSc.Error: error code 75 [0] SNESVISetVariableBounds() at /Users/maths/firedrake/src/petsc/src/snes/impls/vi/vi.c:424 [0] SNESVISetVariableBounds_VI() at /Users/maths/firedrake/src/petsc/src/snes/impls/vi/vi.c:443 [0] Arguments are incompatible [0] Incompatible vector lengths lower bound = 201 solution vector = 302
Hope that one can help me ...
Best regards. Smail.
On Fri, Apr 11, 2025 at 9:11 AM Smail Merabet <merabetsmail@gmail.com> wrote:
Dear firedrake users,
I try to solve a fourth-order 1d nonlinear problem using the Hermite element in firedrake. The functional space for the weak solution is: $$\mathcal{H}\times \mathcal{H}/times H^{2}_{0} $$ where, $$ \mathcal{H} = \left\lbrace f\in H^2([0,\ell]), \; f(0)=0, \; f'(0)=0 \right\rbrace. $$ I can't impose the boundary conditions correctly. When I compile the code, I get the following error: FiniteElement('Hermite', interval, 3), name=None, index=1, component=None), IndexedProxyFunctionSpace(<firedrake.mesh.MeshTopology object at 0x125a2ce80>, FiniteElement('Hermite', interval, 3), name=None, index=2, component=None), name='None_None_None'), Mesh(VectorElement(FiniteElement('Lagrange', interval, 1), dim=1), 0)), 9) defined on incompatible FunctionSpace!
Note the same problem was implemented using FEniCS with C0-interior penalty method, in the reference (sec. 3 eq 3.12-3.13):
Enhanced models for the nonlinear bending of planar rods: localization phenomena and multistability by: Matteo Brunetti , Antonino Favata and Stefano Vidoli.
I attached the codes to this email.
Thanks in advance Smail MERABET
_______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
Thanks Matthew, for the help, in my problem the functions of the mixed space are only the second component is bonded from below by a constant obstacle psi=-0.03 the first component is unconstrained psi = Constant(-0.03) # Create the full function in mixed space obs = Function(Vh, name="Obstacle") upp = Function(Vh, name="Upper") obs_sub1 = obs.sub(1) upp_sub1 = upp.sub(1) # Interpolate into the subfunctions obs_sub1.interpolate(conditional(lt(psi, 0), 0., psi)) upp_sub1.interpolate(Constant(1e10)) the obtained solution (attached figure) is bounded from below but by 0. I attached the code to this email also. Thanks Smail. On Thu, Jul 17, 2025 at 4:11 PM Matthew Knepley <knepley@gmail.com> wrote:
On Thu, Jul 17, 2025 at 11:06 AM Smail Merabet <merabetsmail@gmail.com> wrote:
The lower is a constant obstacle psi =-0.03 I have the plot using freefem++
1. You cannot leave off the Cc if you want an answer from the Firedrake people
2. Did you plot your lower bound vector? The reason I ask is that I am absolutely sure that SNESVI will force the solution to be above the lower bound you input. Therefore I think you are inputting a lower bound different from what you expect.
Thanks,
Matt
On Thu, Jul 17, 2025 at 4:02 PM Matthew Knepley <knepley@gmail.com> wrote:
On Thu, Jul 17, 2025 at 10:49 AM Smail Merabet <merabetsmail@gmail.com> wrote:
Dear Matthew, thanks for reply, I tried this:
psi = Constant(-0.03)
# Create the full function in mixed space
obs = Function(Vh, name="Obstacle")
upp = Function(Vh, name="Upper")
obs_sub1 = obs.sub(1)
upp_sub1 = upp.sub(1)
# Interpolate into the subfunctions
obs_sub1.interpolate(conditional(lt(psi, 0), 0, psi))
I do not understand the Firedrake side, but did you plot your lower bound? I have no idea what conditional gives back.
Thanks,
Matt
upp_sub1.interpolate(Constant(1e10))
But the result doesn't satisfy the constraint $w \geq \psi$
Many Thanks
On Thu, Jul 17, 2025 at 1:22 PM Matthew Knepley <knepley@gmail.com> wrote:
On Thu, Jul 17, 2025 at 7:52 AM Smail Merabet <merabetsmail@gmail.com> wrote:
Dear firedrake users, I try to solve a variational inequality using snes-firedrake. I have two unknown functions say (theta,w) so the functionspace is Vh = MixedFunctionSpace([Th,Wh]) where Th = FunctionSpace(mesh, "Lagrange", 1) Wh = FunctionSpace(mesh, "Lagrange", 2)
the inequality constraint concerns only $w$
I set : obstacle = Function(Vh.sub(1)).interpolate((Constant(-0.03))) upper = Function(Vh.sub(1)).interpolate(Constant(1e10))
The SNES is operating on the full space, but you have interpolated constraints on the subspace. You have to inject those to full space vectors before setting the bounds.
Thanks,
Matt
*solver.solve(bounds=(obstacle, upper)) *
but this line causes the following error:
File "/Users/maths/Desktop/Firedrake/VI/Timoshenko-obs.py", line 64, in <module> solver.solve(bounds=(obstacle, upper)) File "petsc4py/PETSc/Log.pyx", line 188, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func File "petsc4py/PETSc/Log.pyx", line 189, in petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func File "/Users/maths/firedrake/src/firedrake/firedrake/adjoint_utils/variational_solver.py", line 89, in wrapper out = solve(self, **kwargs) File "/Users/maths/firedrake/src/firedrake/firedrake/variational_solver.py", line 312, in solve self.snes.setVariableBounds(lb, ub) File "petsc4py/PETSc/SNES.pyx", line 1978, in petsc4py.PETSc.SNES.setVariableBounds petsc4py.PETSc.Error: error code 75 [0] SNESVISetVariableBounds() at /Users/maths/firedrake/src/petsc/src/snes/impls/vi/vi.c:424 [0] SNESVISetVariableBounds_VI() at /Users/maths/firedrake/src/petsc/src/snes/impls/vi/vi.c:443 [0] Arguments are incompatible [0] Incompatible vector lengths lower bound = 201 solution vector = 302
Hope that one can help me ...
Best regards. Smail.
On Fri, Apr 11, 2025 at 9:11 AM Smail Merabet <merabetsmail@gmail.com> wrote:
> Dear firedrake users, > > I try to solve a fourth-order 1d nonlinear problem using the Hermite > element in firedrake. > The functional space for the weak solution is: > $$\mathcal{H}\times \mathcal{H}/times H^{2}_{0} $$ > where, > $$ > \mathcal{H} = \left\lbrace f\in H^2([0,\ell]), \; f(0)=0, \; f'(0)=0 > \right\rbrace. > $$ > I can't impose the boundary conditions correctly. When I compile the > code, I get the following error: > FiniteElement('Hermite', interval, 3), name=None, index=1, > component=None), IndexedProxyFunctionSpace(<firedrake.mesh.MeshTopology > object at 0x125a2ce80>, FiniteElement('Hermite', interval, 3), name=None, > index=2, component=None), name='None_None_None'), > Mesh(VectorElement(FiniteElement('Lagrange', interval, 1), dim=1), 0)), 9) > defined on incompatible FunctionSpace! > > Note the same problem was implemented using FEniCS with C0-interior > penalty method, in the reference (sec. 3 eq 3.12-3.13): > > Enhanced models for the nonlinear bending of planar rods: > localization phenomena and multistability > by: > Matteo Brunetti , Antonino Favata and Stefano Vidoli. > > I attached the codes to this email. > > Thanks in advance > Smail MERABET > > _______________________________________________ firedrake mailing list firedrake@imperial.ac.uk https://mailman.ic.ac.uk/mailman/listinfo/firedrake
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
I find the mistake, the given obstacle is not correct. Thanks for the help. Best regards Smail. On Thu, Jul 17, 2025 at 4:22 PM Smail Merabet <merabetsmail@gmail.com> wrote:
Thanks Matthew, for the help, in my problem the functions of the mixed space are only the second component is bonded from below by a constant obstacle psi=-0.03 the first component is unconstrained
psi = Constant(-0.03)
# Create the full function in mixed space
obs = Function(Vh, name="Obstacle")
upp = Function(Vh, name="Upper")
obs_sub1 = obs.sub(1)
upp_sub1 = upp.sub(1)
# Interpolate into the subfunctions
obs_sub1.interpolate(conditional(lt(psi, 0), 0., psi))
upp_sub1.interpolate(Constant(1e10))
the obtained solution (attached figure) is bounded from below but by 0.
I attached the code to this email also.
Thanks
Smail.
On Thu, Jul 17, 2025 at 4:11 PM Matthew Knepley <knepley@gmail.com> wrote:
On Thu, Jul 17, 2025 at 11:06 AM Smail Merabet <merabetsmail@gmail.com> wrote:
The lower is a constant obstacle psi =-0.03 I have the plot using freefem++
1. You cannot leave off the Cc if you want an answer from the Firedrake people
2. Did you plot your lower bound vector? The reason I ask is that I am absolutely sure that SNESVI will force the solution to be above the lower bound you input. Therefore I think you are inputting a lower bound different from what you expect.
Thanks,
Matt
On Thu, Jul 17, 2025 at 4:02 PM Matthew Knepley <knepley@gmail.com> wrote:
On Thu, Jul 17, 2025 at 10:49 AM Smail Merabet <merabetsmail@gmail.com> wrote:
Dear Matthew, thanks for reply, I tried this:
psi = Constant(-0.03)
# Create the full function in mixed space
obs = Function(Vh, name="Obstacle")
upp = Function(Vh, name="Upper")
obs_sub1 = obs.sub(1)
upp_sub1 = upp.sub(1)
# Interpolate into the subfunctions
obs_sub1.interpolate(conditional(lt(psi, 0), 0, psi))
I do not understand the Firedrake side, but did you plot your lower bound? I have no idea what conditional gives back.
Thanks,
Matt
upp_sub1.interpolate(Constant(1e10))
But the result doesn't satisfy the constraint $w \geq \psi$
Many Thanks
On Thu, Jul 17, 2025 at 1:22 PM Matthew Knepley <knepley@gmail.com> wrote:
On Thu, Jul 17, 2025 at 7:52 AM Smail Merabet <merabetsmail@gmail.com> wrote:
> Dear firedrake users, I try to solve a variational inequality using > snes-firedrake. > I have two unknown functions say (theta,w) > so the functionspace is Vh = MixedFunctionSpace([Th,Wh]) > where > Th = FunctionSpace(mesh, "Lagrange", 1) > Wh = FunctionSpace(mesh, "Lagrange", 2) > > > the inequality constraint concerns only $w$ > > I set : > obstacle = Function(Vh.sub(1)).interpolate((Constant(-0.03))) > upper = Function(Vh.sub(1)).interpolate(Constant(1e10)) >
The SNES is operating on the full space, but you have interpolated constraints on the subspace. You have to inject those to full space vectors before setting the bounds.
Thanks,
Matt
> *solver.solve(bounds=(obstacle, upper)) * > > but this line causes the following error: > > File "/Users/maths/Desktop/Firedrake/VI/Timoshenko-obs.py", line 64, > in <module> > solver.solve(bounds=(obstacle, upper)) > File "petsc4py/PETSc/Log.pyx", line 188, in > petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func > File "petsc4py/PETSc/Log.pyx", line 189, in > petsc4py.PETSc.Log.EventDecorator.decorator.wrapped_func > File > "/Users/maths/firedrake/src/firedrake/firedrake/adjoint_utils/variational_solver.py", > line 89, in wrapper > out = solve(self, **kwargs) > File > "/Users/maths/firedrake/src/firedrake/firedrake/variational_solver.py", > line 312, in solve > self.snes.setVariableBounds(lb, ub) > File "petsc4py/PETSc/SNES.pyx", line 1978, in > petsc4py.PETSc.SNES.setVariableBounds > petsc4py.PETSc.Error: error code 75 > [0] SNESVISetVariableBounds() at > /Users/maths/firedrake/src/petsc/src/snes/impls/vi/vi.c:424 > [0] SNESVISetVariableBounds_VI() at > /Users/maths/firedrake/src/petsc/src/snes/impls/vi/vi.c:443 > [0] Arguments are incompatible > [0] Incompatible vector lengths lower bound = 201 solution vector = > 302 > > > Hope that one can help me ... > > Best regards. > Smail. > > On Fri, Apr 11, 2025 at 9:11 AM Smail Merabet < > merabetsmail@gmail.com> wrote: > >> Dear firedrake users, >> >> I try to solve a fourth-order 1d nonlinear problem using the >> Hermite element in firedrake. >> The functional space for the weak solution is: >> $$\mathcal{H}\times \mathcal{H}/times H^{2}_{0} $$ >> where, >> $$ >> \mathcal{H} = \left\lbrace f\in H^2([0,\ell]), \; f(0)=0, \; >> f'(0)=0 \right\rbrace. >> $$ >> I can't impose the boundary conditions correctly. When I compile >> the code, I get the following error: >> FiniteElement('Hermite', interval, 3), name=None, index=1, >> component=None), IndexedProxyFunctionSpace(<firedrake.mesh.MeshTopology >> object at 0x125a2ce80>, FiniteElement('Hermite', interval, 3), name=None, >> index=2, component=None), name='None_None_None'), >> Mesh(VectorElement(FiniteElement('Lagrange', interval, 1), dim=1), 0)), 9) >> defined on incompatible FunctionSpace! >> >> Note the same problem was implemented using FEniCS with C0-interior >> penalty method, in the reference (sec. 3 eq 3.12-3.13): >> >> Enhanced models for the nonlinear bending of planar rods: >> localization phenomena and multistability >> by: >> Matteo Brunetti , Antonino Favata and Stefano Vidoli. >> >> I attached the codes to this email. >> >> Thanks in advance >> Smail MERABET >> >> _______________________________________________ > firedrake mailing list > firedrake@imperial.ac.uk > https://mailman.ic.ac.uk/mailman/listinfo/firedrake >
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
-- What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead. -- Norbert Wiener
https://www.cse.buffalo.edu/~knepley/ <http://www.cse.buffalo.edu/~knepley/>
participants (2)
-
Matthew Knepley
-
Smail Merabet