How to calculate second derivatives

[yang]

Hi team

I am using Code_Saturne 7.0.2. Now, I want to know how to use the Code_Saturne source code (c language) to calculate the second derivatives, such as, [Equations will insert Code_Saturne][second-derivates].

I know that we can directly use function 'cs_field_gradient_vector' to calculate first-order gradient, such as (du/dx). do we have similarly method to calculate the second derivatives?

Looking forward to your reply ASAP! Thanks.

Regards,
Yang

[Yvan Fournier]

Hello,

I do not believe we have a simply-usable operator function for this, but I'll check with colleagues in case I missed something.

Best regards,

Yvan

[yang]

Hi Yvan,

Thanks for you reply.

If there is not a simply-usable operator function for second derivatives, I think I can calculate it by finite difference based on the first order derivatives and neighbour cells. As shown below figure, the left and right boundary is symmetry BC, so the second derivatives at cell 'a' can be calculated through 'b' and 'c' gradient based on the difference scheme.

Now, using function 'cs_glob_mesh_adjacencies', I can obtain cell 'b' next to cell a in x direction, but cannot obtain adjacent cell 'c'. How to get adjacent cells following the boundary condition?(i.g. obtain adjacent cell 'c' from cell 'a').

Kind regards,
Yang

[Yvan Fournier]

Hello,

In the global mesh_adjacencies, the cell_cells_idx and cell_cells connectivity should allow you to obtain all cells adjacent through a face. For a given cell (cell_id):

const cs_mesh_adjacencies_t *ma = cs_glob_mesh_adjacencies; // to use "ma" for short notation

Code: Select all

for (cs_lnum_t i = ma->cell_cells_idx[cell_id];  i < ma->cell_cells_idx[cell_id+1]; i++) {
  neighbor_cell_id = ma->cell_cells[i]
  ...
  /do work */
  ...
}

Using cell_cells_e_idx and cell_cells_e, you can also access cells not sharing a face but sharing at least one vertex, of the extended neighborhood is activated (based on gradient computation options). But if you mesh is regular enough you only need the main adjacency.

The adjacency information is directionless, so you need to filter directions on your own.

A simpler solution might be to use the existing operators to compute the gradient of a gradient (as the gradient of a scalar is a vector, you need the vector gradient for a scalar, or the tensor gradient for the gradient of a vector; you would be stuck for the second-order derivatives of a tensor). Also, boundary conditions could be an issue, unless you assume homogeneous Neumann everywhere at least for the second step, but you would have to handle that in a complex manner to do better rewriting your own operator.

Best regards,

Yvan