Dear Daniel
After solving constraints_update.distribute is called. In addition, after assembling the system matrix and residual vector using the first method, constraints_update.condense is called before solving the system. After that the solution is updated and the new residual is calculated based on the first method, so in error calculation the constraints are considered as, BlockVector<double> <https://dealii.org/developer/doxygen/deal.II/classBlockVector.html> error_res(dofs_per_block); for (unsigned int i = 0; i < dof_handler.n_dofs <https://dealii.org/developer/doxygen/deal.II/classDoFHandler.html#a54e7270f2ba6206604f794114b39a2aa>(); ++i) if (!constraints_update.is_constrained <https://dealii.org/developer/doxygen/deal.II/classConstraintMatrix.html#aac0dd11d56f41826d958ade7faa0cfe4> (i)) error_res(i) = system_rhs(i); error_residual.norm = error_res.l2_norm(); constraints update 1 = 0 3 = 0 11 = 0 21 = 0 22 10: 5.00000e-01 22 14: 5.00000e-01 23 15: 5.00000e-01 24 14: 5.00000e-01 24 16: 5.00000e-01 25 15: 5.00000e-01 25 17: 5.00000e-01 26 16: 5.00000e-01 26 40: 5.00000e-01 27 17: 5.00000e-01 27 41: 5.00000e-01 42 = 0 43 = 1.00000e-04 45 = 1.00000e-04 46 = 0 49 = 1.00000e-04 51 = 1.00000e-04 53 = 1.00000e-04 55 = 0 56 40: 5.00000e-01 57 41: 5.00000e-01 57 59: 5.00000e-01 58 = 0 69 = 0 70 = 0 72 = 0 78 = 0 79 = 0 80 = 0 93 87: 5.00000e-01 93 89: 5.00000e-01 94 89: 5.00000e-01 94 90: 5.00000e-01 95 90: 5.00000e-01 95 102: 5.00000e-01 110 102: 5.00000e-01 110 111: 5.00000e-01 constraints hanging nodes 22 10: 5.00000e-01 22 14: 5.00000e-01 23 11: 5.00000e-01 23 15: 5.00000e-01 24 14: 5.00000e-01 24 16: 5.00000e-01 25 15: 5.00000e-01 25 17: 5.00000e-01 26 16: 5.00000e-01 26 40: 5.00000e-01 27 17: 5.00000e-01 27 41: 5.00000e-01 56 40: 5.00000e-01 56 58: 5.00000e-01 57 41: 5.00000e-01 57 59: 5.00000e-01 93 87: 5.00000e-01 93 89: 5.00000e-01 94 89: 5.00000e-01 94 90: 5.00000e-01 95 90: 5.00000e-01 95 102: 5.00000e-01 110 102: 5.00000e-01 110 111: 5.00000e-01 Thank you. Pasha On Thursday, February 2, 2017 at 2:14:34 PM UTC+3:30, Daniel Arndt wrote: > > Pasha, > > >> I checked the the system residual assembly process. For a triangulation >> with 28 cells, the dof index 17 is belong to cells 2, 5, 8 and 27. The >> corresponding cell_rhs values for this index are, >> >> >> cell cell_rhs value >> >> 2 2.88462e+00 >> >> 5 2.88462e+00 >> >> 8 -2.88462e+00 >> >> 27 -1.44231e+00 >> >> >> This dof is not constrained. When the following assembly method is used: >> >> >> *for* (*unsigned* *int* i = 0; i < dofs_per_cell; ++i) >> >> system_rhs(data.local_dof_indices[i]) += data.cell_rhs(i); >> >> >> the system_rhs(17) = 1.4423076, however, using the following method, >> >> >> constraints_update.distribute_local_to_global(data.cell_rhs, >> >> data.local_dof_indices, system_rhs); >> >> >> >> the system_rhs(17) = -1.33e-15. This is the source of the problem because >> I have used the second method in the previous code and the first one >> (similar to step-44) in the new code. I cant understand why the results are >> different. >> > This looks like DoF 17 is constrained, possibly by hanging node > constraints. Can you recheck the constraints stored in `constraints_update` > via `constraints_update.print(std::cout)`? > In any case, you need to use a ConstraintMatrix for assembling the > right-hand side either via calling `distribute_local_to_global` on each > cell or calling `condense` after the loop over all cells. After solving you > would call `constraints.distribute`. > Otherwise, you would not take the hanging node constraints into account. > Are you missing the call to `constraints.condense(system_matrix, > system_rhs);` in the first approach? > > Best, > Daniel > -- 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 [email protected]. For more options, visit https://groups.google.com/d/optout.
