1 /*============================================================================
2 * General functions or variables for the CDO module
3 *============================================================================*/
4
5 /*
6 This file is part of Code_Saturne, a general-purpose CFD tool.
7
8 Copyright (C) 1998-2021 EDF S.A.
9
10 This program is free software; you can redistribute it and/or modify it under
11 the terms of the GNU General Public License as published by the Free Software
12 Foundation; either version 2 of the License, or (at your option) any later
13 version.
14
15 This program is distributed in the hope that it will be useful, but WITHOUT
16 ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
17 FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
18 details.
19
20 You should have received a copy of the GNU General Public License along with
21 this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
22 Street, Fifth Floor, Boston, MA 02110-1301, USA.
23 */
24
25 /*----------------------------------------------------------------------------*/
26
27 /*----------------------------------------------------------------------------
28 * Standard C library headers
29 *----------------------------------------------------------------------------*/
30
31 #include <stdio.h>
32 #include <stdlib.h>
33 #include <string.h>
34 #include <assert.h>
35
36 /*----------------------------------------------------------------------------
37 * Local headers
38 *----------------------------------------------------------------------------*/
39
40 #include "bft_error.h"
41 #include "bft_mem.h"
42
43 #include "cs_log.h"
44
45 /*----------------------------------------------------------------------------
46 * Header for the current file
47 *----------------------------------------------------------------------------*/
48
49 #include "cs_dbg.h"
50
51 /*----------------------------------------------------------------------------*/
52
53 BEGIN_C_DECLS
54
55 /*=============================================================================
56 * Global variables
57 *============================================================================*/
58
59 /*=============================================================================
60 * Local static variables
61 *============================================================================*/
62
63 /*============================================================================
64 * Private function prototypes
65 *============================================================================*/
66
67 /*============================================================================
68 * Public function prototypes
69 *============================================================================*/
70
71 #if defined(DEBUG) && !defined(NDEBUG)
72 /*----------------------------------------------------------------------------*/
73 /*!
74 * \brief Function used to select which element deserves a dump or specific
75 * treatment during a debugging stage
76 *
77 * \param[in] eqp pointer to a cs_equation_param_t structure
78 * \param[in] cm pointer to a cs_cell_mesh_t structure
79 * \param[in] csys pointer to a cs_cell_sys_t structure
80 */
81 /*----------------------------------------------------------------------------*/
82
83 bool
cs_dbg_cw_test(const cs_equation_param_t * eqp,const cs_cell_mesh_t * cm,const cs_cell_sys_t * csys)84 cs_dbg_cw_test(const cs_equation_param_t *eqp,
85 const cs_cell_mesh_t *cm,
86 const cs_cell_sys_t *csys)
87 {
88 CS_UNUSED(eqp);
89
90 #if 1 /* First example: Look for debug information related the cell number 0 */
91 if (cm != NULL)
92 if (cm->c_id == 0)
93 return true;
94 if (csys != NULL)
95 if (csys->c_id == 0)
96 return true;
97 #endif
98
99 #if 0 /* Second example: Only search debug information for a given equation and
100 a requested vertex */
101
102 /* Vertex number 441 belongs to the current cell ? */
103 const cs_lnum_t target_vertex_id = 441;
104 const short int _v = cs_cell_mesh_get_v(target_vertex_id, cm);
105 if (_v < 0)
106 return false;
107
108 if (eqp != NULL)
109 if (strcmp(eqp->name, "Tracer1") == 0)
110 return true;
111 #endif
112
113 #if 0 /* Third example: Look for the cells which have a previous DoF value
114 greater than 0.06 for the face 1 */
115 const cs_lnum_t target_face_id = 1;
116 const short int _f = cs_cell_mesh_get_f(target_face_id, cm);
117
118 if (csys != NULL && _f > -1) {
119 if (csys->val_n[_f] > 0.06)
120 return true;
121 }
122 #endif
123
124 return false; /* Default behavior */
125 }
126
127 /*----------------------------------------------------------------------------*/
128 /*!
129 * \brief Print an array.
130 * Print into the file f if given otherwise open a new file named
131 * fname if given otherwise print into the standard output
132 * The usage of threshold allows one to compare more easier arrays
133 * without taking into account numerical roundoff.
134 *
135 * \param[in] fp pointer to a file structure or NULL
136 * \param[in] fname filename or NULL
137 * \param[in] thd threshold (below this value --> set 0)
138 * \param[in] n_elts size of the array
139 * \param[in] array list of values to dump
140 * \param[in] n_cols print array with n_cols columns
141 */
142 /*----------------------------------------------------------------------------*/
143
144 void
cs_dbg_array_fprintf(FILE * fp,const char * fname,cs_real_t thd,cs_lnum_t n_elts,const cs_real_t array[],int n_cols)145 cs_dbg_array_fprintf(FILE *fp,
146 const char *fname,
147 cs_real_t thd,
148 cs_lnum_t n_elts,
149 const cs_real_t array[],
150 int n_cols)
151 {
152 FILE *fout = stdout;
153 if (fp != NULL)
154 fout = fp;
155 else if (fname != NULL) {
156 fout = fopen(fname, "w");
157 }
158
159 fprintf(fout, "array %p\n", (const void *)array);
160
161 if (array == NULL)
162 return;
163
164 if (n_cols < 1) n_cols = 1;
165 int n_rows = n_elts/n_cols;
166
167 for (cs_lnum_t i = 0; i < n_rows; i++) {
168 for (cs_lnum_t j = i*n_cols; j < (i+1)*n_cols; j++) {
169 if (fabs(array[j]) < thd)
170 fprintf(fout, " % -8.5e", 0.);
171 else
172 fprintf(fout, " % -8.5e", array[j]);
173 }
174 fprintf(fout, "\n");
175 }
176
177 if (n_rows*n_cols < n_elts) {
178 for (cs_lnum_t j = n_rows*n_cols; j < n_elts; j++) {
179 if (fabs(array[j]) < thd)
180 fprintf(fout, " % -8.5e", 0.);
181 else
182 fprintf(fout, " % -8.5e", array[j]);
183 }
184 fprintf(fout, "\n");
185 }
186
187 if (fout != stdout && fout != fp)
188 fclose(fout);
189 }
190
191 /*----------------------------------------------------------------------------*/
192 /*!
193 * \brief In debug mode, print into a file the solution and its right-hand
194 * side
195 *
196 * \param[in] eqname name of the related equation
197 * \param[in] nt number of time step
198 * \param[in] level level of debug
199 * \param[in] sol solution array
200 * \param[in] rhs rhs array
201 * \param[in] size size of the array to print
202 */
203 /*----------------------------------------------------------------------------*/
204
205 void
cs_dbg_fprintf_system(const char * eqname,int nt,int level,const cs_real_t * sol,const cs_real_t * rhs,cs_lnum_t size)206 cs_dbg_fprintf_system(const char *eqname,
207 int nt,
208 int level,
209 const cs_real_t *sol,
210 const cs_real_t *rhs,
211 cs_lnum_t size)
212 {
213 int len = strlen(eqname) + strlen("-sol-.log") + 4 + 1;
214 char *filename = NULL;
215
216 BFT_MALLOC(filename, len, char);
217
218 sprintf(filename, "%s-sol-%04d.log", eqname, nt);
219 if (sol != NULL && level > 4)
220 cs_dbg_array_fprintf(NULL, filename, 1e-16, size, sol, 6);
221
222 sprintf(filename, "%s-rhs-%04d.log", eqname, nt);
223 if (rhs != NULL && level > 5)
224 cs_dbg_array_fprintf(NULL, filename, 1e-16, size, rhs, 6);
225
226 BFT_FREE(filename);
227 }
228
229 /*----------------------------------------------------------------------------*/
230 /*!
231 * \brief In debug mode, dump an array of double into the log
232 *
233 * \param[in] header header message to write
234 * \param[in] size number of elements in array
235 * \param[in] array pointer to the array of values
236 * \param[in] n_cols print array with n_cols columns
237 */
238 /*----------------------------------------------------------------------------*/
239
240 void
cs_dbg_darray_to_listing(const char * header,const cs_lnum_t size,const cs_real_t array[],int n_cols)241 cs_dbg_darray_to_listing(const char *header,
242 const cs_lnum_t size,
243 const cs_real_t array[],
244 int n_cols)
245 {
246 cs_log_printf(CS_LOG_DEFAULT, "\nDUMP>> %s\n", header);
247
248 if (n_cols < 1) n_cols = 1;
249 int n_rows = size/n_cols;
250
251 for (cs_lnum_t i = 0; i < n_rows; i++) {
252 for (cs_lnum_t j = i*n_cols; j < (i+1)*n_cols; j++)
253 cs_log_printf(CS_LOG_DEFAULT, " (%04ld) % 6.4e", (long)j, array[j]);
254 cs_log_printf(CS_LOG_DEFAULT, "\n");
255 }
256
257 if (n_rows*n_cols < size) {
258 for (cs_lnum_t j = n_rows*n_cols; j < size; j++)
259 cs_log_printf(CS_LOG_DEFAULT, " (%04ld) % 6.4e", (long)j, array[j]);
260 cs_log_printf(CS_LOG_DEFAULT, "\n");
261 }
262
263 }
264
265 /*----------------------------------------------------------------------------*/
266 /*!
267 * \brief In debug mode, dump an array of integer into the log
268 *
269 * \param[in] header header message to write
270 * \param[in] size number of elements in array
271 * \param[in] array pointer to the array of values
272 * \param[in] n_cols print array with n_cols columns
273 */
274 /*----------------------------------------------------------------------------*/
275
276 void
cs_dbg_iarray_to_listing(const char * header,const cs_lnum_t size,const cs_lnum_t array[],int n_cols)277 cs_dbg_iarray_to_listing(const char *header,
278 const cs_lnum_t size,
279 const cs_lnum_t array[],
280 int n_cols)
281 {
282 cs_log_printf(CS_LOG_DEFAULT, "\nDUMP>> %s\n", header);
283
284 if (n_cols < 1) n_cols = 1;
285 int n_rows = size/n_cols;
286
287 for (cs_lnum_t i = 0; i < n_rows; i++) {
288 for (cs_lnum_t j = i*n_cols; j < (i+1)*n_cols; j++)
289 cs_log_printf(CS_LOG_DEFAULT, " (%04ld) % 6ld", (long)j, (long)array[j]);
290 cs_log_printf(CS_LOG_DEFAULT, "\n");
291 }
292
293 if (n_rows*n_cols < size) {
294 for (cs_lnum_t j = n_rows*n_cols; j < size; j++)
295 cs_log_printf(CS_LOG_DEFAULT, " (%04ld) % 6ld", (long)j, (long)array[j]);
296 cs_log_printf(CS_LOG_DEFAULT, "\n");
297 }
298
299 }
300
301 /*----------------------------------------------------------------------------*/
302 /*!
303 * \brief In debug mode, dump a linear system. Case of scalar-valued entries.
304 *
305 * \param[in] eqname name of the equation related to the current system
306 * \param[in] matrix pointer to the matrix to dump
307 */
308 /*----------------------------------------------------------------------------*/
309
310 void
cs_dbg_dump_local_scalar_msr_matrix(const char * name,const cs_matrix_t * matrix)311 cs_dbg_dump_local_scalar_msr_matrix(const char *name,
312 const cs_matrix_t *matrix)
313 {
314 if (cs_matrix_get_type(matrix) != CS_MATRIX_MSR)
315 return;
316
317 cs_log_printf(CS_LOG_DEFAULT, "\nDUMP MSR MATRIX FOR THE EQUATION >> %s\n",
318 name);
319
320 const cs_lnum_t size = cs_matrix_get_n_rows(matrix);
321 const cs_lnum_t *row_index, *col_id;
322 const cs_real_t *d_val, *x_val;
323
324 cs_matrix_get_msr_arrays(matrix, &row_index, &col_id, &d_val, &x_val);
325
326 for (cs_lnum_t i = 0; i < size; i++) {
327
328 const cs_lnum_t *idx = row_index + i;
329
330 cs_log_printf(CS_LOG_DEFAULT, "%4ld |D|% -6.4e |E", (long)i, d_val[i]);
331 for (cs_lnum_t j = idx[0]; j < idx[1]; j++)
332 cs_log_printf(CS_LOG_DEFAULT, "|% -6.4e c%4ld", x_val[j], (long)col_id[j]);
333 cs_log_printf(CS_LOG_DEFAULT, "\n");
334
335 } /* Loop on rows */
336 }
337
338 /*----------------------------------------------------------------------------*/
339 /*!
340 * \brief In debug mode, dump a linear system
341 *
342 * \param[in] eqname name of the equation related to the current system
343 * \param[in] size number of elements in array
344 * \param[in] verbosity level of display
345 * \param[in] x solution array
346 * \param[in] b right-hand side
347 * \param[in] row_index index on row entries (column id and extra-diag values)
348 * \param[in] col_id list of column id
349 * \param[in] xval array of extra-diagonal values
350 * \param[in] dval array of diagonal values
351 */
352 /*----------------------------------------------------------------------------*/
353
354 void
cs_dbg_dump_linear_system(const char * eqname,cs_lnum_t size,int verbosity,const cs_real_t x[],const cs_real_t b[],const cs_lnum_t row_index[],const cs_lnum_t col_id[],const cs_real_t xval[],const cs_real_t dval[])355 cs_dbg_dump_linear_system(const char *eqname,
356 cs_lnum_t size,
357 int verbosity,
358 const cs_real_t x[],
359 const cs_real_t b[],
360 const cs_lnum_t row_index[],
361 const cs_lnum_t col_id[],
362 const cs_real_t xval[],
363 const cs_real_t dval[])
364 {
365 cs_log_printf(CS_LOG_DEFAULT, "\nDUMP LINEAR SYSTEM FOR THE EQUATION >> %s\n",
366 eqname);
367
368 int n_dump_cols = 8;
369 int n_dump_rows = size/n_dump_cols;
370
371 cs_log_printf(CS_LOG_DEFAULT, " >> SOLUTION\n");
372 for (cs_lnum_t i = 0; i < n_dump_rows; i++) {
373 for (cs_lnum_t j = i*n_dump_cols; j < (i+1)*n_dump_cols; j++)
374 cs_log_printf(CS_LOG_DEFAULT, "%4ld % -6.4e |", (long)j, x[j]);
375 cs_log_printf(CS_LOG_DEFAULT, "\n");
376 }
377 if (n_dump_rows*n_dump_cols < size) {
378 for (cs_lnum_t j = n_dump_rows*n_dump_cols; j < size; j++)
379 cs_log_printf(CS_LOG_DEFAULT, "%4ld % -6.4e |", (long)j, x[j]);
380 cs_log_printf(CS_LOG_DEFAULT, "\n");
381 }
382
383 cs_log_printf(CS_LOG_DEFAULT, " >> RIGHT-HAND SIDE\n");
384 for (cs_lnum_t i = 0; i < n_dump_rows; i++) {
385 for (cs_lnum_t j = i*n_dump_cols; j < (i+1)*n_dump_cols; j++)
386 cs_log_printf(CS_LOG_DEFAULT, "%4ld % -6.4e |", (long)j, b[j]);
387 cs_log_printf(CS_LOG_DEFAULT, "\n");
388 }
389 if (n_dump_rows*n_dump_cols < size) {
390 for (cs_lnum_t j = n_dump_rows*n_dump_cols; j < size; j++)
391 cs_log_printf(CS_LOG_DEFAULT, "%4ld % -6.4e |", (long)j, b[j]);
392 cs_log_printf(CS_LOG_DEFAULT, "\n");
393 }
394
395 if (verbosity > 2) {
396
397 cs_log_printf(CS_LOG_DEFAULT, " >> DIAGONAL ENTRIES\n");
398 for (cs_lnum_t i = 0; i < n_dump_rows; i++) {
399 for (cs_lnum_t j = i*n_dump_cols; j < (i+1)*n_dump_cols; j++)
400 cs_log_printf(CS_LOG_DEFAULT, "%4ld % -6.4e |", (long)j, dval[j]);
401 cs_log_printf(CS_LOG_DEFAULT, "\n");
402 }
403 if (n_dump_rows*n_dump_cols < size) {
404 for (cs_lnum_t j = n_dump_rows*n_dump_cols; j < size; j++)
405 cs_log_printf(CS_LOG_DEFAULT, "%4ld % -6.4e |", (long)j, dval[j]);
406 cs_log_printf(CS_LOG_DEFAULT, "\n");
407 }
408
409 cs_log_printf(CS_LOG_DEFAULT, " >> EXTRA-DIAGONAL ENTRIES\n");
410 for (cs_lnum_t i = 0; i < size; i++) {
411
412 const cs_lnum_t *idx = row_index + i;
413 const cs_lnum_t *_col = col_id + idx[0];
414 const cs_real_t *_val = xval + idx[0];
415 const int n_entries = idx[1] - idx[0];
416
417 int _n_cols = 6;
418 int _n_rows = n_entries/_n_cols;
419
420 for (cs_lnum_t ii = 0; ii < _n_rows; ii++) {
421 cs_log_printf(CS_LOG_DEFAULT, "ROW%4ld >> ", (long)i);
422 for (cs_lnum_t jj = ii*_n_cols; jj < (ii+1)*_n_cols; jj++)
423 cs_log_printf(CS_LOG_DEFAULT, "%4ld: % -6.4e |",
424 (long)_col[jj], _val[jj]);
425 cs_log_printf(CS_LOG_DEFAULT, "\n");
426 }
427 if (_n_rows*_n_cols < n_entries) {
428 cs_log_printf(CS_LOG_DEFAULT, "ROW%4ld >> ", (long)i);
429 for (cs_lnum_t jj = _n_rows*_n_cols; jj < n_entries; jj++)
430 cs_log_printf(CS_LOG_DEFAULT, "%4ld: % -6.4e |",
431 (long)_col[jj], _val[jj]);
432 cs_log_printf(CS_LOG_DEFAULT, "\n");
433 }
434
435 }
436
437 } /* verbosity */
438 }
439 #endif /* Only in debug mode */
440
441 /*----------------------------------------------------------------------------*/
442
443 END_C_DECLS
444