Thanks Markus, I have all different codes using "fe_values.shape_grad()" and "fe_values.shape_grad()"..., not using "FEValues::Extractors class". This part of the code is based on version 6.1.0 of DealII.
It is because of that I wanted to know if I had programmed the correct code (bellow this email) in order to get the "full" diffusion matrix. I would be very grateful if someone someone told me if my code is right or what I would have to change? Thanks again, best regards. Isa > Hello Isa, > > I would propose to use the FEValues::Extractors class to access the > single solution components U and V. You can find more about this > approach in the module on vector-valued problems. > > Best Regards, > Markus > > > > Am 18.02.11 18:52, schrieb [email protected]: >> Hello, >> >> I'd like to assemble the "complete" diffusion matrix in the 2D >> Navier-Stokes equation: >> >> for component U: d/dx(nu*d/dx(U))+d/dy(nu*d/dy(U)) + >> d/dx(nu*d/dx(U))+d/dy(nu*d/dx(V)) >> >> for component V: d/dx(nu*d/dx(V))+d/dy(nu*d/dy(V)) + >> d/dy(nu*d/dy(V))+d/dx(nu*d/dy(U)) >> >> nu is the viscosity and it isn't constant, so I do: >> >> So, I do: >> >> ///////////////////////////////////////////////////// >> for (; cell!=endc; ++cell) >> { >> >> cell_matrix=0.0; >> cell_rhs=0.0; >> fe_values.reinit (cell); >> . >> . >> . >> >> for (unsigned int i=0; i<dofs_per_cell; ++i) >> { >> //0=U,1=V,2=P >> const unsigned int comp_i=fe_total.system_to_component_index(i).first; >> for (unsigned int j=0; j<dofs_per_cell; ++j) >> { >> const unsigned int >> comp_j=fe_total.system_to_component_index(j).first; >> >> for (unsigned int q_point=0; q_point<n_q_points; ++q_point) >> { >> if (comp_i == comp_j) >> { >> >> //To assemble d/dx(nu*d/dx(U))+d/dy(nu*d/dy(U)) for component U >> //To assemble d/dx(nu*d/dx(V))+d/dy(nu*d/dy(V)) for component V: >> >> if (comp_i != DIMENSION) >> cell_matrix(i,j)+=(fe_values.JxW(q_point)* >> ( >> (nu*fe_values.shape_grad(i,q_point) >> * >> fe_values.shape_grad(j,q_point)) >> ) >> ); >> >> //To assemble d/dx(nu*d/dx(U)) for component U: >> >> else if (comp_i == 0) >> cell_matrix(i,j)+=(fe_values.JxW(q_point)* >> ( >> (nu*fe_values.shape_grad(i,q_point)[comp_i] >> * >> fe_values.shape_grad(j,q_point))[comp_i] >> ) >> ); >> //To assemble d/dy(nu*d/dy(V)) for component V: >> >> else if (comp_i == 1) >> cell_matrix(i,j)+=(fe_values.JxW(q_point)* >> ( >> (nu*fe_values.shape_grad(i,q_point[comp_i] >> * >> fe_values.shape_grad(j,q_point))[comp_i] >> ) >> ); >> >> >> } >> if ((comp_i != comp_j)&& (comp_i != >> DIMENSION)&&(comp_j!=DIMENSION)) >> { >> //To assemble d/dy(nu*d/dx(V)) for component U >> //To assemble d/dx(nu*d/dy(U)) for component V >> cell_matrix(i,j)+=(fe_values.JxW(q_point)* >> ( >> (nu*fe_values.shape_grad(i,q_point)[comp_j] >> * >> fe_values.shape_grad(j,q_point))[comp_i] >> ) >> ); >> >> } >> >> } >> } >> } >> . >> . >> . >> } >> ///////////////////////////////////////////////////// >> >> Is it the way to create this complete diffusion matrix? >> >> Thanks in advance. >> Best regards. >> Isa >> >> >> >> >> _______________________________________________ >> 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
