Kyle, what you are asking is a very difficult question. If you take a look at the (new) step-86 branch (https://github.com/dealii/dealii/pull/15540), you’ll find a similar approach to what you describe.
Everything you said is correct, with possibly one problem (that I have encountered very recently, and may be related to what you are observing). Boundary conditions and hanging nodes have a very delicate interplay. In particular, when you *build* the constraints, if one of the hanging nodes (dependent) has one of its dependencies on the boundary, then when you add non-homogeneous boundary conditions to the constraints, the hanging node constraints are also modified to take into account the values of the boundary conditions. This computation is done *when you compute the constraints for BC*, and it holds the same values for the rest of the life of the affine constraints. When you `distribute`, the constraints are computed with the BC that you had at the beginning. This information is retained even if you change later your (time dependent) boundary conditions, and may result in the hanging nodes being compute in a wrong manner (in your solution, you may see discontinuous solutions instead of continuous ones where hanging nodes interact with boundary nodes). This affects the final call to “distribute”, not the “distribute_local_to_global” (with false as last argument), since that call condensates away the constraints. So your process is correct, but `distribute` does not do the right job, because the constraints on the hanging nodes have not been updated to take into account the new boundary conditions. To verify that this is the case, try doing the following when you update the boundary conditions: clear the constraints, and recreate them from scratch, and see if this solves your problem. The `distribute_local_to_global` will behave exactly in the same way, but `distribute` should now do the right thing. L. > On 20 Jul 2023, at 22:28, Kyle Schwiebert <kjsch...@mtu.edu> wrote: > > Hello all, > > Thank you for taking the time to look at my question. First, I'll ask a > couple of basic questions about the built-in functions, and then I'll give a > few details of why I ask. > > Does the inhomogeneity-handling call distribute_local_to_global(...,false) do > the following: > > Say that we have a constraint x_i = a. I would guess that the function takes > the following steps: (matrix is M, RHS is right hand side) > 1) Set M(i,i) = L and set RHS(i) = 0, where L is somehow "nice" for the > matrix's spectrum. > 2) Set the rest of the entries of the ith row and column of M = 0. Perhaps > technically we don't bother to zero the column since we plan to call > AffineConstraints::set_zero. > 3) For each DoF j which is coupled to DoF i and is not itself constrained, we > do RHS(j) -= a*M(j,i). > > Now considering the same constraint with instead a call to > distribute_local_to_global(local_rhs,dof_inds,global_rhs) I would guess that > we just set RHS(i) = 0 and do nothing else special. > > Most importantly, is the above correct? > > If so, I am having trouble understanding where I might be going wrong. I'm > solving a time dependent coupled problem with a decoupled solver. By > decoupled solver I mean that information from the second domain only appears > in the right hand side. To make matters one step more complex, one term in my > matrix is time-dependent and so I do not recompute the whole matrix. My > procedure is this: > > 1) Assemble the static part of the matrix with a zero right hand side. This > matrix will never be touched again. The vector, RHS_bc, which will also never > be touched again contains no nonzero entries except for those created by the > call to distribute_local_to_global(...,false). > > 2) I then assemble the time-dependent terms into a separate matrix and > vector, RHS_td, using distribute_local_to_global(...,false) and add the > RHS_bc to RHS_td. > > 3) I then assemble the coupled terms from the other solver into RHS_c with > distribute_local_to_global(local_rhs,dof_inds,RHS_c). > > 4) I then add RHS_c to RHS_td and solve (being careful to call > constraints.set_zero(solution) before the solve and > constraints.distribute(solution) after the solve). I am of course using a > custom class similar to the linear operator classes to bundle the vmult() > outputs of the static and time-dependent matrices. > > However, I am seeing blowup in my solution and I'm seeing some indication > that the inhomogeneous boundary conditions are to blame. Note that I also > tried calling constraints.set_zero(RHS_c) in step 4 above and as I expected > that had no change (since presumably constrained entries are already zero). > > Thank you again for taking the time to engage with my question and in advance > for any tips you are able to offer. > > -Kyle > > -- > The deal.II project is located at http://www.dealii.org/ > For mailing list/forum options, see > https://groups.google.com/d/forum/dealii?hl=en > --- > You received this message because you are subscribed to the Google Groups > "deal.II User Group" group. > To unsubscribe from this group and stop receiving emails from it, send an > email to dealii+unsubscr...@googlegroups.com. > To view this discussion on the web visit > https://groups.google.com/d/msgid/dealii/0589f784-e1ed-4699-bf30-a7358b3acee1n%40googlegroups.com. -- The deal.II project is located at http://www.dealii.org/ For mailing list/forum options, see https://groups.google.com/d/forum/dealii?hl=en --- You received this message because you are subscribed to the Google Groups "deal.II User Group" group. To unsubscribe from this group and stop receiving emails from it, send an email to dealii+unsubscr...@googlegroups.com. To view this discussion on the web visit https://groups.google.com/d/msgid/dealii/38123EFB-6050-4AFD-AC45-E623412CD3BA%40gmail.com.
