## Re: [Libmesh-devel] User specified adjoint boundary conditions

 Re: [Libmesh-devel] User specified adjoint boundary conditions From: Roy Stogner - 2011-02-10 20:12:53 ``` >> Dirichlet(penalty)->non-Dirichlet is even worse.  Even if you edit the >> matrix again, floating point error makes it impossible to just >> subtract off the penalty terms; x + O(1/epsilon) - O(1/epsilon) isn't >> x, it's garbage. > > Cant we just reinit those edges ? No, thanks to our sparse-matrix-based data structure and element+side based assembly loops. If you throw out the coefficient connecting two side DoFs, you can't just reinit the side, you've got to reinit the element (once the BC becomes non-Dirichlet, the interior shape function based coupling of side DoFs becomes important again). But that would invalidate the coefficients connecting interior dofs, which also have support on other elements, etc. --- Roy```

 [Libmesh-devel] User specified adjoint boundary conditions From: Vikram Garg - 2011-02-09 17:30:50 ```All, Currently our adjoints framework assembles the discrete adjoint problem by transposing the matrix representing the discrete forward problem. This is simple and works as long as we have (or only desire) the discrete adjoint for an adjoint consistent formulation. Of course, this might not be the case for all problems. In particular, while computing boundary QoIs that are ill-posed due to regularity issues, one might want to follow the approach described in a paper by Giles, Larson et al. Here one imposes a Dirichlet bc for the adjoint problem on the boundary in question, modifying the trial space. This ensures that we are solving a well posed adjoint problem. One way of extending the current framework to include such adjoint formulations would be to have the user specify a side_qoi_constraint along the lines of side_qoi_derivative and side_constraint. The user would specify the stiffness matrix entries (penalty terms) for the Dirichlet bc in side_qoi_residual and the corresponding residual term in side_qoi_derivative. We can then include side_qoi_constraint in the assembly before solving the adjoint. Does this sound reasonable ? Of course any other ideas are welcome. Thanks. -- Vikram Garg PhD Candidate Institute for Computational and Engineering Sciences The University of Texas at Austin http://users.ices.utexas.edu/~vikram/ http://www.runforindia.org/runners/vikramg ```
 Re: [Libmesh-devel] User specified adjoint boundary conditions From: Roy Stogner - 2011-02-10 17:25:33 ```On Wed, 9 Feb 2011, Vikram Garg wrote: > Currently our adjoints framework assembles the discrete adjoint > problem by transposing the matrix representing the discrete forward > problem. This is simple and works as long as we have (or only desire) > the discrete adjoint for an adjoint consistent formulation. Of course, > this might not be the case for all problems. > > In particular, while computing boundary QoIs that are ill-posed due to > regularity issues, one might want to follow the approach described in > a paper by Giles, Larson et al. Here one imposes a Dirichlet bc for > the adjoint problem on the boundary in question, modifying the trial > space. This ensures that we are solving a well posed adjoint problem. The problem class I'm most familiar with is the one we looked at when first starting with the automatic adjoints code, where the discrete adjoint solution looked ill-posed but surprisingly worked, the more traditional problem would have replaced the automatic homogeneous Dirichlet boundary condition with a heterogeneous Dirichlet B.C. Worst case I know of in these problems, the automatic adjoint does fail but the homogeneous->heterogeneous Dirichlet tricks like in that Giles/Larson paper look similar. That's (in situations where we need it) doable within the current framework, for anyone using penalty Dirichlet BCs, no? You've already got the correct penalty matrix, and you just need to add the proper penalty forcing term to your qoi_derivative. We'll have to worry again when we extend the library: when we add an official in-library method for strict Dirichlet BCs we'll need to remember to allow optional separate Dirichlet BCs to be used on adjoint problems, and when we do transient adjoints then the penalty really will need to be in a constraint-type rather than a time-derivative-type term. But what I'm worried about NOW: are there any cases where we'd want more complicated boundary condition changes? Any Neumann/mixed primal -> Dirichlet adjoint, or Dirichlet primal -> Neumann/mixed adjoint is pretty much impossible right now, for example; any chance we'd be painting ourselves into a corner if a feature designed now wasn't flexible enough to handle those in the future? > One way of extending the current framework to include such adjoint > formulations would be to have the user specify a side_qoi_constraint > along the lines of side_qoi_derivative and side_constraint. The user > would specify the stiffness matrix entries (penalty terms) for the > Dirichlet bc in side_qoi_residual and the corresponding residual term > in side_qoi_derivative. We can then include side_qoi_constraint in the > assembly before solving the adjoint. > > Does this sound reasonable ? Of course any other ideas are welcome. The only additional issue with side_qoi_constraint (I'd name it side_adjoint_constraint, since it wouldn't end up being called in qoi evaluations) only makes sense for the DiffSystem-based adjoint. The rest of the adjoint features work for LinearImplicitSystem and NonlinearImplicitSystem based codes as well, and I'd rather any new feature work as widely, if possible. On the other hand, side_qoi_constraint might be the proper way to do this for transient adjoint problems, which are going to be difficult-to-impossible to do in TransientSystem anyway, so maybe the thing to do is to get the DiffSystem API sorted out now and leave the other systems stuck with qoi_derivative hacks. --- Roy ```
 Re: [Libmesh-devel] User specified adjoint boundary conditions From: Vikram Garg - 2011-02-10 18:03:44 ```On Thu, Feb 10, 2011 at 11:25 AM, Roy Stogner wrote: > > On Wed, 9 Feb 2011, Vikram Garg wrote: > >>    Currently our adjoints framework assembles the discrete adjoint >> problem by transposing the matrix representing the discrete forward >> problem. This is simple and works as long as we have (or only desire) >> the discrete adjoint for an adjoint consistent formulation. Of course, >> this might not be the case for all problems. >> >> In particular, while computing boundary QoIs that are ill-posed due to >> regularity issues, one might want to follow the approach described in >> a paper by Giles, Larson et al. Here one imposes a Dirichlet bc for >> the adjoint problem on the boundary in question, modifying the trial >> space. This ensures that we are solving a well posed adjoint problem. > > The problem class I'm most familiar with is the one we looked at when > first starting with the automatic adjoints code, where the discrete > adjoint solution looked ill-posed but surprisingly worked, the more > traditional problem would have replaced the automatic homogeneous > Dirichlet boundary condition with a heterogeneous Dirichlet B.C. > Worst case I know of in these problems, the automatic adjoint does > fail but the homogeneous->heterogeneous Dirichlet tricks like in that > Giles/Larson paper look similar. > > That's (in situations where we need it) doable within the current > framework, for anyone using penalty Dirichlet BCs, no?  You've already > got the correct penalty matrix, and you just need to add the proper > penalty forcing term to your qoi_derivative. > Not always. Think about Navier Stokes in a channel. There is no Dirichlet bc on the pressure on the channel walls. But if our QoI involved the pressure on that boundary, then the approach of Giles et al would require for us to apply a Dirichlet condition on the adjoint pressure. But with the current setup we wouldnt have the penalty matrix for the pressure on that boundary. > We'll have to worry again when we extend the library: when we add an > official in-library method for strict Dirichlet BCs we'll need to > remember to allow optional separate Dirichlet BCs to be used on > adjoint problems, and when we do transient adjoints then the penalty > really will need to be in a constraint-type rather than a > time-derivative-type term. > > But what I'm worried about NOW: are there any cases where we'd want > more complicated boundary condition changes?  Any Neumann/mixed primal > -> Dirichlet adjoint, or Dirichlet primal -> Neumann/mixed adjoint is > pretty much impossible right now, for example; any chance we'd be > painting ourselves into a corner if a feature designed now wasn't > flexible enough to handle those in the future? > Yeah, I thought about this as well. But wouldnt the side_adjoint_constraint and side_adjoint_derivative be able to handle those cases ? >> One way of extending the current framework to include such adjoint >> formulations would be to have the user specify a side_qoi_constraint >> along the lines of side_qoi_derivative and side_constraint. The user >> would specify the stiffness matrix entries (penalty terms) for the >> Dirichlet bc in side_qoi_residual and the corresponding residual term >> in side_qoi_derivative. We can then include side_qoi_constraint in the >> assembly before solving the adjoint. >> >> Does this sound reasonable ? Of course any other ideas are welcome. > > The only additional issue with side_qoi_constraint (I'd name it > side_adjoint_constraint, since it wouldn't end up being called in qoi > evaluations) only makes sense for the DiffSystem-based adjoint.  The > rest of the adjoint features work for LinearImplicitSystem and > NonlinearImplicitSystem based codes as well, and I'd rather any new > feature work as widely, if possible. > > On the other hand, side_qoi_constraint might be the proper way to do > this for transient adjoint problems, which are going to be > difficult-to-impossible to do in TransientSystem > anyway, so maybe the thing to do is to get the DiffSystem API sorted > out now and leave the other systems stuck with qoi_derivative hacks. I would say lets get the DiffSystem API sorted out first. Unless there are lots of people who are wanting to use adjoints (with user specified boundary conditions) and are using non DiffSystem frameworks. > --- > Roy Thanks. -- Vikram Garg PhD Candidate Institute for Computational and Engineering Sciences The University of Texas at Austin http://users.ices.utexas.edu/~vikram/ http://www.runforindia.org/runners/vikramg ```
 Re: [Libmesh-devel] User specified adjoint boundary conditions From: Roy Stogner - 2011-02-10 18:31:17 ```On Thu, 10 Feb 2011, Vikram Garg wrote: > Not always. Think about Navier Stokes in a channel. There is no > Dirichlet bc on the pressure on the channel walls. But if our QoI > involved the pressure on that boundary, then the approach of Giles et > al would require for us to apply a Dirichlet condition on the adjoint > pressure. But with the current setup we wouldnt have the penalty > matrix for the pressure on that boundary. Yeah, that would fall into one of my "what I'm worried about NOW" categories... except that one's especially conceptually tricky, isn't it, due to the different function spaces involved? Do the Giles et al tricks work the same way? Going from H^1 to an affine transform of H^1_0 is a smaller jump than going there from L_2. >> But what I'm worried about NOW: are there any cases where we'd want >> more complicated boundary condition changes?  Any Neumann/mixed primal >> -> Dirichlet adjoint, or Dirichlet primal -> Neumann/mixed adjoint is >> pretty much impossible right now, for example; any chance we'd be >> painting ourselves into a corner if a feature designed now wasn't >> flexible enough to handle those in the future? > > Yeah, I thought about this as well. But wouldnt the > side_adjoint_constraint and side_adjoint_derivative be able to handle > those cases ? non-Dirichlet->Dirichlet we might be able to handle with penalty terms using adjoint_constraint equations that, unlike the qoi_derivative stuff, add to both the system rhs and matrix.... but then this trashes the matrix, for the more complicated adjoint-based Hessian and Hessian-vector tricks that rely on multiple adjoint + forward sensitivity solves. Dirichlet(penalty)->non-Dirichlet is even worse. Even if you edit the matrix again, floating point error makes it impossible to just subtract off the penalty terms; x + O(1/epsilon) - O(1/epsilon) isn't x, it's garbage. > I would say lets get the DiffSystem API sorted out first. "First" is probably appropriate; that's basically the way we did it with all the other adjoint APIs. But we don't want to stop thinking ahead to how they can be generalized. > Unless there are lots of people who are wanting to use adjoints > (with user specified boundary conditions) and are using non > DiffSystem frameworks. Not sure about the user-specified boundary conditions, but "wanting to use adjoints while using non-DiffSystem frameworks" is true for all the INL-based folks; loosen that very slightly to "wanting to use adjoints even when stuck with non-DiffSystem frameworks" and it's emphatically true for me too. ;-) --- Roy```
 Re: [Libmesh-devel] User specified adjoint boundary conditions From: Vikram Garg - 2011-02-10 19:47:39 ```On Thu, Feb 10, 2011 at 12:31 PM, Roy Stogner wrote: > > On Thu, 10 Feb 2011, Vikram Garg wrote: > >> Not always. Think about Navier Stokes in a channel. There is no >> Dirichlet bc on the pressure on the channel walls. But if our QoI >> involved the pressure on that boundary, then the approach of Giles et >> al would require for us to apply a Dirichlet condition on the adjoint >> pressure. But with the current setup we wouldnt have the penalty >> matrix for the pressure on that boundary. > > Yeah, that would fall into one of my "what I'm worried about NOW" > categories... except that one's especially conceptually tricky, isn't > it, due to the different function spaces involved?  Do the Giles et al > tricks work the same way?  Going from H^1 to an affine transform of > H^1_0 is a smaller jump than going there from L_2. > Well the trick works for grad u which is in the same function space as p. So there should be a way to make this work. The QoIs in the paper are lift and drag, you write out the stress tensor for the Stokes equation which of course includes both the velocity gradient and the pressure tensors, and end up with an analogous case to the earlier scalar example. I think we just need to pick the right function spaces to isolate the pressure QoI of our interest. >>> But what I'm worried about NOW: are there any cases where we'd want >>> more complicated boundary condition changes?  Any Neumann/mixed primal >>> -> Dirichlet adjoint, or Dirichlet primal -> Neumann/mixed adjoint is >>> pretty much impossible right now, for example; any chance we'd be >>> painting ourselves into a corner if a feature designed now wasn't >>> flexible enough to handle those in the future? >> >> Yeah, I thought about this as well. But wouldnt the >> side_adjoint_constraint and side_adjoint_derivative be able to handle >> those cases ? > > non-Dirichlet->Dirichlet we might be able to handle with penalty terms > using adjoint_constraint equations that, unlike the qoi_derivative > stuff, add to both the system rhs and matrix.... but then this trashes > the matrix, for the more complicated adjoint-based Hessian and > Hessian-vector tricks that rely on multiple adjoint + forward > sensitivity solves. > > Dirichlet(penalty)->non-Dirichlet is even worse.  Even if you edit the > matrix again, floating point error makes it impossible to just > subtract off the penalty terms; x + O(1/epsilon) - O(1/epsilon) isn't > x, it's garbage. > Cant we just reinit those edges ? Even if we have to reassemble from the sides all over again it shouldnt cost us much. I believe all we'll have to do is to add a system.assemble_sides for the forward problem ? >> I would say lets get the DiffSystem API sorted out first. > > "First" is probably appropriate; that's basically the way we did it > with all the other adjoint APIs.  But we don't want to stop thinking > ahead to how they can be generalized. > I havent really used non DiffSystem based apps for a while now. I'll need to talk to you a little bit before I can contribute meaningfully to generalizing to non-DiffSystem apps. >> Unless there are lots of people who are wanting to use adjoints >> (with user specified boundary conditions) and are using non >> DiffSystem frameworks. > > Not sure about the user-specified boundary conditions, but "wanting to > use adjoints while using non-DiffSystem frameworks" is true for all > the INL-based folks; loosen that very slightly to "wanting to use > adjoints even when stuck with non-DiffSystem frameworks" and it's > emphatically true for me too.  ;-) > --- > Roy Thanks. -- Vikram Garg PhD Candidate Institute for Computational and Engineering Sciences The University of Texas at Austin http://users.ices.utexas.edu/~vikram/ http://www.runforindia.org/runners/vikramg ```
 Re: [Libmesh-devel] User specified adjoint boundary conditions From: Roy Stogner - 2011-02-10 20:12:53 ``` >> Dirichlet(penalty)->non-Dirichlet is even worse.  Even if you edit the >> matrix again, floating point error makes it impossible to just >> subtract off the penalty terms; x + O(1/epsilon) - O(1/epsilon) isn't >> x, it's garbage. > > Cant we just reinit those edges ? No, thanks to our sparse-matrix-based data structure and element+side based assembly loops. If you throw out the coefficient connecting two side DoFs, you can't just reinit the side, you've got to reinit the element (once the BC becomes non-Dirichlet, the interior shape function based coupling of side DoFs becomes important again). But that would invalidate the coefficients connecting interior dofs, which also have support on other elements, etc. --- Roy```