Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Should d2udt2 be set to 0 when Dirichlet boundary condition is applied? #4181

Open
Wayne901 opened this issue Mar 8, 2024 · 4 comments
Open

Comments

@Wayne901
Copy link

Wayne901 commented Mar 8, 2024

I have a little suggestion/question about ex23. When the dirichlet boundary conditions are applied, dudt and d2udt2 should be 0. However when I run the wave euqation example ex23 and output the d2udt2 vector in ImplicitSolver function, d2udt2 is not zero at the boundary DoFs. I set d2udt2's boundary entries to 0, and the final result is slightly different to the original code (./ex23 -m ../../data/disc-nurbs.mesh -r 3 -o 4 -tf 5 --dirichlet)
Result from original code:
Wave_d2dut2_notset
Result when boundary entris of d2udt2 set to 0:
Wave_d2dut2_setzero

I think the nonzero boundary DoFs in d2udt2 comes from solving the linear system
$[M+\beta (\Delta t)^2K]\frac{d^2u(t+\Delta t)}{dt^2}=-Ku_{ImplicitSolve}$? But why would the solution on boundary entries be nonzero when $u_{ImplicitSolve}$'s are 0 there?

Please let me know if my choice of setting d2udt2's boundary DoFs 0 is correct.

Thank

@Wayne901
Copy link
Author

May someone give some hint on this problem? I'm still confused about it. The different result is not observed when t_final <= 4. I set d2udt2 to 0 after T_solver with

T_solver.Mult(z, d2udt2);
for (int i = 0; i < ess_tdof_list.Size(); i++)
  d2udt2[ess_tdof_list[i]] = 0.0;

@Wayne901
Copy link
Author

I wrote another wave equation code based on the coupled first-order equations, like ex10.
The result at t=5 is same as the ex23 with d2udt2 set to 0.

@IdoAkkerman
Copy link
Contributor

I have tried to address the issue in a small PR #4233. I added as step after solving for the acceleration where the acceleration is hardwired set to zero.

I think, this should technically not have been necessary, as the solve should have resulted in zero accelerations on those places. But apparently round-off and/or the approximate nature of the linear solves create for excessive slip of the solution.

@IdoAkkerman
Copy link
Contributor

I have checked the BC values and after zeroing the accelerations the solution values at the BC were of the order e-105 which is effectively zero.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

2 participants