Hi Wolfgang, right so in the Stokes subsystem of my equations, in simplified form I am trying to solve:
div tau - grad p = rhs1 div v = rhs2 where my tau is 2 epsilon as in step-22 The boundary conditions I now want to implement to solve the 'real' problem is: on the sides (boundary 0): zero tangential stress and no normal flux, so we have v dot n = 0 and n.(-pI + tau)t = 0. at the top (boundary 1): zero tangential stress and prescribed normal component of normal stress, so that n.(-pI+tau)n= top stress value, and similarly, on the bottom (boundary 2): zero tangential stress and prescribed normal component of normal stress, so that n.(-pI+tau)n= bottom stress value. n and t here are the normal and tangential vectors relative to the surface at the respective boundaries, respectively. The test code works with Dirichlet conditions around the boundary it seems, but i need to implement the problem as above. Since the last message, I realised that I can just expand n.(-pI+tau) into normal and tangential components, which means that for the above real problem, I can just implement Neumann conditions with top stress value*normal vector to the surface, as the tangential stress are zero and therefore in the expansions, just simply equal zero. ie. I have n.(-pI+tau) = stress_value*normal vector. This is for the top and bottom boundaries. For the sides, I use no_normal_flux_constraints in vectortools. This means that I can indeed implement the boundary condition in the weak form of my equation. I am now having a small issue with this, and tried to look up more information about fe_face_values, but I am a little stuck. I am doing (eg. for boundary id 1, and i have the same for 2): for (unsigned int face_number=0; face_number<GeometryInfo<dim>::faces_per_cell; ++face_number) if (cell->face(face_number)->at_boundary() && (cell->face(face_number)->boundary_id() == 1)) { fe_face_values.reinit (cell, face_number); for (unsigned int q_point=0; q_point<n_face_q_points; ++q_point) { for (unsigned int i=0; i<dofs_per_cell; ++i) local_rhs(i) += data::topstress*(fe_face_values.normal_vector(q_point) * fe_face_values.shape_value(i,q_point) * fe_face_values.JxW(q_point)); } } I can't seem to find info about the fe_face_values.normal_vector, but it seems I am struggling with the types (tensor/double/vector) in the calculation - I am getting the error no match for óperator+= '(operand types are 'double'' and 'dealii::Tensor<1, 2, double>') I understand that the condition is on all velocities and pressure, not just the velocities, so I didn't use a component mask or similar. In the case of the mixed formulation like this one, how am I supposed to be extracting the normal vector correctly in a way I am applying it to the whole system? I hope that makes it clearer! Thanks for your help again. On Monday, December 4, 2017 at 9:17:20 PM UTC+1, Wolfgang Bangerth wrote: > > > Jane, > > > I get that the stress boundary can be put into the weak form, but only > > really when you have Neumann conditions. > > Correct, but prescribing the stress in an elasticity problem is exactly > a Neumann boundary condition. > > > > And why would I have to use a > > compute_nonzero_normal_flux for a zero one? > > I thought you had a nonzero displacement/velocity on parts of the > boundary? If it is zero, then of course you would want to use a > different function. > > > > I believe it my have been my fault not being clearer. > > > > On the 'side' boundaries, I have no normal flux. > > On the other boundaries, I have a prescribed nonzero inhomogeneous > > normal component of the normal stress (so this is only one of the > > components). > > On all boundaries, I have zero tangential stresses (which sorts out the > > other 2 components in 2D in additions to the normal component conditions > > from the above). > > Can you just put that into formulas, to make things really clear? State > the equations you want to solve and for each type of boundary what > conditions hold there. > > Maybe that would make it easier for everyone to understand what the > other side is saying :-) > > > > Pardon my potential ignorance, but I'm unsure how I can put these in the > > weak form when i have component-wise conditions (partial conditions). > > I had thought to use compute_no_normal_flux_constraints for the side > > boundaries and compute_nonzero_normal_flux_constraints for the other two > > as constraints, but this doesn't seem to work. > > Be specific here as well: what doesn't work? Do you get an assertion, or > a wrong solution, or a compiler error, ...? > > Best > W. > > -- > ------------------------------------------------------------------------ > Wolfgang Bangerth email: bang...@colostate.edu > <javascript:> > www: http://www.math.colostate.edu/~bangerth/ > -- 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 dealii+unsubscr...@googlegroups.com. For more options, visit https://groups.google.com/d/optout.