Hi Michael, Thank you for the suggestions. I will look at the filtered matrix class. I understand that distribute_local_to_global() usually requires matrix and righthand side at the same time. That's why I need an auxiliary matrix matrix_for_bc. I was basically following tutorial-34 and my Poisson test is much simpler than the hydrodynamic system solved there. So I guess I am using ConstraintMatrix in a legitimate way ... somehow something goes wrong.
Also I have one more question on the inhomogeneous boundary condition: in this case is VectorTools::project() forbidden? Thanks, Huan On Tue, 2010-05-18 at 04:12 +0200, Michael Rapson wrote: > Hi Huan, > > I also sometimes separate the matrix and right hand side assembly. > However, care must be taken because the method > distribute_local_to_global uses both the right hand side matrix and > the system matrix. Basically, the algorithm scales the Dirichlet > boundary condition relative to the size of the terms in the matrix. In > generally it will not work to apply the method to matrix and right > hand side separately (hence the method which takes both) or even to > apply the method again on a matrix that has had some of the terms > assembled already. > > One option in the library is to use the FilteredMatrix class. (Which > should work as long as you don't need to add terms to the matrix.) > Another is to keep a copy of the matrix without the boundary > conditions applied which you keep initializing from, completing with > any extra terms and then use to apply the boundary conditions to. > > Regards, > Michael > > On Mon, May 17, 2010 at 8:18 PM, Huan Sun <[email protected]> wrote: > > Hi Markus, > > > > Yes, it works to use the original local stiffness matrix. But I want to > > later separate the assemblage of the matrix and the right-hand side, > > similar to what tutorial 32 does. So I was testing this easier case. > > Somehow it doesn't work with this auxiliary matrix. I am just wondering > > if there is anything wrong or missing in the code. > > > > Thanks, > > > > Huan > > > > > > On Mon, 2010-05-17 at 08:29 +0200, Markus Bürg wrote: > >> Hello Huan Sun, > >> > >> why do you create another matrix for the boundary conditions? You have > >> already created the problem's matrix for this cell right above. Thus > >> you can hand over just this matrix: > >> constraints.distribute_local_to_global (local_rhs, local_dof_indices, > >> rhs, local_matrix); > >> > >> Best Regards, > >> Markus > >> > >> > >> > >> Am 17.05.10 08:12, schrieb Huan Sun: > >> > Hi all, > >> > > >> > I was following tutorial 22 and 34 to implement inhomogeneous Dirichlet > >> > boundary conditions using the ConsraintMatrix class but couldn't get it > >> > right. I think the problem mainly lies in the assemblage part (i was > >> > solving the basic Poisson's eq.) > >> > > >> > if(constraints.is_inhomogeneously_constrained > >> > (local_dof_indices[i])){ > >> > for(unsigned int j=0; j <dofs_per_cell; ++j){ > >> > matrix_for_bc(j,i) += > >> > grad_basis_phi[i] * grad_basis_phi[j] * > >> > fe_values.JxW(q); > >> > }//j-loop > >> > }//if > >> > ... > >> > constraints.distribute_local_to_global(local_rhs, > >> > local_dof_indices, > >> > rhs, > >> > matrix_for_bc); > >> > > >> > However, when I comment out the if statement, the program seems to work > >> > correctly. The code is attached in this email. > >> > Your help will be greatly appreciated. > >> > > >> > Thanks! > >> > > >> > Huan > >> > > >> > > >> > _______________________________________________ > >> > dealii mailing list http://poisson.dealii.org/mailman/listinfo/dealii > >> > > >> _______________________________________________ > >> dealii mailing list http://poisson.dealii.org/mailman/listinfo/dealii > > > > > > _______________________________________________ > > dealii mailing list http://poisson.dealii.org/mailman/listinfo/dealii > > _______________________________________________ dealii mailing list http://poisson.dealii.org/mailman/listinfo/dealii
