Re: [petsc-users] Pre-check before each line search function evaluation
On Thu, Nov 30, 2023 at 5:08 PM Alexander Lindsay
wrote:
> Hi Matt, your derivation is spot on. However, the problem is not linear,
> which is why I am using SNES. So you need to replace
>
> U = A^{-1} f - A^{-1} B L
>
> with
>
> dU = A^{-1} f - A^{-1} B dL
>
I see two cases:
1) There is an easy nonlinear elimination for U. In this case, you do
this to get U_1.
2) There is only a linear elimination. In this case, I don't think the
nonlinear system should be phrased
only on L, but rather on (U, L) itself. The linear elimination can be
used as an excellent preconditioner
for the Newton system.
Thanks,
Matt
> On Thu, Nov 30, 2023 at 1:47 PM Matthew Knepley wrote:
>
>> On Thu, Nov 30, 2023 at 4:23 PM Alexander Lindsay <
>> alexlindsay...@gmail.com> wrote:
>>
>>> If someone passes me just L, where L represents the "global" degrees of
>>> freedom, in this case they represent unknowns on the trace of the mesh,
>>> this is insufficient information for me to evaluate my function. Because in
>>> truth my degrees of freedom are the sum of the trace unknowns (the unknowns
>>> in the global solution vector) and the eliminated unknowns which are
>>> entirely local to each element. So I will say my dofs are L + U.
>>>
>>
>> I want to try and reduce this to the simplest possible thing so that I
>> can understand. We have some system which has two parts to the solution, L
>> and U. If this problem is linear, we have
>>
>> / A B \ / U \ = / f \
>> \ C D / \ L / \ g /
>>
>> and we assume that A is easily invertible, so that
>>
>> U + A^{-1} B L = f
>> U = f - A^{-1} B L
>>
>> C U + D L = g
>> C (f - A^{-1} B L) + D L = g
>> (D - C A^{-1} B) L = g - C f
>>
>> where I have reproduced the Schur complement derivation. Here, given any
>> L, I can construct the corresponding U by inverting A. I know your system
>> may be different, but if you are only solving for L,
>> it should have this property I think.
>>
>> Thus, if the line search generates a new L, say L_1, I should be able to
>> get U_1 by just plugging in. If this is not so, can you write out the
>> equations so we can see why this is not true?
>>
>> Thanks,
>>
>> Matt
>>
>>
>>> I start with some initial guess L0 and U0. I perform a finite element
>>> assembly procedure on each element which gives me things like K_LL, K_UL,
>>> K_LU, K_UU, F_U, and F_L. I can do some math:
>>>
>>> K_LL = -K_LU * K_UU^-1 * K_UL + K_LL
>>> F_L = -K_LU * K_UU^-1 * F_U + F_L
>>>
>>> And then I feed K_LL and F_L into the global system matrix and vector
>>> respectively. I do something (like a linear solve) which gives me an
>>> increment to L, I'll call it dL. I loop back through and do a finite
>>> element assembly again using **L0 and U0** (or one could in theory save off
>>> the previous assemblies) to once again obtain the same K_LL, K_UL, K_LU,
>>> K_UU, F_U, F_L. And now I can compute the increment for U, dU, according to
>>>
>>> dU = K_UU^-1 * (-F_U - K_UL * dL)
>>>
>>> Armed now with both dL and dU, I am ready to perform a new residual
>>> evaluation with (L0 + dL, U0 + dU) = (L1, U1).
>>>
>>> The key part is that I cannot get U1 (or more generally an arbitrary U)
>>> just given L1 (or more generally an arbitrary L). In order to get U1, I
>>> must know both L0 and dL (and U0 of course). This is because at its core U
>>> is not some auxiliary vector; it represents true degrees of freedom.
>>>
>>> On Thu, Nov 30, 2023 at 12:32 PM Barry Smith wrote:
>>>
Why is this all not part of the function evaluation?
> On Nov 30, 2023, at 3:25 PM, Alexander Lindsay <
alexlindsay...@gmail.com> wrote:
>
> Hi I'm looking at the sources, and I believe the answer is no, but is
there a dedicated callback that is akin to SNESLineSearchPrecheck but is
called before *each* function evaluation in a line search method? I am
using a Hybridized Discontinuous Galerkin method in which most of the
degrees of freedom are eliminated from the global system. However, an
accurate function evaluation requires that an update to the "global" dofs
also trigger an update to the eliminated dofs.
>
> I can almost certainly do this manually but I believe it would be
more prone to error than a dedicated callback.
>>
>> --
>> 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/
>>
Re: [petsc-users] Pre-check before each line search function evaluation
Hi Matt, your derivation is spot on. However, the problem is not linear,
which is why I am using SNES. So you need to replace
U = A^{-1} f - A^{-1} B L
with
dU = A^{-1} f - A^{-1} B dL
On Thu, Nov 30, 2023 at 1:47 PM Matthew Knepley wrote:
> On Thu, Nov 30, 2023 at 4:23 PM Alexander Lindsay <
> alexlindsay...@gmail.com> wrote:
>
>> If someone passes me just L, where L represents the "global" degrees of
>> freedom, in this case they represent unknowns on the trace of the mesh,
>> this is insufficient information for me to evaluate my function. Because in
>> truth my degrees of freedom are the sum of the trace unknowns (the unknowns
>> in the global solution vector) and the eliminated unknowns which are
>> entirely local to each element. So I will say my dofs are L + U.
>>
>
> I want to try and reduce this to the simplest possible thing so that I can
> understand. We have some system which has two parts to the solution, L and
> U. If this problem is linear, we have
>
> / A B \ / U \ = / f \
> \ C D / \ L / \ g /
>
> and we assume that A is easily invertible, so that
>
> U + A^{-1} B L = f
> U = f - A^{-1} B L
>
> C U + D L = g
> C (f - A^{-1} B L) + D L = g
> (D - C A^{-1} B) L = g - C f
>
> where I have reproduced the Schur complement derivation. Here, given any
> L, I can construct the corresponding U by inverting A. I know your system
> may be different, but if you are only solving for L,
> it should have this property I think.
>
> Thus, if the line search generates a new L, say L_1, I should be able to
> get U_1 by just plugging in. If this is not so, can you write out the
> equations so we can see why this is not true?
>
> Thanks,
>
> Matt
>
>
>> I start with some initial guess L0 and U0. I perform a finite element
>> assembly procedure on each element which gives me things like K_LL, K_UL,
>> K_LU, K_UU, F_U, and F_L. I can do some math:
>>
>> K_LL = -K_LU * K_UU^-1 * K_UL + K_LL
>> F_L = -K_LU * K_UU^-1 * F_U + F_L
>>
>> And then I feed K_LL and F_L into the global system matrix and vector
>> respectively. I do something (like a linear solve) which gives me an
>> increment to L, I'll call it dL. I loop back through and do a finite
>> element assembly again using **L0 and U0** (or one could in theory save off
>> the previous assemblies) to once again obtain the same K_LL, K_UL, K_LU,
>> K_UU, F_U, F_L. And now I can compute the increment for U, dU, according to
>>
>> dU = K_UU^-1 * (-F_U - K_UL * dL)
>>
>> Armed now with both dL and dU, I am ready to perform a new residual
>> evaluation with (L0 + dL, U0 + dU) = (L1, U1).
>>
>> The key part is that I cannot get U1 (or more generally an arbitrary U)
>> just given L1 (or more generally an arbitrary L). In order to get U1, I
>> must know both L0 and dL (and U0 of course). This is because at its core U
>> is not some auxiliary vector; it represents true degrees of freedom.
>>
>> On Thu, Nov 30, 2023 at 12:32 PM Barry Smith wrote:
>>
>>>
>>> Why is this all not part of the function evaluation?
>>>
>>>
>>> > On Nov 30, 2023, at 3:25 PM, Alexander Lindsay <
>>> alexlindsay...@gmail.com> wrote:
>>> >
>>> > Hi I'm looking at the sources, and I believe the answer is no, but is
>>> there a dedicated callback that is akin to SNESLineSearchPrecheck but is
>>> called before *each* function evaluation in a line search method? I am
>>> using a Hybridized Discontinuous Galerkin method in which most of the
>>> degrees of freedom are eliminated from the global system. However, an
>>> accurate function evaluation requires that an update to the "global" dofs
>>> also trigger an update to the eliminated dofs.
>>> >
>>> > I can almost certainly do this manually but I believe it would be more
>>> prone to error than a dedicated callback.
>>>
>>>
>
> --
> 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/
>
Re: [petsc-users] Pre-check before each line search function evaluation
On Thu, Nov 30, 2023 at 4:23 PM Alexander Lindsay
wrote:
> If someone passes me just L, where L represents the "global" degrees of
> freedom, in this case they represent unknowns on the trace of the mesh,
> this is insufficient information for me to evaluate my function. Because in
> truth my degrees of freedom are the sum of the trace unknowns (the unknowns
> in the global solution vector) and the eliminated unknowns which are
> entirely local to each element. So I will say my dofs are L + U.
>
I want to try and reduce this to the simplest possible thing so that I can
understand. We have some system which has two parts to the solution, L and
U. If this problem is linear, we have
/ A B \ / U \ = / f \
\ C D / \ L / \ g /
and we assume that A is easily invertible, so that
U + A^{-1} B L = f
U = f - A^{-1} B L
C U + D L = g
C (f - A^{-1} B L) + D L = g
(D - C A^{-1} B) L = g - C f
where I have reproduced the Schur complement derivation. Here, given any L,
I can construct the corresponding U by inverting A. I know your system may
be different, but if you are only solving for L,
it should have this property I think.
Thus, if the line search generates a new L, say L_1, I should be able to
get U_1 by just plugging in. If this is not so, can you write out the
equations so we can see why this is not true?
Thanks,
Matt
> I start with some initial guess L0 and U0. I perform a finite element
> assembly procedure on each element which gives me things like K_LL, K_UL,
> K_LU, K_UU, F_U, and F_L. I can do some math:
>
> K_LL = -K_LU * K_UU^-1 * K_UL + K_LL
> F_L = -K_LU * K_UU^-1 * F_U + F_L
>
> And then I feed K_LL and F_L into the global system matrix and vector
> respectively. I do something (like a linear solve) which gives me an
> increment to L, I'll call it dL. I loop back through and do a finite
> element assembly again using **L0 and U0** (or one could in theory save off
> the previous assemblies) to once again obtain the same K_LL, K_UL, K_LU,
> K_UU, F_U, F_L. And now I can compute the increment for U, dU, according to
>
> dU = K_UU^-1 * (-F_U - K_UL * dL)
>
> Armed now with both dL and dU, I am ready to perform a new residual
> evaluation with (L0 + dL, U0 + dU) = (L1, U1).
>
> The key part is that I cannot get U1 (or more generally an arbitrary U)
> just given L1 (or more generally an arbitrary L). In order to get U1, I
> must know both L0 and dL (and U0 of course). This is because at its core U
> is not some auxiliary vector; it represents true degrees of freedom.
>
> On Thu, Nov 30, 2023 at 12:32 PM Barry Smith wrote:
>
>>
>> Why is this all not part of the function evaluation?
>>
>>
>> > On Nov 30, 2023, at 3:25 PM, Alexander Lindsay <
>> alexlindsay...@gmail.com> wrote:
>> >
>> > Hi I'm looking at the sources, and I believe the answer is no, but is
>> there a dedicated callback that is akin to SNESLineSearchPrecheck but is
>> called before *each* function evaluation in a line search method? I am
>> using a Hybridized Discontinuous Galerkin method in which most of the
>> degrees of freedom are eliminated from the global system. However, an
>> accurate function evaluation requires that an update to the "global" dofs
>> also trigger an update to the eliminated dofs.
>> >
>> > I can almost certainly do this manually but I believe it would be more
>> prone to error than a dedicated callback.
>>
>>
--
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/
Re: [petsc-users] Pre-check before each line search function evaluation
I can do exactly what I want using SNESLineSearchPrecheck and -snes_linesearch_type basic ... I just can't use any more exotic line searches On Thu, Nov 30, 2023 at 1:22 PM Alexander Lindsay wrote: > If someone passes me just L, where L represents the "global" degrees of > freedom, in this case they represent unknowns on the trace of the mesh, > this is insufficient information for me to evaluate my function. Because in > truth my degrees of freedom are the sum of the trace unknowns (the unknowns > in the global solution vector) and the eliminated unknowns which are > entirely local to each element. So I will say my dofs are L + U. > > I start with some initial guess L0 and U0. I perform a finite element > assembly procedure on each element which gives me things like K_LL, K_UL, > K_LU, K_UU, F_U, and F_L. I can do some math: > > K_LL = -K_LU * K_UU^-1 * K_UL + K_LL > F_L = -K_LU * K_UU^-1 * F_U + F_L > > And then I feed K_LL and F_L into the global system matrix and vector > respectively. I do something (like a linear solve) which gives me an > increment to L, I'll call it dL. I loop back through and do a finite > element assembly again using **L0 and U0** (or one could in theory save off > the previous assemblies) to once again obtain the same K_LL, K_UL, K_LU, > K_UU, F_U, F_L. And now I can compute the increment for U, dU, according to > > dU = K_UU^-1 * (-F_U - K_UL * dL) > > Armed now with both dL and dU, I am ready to perform a new residual > evaluation with (L0 + dL, U0 + dU) = (L1, U1). > > The key part is that I cannot get U1 (or more generally an arbitrary U) > just given L1 (or more generally an arbitrary L). In order to get U1, I > must know both L0 and dL (and U0 of course). This is because at its core U > is not some auxiliary vector; it represents true degrees of freedom. > > On Thu, Nov 30, 2023 at 12:32 PM Barry Smith wrote: > >> >> Why is this all not part of the function evaluation? >> >> >> > On Nov 30, 2023, at 3:25 PM, Alexander Lindsay < >> alexlindsay...@gmail.com> wrote: >> > >> > Hi I'm looking at the sources, and I believe the answer is no, but is >> there a dedicated callback that is akin to SNESLineSearchPrecheck but is >> called before *each* function evaluation in a line search method? I am >> using a Hybridized Discontinuous Galerkin method in which most of the >> degrees of freedom are eliminated from the global system. However, an >> accurate function evaluation requires that an update to the "global" dofs >> also trigger an update to the eliminated dofs. >> > >> > I can almost certainly do this manually but I believe it would be more >> prone to error than a dedicated callback. >> >>
Re: [petsc-users] Pre-check before each line search function evaluation
If someone passes me just L, where L represents the "global" degrees of freedom, in this case they represent unknowns on the trace of the mesh, this is insufficient information for me to evaluate my function. Because in truth my degrees of freedom are the sum of the trace unknowns (the unknowns in the global solution vector) and the eliminated unknowns which are entirely local to each element. So I will say my dofs are L + U. I start with some initial guess L0 and U0. I perform a finite element assembly procedure on each element which gives me things like K_LL, K_UL, K_LU, K_UU, F_U, and F_L. I can do some math: K_LL = -K_LU * K_UU^-1 * K_UL + K_LL F_L = -K_LU * K_UU^-1 * F_U + F_L And then I feed K_LL and F_L into the global system matrix and vector respectively. I do something (like a linear solve) which gives me an increment to L, I'll call it dL. I loop back through and do a finite element assembly again using **L0 and U0** (or one could in theory save off the previous assemblies) to once again obtain the same K_LL, K_UL, K_LU, K_UU, F_U, F_L. And now I can compute the increment for U, dU, according to dU = K_UU^-1 * (-F_U - K_UL * dL) Armed now with both dL and dU, I am ready to perform a new residual evaluation with (L0 + dL, U0 + dU) = (L1, U1). The key part is that I cannot get U1 (or more generally an arbitrary U) just given L1 (or more generally an arbitrary L). In order to get U1, I must know both L0 and dL (and U0 of course). This is because at its core U is not some auxiliary vector; it represents true degrees of freedom. On Thu, Nov 30, 2023 at 12:32 PM Barry Smith wrote: > > Why is this all not part of the function evaluation? > > > > On Nov 30, 2023, at 3:25 PM, Alexander Lindsay > wrote: > > > > Hi I'm looking at the sources, and I believe the answer is no, but is > there a dedicated callback that is akin to SNESLineSearchPrecheck but is > called before *each* function evaluation in a line search method? I am > using a Hybridized Discontinuous Galerkin method in which most of the > degrees of freedom are eliminated from the global system. However, an > accurate function evaluation requires that an update to the "global" dofs > also trigger an update to the eliminated dofs. > > > > I can almost certainly do this manually but I believe it would be more > prone to error than a dedicated callback. > >
Re: [petsc-users] Pre-check before each line search function evaluation
Why is this all not part of the function evaluation? > On Nov 30, 2023, at 3:25 PM, Alexander Lindsay > wrote: > > Hi I'm looking at the sources, and I believe the answer is no, but is there a > dedicated callback that is akin to SNESLineSearchPrecheck but is called > before *each* function evaluation in a line search method? I am using a > Hybridized Discontinuous Galerkin method in which most of the degrees of > freedom are eliminated from the global system. However, an accurate function > evaluation requires that an update to the "global" dofs also trigger an > update to the eliminated dofs. > > I can almost certainly do this manually but I believe it would be more prone > to error than a dedicated callback.
[petsc-users] Pre-check before each line search function evaluation
Hi I'm looking at the sources, and I believe the answer is no, but is there a dedicated callback that is akin to SNESLineSearchPrecheck but is called before *each* function evaluation in a line search method? I am using a Hybridized Discontinuous Galerkin method in which most of the degrees of freedom are eliminated from the global system. However, an accurate function evaluation requires that an update to the "global" dofs also trigger an update to the eliminated dofs. I can almost certainly do this manually but I believe it would be more prone to error than a dedicated callback.
Re: [petsc-users] Reading VTK files in PETSc
I assume you're working with a DA, in which case you can write in HDF5 format and add an Xdmf header so Paraview can read it. The utility lib/petsc/bin/petsc_gen_xdmf.py should be able to handle it. I haven't written support for it (my problems are unstructured so I've focused on DMPlex), but the CGNS format supports structured meshes and that would be an efficient parallel format that doesn't need a header. "Kevin G. Wang" writes: > Hi Jed, > > Thanks for your help! It does not have to be VTK (.vtr in my case). But I > am trying to have a "seamless" workflow between reading, writing, and > visualization, without the need of format conversions. I opted for VTK > simply because it is easy to write, and can be directly visualized (using > Paraview). > > Could you please advise as to what is the best practice in my scenario? My > meshes are Cartesian, but non-uniform. > > Thanks, > Kevin > > On Thu, Nov 30, 2023 at 1:02 AM Jed Brown wrote: > >> Is it necessary that it be VTK format or can it be PETSc's binary format >> or a different mesh format? VTK (be it legacy .vtk or the XML-based .vtu, >> etc.) is a bad format for parallel reading, no matter how much effort might >> go into an implementation. >> >> "Kevin G. Wang" writes: >> >> > Good morning everyone. >> > >> > I use the following functions to output parallel vectors --- "globalVec" >> in >> > this example --- to VTK files. It works well, and is quite convenient. >> > >> > ~ >> > PetscViewer viewer; >> > PetscViewerVTKOpen(PetscObjectComm((PetscObject)*dm), filename, >> > FILE_MODE_WRITE, ); >> > VecView(globalVec, viewer); >> > PetscViewerDestroy(); >> > ~ >> > >> > Now, I am trying to do the opposite. I would like to read the VTK files >> > generated by PETSc back into memory, and assign each one to a Vec. Could >> > someone let me know how this can be done? >> > >> > Thanks! >> > Kevin >> > >> > >> > -- >> > Kevin G. Wang, Ph.D. >> > Associate Professor >> > Kevin T. Crofton Department of Aerospace and Ocean Engineering >> > Virginia Tech >> > 1600 Innovation Dr., VTSS Rm 224H, Blacksburg, VA 24061 >> > Office: (540) 231-7547 | Mobile: (650) 862-2663 >> > URL: https://www.aoe.vt.edu/people/faculty/wang.html >> > Codes: https://github.com/kevinwgy >> > > > -- > Kevin G. Wang, Ph.D. > Associate Professor > Kevin T. Crofton Department of Aerospace and Ocean Engineering > Virginia Tech > 1600 Innovation Dr., VTSS Rm 224H, Blacksburg, VA 24061 > Office: (540) 231-7547 | Mobile: (650) 862-2663 > URL: https://www.aoe.vt.edu/people/faculty/wang.html > Codes: https://github.com/kevinwgy
