Isaac:

On 10/4/23 16:05, Isaac Bannerman wrote:
I am using dealii for a CFD solver and the solution vector which is passed to the linear algebra solver is a BlockVector where block zero is the velocity and block 1 is the pressure. I would like to obtain the global DOF numbers of all the entries in the block vector; I would therefore like to know if the  the indices of the solution vector (especially for block zero) corresponds to the global DOF  numbers.   The following lines (ignoring other details)  illustrates the question I am asking:

//solution vector
BlockVector<double> newton_update;

//the solution steps
SolverFGMRES<BlockVector<double>> gmres(solver_control, vector_memory);
gmres.solve(system_matrix, newton_update, system_rhs, *preconditioner);

I want to know for instance if index "44" used in the snippet below will correspond to global dof number 44.

  newton_update.block(0)[44]

In the case where the answer to the above question is no, could someone please suggest a way to obtain the global DOF number of each entry in the BlockVector?

At the end of the day, the connection between DoF numbers and vector elements is what you make it. If you take a look at what step-3 does, for example, you see in the assembly this:

      // get DoF indices
      cell->get_dof_indices(local_dof_indices);

      for (const unsigned int i : ...)
        // Use DoF indices as index into a vector
        system_rhs(local_dof_indices[i]) += cell_rhs(i);

So yes, generally the DoF indices correspond 1:1 to vector entries.

Of course, when you split a vector into blocks, the indices *within these blocks* do not necessarily correspond to DoF indices. So in your example above,
  newton_update[44]
corresponds to DoF index 44. The use of
  newton_update.block(0)[44]
also corresponds to DoF index 44, but only because block 0 starts at index 0. On the other hand,
  newton_update.block(1)[44]
does *not* correspond to the DoF with index 44, but to the DoF with index
  newton_update.block(0).size() + 44


My ultimate goal is the obtain the DOF and real coordinates (x,y,z) of each entry in the block vector. Fortunately, I already obtained the DOFs and their corresponding coordinates in the real space using DoFTools::map_dofs_to_support_points(...) ;  I want to get the DOFs for the BlockVector entries so I can match the the coordinates to the entries in the BlockVector.  Any suggestion on a better way to do it will also be much appreciated.

As shown above, if you index into a block, you need to take into account the offset of that block within the vector as a whole.

Best
 W.


--
------------------------------------------------------------------------
Wolfgang Bangerth          email:                 bange...@colostate.edu
                           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.
To view this discussion on the web visit 
https://groups.google.com/d/msgid/dealii/0c8cdbc3-b0b7-b918-0744-c72b7e6787b7%40colostate.edu.

Reply via email to