<div dir="ltr"><div class="gmail_extra"><div class="gmail_quote">On Tue, Nov 25, 2014 at 4:08 AM, Olivier Bonnefon <span dir="ltr"><<a href="mailto:olivier.bonnefon@avignon.inra.fr" target="_blank">olivier.bonnefon@avignon.inra.fr</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
<div text="#000000" bgcolor="#FFFFFF">
<div>Hello,<br>
<br>
Thank you for your answer, I agree with you.<br>
<br>
I have adapted the ex12 example for the system:<br>
<br>
-\nabla \nabla u + u + f = 0 in \Omega (1)<br>
<br>
with f = 4-x^2-y^2<br>
The solution is still x^2-y^2.<br>
<br>
It consists in adding the following gradient function:<br>
<br>
void g0_uu(const PetscScalar u[], const PetscScalar u_t[], const
PetscScalar u_x[], const PetscScalar a[], const PetscScalar a_t[],
const PetscScalar a_x[], const PetscReal x[], PetscScalar g0[])<br>
{<br>
g0[0] = 1.0;<br>
}<br>
...<br>
//and set gradient<br>
ierr = PetscDSSetJacobian(prob, 0, 0, g0_uu, NULL, NULL,
g3_uu);CHKERRQ(ierr);<br>
...<br>
Of course, I have modified the f0_u function.<br>
<br>
<br>
With Dirichlet BC, I get the solution.<br>
<br>
With Neumann BC, the solution is still shifted down (-2/3). The
problem (1) with Neumann BC has a unique solution.<br>
How remove the condition \int_\Omega u = 0 ?<br></div></div></blockquote><div><br></div><div>That condition is one way of imposing boundary conditions for pressure. When I call MatSetNullSpace(), I am implying</div><div>this condition. You can take out that call, but the system will be singular. You would have to add a condition on the</div><div>pressure somehow, like fixing it at a point. This tends to make the system more ill-conditioned, but if that is appropriate</div><div>for you, then you make a similar call to</div><div><br></div><div> ierr = DMPlexAddBoundary(cdm, user->bcType == DIRICHLET, "wall", user->bcType == NEUMANN ? "boundary" : "marker", 0, user->exactFuncs[0], 1, &id, user);CHKERRQ(ierr);</div><div><br></div><div>but for field 1. To test , you can jsut change 0 to 1 in this call. However, this sets the entire boundary, and you really only</div><div>need to set the pressure at one point, so you could make another label that only has a single point and use that instead.</div><div><br></div><div> Thanks,</div><div><br></div><div> Matt</div><div> </div><blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex"><div text="#000000" bgcolor="#FFFFFF"><div>
Thanks<span class=""><font color="#888888"><br>
Olivier Bonnefon</font></span><div><div class="h5"><br>
<br>
<br>
On 11/06/2014 07:14 AM, Matthew Knepley wrote:<br>
</div></div></div><div><div class="h5">
<blockquote type="cite">
<div dir="ltr">
<div class="gmail_extra">
<div class="gmail_quote">On Mon, Nov 3, 2014 at 9:46 AM,
Olivier Bonnefon <span dir="ltr"><<a href="mailto:olivier.bonnefon@avignon.inra.fr" target="_blank">olivier.bonnefon@avignon.inra.fr</a>></span>
wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
<div text="#000000" bgcolor="#FFFFFF">
<div>Hello,<br>
<br>
Thank for your answer, I'll explain my trouble:<br>
<br>
My problem is that the BC Neumann leads to wrong
result.<br>
<br>
With Dirichlet BC, I get:<br>
-------------------------------<br>
<br>
> ./ex12 -run_type full -refinement_limit 0.0
-bc_type dirichlet -interpolate 0 -petscspace_order 1
-show_initial -show_solution -dm_plex_print_fem 1<br>
...<br>
...<br>
Solution<br>
Vec Object:potential 1 MPI processes<br>
type: seq<br>
0.5<br>
<br>
This result is correct.<br>
<br>
With Neuman BC, I get:<br>
--------------------------------<br>
<br>
> ./ex12 -run_type full -refinement_limit 0.0
-bc_type neumann -interpolate 1 -petscspace_order 1
-show_initial -dm_plex_print_fem 1 -show_solution
-bd_petscspace_order 1 -snes_linesearch_monitor
-snes_monitor -ksp_monitor_true_residual
-snes_converged_reason -ksp_converged_reason<br>
....<br>
....<br>
<br>
<br>
Solution<br>
Vec Object:potential 1 MPI processes<br>
type: seq<br>
-0.75<br>
-0.583333<br>
0.0833333<br>
-0.583333<br>
-0.333333<br>
0.416667<br>
0.0833333<br>
0.416667<br>
1.25<br>
<br>
<br>
That is not the values of the solution x*x+y*y.<br>
</div>
</div>
</blockquote>
<div><br>
</div>
<div>Sorry, this is poor documentation on my part. I will
fix it in the example. The Naumann solution</div>
<div>explicitly discards the constant part, meaning</div>
<div><br>
</div>
<div> \int_\Omega u = 0</div>
<div><br>
</div>
<div>If we look at my solution</div>
<div><br>
</div>
<div> \int^1_0 dx \int^1_0 dy x^2 + y^2 = 2/3</div>
<div><br>
</div>
<div>However, this is for the continuum result, whereas we
are enforcing this in the discrete case.</div>
<div>Thus, what we really get is 0.75, since we are
integrating linear patches instead of quadratic</div>
<div>patches. After shifting by that, you still do not get
exactly x^2 + y^2 since there is some</div>
<div>discretization error. If you run with P2 you should get
the exact answer, shifted down to eliminate</div>
<div>the DC component.</div>
<div><br>
</div>
<div> Thanks,</div>
<div><br>
</div>
<div> Matt</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">
<div text="#000000" bgcolor="#FFFFFF">
<div> <br>
I tried many ksp options. <br>
Moreover, the neumann BC with "-run_type full" is not
cover in the list <a href="https://bitbucket.org/petsc/petsc/src/fced3c3f9e703542693913793d15321603e40fe6/config/builder.py?at=master#cl-257" target="_blank">https://bitbucket.org/petsc/petsc/src/fced3c3f9e703542693913793d15321603e40fe6/config/builder.py?at=master#cl-257</a><br>
<br>
Do you know what is wrong ?<br>
<br>
Thanks,<br>
<br>
Olivier Bonnefon<br>
<br>
On 10/31/2014 06:50 PM, Matthew Knepley wrote:<br>
</div>
<blockquote type="cite">
<div dir="ltr">
<div class="gmail_extra">
<div class="gmail_quote">On Fri, Oct 31, 2014 at
10:43 AM, Olivier Bonnefon <span dir="ltr"><<a href="mailto:olivier.bonnefon@avignon.inra.fr" target="_blank">olivier.bonnefon@avignon.inra.fr</a>></span>
wrote:<br>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">Hello,<br>
<br>
I'm working on the snes/examples/tutorial/ex12
version 3.5.2.<br>
<br>
I didn't succed to run the simplest case with
Neumann BC:<br>
<br>
./ex12 -run_type full -refinement_limit 0.0
-bc_type neumann -interpolate 1
-petscspace_order 1 -show_initial
-dm_plex_print_fem 1 -show_solution
-bd_petscspace_order 1
-snes_linesearch_monitor -snes_monitor
-ksp_monitor_true_residual
-snes_converged_reason -ksp_converged_reason<br>
<br>
This leads to dofs negatives values.<br>
</blockquote>
<div><br>
</div>
<div>I do not understand what you mean. Please
always mail the full error message.</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">Do
you know the options to get a correct result
with Neumann BC ?<br>
</blockquote>
<div><br>
</div>
<div>There are some tests here:</div>
<div><br>
</div>
<div><a href="https://bitbucket.org/petsc/petsc/src/fced3c3f9e703542693913793d15321603e40fe6/config/builder.py?at=master#cl-257" target="_blank">https://bitbucket.org/petsc/petsc/src/fced3c3f9e703542693913793d15321603e40fe6/config/builder.py?at=master#cl-257</a><br>
</div>
<div><br>
</div>
<div> Matt</div>
<div> </div>
<blockquote class="gmail_quote" style="margin:0px 0px 0px 0.8ex;border-left-width:1px;border-left-color:rgb(204,204,204);border-left-style:solid;padding-left:1ex">Regards,<br>
Olivier Bonnefon<br>
<br>
</blockquote>
</div>
<br>
<br clear="all">
<span><font color="#888888">
<div><br>
</div>
-- <br>
What most experimenters take for granted
before they begin their experiments is
infinitely more interesting than any results
to which their experiments lead.<br>
-- Norbert Wiener </font></span></div>
<span><font color="#888888"> </font></span></div>
<span><font color="#888888"> </font></span></blockquote>
<span><font color="#888888"> <br>
<br>
<pre cols="72">--
Olivier Bonnefon
INRA PACA-Avignon, Unité BioSP
Tel: <a href="tel:%2B33%20%280%294%2032%2072%2021%2058" value="+33432722158" target="_blank">+33 (0)4 32 72 21 58</a></pre>
</font></span></div>
</blockquote>
</div>
<br>
<br clear="all">
<div><br>
</div>
-- <br>
<div>What most experimenters take for
granted before they begin their experiments is infinitely
more interesting than any results to which their experiments
lead.<br>
-- Norbert Wiener</div>
</div>
</div>
</blockquote>
<br>
<br>
<pre cols="72">--
Olivier Bonnefon
INRA PACA-Avignon, Unité BioSP
Tel: <a href="tel:%2B33%20%280%294%2032%2072%2021%2058" value="+33432722158" target="_blank">+33 (0)4 32 72 21 58</a></pre>
</div></div></div>
</blockquote></div><br><br clear="all"><div><br></div>-- <br><div class="gmail_signature">What most experimenters take for granted before they begin their experiments is infinitely more interesting than any results to which their experiments lead.<br>-- Norbert Wiener</div>
</div></div>