Hi, !697 triggered a discussion of how to handle boundary conditions for vector-valued functions, in this case for the `CSFVectorFunction` family. The primary questions are: 1. Which of the following two options do we choose: a. The boundary condition for a `CSFVectorFunction` is a property of that function alone, i.e. we enforce that the individual component functions always have the same attribute. b. We allow that the individual component functions carry different attributes. 2. If we go for the second option (b), this induces the question what the method `CSFVectorFunction::getBoundaryCondition()` should return. The suggestion was that in this case we return an array of size dimension with the boundary condition attributes of the scalar sub-functions. **Note:** The current implementation is (blame on me) inconsistent in this respect. It started out assuming option (a), but does not enforce this and even offers a method to set the attribute on individual components. In order to better understand the pros and cons of the approaches let us look at the *free-slip* boundary condition as an example: ```math \mathbf{u} \cdot \mathbf{n} = 0 \enspace,\quad \frac{\partial (\mathbf{u} \cdot \mathbf{t})}{\partial \mathbf{n} } = 0 ``` with normal and tangential vectors $\mathbf{n}$ and $\mathbf{t}$. **PRO for (b):** - In the case of a cartesian box, the free-slip condition on a face simply translates into having a homogeneous Dirichlet boundary condition for the scalar sub-function of the normal component and homogeneous Neumann boundary conditions for the two other sub-functions. If we go for (b), this can nicely be used to handle the case. - Always setting a well-defined boundary condition on each sub-functions is IMHO good. And, without the need to enforce uniqueness over the sub-functions, is simple in the current implementation. **CONS for (b):** - In any case where the domain boundaries are not axis-aligned, the splitting cannot be exploited. Especially in the case of an annulus or a thick spherical shell, we then need to set the boundary condition of the scalar sub-functions to `FreeslipBoundary`. Mathematically this feels a little unnatural for a scalar function. - We would handle the `CSFVectorFunction` case differently from any other vector-valued FEM function space. An *outflow* conditions of the form ```math \frac{\partial (\mathbf{u} \cdot \mathbf{n})}{\partial \mathbf{n} } = 0 \enspace,\quad \frac{\partial (\mathbf{u} \cdot \mathbf{t})}{\partial \mathbf{n} } = 0 ``` could be treated in either scenario: - in (a) we could set the attribute of the `CSFVectorFunction` to `NeumannBoundary` - in (b) we could do the same for the scalar sub-functions Any input on the question is welcome. Cheers Marcus