programmer's documentation
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Groups Pages
cs_time_plot.h
Go to the documentation of this file.
1 #ifndef __CS_TIME_PLOT_H__
2 #define __CS_TIME_PLOT_H__
3 
4 /*============================================================================
5  * Time_Plot helper structures
6  *============================================================================*/
7 
8 /*
9  This file is part of Code_Saturne, a general-purpose CFD tool.
10 
11  Copyright (C) 1998-2013 EDF S.A.
12 
13  This program is free software; you can redistribute it and/or modify it under
14  the terms of the GNU General Public License as published by the Free Software
15  Foundation; either version 2 of the License, or (at your option) any later
16  version.
17 
18  This program is distributed in the hope that it will be useful, but WITHOUT
19  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20  FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
21  details.
22 
23  You should have received a copy of the GNU General Public License along with
24  this program; if not, write to the Free Software Foundation, Inc., 51 Franklin
25  Street, Fifth Floor, Boston, MA 02110-1301, USA.
26 */
27 
28 /*----------------------------------------------------------------------------*/
29 
30 /*----------------------------------------------------------------------------
31  * Local headers
32  *----------------------------------------------------------------------------*/
33 
34 #include "cs_base.h"
35 
36 /*----------------------------------------------------------------------------*/
37 
39 
40 /*============================================================================
41  * Macro definitions
42  *============================================================================*/
43 
44 /*============================================================================
45  * Type definitions
46  *============================================================================*/
47 
48 typedef struct _cs_time_plot_t cs_time_plot_t;
49 
50 /*============================================================================
51  * Local type definitions
52  *============================================================================*/
53 
54 /* Type of 1D plot file format */
55 
56 typedef enum {
57  CS_TIME_PLOT_DAT, /* .dat file (usable by Qtplot or Grace) */
58  CS_TIME_PLOT_CSV /* .csv file (readable by ParaView or spreadsheat) */
60 
61 /*============================================================================
62  * Global variables
63  *============================================================================*/
64 
65 /*=============================================================================
66  * Public function prototypes for Fortran API
67  *============================================================================*/
68 
69 /*----------------------------------------------------------------------------
70  * Create a writer for time plot probe-type data.
71  *
72  * This subroutine should only be called by one rank for a given data series.
73  *
74  * subroutine tppini (tplnum, tplnam, tplpre, tplfmt, idtvar,
75  * *****************
76  * ntflsh, wtflsh, nprb, lstprb, xyzprb,
77  * lnam, lpre)
78  *
79  * integer tplnum : <-- : number of plot to create (> 0)
80  * character tplnam : <-- : name of associated plot
81  * character tplpre : <-- : prefix for associated file
82  * integer tplfmt : <-- : associated format
83  * (1: dat, 2: csv, 3: both)
84  * integer idtvar : <-- : calculation time dependency
85  * integer ntflsh : <-- : file write every ntflsh output
86  * time steps if > 0 (file kept
87  * open otherwise)
88  * integer wtflsh : <-- : file flush forced every wtflsh
89  * elapsed seconds if > 0
90  * integer nprb : <-- : number of probes
91  * integer lstprb : <-- : list of probes (1 to n)
92  * double precision xyzprb : <-- : probe coordinates
93  * integer lnam : <-- : name length
94  * integer lpre : <-- : prefix length
95  *----------------------------------------------------------------------------*/
96 
97 void CS_PROCF (tppini, TPPINI)
98 (
99  const cs_int_t *tplnum,
100  const char *tplnam,
101  const char *tplpre,
102  const cs_int_t *tplfmt,
103  const cs_int_t *idtvar,
104  const cs_int_t *ntflsh,
105  const cs_real_t *wtflsh,
106  const cs_int_t *nprb,
107  const cs_int_t *lstprb,
108  const cs_real_t *xyzprb,
109  const cs_int_t *lnam,
110  const cs_int_t *lpre
111  CS_ARGF_SUPP_CHAINE /* (possible 'length' arguments added
112  by many Fortran compilers) */
113 );
114 
115 /*----------------------------------------------------------------------------
116  * Create a writer for time plot structure-type data.
117  *
118  * This subroutine should only be called by one rank for a given data series.
119  *
120  * subroutine tpsini (tplnum, tplnam, tplpre, tplfmt, idtvar,
121  * *****************
122  * ntflsh, wtflsh, nprb, lstprb, xyzprb,
123  * lnam, lpre)
124  *
125  * integer tplnum : <-- : number of plot to create (> 0)
126  * character tplnam : <-- : name of associated plot
127  * character tplpre : <-- : prefix for associated file
128  * integer tplfmt : <-- : associated format
129  * (1: dat, 2: csv, 3: both)
130  * integer idtvar : <-- : calculation time dependency
131  * integer ntflsh : <-- : file write every ntflsh output
132  * time steps if > 0 (file kept
133  * open otherwise)
134  * integer wtflsh : <-- : file flush forced every wtflsh
135  * elapsed seconds if > 0
136  * integer nstru : <-- : number of structures
137  * double precision xmstru : <-- : mass matrixes
138  * double precision xcstru : <-- : damping matrixes
139  * double precision xkstru : <-- : stiffness matrixes
140  * integer lnam : <-- : name length
141  * integer lpre : <-- : prefix length
142  *----------------------------------------------------------------------------*/
143 
144 void CS_PROCF (tpsini, TPPINI)
145 (
146  const cs_int_t *tplnum,
147  const char *tplnam,
148  const char *tplpre,
149  const cs_int_t *tplfmt,
150  const cs_int_t *idtvar,
151  const cs_int_t *ntflsh,
152  const cs_real_t *wtflsh,
153  const cs_int_t *nstru,
154  const cs_real_t *xmstru,
155  const cs_real_t *xcstru,
156  const cs_real_t *xkstru,
157  const cs_int_t *lnam,
158  const cs_int_t *lpre
159  CS_ARGF_SUPP_CHAINE /* (possible 'length' arguments added
160  by many Fortran compilers) */
161 );
162 
163 /*----------------------------------------------------------------------------
164  * Finalize a writer for time plot data.
165  *
166  * This subroutine should only be called by one rank for a given data series.
167  *
168  * subroutine tplend (tplnum)
169  * *****************
170  *
171  * integer tplnum : <-- : number of plot to create (> 0)
172  * integer tplfmt : <-- : associated format
173  * (1: dat, 2: csv, 3: both)
174  *----------------------------------------------------------------------------*/
175 
176 void CS_PROCF (tplend, TPLEND)
177 (
178  const cs_int_t *tplnum,
179  const cs_int_t *tplfmt
180 );
181 
182 /*----------------------------------------------------------------------------
183  * Write time plot values.
184  *
185  * subroutine tplwri (tplnum, tplfmt, nprb, ntcabs, ttcabs, valprb)
186  * *****************
187  *
188  * integer tplnum : <-- : number of associated plot (> 0)
189  * integer tplfmt : <-- : associated format
190  * (1: dat, 2: csv, 3: both)
191  * integer nprb : <-- : number of probes
192  * integer ntcabs : <-- : current time step number
193  * double precision ttcabs : <-- : current time value
194  * double precision valprb : <-- : probe values
195  *----------------------------------------------------------------------------*/
196 
197 void CS_PROCF (tplwri, TPLWRI)
198 (
199  const cs_int_t *tplnum,
200  const cs_int_t *tplfmt,
201  const cs_int_t *nprb,
202  const cs_int_t *ntcabs,
203  const cs_real_t *ttcabs,
204  const cs_real_t *valprb
205 );
206 
207 /*----------------------------------------------------------------------------
208  * Return the number of time plots accessible through the Fortran API
209  *
210  * This subroutine will only return the number of time plots defined by the
211  * local rank
212  *
213  * subroutine tplnbr (ntpl)
214  * *****************
215  *
216  * integer ntpl : --> : number of time plots defined
217  *----------------------------------------------------------------------------*/
218 
219 void CS_PROCF (tplnbr, TPLNBR)
220 (
221  cs_int_t *ntpl
222 );
223 
224 /*=============================================================================
225  * Public function prototypes
226  *============================================================================*/
227 
228 /*----------------------------------------------------------------------------
229  * Initialize a plot file writer for probe-type plots
230  *
231  * This function should only be called by one rank for a given data series.
232  *
233  * parameters:
234  * plot_name <-- plot (variable) name
235  * file_prefix <-- file name prefix
236  * format <-- associated file format
237  * use_iteration <-- should we use the iteration number instead of the
238  * physical time ?
239  * flush_wtime <-- elapsed time interval between file flushes
240  * (if < 0, no forced flush)
241  * n_buffer_steps <-- number of time steps in output buffer if
242  * file is not to be kept open
243  * n_probes <-- number of probes associated with this plot
244  * probe_list <-- numbers (1 to n) of probes if filtered, or NULL
245  * probe_coords <-- probe coordinates, or NULL
246  *
247  * returns:
248  * pointer to new time plot writer
249  *----------------------------------------------------------------------------*/
250 
251 cs_time_plot_t *
252 cs_time_plot_init_probe(const char *plot_name,
253  const char *file_prefix,
254  cs_time_plot_format_t format,
255  bool use_iteration,
256  double flush_wtime,
257  int n_buffer_steps,
258  int n_probes,
259  const int *probe_list,
260  const cs_real_t probe_coords[]);
261 
262 /*----------------------------------------------------------------------------
263  * Initialize a plot file writer for structure-type plots
264  *
265  * This function should only be called by one rank for a given data series.
266  *
267  * parameters:
268  * plot_name <-- plot (variable) name
269  * file_prefix <-- file name prefix
270  * format <-- associated file format
271  * use_iteration <-- should we use the iteration number instead of the
272  * physical time ?
273  * flush_wtime <-- elapsed time interval between file flushes
274  * (if < 0, no forced flush)
275  * n_buffer_steps <-- number of time steps in output buffer if
276  * file is not to be kept open
277  * n_structures <-- number of structures associated with this plot
278  * mass_matrixes <-- mass matrix coefficients (3x3 blocks)
279  * damping_matrixes <-- damping matrix coefficients (3x3 blocks)
280  * stiffness_matrixes <-- stiffness matrix coefficients (3x3 blocks)
281  *
282  * returns:
283  * pointer to new time plot writer
284  *----------------------------------------------------------------------------*/
285 
286 cs_time_plot_t *
287 cs_time_plot_init_struct(const char *plot_name,
288  const char *file_prefix,
289  cs_time_plot_format_t format,
290  bool use_iteration,
291  double flush_wtime,
292  int n_buffer_steps,
293  int n_structures,
294  const cs_real_t mass_matrixes[],
295  const cs_real_t damping_matrixes[],
296  const cs_real_t stiffness_matrixes[]);
297 
298 /*----------------------------------------------------------------------------
299  * Finalize time plot writer for a given variable
300  *
301  * This function should only be called by one rank for a given data series.
302  *
303  * parameters:
304  * p <-> time plot values file handler
305  *----------------------------------------------------------------------------*/
306 
307 void
308 cs_time_plot_finalize(cs_time_plot_t **p);
309 
310 /*----------------------------------------------------------------------------
311  * Write time plot values
312  *
313  * This function should only be called by one rank for a given data series.
314  *
315  * parameters:
316  * p <-- pointer to associated plot structure
317  * tn <-- associated time step number
318  * t <-- associated time value
319  * n_vals <-- number of associated time values
320  * vals <-- associated time values
321  *----------------------------------------------------------------------------*/
322 
323 void
324 cs_time_plot_vals_write(cs_time_plot_t *p,
325  int tn,
326  double t,
327  int n_vals,
328  const cs_real_t vals[]);
329 
330 /*----------------------------------------------------------------------------*/
331 
333 
334 #endif /* __CS_PROBE_H__ */