## Re: [Libmesh-users] how does algebraic method constrain hanging nodes?

 Re: [Libmesh-users] how does algebraic method constrain hanging nodes? From: Benjamin Kirk - 2008-02-13 02:14:08 ```> Now, I am wondering how algebraic method constrains hanging nodes. Whether > does it set the values at hanging nodes to zero? I output the matrix > assembled in libmesh using PetSC function, I find lots of zero values in > matrix. Because the matrix in Petsc is stored using sparse compressed strage > format, zero values should not appear. Could you give me some hints? The constraint u_f = 1/2 u_c1 + 1/2 u_c2 is a typical constraint - it says the fine value is the average of the coarse values on an edge. It can be imposed through the linear system by zeroing the entire row for u_f and the right-hand-side entry and replacing it with [u_c1] [ 0 ... -1/2 1 -1/2 ...0 ] [u_f ] = [0] [u_c2] Strictly speaking, you could eliminate the u_f row & column from the sparse matrix and the resulting linear system would be 1 dof smaller. In practice though that is more trouble than it is worth. It would take a pathological case (like every-other element in a checkerboard pattern is refined) for the savings gained by reducing the size of the linear system to be worth the effort. Now, the reason there are a number of zeros in that row is simple... In 2D for linear quadrilateral elements a typical row will have 9 nonzero entries, and this is what gets allocated for that row. A typical linear constraint only involves 3 dofs, though, so you are left with some zeros. A quadratic Lagrange quadrilateral typically has 25 nonzeros on a row, but a constraint involves 4 dofs. So, there are 21 stored zeros in that row. This functionality is implemented in the DofMap::constrain_element_matrix_an_vector() function -- see http://libmesh.sourceforge.net/doxygen/classDofMap.php It takes the element matrix and vector written in terms of unconstrained dofs and alters them so that they are written entirely in terms of unconstrained dofs. Finally, Is PETSc outputting a compressed matrix, or is it giving you the square matrix and padding with zeros itself? -Ben ```

 [Libmesh-users] how does algebraic method constrain hanging nodes? From: Yujie - 2008-02-13 01:47:46 Attachments: Message as HTML ```Hi, everyone Now, I am wondering how algebraic method constrains hanging nodes. Whether does it set the values at hanging nodes to zero? I output the matrix assembled in libmesh using PetSC function, I find lots of zero values in matrix. Because the matrix in Petsc is stored using sparse compressed strage format, zero values should not appear. Could you give me some hints? thanks a lot. Regards, Yujie ```
 Re: [Libmesh-users] how does algebraic method constrain hanging nodes? From: Roy Stogner - 2008-02-13 02:04:05 ```On Tue, 12 Feb 2008, Yujie wrote: > Now, I am wondering how algebraic method constrains hanging nodes. By default, the matrix line corresponding to each hanging node is the constraint equation which sets its value. This constrains each hanging node to whatever precision your linear solver gives; the hanging node degree of freedom coefficient is then not zero, it's whatever value is correct for the degree of freedom equation on the fine elements. > Whether does it set the values at hanging nodes to zero? With non-default options set, the matrix line for each hanging node is (like the matrix column) just a "1" on the diagonal, with a corresponding 0 on the residual vector, so that the node coefficient is zero immediately after the linear solve. In this case, user code (or higher level library code like NewtonSolver::solve()) must then use the constraint equations to set the hanging node coefficients to their correct values. > I output the matrix assembled in libmesh using PetSC function, I > find lots of zero values in matrix. Because the matrix in Petsc is > stored using sparse compressed strage format, zero values should not > appear. That is incorrect. The sparsity pattern we tell PETSc to use comes from the connectivity of your mesh, which is much more efficient than building it on the fly as you assemble the matrix with your particular equation. Some multi-variable equations (Navier-Stokes, in particular) and/or some choices of basis functions can then leave you with zero entries. If your equation is of the former type you can use DofObject::_dof_coupling to tell libMesh not to bother allocating matrix entries that you know will always be zero. However, it's also possible that we're allocating too many entries in the sparsity pattern. Can you come up with a specific example of two DoFs that you know shouldn't be coupled? --- Roy ```
 Re: [Libmesh-users] how does algebraic method constrain hanging nodes? From: Benjamin Kirk - 2008-02-13 02:14:08 ```> Now, I am wondering how algebraic method constrains hanging nodes. Whether > does it set the values at hanging nodes to zero? I output the matrix > assembled in libmesh using PetSC function, I find lots of zero values in > matrix. Because the matrix in Petsc is stored using sparse compressed strage > format, zero values should not appear. Could you give me some hints? The constraint u_f = 1/2 u_c1 + 1/2 u_c2 is a typical constraint - it says the fine value is the average of the coarse values on an edge. It can be imposed through the linear system by zeroing the entire row for u_f and the right-hand-side entry and replacing it with [u_c1] [ 0 ... -1/2 1 -1/2 ...0 ] [u_f ] = [0] [u_c2] Strictly speaking, you could eliminate the u_f row & column from the sparse matrix and the resulting linear system would be 1 dof smaller. In practice though that is more trouble than it is worth. It would take a pathological case (like every-other element in a checkerboard pattern is refined) for the savings gained by reducing the size of the linear system to be worth the effort. Now, the reason there are a number of zeros in that row is simple... In 2D for linear quadrilateral elements a typical row will have 9 nonzero entries, and this is what gets allocated for that row. A typical linear constraint only involves 3 dofs, though, so you are left with some zeros. A quadratic Lagrange quadrilateral typically has 25 nonzeros on a row, but a constraint involves 4 dofs. So, there are 21 stored zeros in that row. This functionality is implemented in the DofMap::constrain_element_matrix_an_vector() function -- see http://libmesh.sourceforge.net/doxygen/classDofMap.php It takes the element matrix and vector written in terms of unconstrained dofs and alters them so that they are written entirely in terms of unconstrained dofs. Finally, Is PETSc outputting a compressed matrix, or is it giving you the square matrix and padding with zeros itself? -Ben ```
 Re: [Libmesh-users] how does algebraic method constrain hanging nodes? From: Lorenzo Botti - 2008-02-14 00:28:06 ```Il giorno 13/feb/08, alle ore 03:03, Roy Stogner ha scritto: >> I output the matrix assembled in libmesh using PetSC function, I >> find lots of zero values in matrix. Because the matrix in Petsc is >> stored using sparse compressed strage format, zero values should not >> appear. > > That is incorrect. The sparsity pattern we tell PETSc to use comes > from the connectivity of your mesh, which is much more efficient than > building it on the fly as you assemble the matrix with your particular > equation. Some multi-variable equations (Navier-Stokes, in > particular) and/or some choices of basis functions can then leave you > with zero entries. If your equation is of the former type you can use > DofObject::_dof_coupling to tell libMesh not to bother allocating > matrix entries that you know will always be zero. Hi Roy, thanks for the explanations but i have a doubt. I have tried to use DofMap::_dof_coupling but it seems that an issue with amr arises... maybe I miss something obvious... For example if I solve an advection diffusion equation in 3D with a semi-implicit scheme I want to save the memory required to store the Kuv, Kuw, Kvw, (and also Kvu, Kvw, Kwu,) blocks of each of my element matrices. With _dof_coupling I can obtain the right sparsity pattern but then at assembly time I'll need a sparse element matrix. I can instead decide to assemble only Kuu, Kvv, Kww as tree different dense matrices. This works without amr, It saves a lot of memory. If I want to use amr the problem is that the method DofMap::constrain_element_matrix_and_vector wants to constrain all my system variables at the same time preventing me to assemble Kuu, Kvv, Kww as different matrices. Is there a simple solution? Lorenzo ```
 Re: [Libmesh-users] how does algebraic method constrain hanging nodes? From: Roy Stogner - 2008-02-14 01:23:23 ``` On Thu, 14 Feb 2008, Lorenzo Botti wrote: > I have tried to use DofMap::_dof_coupling but it seems that an issue > with amr arises... maybe I miss something obvious... > For example if I solve an advection diffusion equation in 3D with a > semi-implicit scheme I want to save the memory required to store the Kuv, > Kuw, Kvw, (and also Kvu, Kvw, Kwu,) blocks of each of my element matrices. > With _dof_coupling I can obtain the right sparsity pattern but then at > assembly time I'll need a sparse element matrix. I can instead decide to > assemble only Kuu, Kvv, Kww as tree different dense matrices. This works > without amr, It saves a lot of memory. If this saves a lot of memory, you may be doing something wrong. A triquadratic element has 27 degrees of freedom per component. If you're using the default 8-byte floating point numbers, then using three 27x27 matrices would take up about 17kB, whereas using a (27*3)*(27*3) matrix would take up about 51kB. In any significant 3D problem, 51kB should be lost in the noise. But you may be right that the constrain_element* functions are ignoring _dof_coupling; I'll leave it to Ben to check on that. --- Roy ```