Re: [petsc-users] Doubt about BT and BASIC NEWTONLS
By "scaling" I mean you may need to rescale the residual terms from the
Dirichlet boundary conditions to be the same order (in terms of the
characteristic size of the elements) as they are for the other residuals.
It seems in your formulation the Jacobian entries of the residual for
Dirichlet boundary conditions are 1 on the diagonal while the Jacobian entries
of the residual for the PDE entries will have a scaling related to the finite
element method you are using, for example in 3d with an operator grad dot grad
the Jacobian entry will be order dx = (dx)^3 * (1/dx)^2. Thus you should scale
the Dirichlet boundary condition residuals by dx to get the same scaling.
Barry
> On Nov 16, 2021, at 1:37 PM, Francesc Levrero-Florencio
> wrote:
>
> Hi Barry,
>
> Thanks for the quick reply. I am not sure what you meant by "scaling". In our
> code, every time our residual is called we do the following in this order:
> - Apply updated Dirichlet BCs in order to assemble the internal forces, and
> add them to the residual.
> - Assemble the external forces, and subtract them from the residual.
> - Update the residual to take into account Dirichet BCs, so we set the
> corresponding indices of the residual to - (dirichlet_bc_value -
> solution_value), the latter obtained directly from the solution vector of
> PETSc. However, note that in this case all Dirichlet BCs are homogeneous, so
> this value would be just "solution_value", or 0 throughout the simulation.
> This ensures that the solution vector coming from PETSc has the Dirichlet BCs
> applied.
>
> This is basically what we do in the residual. If I am not mistaken, the
> lambda scaling factor from the line-search would scale the factor
> "jacobian^{-1} * residual", so it would scale equally both internal and
> external forces.
>
> Regards,
> Francesc.
>
> On Tue, Nov 16, 2021 at 5:38 PM Barry Smith > wrote:
>
> Perhaps this behavior is the result of a "scaling" issue in how various
> terms affect the residual? In particular perhaps the terms for enforcing
> boundary conditions are scaled differently than terms for the PDE enforcement?
>
>
>
> > On Nov 16, 2021, at 11:19 AM, Francesc Levrero-Florencio
> > mailto:f.levrero-floren...@onscale.com>>
> > wrote:
> >
> > Dear PETSc team and users,
> >
> > We are running a simple cantilever beam bending, where the profile of the
> > beam is I-shaped, where we apply a bending force on one end and fully
> > constrained displacements on the other end. The formulation is a large
> > strain formulation in Total Lagrangian form, where the material of the beam
> > is a Saint Venant-Kirchhoff hyperelastic material that uses the same
> > constants as steel (200E9 GPa Young’s modulus and 0.3 Poisson’s ratio).
> >
> > The simulation finishes in the requested number of time steps by using the
> > “basic” type of line-search in the SNES (i.e. traditional Newton method
> > without line-search) in a reasonable number of Newton iterations per time
> > step (3 or 4 iterations). However, when using the “bt” (or “l2”, and in
> > fact no type of line-search ends up yielding convergence) line-search type,
> > the convergence never happens within the SNES default maximum number of
> > iterations of 50.
> >
> > During solving with traditional Newton, the general behaviour of each time
> > step is that the norm of the residual increases on the second call to the
> > residual function, but then hugely decreases in the following one or two,
> > up to the point where convergence is achieved. Using “bt” line-search, the
> > line-search discards the step at lambda=1 immediately because the norm of
> > the residual is higher than that produced in the first call to the residual
> > function, cutting down the value of lambda to a value significantly lower
> > than 1. The simulation then progresses in following Newton iterations in a
> > similar fashion, the line-search step at lambda=1 is always discarded, and
> > then smaller steps are taken but convergence never occurs, for even the
> > first time step.
> >
> > Here are a few values of the relevant norms (using traditional Newton) in
> > the first time step:
> >
> > BASIC NEWTON LS
> > Norm of the internal forces is 0
> > Norm of the external forces is 1374.49
> > Norm of the residual is 1374.49
> > Norm of the solution with Dirichlet BCs is 0
> > Number of SNES iteration is 0
> > -
> > Norm of the internal forces is 113498
> > Norm of the external forces is 1374.49
> > Norm of the residual is 105053
> > Norm of the solution with Dirichlet BCs is 0.441466
> > Number of SNES iteration is 0
> > -
> > Norm of the internal forces is 42953.5
> > Norm of the external forces is 1374.49
> > Norm of the residual is 11.3734
> > Norm of the solution with Dirichlet BCs is 0.441438
>
Re: [petsc-users] Doubt about BT and BASIC NEWTONLS
Hi Barry,
Thanks for the quick reply. I am not sure what you meant by "scaling". In
our code, every time our residual is called we do the following in this
order:
- Apply updated Dirichlet BCs in order to assemble the internal forces, and
add them to the residual.
- Assemble the external forces, and subtract them from the residual.
- Update the residual to take into account Dirichet BCs, so we set the
corresponding indices of the residual to - (dirichlet_bc_value -
solution_value), the latter obtained directly from the solution vector of
PETSc. However, note that in this case all Dirichlet BCs are homogeneous,
so this value would be just "solution_value", or 0 throughout the
simulation. This ensures that the solution vector coming from PETSc has the
Dirichlet BCs applied.
This is basically what we do in the residual. If I am not mistaken, the
lambda scaling factor from the line-search would scale the factor
"jacobian^{-1} * residual", so it would scale equally both internal and
external forces.
Regards,
Francesc.
On Tue, Nov 16, 2021 at 5:38 PM Barry Smith wrote:
>
> Perhaps this behavior is the result of a "scaling" issue in how various
> terms affect the residual? In particular perhaps the terms for enforcing
> boundary conditions are scaled differently than terms for the PDE
> enforcement?
>
>
>
> > On Nov 16, 2021, at 11:19 AM, Francesc Levrero-Florencio <
> f.levrero-floren...@onscale.com> wrote:
> >
> > Dear PETSc team and users,
> >
> > We are running a simple cantilever beam bending, where the profile of
> the beam is I-shaped, where we apply a bending force on one end and fully
> constrained displacements on the other end. The formulation is a large
> strain formulation in Total Lagrangian form, where the material of the beam
> is a Saint Venant-Kirchhoff hyperelastic material that uses the same
> constants as steel (200E9 GPa Young’s modulus and 0.3 Poisson’s ratio).
> >
> > The simulation finishes in the requested number of time steps by using
> the “basic” type of line-search in the SNES (i.e. traditional Newton method
> without line-search) in a reasonable number of Newton iterations per time
> step (3 or 4 iterations). However, when using the “bt” (or “l2”, and in
> fact no type of line-search ends up yielding convergence) line-search type,
> the convergence never happens within the SNES default maximum number of
> iterations of 50.
> >
> > During solving with traditional Newton, the general behaviour of each
> time step is that the norm of the residual increases on the second call to
> the residual function, but then hugely decreases in the following one or
> two, up to the point where convergence is achieved. Using “bt” line-search,
> the line-search discards the step at lambda=1 immediately because the norm
> of the residual is higher than that produced in the first call to the
> residual function, cutting down the value of lambda to a value
> significantly lower than 1. The simulation then progresses in following
> Newton iterations in a similar fashion, the line-search step at lambda=1 is
> always discarded, and then smaller steps are taken but convergence never
> occurs, for even the first time step.
> >
> > Here are a few values of the relevant norms (using traditional Newton)
> in the first time step:
> >
> > BASIC NEWTON LS
> > Norm of the internal forces is 0
> > Norm of the external forces is 1374.49
> > Norm of the residual is 1374.49
> > Norm of the solution with Dirichlet BCs is 0
> > Number of SNES iteration is 0
> > -
> > Norm of the internal forces is 113498
> > Norm of the external forces is 1374.49
> > Norm of the residual is 105053
> > Norm of the solution with Dirichlet BCs is 0.441466
> > Number of SNES iteration is 0
> > -
> > Norm of the internal forces is 42953.5
> > Norm of the external forces is 1374.49
> > Norm of the residual is 11.3734
> > Norm of the solution with Dirichlet BCs is 0.441438
> > Number of SNES iteration is 1
> >
> > Here are a few values of the relevant norms (using “bt” line-search) in
> the first time step:
> >
> > BT NEWTONLS
> > Norm of the internal forces is 0
> > Norm of the external forces is 1374.49
> > Norm of the residual is 1374.49
> > Norm of the solution with Dirichlet BCs is 0
> > Number of SNES iteration is 0
> > -
> > Norm of the internal forces is 113498
> > Norm of the external forces is 1374.49
> > Norm of the residual is 105053
> > Norm of the solution with Dirichlet BCs is 0.441466
> > Number of SNES iteration is 0
> > (I assume this is the first try at lambda=1)
> > -
> > Norm of the internal forces is 4422.12
> > Norm of the external forces is 1374.49
> > Norm of the residual is 1622.74
> > Norm of the solution with Dirichlet BCs is 0.0441466
> > Num
Re: [petsc-users] Doubt about BT and BASIC NEWTONLS
Perhaps this behavior is the result of a "scaling" issue in how various terms affect the residual? In particular perhaps the terms for enforcing boundary conditions are scaled differently than terms for the PDE enforcement? > On Nov 16, 2021, at 11:19 AM, Francesc Levrero-Florencio > wrote: > > Dear PETSc team and users, > > We are running a simple cantilever beam bending, where the profile of the > beam is I-shaped, where we apply a bending force on one end and fully > constrained displacements on the other end. The formulation is a large strain > formulation in Total Lagrangian form, where the material of the beam is a > Saint Venant-Kirchhoff hyperelastic material that uses the same constants as > steel (200E9 GPa Young’s modulus and 0.3 Poisson’s ratio). > > The simulation finishes in the requested number of time steps by using the > “basic” type of line-search in the SNES (i.e. traditional Newton method > without line-search) in a reasonable number of Newton iterations per time > step (3 or 4 iterations). However, when using the “bt” (or “l2”, and in fact > no type of line-search ends up yielding convergence) line-search type, the > convergence never happens within the SNES default maximum number of > iterations of 50. > > During solving with traditional Newton, the general behaviour of each time > step is that the norm of the residual increases on the second call to the > residual function, but then hugely decreases in the following one or two, up > to the point where convergence is achieved. Using “bt” line-search, the > line-search discards the step at lambda=1 immediately because the norm of the > residual is higher than that produced in the first call to the residual > function, cutting down the value of lambda to a value significantly lower > than 1. The simulation then progresses in following Newton iterations in a > similar fashion, the line-search step at lambda=1 is always discarded, and > then smaller steps are taken but convergence never occurs, for even the first > time step. > > Here are a few values of the relevant norms (using traditional Newton) in the > first time step: > > BASIC NEWTON LS > Norm of the internal forces is 0 > Norm of the external forces is 1374.49 > Norm of the residual is 1374.49 > Norm of the solution with Dirichlet BCs is 0 > Number of SNES iteration is 0 > - > Norm of the internal forces is 113498 > Norm of the external forces is 1374.49 > Norm of the residual is 105053 > Norm of the solution with Dirichlet BCs is 0.441466 > Number of SNES iteration is 0 > - > Norm of the internal forces is 42953.5 > Norm of the external forces is 1374.49 > Norm of the residual is 11.3734 > Norm of the solution with Dirichlet BCs is 0.441438 > Number of SNES iteration is 1 > > Here are a few values of the relevant norms (using “bt” line-search) in the > first time step: > > BT NEWTONLS > Norm of the internal forces is 0 > Norm of the external forces is 1374.49 > Norm of the residual is 1374.49 > Norm of the solution with Dirichlet BCs is 0 > Number of SNES iteration is 0 > - > Norm of the internal forces is 113498 > Norm of the external forces is 1374.49 > Norm of the residual is 105053 > Norm of the solution with Dirichlet BCs is 0.441466 > Number of SNES iteration is 0 > (I assume this is the first try at lambda=1) > - > Norm of the internal forces is 4422.12 > Norm of the external forces is 1374.49 > Norm of the residual is 1622.74 > Norm of the solution with Dirichlet BCs is 0.0441466 > Number of SNES iteration is 0 > Line search: gnorm after quadratic fit 1.622742343614e+03 > (I assume that in this line-search iteration 0.05 < lambda < 1, but the > corresponding residual is not smaller than the one in the first call, 1622 > > 1374) > - > Norm of the internal forces is 2163.76 > Norm of the external forces is 1374.49 > Norm of the residual is 1331.88 > Norm of the solution with Dirichlet BCs is 0.0220733 > Number of SNES iteration is 0 > Line search: Cubically determined step, current gnorm 1.331884625811e+03 > lambda=5.0003e-02 > (This is the accepted lambda for the current Newton iteration) > - > Norm of the internal forces is 104020 > Norm of the external forces is 1374.49 > Norm of the residual is 94739 > Norm of the solution with Dirichlet BCs is 0.441323 > Number of SNES iteration is 1 > > My question would be, any idea on how to deal with this situation? I would > imagine a “hack” would be to bypass the first residual norm, and have the > line-search use the following one as the “base norm” to check its reduct
[petsc-users] Doubt about BT and BASIC NEWTONLS
Dear PETSc team and users, We are running a simple cantilever beam bending, where the profile of the beam is I-shaped, where we apply a bending force on one end and fully constrained displacements on the other end. The formulation is a large strain formulation in Total Lagrangian form, where the material of the beam is a Saint Venant-Kirchhoff hyperelastic material that uses the same constants as steel (200E9 GPa Young’s modulus and 0.3 Poisson’s ratio). The simulation finishes in the requested number of time steps by using the “basic” type of line-search in the SNES (i.e. traditional Newton method without line-search) in a reasonable number of Newton iterations per time step (3 or 4 iterations). However, when using the “bt” (or “l2”, and in fact no type of line-search ends up yielding convergence) line-search type, the convergence never happens within the SNES default maximum number of iterations of 50. During solving with traditional Newton, the general behaviour of each time step is that the norm of the residual increases on the second call to the residual function, but then hugely decreases in the following one or two, up to the point where convergence is achieved. Using “bt” line-search, the line-search discards the step at lambda=1 immediately because the norm of the residual is higher than that produced in the first call to the residual function, cutting down the value of lambda to a value significantly lower than 1. The simulation then progresses in following Newton iterations in a similar fashion, the line-search step at lambda=1 is always discarded, and then smaller steps are taken but convergence never occurs, for even the first time step. Here are a few values of the relevant norms (using traditional Newton) in the first time step: BASIC NEWTON LS Norm of the internal forces is 0 Norm of the external forces is 1374.49 Norm of the residual is 1374.49 Norm of the solution with Dirichlet BCs is 0 Number of SNES iteration is 0 - Norm of the internal forces is 113498 Norm of the external forces is 1374.49 Norm of the residual is 105053 Norm of the solution with Dirichlet BCs is 0.441466 Number of SNES iteration is 0 - Norm of the internal forces is 42953.5 Norm of the external forces is 1374.49 Norm of the residual is 11.3734 Norm of the solution with Dirichlet BCs is 0.441438 Number of SNES iteration is 1 Here are a few values of the relevant norms (using “bt” line-search) in the first time step: BT NEWTONLS Norm of the internal forces is 0 Norm of the external forces is 1374.49 Norm of the residual is 1374.49 Norm of the solution with Dirichlet BCs is 0 Number of SNES iteration is 0 - Norm of the internal forces is 113498 Norm of the external forces is 1374.49 Norm of the residual is 105053 Norm of the solution with Dirichlet BCs is 0.441466 Number of SNES iteration is 0 (I assume this is the first try at lambda=1) - Norm of the internal forces is 4422.12 Norm of the external forces is 1374.49 Norm of the residual is 1622.74 Norm of the solution with Dirichlet BCs is 0.0441466 Number of SNES iteration is 0 Line search: gnorm after quadratic fit 1.622742343614e+03 (I assume that in this line-search iteration 0.05 < lambda < 1, but the corresponding residual is not smaller than the one in the first call, 1622 > 1374) - Norm of the internal forces is 2163.76 Norm of the external forces is 1374.49 Norm of the residual is 1331.88 Norm of the solution with Dirichlet BCs is 0.0220733 Number of SNES iteration is 0 Line search: Cubically determined step, current gnorm 1.331884625811e+03 lambda=5.0003e-02 (This is the accepted lambda for the current Newton iteration) - Norm of the internal forces is 104020 Norm of the external forces is 1374.49 Norm of the residual is 94739 Norm of the solution with Dirichlet BCs is 0.441323 Number of SNES iteration is 1 My question would be, any idea on how to deal with this situation? I would imagine a “hack” would be to bypass the first residual norm, and have the line-search use the following one as the “base norm” to check its reduction in further iterations, but we are open to any ideas. Thanks for your help in advance and please keep up the good work! Regards, Francesc.
Re: [petsc-users] How to apply boundary conditions when using TS?
On Mon, Nov 15, 2021 at 10:24 PM zhfreewill wrote:
> Hi,
> I'm new to PETSc. I have been confused about how to apply boundary
> conditions (BCs) (e.g., the Dirichlet and the Neumann type) when using the
> time-stepping (TS) object of PETSc.
>
> I have read the example codes (
> https://petsc.org/release/src/ts/tutorials/index.html) and vaguely found
> that BCs seems to be set as follows:
>
> 1. when evaluating a matrix of the right hand side function g(x) (e.g.,
> using a function like RHSMatrixHeat)
> 2. when evaluating a matrix of the Jacobian matrix of g'(x) (e.g., using a
> function like FormJacobian)
>
> For instance, in the ex1.4 (
> https://petsc.org/release/src/ts/tutorials/ex4.c.html), with comments at
> the beginning states
>
> 17: /*
>
> 19:This program solves the one-dimensional heat equation (also called
> the
> 20:diffusion equation),
> 21:u_t = u_xx,
> 22:on the domain 0 <= x <= 1, with the boundary conditions
> 23:u(t,0) = 0, u(t,1) = 0,
> 24:and the initial condition
> 25:u(0,x) = sin(6*pi*x) + 3*sin(2*pi*x).
> 26:This is a linear, second-order, parabolic equation.
>
> 28:We discretize the right-hand side using finite differences with
> 29:uniform grid spacing h:
> 30:u_xx = (u_{i+1} - 2u_{i} + u_{i-1})/(h^2)
> ...
> 47:
> - */
>
> and Dirichlet BCs, u(t,0) = 0 and u(t,1) = 0, are to set on the both side
> of the discretized domain, the corresponding codes can be found on the
> following lines in a function RHSMatrixHeat:
>
> /*
> - */
> 463: RHSMatrixHeat - User-provided routine to compute the right-hand-side
> matrix for the heat equation. */
> 491: PetscErrorCode RHSMatrixHeat(TS ts,PetscReal t,Vec X,Mat AA,Mat
> BB,void *ctx)
> 492: {
> ...
> 505: /* Set matrix rows corresponding to boundary data */
> 509: if (mstart == 0) { /* first processor only */
> 510: v[0] = 1.0;
> 511: MatSetValues(A,1,&mstart,1,&mstart,v,INSERT_VALUES);
> 512: mstart++;
> 513: }
> 515: if (mend == appctx->m) { /* last processor only */
> 516: mend--;
> 517: v[0] = 1.0;
> 518: MatSetValues(A,1,&mend,1,&mend,v,INSERT_VALUES);
> 519: }
>
> 521: /* Set matrix rows corresponding to interior data. We construct
> the matrix one row at a time. */
> 525: v[0] = sone; v[1] = stwo; v[2] = sone;
> 526: for (i=mstart; i 527: idx[0] = i-1; idx[1] = i; idx[2] = i+1;
> 528: MatSetValues(A,1,&i,3,idx,v,INSERT_VALUES);
> 529: }
> ...
> 550: }
> /*
> - */
>
>
> My questions are:
>
> 1. How do these lines of code embody the Dirichlet BCs u(t,0) = 0 and
> u(t,1) = 0.
>
This is the finite difference version of BC. It is just moving the term for
a known value to the forcing. However, since
the BC is 0, then we just do not compute the term. If the BC were nonzero,
you would also see the term in the
residual. For example,
https://gitlab.com/petsc/petsc/-/blob/main/src/snes/tutorials/ex19.c#L271
> 2. Similarly, if I want to apply a Neumann type BCs in these lines, How
> can I do?
>
For finite differences, you usually just approximate the derivative as the
one-sided difference of the last two values,
and set this equal to the known value. For finite elements, Neumann
conditions are just an addition weak form
integrated over the boundary.
Thanks,
Matt
> Any suggestions or documents?
>
> Best,
> Qw
>
--
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/
