Report Typos and Errors    
Semi-Lagrangian Library
Modular library for kinetic and gyrokinetic simulations of plasmas in fusion energy devices.
sll_m_xdmf.F90
Go to the documentation of this file.
1 !**************************************************************
2 ! Copyright INRIA
3 !
4 ! This code SeLaLib (for Semi-Lagrangian-Library)
5 ! is a parallel library for simulating the plasma turbulence
6 ! in a tokamak.
7 !
8 ! This software is governed by the CeCILL-B license
9 ! under French law and abiding by the rules of distribution
10 ! of free software. You can use, modify and redistribute
11 ! the software under the terms of the CeCILL-B license as
12 ! circulated by CEA, CNRS and INRIA at the following URL
13 ! "http://www.cecill.info".
14 !**************************************************************
15 
24 module sll_m_xdmf
25 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
26 #include "sll_assert.h"
27 #include "sll_errors.h"
28 #include "sll_memory.h"
29 #include "sll_working_precision.h"
30 
31  use sll_m_ascii_io, only: &
33 
34  use sll_m_utilities, only: &
36 
37  use sll_m_xml_io, only: &
42 
43 #ifdef NOHDF5
44  use sll_m_binary_io, only: &
48 
49 #else
50  use sll_m_hdf5_io_serial, only: &
55 
56 #endif
57  implicit none
58 
59  public :: &
60  sll_s_plot_f, &
71 
72  private
73 !+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
74 
76  interface sll_o_xdmf_open
77  module procedure sll_xdmf_open_2d
78  module procedure sll_xdmf_open_3d
79  end interface
80 
83  module procedure sll_xdmf_array_2d
84  module procedure sll_xdmf_array_3d
85  end interface
86 
87 contains
88 
90  subroutine sll_xdmf_set_time(file_id, time)
91  sll_real64, intent(in) :: time
92  sll_int32, intent(in) :: file_id
93  logical :: i_opened
94 
95  inquire (file_id, opened=i_opened)
96 
97  if (i_opened) then
98  write (file_id, "(a13,g15.3,a3)") "<Time Value='", time, "'/>"
99  else
100  sll_error("sll_xdmf_set_time", "This xdmf file is not not open.")
101  end if
102 
103  end subroutine sll_xdmf_set_time
104 
106  subroutine sll_xdmf_open_2d(file_name, mesh_name, nnodes_x1, nnodes_x2, file_id, error)
107 
108  character(len=*), intent(in) :: file_name
109  character(len=*), intent(in) :: mesh_name
110  sll_int32, intent(in) :: nnodes_x1
111  sll_int32, intent(in) :: nnodes_x2
112  sll_int32, intent(out) :: file_id
113  sll_int32, intent(out) :: error
114 
115  call sll_s_xml_file_create(trim(file_name), file_id, error)
116  call sll_o_xml_grid_geometry(file_id, trim(mesh_name), &
117  nnodes_x1, nnodes_x2)
118 
119  end subroutine sll_xdmf_open_2d
120 
122  subroutine sll_xdmf_open_3d(file_name, &
123  mesh_name, &
124  nnodes_x1, &
125  nnodes_x2, &
126  nnodes_x3, &
127  file_id, &
128  error)
129 
130  character(len=*), intent(in) :: file_name
131  character(len=*), intent(in) :: mesh_name
132  sll_int32, intent(out) :: file_id
133  sll_int32, intent(out) :: error
134  sll_int32, intent(in) :: nnodes_x1
135  sll_int32, intent(in) :: nnodes_x2
136  sll_int32, intent(in) :: nnodes_x3
137 
138  call sll_s_xml_file_create(trim(file_name), file_id, error)
139  call sll_o_xml_grid_geometry(file_id, trim(mesh_name), &
140  nnodes_x1, nnodes_x2, nnodes_x3)
141 
142  end subroutine sll_xdmf_open_3d
143 
145  subroutine sll_xdmf_array_2d(mesh_name, array, array_name, error, xmffile_id, center)
146 
147  character(len=*), intent(in) :: mesh_name
148  sll_real64, intent(in) :: array(:, :)
149  character(len=*), intent(in) :: array_name
150  sll_int32, intent(out) :: error
151  sll_int32 :: npoints_x1
152  sll_int32 :: npoints_x2
153  sll_int32, intent(in), optional :: xmffile_id
154  character(len=4), optional :: center
155 
156 #ifndef NOHDF5
157  type(sll_t_hdf5_ser_handle) :: hfile_id
158 #else
159  sll_int32 :: file_id
160 #endif
161 
162  npoints_x1 = size(array, 1)
163  npoints_x2 = size(array, 2)
164 
165 #ifdef NOHDF5
166  call sll_s_binary_file_create(trim(mesh_name)//"-"//trim(array_name)//".bin", &
167  file_id, error)
168  call sll_o_binary_write_array(file_id, array, error)
169  call sll_s_binary_file_close(file_id, error)
170 #else
171  call sll_s_hdf5_ser_file_create(trim(mesh_name)//"-"//trim(array_name)//".h5", &
172  hfile_id, error)
173  call sll_o_hdf5_ser_write_array(hfile_id, array, "/"//trim(array_name), error)
174  call sll_s_hdf5_ser_file_close(hfile_id, error)
175 #endif
176 
177  if (present(xmffile_id) .and. present(center)) then
178 #ifdef NOHDF5
179  call sll_o_xml_field(xmffile_id, trim(array_name), &
180  trim(mesh_name)//"-"//trim(array_name)//".bin", &
181  npoints_x1, npoints_x2, 'Binary', center)
182 #else
183  call sll_o_xml_field( &
184  xmffile_id, &
185  trim(array_name), &
186  trim(mesh_name)//"-"//trim(array_name)//".h5:/"//trim(array_name), &
187  npoints_x1, &
188  npoints_x2, &
189  'HDF', &
190  center)
191 #endif
192  end if
193  end subroutine sll_xdmf_array_2d
194 
196  subroutine sll_xdmf_array_3d(mesh_name, array, array_name, error, xmffile_id, center)
197 
198  character(len=*), intent(in) :: mesh_name
199  sll_real64, intent(in) :: array(:, :, :)
200  character(len=*), intent(in) :: array_name
201  sll_int32, intent(out) :: error
202  sll_int32, intent(in), optional :: xmffile_id
203  character(len=4), optional :: center
204  sll_int32 :: npoints_x1
205  sll_int32 :: npoints_x2
206  sll_int32 :: npoints_x3
207 
208 #ifndef NOHDF5
209  type(sll_t_hdf5_ser_handle) :: hfile_id
210 #else
211  sll_int32 :: file_id
212 #endif
213 
214  npoints_x1 = size(array, 1)
215  npoints_x2 = size(array, 2)
216  npoints_x3 = size(array, 3)
217 
218 #ifdef NOHDF5
219  call sll_s_binary_file_create(trim(mesh_name)//"-"//trim(array_name)//".bin", &
220  file_id, error)
221  call sll_o_binary_write_array(file_id, array, error)
222  call sll_s_binary_file_close(file_id, error)
223 #else
224  call sll_s_hdf5_ser_file_create(trim(mesh_name)//"-"//trim(array_name)//".h5", &
225  hfile_id, error)
226  call sll_o_hdf5_ser_write_array(hfile_id, array, "/"//trim(array_name), error)
227  call sll_s_hdf5_ser_file_close(hfile_id, error)
228 #endif
229 
230  if (present(xmffile_id) .and. present(center)) then
231 
232 #ifdef NOHDF5
233  call sll_o_xml_field(xmffile_id, trim(array_name), &
234  trim(mesh_name)//"-"//trim(array_name)//".bin", &
235  npoints_x1, npoints_x2, npoints_x3, 'Binary', center)
236 #else
237  call sll_o_xml_field( &
238  xmffile_id, &
239  trim(array_name), &
240  trim(mesh_name)//"-"//trim(array_name)//".h5:/"//trim(array_name), &
241  npoints_x1, &
242  npoints_x2, &
243  npoints_x3, &
244  'HDF', &
245  center)
246 #endif
247 
248  end if
249 
250  end subroutine sll_xdmf_array_3d
251 
255  subroutine sll_s_xdmf_corect2d_nodes(file_name, &
256  array, &
257  array_name, &
258  eta1_min, &
259  delta_eta1, &
260  eta2_min, &
261  delta_eta2, &
262  file_format, &
263  iplot, &
264  time)
265 
266  sll_real64, intent(in) :: array(:, :)
267  character(len=*), intent(in) :: file_name
268  character(len=*), intent(in) :: array_name
269  sll_int32 :: error
270  sll_real64 :: eta1_min
271  sll_real64 :: eta2_min
272  sll_real64 :: delta_eta1
273  sll_real64 :: delta_eta2
274  sll_int32 :: file_id
275  sll_int32 :: nx1
276  sll_int32 :: nx2
277  character(len=4), optional :: file_format
278  sll_int32, optional :: iplot
279  sll_real64, optional :: time
280 
281  character(len=4) :: cplot
282 #ifndef NOHDF5
283  type(sll_t_hdf5_ser_handle) :: hfile_id
284 #endif
285 
286  nx1 = size(array, 1)
287  nx2 = size(array, 2)
288 
289  if (present(iplot)) then
290  call sll_s_int2string(iplot, cplot)
291  call sll_s_xml_file_create(file_name//cplot//".xmf", file_id, error)
292  else
293  call sll_s_xml_file_create(file_name//".xmf", file_id, error)
294  end if
295 
296  write (file_id, "(a)") "<Grid Name='mesh' GridType='Uniform'>"
297  if (present(time)) then
298  write (file_id, "(a,f12.6,a)") "<Time Value='", time, "'/>"
299  end if
300  write (file_id, "(a,2i5,a)") "<Topology TopologyType='2DCoRectMesh' NumberOfElements='", &
301  nx2, nx1, "'/>"
302  write (file_id, "(a)") "<Geometry GeometryType='ORIGIN_DXDY'>"
303  write (file_id, "(a)") "<DataItem Dimensions='2' NumberType='Float' Format='XML'>"
304  write (file_id, "(2f12.5)") eta1_min, eta2_min
305  write (file_id, "(a)") "</DataItem>"
306  write (file_id, "(a)") "<DataItem Dimensions='2' NumberType='Float' Format='XML'>"
307  write (file_id, "(2f12.5)") delta_eta1, delta_eta2
308  write (file_id, "(a)") "</DataItem>"
309  write (file_id, "(a)") "</Geometry>"
310  write (file_id, "(a)") "<Attribute Name='"//array_name//"' AttributeType='Scalar' Center='Node'>"
311 
312  if (present(file_format)) then
313  if (file_format == "HDF5") then
314  write (file_id, "(a,2i5,a)") "<DataItem Dimensions='", nx2, nx1, &
315  "' NumberType='Float' Precision='8' Format='HDF'>"
316  if (present(iplot)) then
317  write (file_id, "(a)") array_name//cplot//".h5:/node_values"
318  else
319  write (file_id, "(a)") array_name//".h5:/node_values"
320  end if
321 #ifndef NOHDF5
322  if (present(iplot)) then
323  call sll_s_hdf5_ser_file_create(array_name//cplot//".h5", hfile_id, error)
324  else
325  call sll_s_hdf5_ser_file_create(array_name//".h5", hfile_id, error)
326  end if
327  call sll_o_hdf5_ser_write_array(hfile_id, array, "/node_values", error)
328  call sll_s_hdf5_ser_file_close(hfile_id, error)
329 #endif
330  end if
331  else
332  write (file_id, "(a,2i5,a)") "<DataItem Dimensions='", nx2, nx1, &
333  "' NumberType='Float' Precision='8' Format='XML'>"
334  call sll_o_ascii_write_array(file_id, array, error)
335  end if
336  write (file_id, "(a)") "</DataItem>"
337  write (file_id, "(a)") "</Attribute>"
338  call sll_s_xml_file_close(file_id, error)
339 
340  end subroutine sll_s_xdmf_corect2d_nodes
341 
345  subroutine sll_s_xdmf_corect3d_nodes(file_name, &
346  array, &
347  array_name, &
348  eta1_min, &
349  delta_eta1, &
350  eta2_min, &
351  delta_eta2, &
352  eta3_min, &
353  delta_eta3, &
354  file_format, &
355  iplot)
356 
357  sll_real64, intent(in) :: array(:, :, :)
358  character(len=*), intent(in) :: file_name
359  character(len=*), intent(in) :: array_name
360  sll_int32 :: error
361  sll_real64 :: eta1_min
362  sll_real64 :: eta2_min
363  sll_real64 :: eta3_min
364  sll_real64 :: delta_eta1
365  sll_real64 :: delta_eta2
366  sll_real64 :: delta_eta3
367  sll_int32 :: file_id
368  sll_int32 :: nx1
369  sll_int32 :: nx2
370  sll_int32 :: nx3
371  character(len=4), optional :: file_format
372  sll_int32, optional :: iplot
373 
374  character(len=4) :: cplot
375  logical :: hdf5_file_format
376 #ifndef NOHDF5
377  type(sll_t_hdf5_ser_handle) :: hfile_id
378 #endif
379 
380  nx1 = size(array, 1)
381  nx2 = size(array, 2)
382  nx3 = size(array, 3)
383 
384  if (present(iplot)) then
385  call sll_s_int2string(iplot, cplot)
386  call sll_s_xml_file_create(file_name//cplot//".xmf", file_id, error)
387  else
388  call sll_s_xml_file_create(file_name//".xmf", file_id, error)
389  end if
390 
391  write (file_id, "(a)") "<Grid Name='mesh' GridType='Uniform'>"
392  write (file_id, "(a,3i5,a)") "<Topology TopologyType='3DCoRectMesh' NumberOfElements='", &
393  nx3, nx2, nx1, "'/>"
394  write (file_id, "(a)") "<Geometry GeometryType='ORIGIN_DXDYDZ'>"
395  write (file_id, "(a)") "<DataItem Dimensions='3' NumberType='Float' Format='XML'>"
396  write (file_id, "(3f12.5)") eta1_min, eta2_min, eta3_min
397  write (file_id, "(a)") "</DataItem>"
398  write (file_id, "(a)") "<DataItem Dimensions='3' NumberType='Float' Format='XML'>"
399  write (file_id, "(3f12.5)") delta_eta1, delta_eta2, delta_eta3
400  write (file_id, "(a)") "</DataItem>"
401  write (file_id, "(a)") "</Geometry>"
402  write (file_id, "(a)") "<Attribute Name='"//array_name//"' AttributeType='Scalar' Center='Node'>"
403 
404  hdf5_file_format = .false.
405  if (present(file_format)) then
406  if (file_format == "HDF5") hdf5_file_format = .true.
407  end if
408  if (hdf5_file_format) then
409  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
410  "' NumberType='Float' Precision='8' Format='HDF'>"
411  if (present(iplot)) then
412  write (file_id, "(a)") array_name//cplot//".h5:/node_values"
413  else
414  write (file_id, "(a)") array_name//".h5:/node_values"
415  end if
416 #ifndef NOHDF5
417  if (present(iplot)) then
418  call sll_s_hdf5_ser_file_create(array_name//cplot//".h5", hfile_id, error)
419  else
420  call sll_s_hdf5_ser_file_create(array_name//".h5", hfile_id, error)
421  end if
422  call sll_o_hdf5_ser_write_array(hfile_id, array, "/node_values", error)
423  call sll_s_hdf5_ser_file_close(hfile_id, error)
424 #endif
425  else
426  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
427  "' NumberType='Float' Precision='8' Format='XML'>"
428  call sll_o_ascii_write_array(file_id, array, error)
429  end if
430  write (file_id, "(a)") "</DataItem>"
431  write (file_id, "(a)") "</Attribute>"
432  call sll_s_xml_file_close(file_id, error)
433 
434  end subroutine sll_s_xdmf_corect3d_nodes
435 
439  subroutine sll_s_xdmf_rect2d_nodes(file_name, &
440  array, &
441  array_name, &
442  eta1, &
443  eta2, &
444  file_format, &
445  iplot, &
446  time)
447 
448  sll_real64, intent(in) :: array(:, :)
449  sll_real64, intent(in) :: eta1(:)
450  sll_real64, intent(in) :: eta2(:)
451  character(len=*), intent(in) :: file_name
452  character(len=*), intent(in) :: array_name
453  sll_int32 :: error
454  sll_int32 :: file_id
455  sll_int32 :: nx1
456  sll_int32 :: nx2
457  character(len=4), optional :: file_format
458  sll_int32, optional :: iplot
459  sll_real64, optional :: time
460  sll_int32 :: i, j
461  character(len=4) :: cplot
462  logical :: hdf5_file_format
463 #ifndef NOHDF5
464  type(sll_t_hdf5_ser_handle) :: hfile_id
465 #endif
466 
467  nx1 = size(array, 1)
468  nx2 = size(array, 2)
469 
470  sll_assert_always(nx1 == size(eta1))
471  sll_assert_always(nx2 == size(eta2))
472 
473  if (present(iplot)) then
474  call sll_s_int2string(iplot, cplot)
475  call sll_s_xml_file_create(file_name//cplot//".xmf", file_id, error)
476  else
477  call sll_s_xml_file_create(file_name//".xmf", file_id, error)
478  end if
479  write (file_id, "(a)") "<Grid Name='mesh' GridType='Uniform'>"
480  if (present(time)) then
481  write (file_id, "(a,f12.6,a)") "<Time Value='", time, "'/>"
482  end if
483  write (file_id, "(a,2i5,a)") "<Topology TopologyType='2DRectMesh' NumberOfElements='", &
484  nx2, nx1, "'/>"
485  write (file_id, "(a)") "<Geometry GeometryType='VXVY'>"
486  write (file_id, "(a,i5,a)") "<DataItem Dimensions='", nx1, &
487  "' NumberType='Float' Format='XML'>"
488  write (file_id, *) (eta1(i), i=1, nx1)
489  write (file_id, "(a)") "</DataItem>"
490  write (file_id, "(a,i5,a)") "<DataItem Dimensions='", nx2, &
491  "' NumberType='Float' Format='XML'>"
492  write (file_id, *) (eta2(j), j=1, nx2)
493  write (file_id, "(a)") "</DataItem>"
494  write (file_id, "(a)") "</Geometry>"
495  write (file_id, "(a)") "<Attribute Name='"//array_name//"' AttributeType='Scalar' Center='Node'>"
496  hdf5_file_format = .false.
497  if (present(file_format)) then
498  if (file_format == "HDF5") hdf5_file_format = .true.
499  end if
500  if (hdf5_file_format) then
501  write (file_id, "(a,2i5,a)") "<DataItem Dimensions='", nx2, nx1, &
502  "' NumberType='Float' Precision='8' Format='HDF'>"
503  if (present(iplot)) then
504  write (file_id, "(a)") array_name//cplot//".h5:/node_values"
505  else
506  write (file_id, "(a)") array_name//".h5:/node_values"
507  end if
508 #ifndef NOHDF5
509  if (present(iplot)) then
510  call sll_s_hdf5_ser_file_create(array_name//cplot//".h5", hfile_id, error)
511  else
512  call sll_s_hdf5_ser_file_create(array_name//".h5", hfile_id, error)
513  end if
514  call sll_o_hdf5_ser_write_array(hfile_id, array, "/node_values", error)
515  call sll_s_hdf5_ser_file_close(hfile_id, error)
516 #endif
517  else
518  write (file_id, "(a,2i5,a)") "<DataItem Dimensions='", nx2, nx1, &
519  "' NumberType='Float' Precision='4' Format='XML'>"
520  call sll_o_ascii_write_array(file_id, array, error)
521  end if
522  write (file_id, "(a)") "</DataItem>"
523  write (file_id, "(a)") "</Attribute>"
524  call sll_s_xml_file_close(file_id, error)
525 
526  end subroutine sll_s_xdmf_rect2d_nodes
527 
531  subroutine sll_s_xdmf_rect3d_nodes(file_name, &
532  array, &
533  array_name, &
534  eta1, &
535  eta2, &
536  eta3, &
537  file_format, &
538  iplot)
539 
540  sll_real64, intent(in) :: array(:, :, :)
541  sll_real64, intent(in) :: eta1(:)
542  sll_real64, intent(in) :: eta2(:)
543  sll_real64, intent(in) :: eta3(:)
544  character(len=*), intent(in) :: file_name
545  character(len=*), intent(in) :: array_name
546  sll_int32 :: error
547  sll_int32 :: file_id
548  sll_int32 :: nx1
549  sll_int32 :: nx2
550  sll_int32 :: nx3
551  character(len=4), optional :: file_format
552  sll_int32, optional :: iplot
553  sll_int32 :: i, j, k
554 #ifndef NOHDF5
555  type(sll_t_hdf5_ser_handle) :: hfile_id
556 #endif
557  character(len=4) :: cplot
558  logical :: hdf5_file_format
559 
560  nx1 = size(array, 1)
561  nx2 = size(array, 2)
562  nx3 = size(array, 3)
563 
564  sll_assert_always(nx1 == size(eta1))
565  sll_assert_always(nx2 == size(eta2))
566  sll_assert_always(nx3 == size(eta3))
567 
568  if (present(iplot)) then
569  call sll_s_int2string(iplot, cplot)
570  call sll_s_xml_file_create(file_name//cplot//".xmf", file_id, error)
571  else
572  call sll_s_xml_file_create(file_name//".xmf", file_id, error)
573  end if
574 
575  write (file_id, "(a)") "<Grid Name='mesh' GridType='Uniform'>"
576  write (file_id, "(a,3i5,a)") "<Topology TopologyType='3DRectMesh' NumberOfElements='", &
577  nx3, nx2, nx1, "'/>"
578  write (file_id, "(a)") "<Geometry GeometryType='VXVYVZ'>"
579  write (file_id, "(a,i5,a)") "<DataItem Dimensions='", nx1, &
580  "' NumberType='Float' Format='XML'>"
581  write (file_id, *) (eta1(i), i=1, nx1)
582  write (file_id, "(a)") "</DataItem>"
583  write (file_id, "(a,i5,a)") "<DataItem Dimensions='", nx2, &
584  "' NumberType='Float' Format='XML'>"
585  write (file_id, *) (eta2(j), j=1, nx2)
586  write (file_id, "(a)") "</DataItem>"
587  write (file_id, "(a,i5,a)") "<DataItem Dimensions='", nx3, &
588  "' NumberType='Float' Format='XML'>"
589  write (file_id, *) (eta3(k), k=1, nx3)
590  write (file_id, "(a)") "</DataItem>"
591  write (file_id, "(a)") "</Geometry>"
592  write (file_id, "(a)") "<Attribute Name='"//array_name//"' AttributeType='Scalar' Center='Node'>"
593  hdf5_file_format = .false.
594  if (present(file_format)) then
595  if (file_format == "HDF5") hdf5_file_format = .true.
596  end if
597  if (hdf5_file_format) then
598  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
599  "' NumberType='Float' Precision='8' Format='HDF'>"
600  if (present(iplot)) then
601  write (file_id, "(a)") array_name//cplot//".h5:/node_values"
602  else
603  write (file_id, "(a)") array_name//".h5:/node_values"
604  end if
605 #ifndef NOHDF5
606  if (present(iplot)) then
607  call sll_s_hdf5_ser_file_create(array_name//cplot//".h5", hfile_id, error)
608  else
609  call sll_s_hdf5_ser_file_create(array_name//".h5", hfile_id, error)
610  end if
611  call sll_o_hdf5_ser_write_array(hfile_id, array, "/node_values", error)
612  call sll_s_hdf5_ser_file_close(hfile_id, error)
613 #endif
614  else
615  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
616  "' NumberType='Float' Precision='4' Format='XML'>"
617  call sll_o_ascii_write_array(file_id, array, error)
618  end if
619  write (file_id, "(a)") "</DataItem>"
620  write (file_id, "(a)") "</Attribute>"
621  call sll_s_xml_file_close(file_id, error)
622 
623  end subroutine sll_s_xdmf_rect3d_nodes
624 
628  subroutine sll_s_xdmf_curv2d_nodes(file_name, &
629  array, &
630  array_name, &
631  eta1, &
632  eta2, &
633  file_format, &
634  iplot)
635 
636  sll_real64, intent(in) :: array(:, :)
637  sll_real64, intent(in) :: eta1(:, :)
638  sll_real64, intent(in) :: eta2(:, :)
639  character(len=*), intent(in) :: file_name
640  character(len=*), intent(in) :: array_name
641  sll_int32 :: error
642  sll_int32 :: file_id
643  sll_int32 :: nx1
644  sll_int32 :: nx2
645  character(len=4), optional :: file_format
646  sll_int32, optional :: iplot
647  character(len=4) :: cplot
648  logical :: hdf5_file_format
649 #ifndef NOHDF5
650  type(sll_t_hdf5_ser_handle) :: hfile_id
651 #endif
652 
653  nx1 = size(array, 1)
654  nx2 = size(array, 2)
655 
656  sll_assert_always(nx1 == size(eta1, 1))
657  sll_assert_always(nx2 == size(eta2, 2))
658 
659  if (present(iplot)) then
660  call sll_s_int2string(iplot, cplot)
661  call sll_s_xml_file_create(file_name//cplot//".xmf", file_id, error)
662  else
663  call sll_s_xml_file_create(file_name//".xmf", file_id, error)
664  end if
665 
666  write (file_id, "(a)") "<Grid Name='mesh' GridType='Uniform'>"
667  write (file_id, "(a,2i5,a)") "<Topology TopologyType='2DSMesh' NumberOfElements='", &
668  nx2, nx1, "'/>"
669  write (file_id, "(a)") "<Geometry GeometryType='X_Y'>"
670 
671  hdf5_file_format = .false.
672  if (present(file_format)) then
673  if (file_format == "HDF5") hdf5_file_format = .true.
674  end if
675  if (hdf5_file_format) then
676 
677  write (file_id, "(a,2i5,a)") "<DataItem Dimensions='", nx2, nx1, &
678  "' NumberType='Float' Precision='8' Format='HDF'>"
679  write (file_id, "(a)") array_name//".h5:/x1_values"
680  write (file_id, "(a)") "</DataItem>"
681  write (file_id, "(a,2i5,a)") "<DataItem Dimensions='", nx2, nx1, &
682  "' NumberType='Float' Precision='8' Format='HDF'>"
683  if (present(iplot)) then
684  write (file_id, "(a)") array_name//cplot//".h5:/x2_values"
685  else
686  write (file_id, "(a)") array_name//".h5:/x2_values"
687  end if
688  write (file_id, "(a)") "</DataItem>"
689 
690 #ifndef NOHDF5
691  call sll_s_hdf5_ser_file_create(file_name//".h5", hfile_id, error)
692  call sll_o_hdf5_ser_write_array(hfile_id, eta1, "/x1_values", error)
693  call sll_o_hdf5_ser_write_array(hfile_id, eta2, "/x2_values", error)
694 #endif
695 
696  else
697 
698  write (file_id, "(a,2i5,a)") "<DataItem Dimensions='", nx2, nx1, &
699  "' NumberType='Float' Precision='4' Format='XML'>"
700  call sll_o_ascii_write_array(file_id, eta1, error)
701  write (file_id, "(a)") "</DataItem>"
702  write (file_id, "(a,2i5,a)") "<DataItem Dimensions='", nx2, nx1, &
703  "' NumberType='Float' Precision='4' Format='XML'>"
704  call sll_o_ascii_write_array(file_id, eta2, error)
705  write (file_id, "(a)") "</DataItem>"
706 
707  end if
708 
709  write (file_id, "(a)") "</Geometry>"
710  write (file_id, "(a)") &
711  "<Attribute Name='"//array_name//"' AttributeType='Scalar' Center='Node'>"
712 
713 #ifndef NOHDF5
714  if (hdf5_file_format) then
715 
716  write (file_id, "(a,2i5,a)") "<DataItem Dimensions='", nx2, nx1, &
717  "' NumberType='Float' Precision='8' Format='HDF'>"
718  if (present(iplot)) then
719  write (file_id, "(a)") array_name//cplot//".h5:/node_values"
720  else
721  write (file_id, "(a)") array_name//".h5:/node_values"
722  end if
723 
724  call sll_o_hdf5_ser_write_array(hfile_id, array, "/node_values", error)
725  call sll_s_hdf5_ser_file_close(hfile_id, error)
726 
727  else
728 #endif
729 
730  write (file_id, "(a,2i5,a)") "<DataItem Dimensions='", nx2, nx1, &
731  "' NumberType='Float' Precision='4' Format='XML'>"
732  call sll_o_ascii_write_array(file_id, array, error)
733 
734 #ifndef NOHDF5
735  end if
736 #endif
737 
738  write (file_id, "(a)") "</DataItem>"
739  write (file_id, "(a)") "</Attribute>"
740  call sll_s_xml_file_close(file_id, error)
741 
742  end subroutine sll_s_xdmf_curv2d_nodes
743 
760 
761  subroutine sll_s_xdmf_curv3d_nodes(file_name, &
762  array, &
763  array_name, &
764  eta1, &
765  eta2, &
766  eta3, &
767  file_format, &
768  iplot)
769 
770  sll_real64, intent(in) :: array(:, :, :)
771  sll_real64, intent(in) :: eta1(:, :, :)
772  sll_real64, intent(in) :: eta2(:, :, :)
773  sll_real64, intent(in) :: eta3(:, :, :)
774  character(len=*), intent(in) :: file_name
775  character(len=*), intent(in) :: array_name
776  sll_int32 :: error
777  sll_int32 :: file_id
778  sll_int32 :: nx1
779  sll_int32 :: nx2
780  sll_int32 :: nx3
781  character(len=4), optional :: file_format
782  sll_int32, optional :: iplot
783 
784  character(len=4) :: cplot
785  logical :: hdf5_file_format
786 #ifndef NOHDF5
787  type(sll_t_hdf5_ser_handle) :: hfile_id
788 #endif
789 
790  nx1 = size(array, 1)
791  nx2 = size(array, 2)
792  nx3 = size(array, 3)
793 
794  sll_assert_always(nx1 == size(eta1, 1))
795  sll_assert_always(nx2 == size(eta2, 2))
796  sll_assert_always(nx3 == size(eta3, 3))
797 
798  if (present(iplot)) then
799  call sll_s_int2string(iplot, cplot)
800  call sll_s_xml_file_create(file_name//cplot//".xmf", file_id, error)
801  else
802  call sll_s_xml_file_create(file_name//".xmf", file_id, error)
803  end if
804 
805  write (file_id, "(a)") "<Grid Name='mesh' GridType='Uniform'>"
806  write (file_id, "(a,3i5,a)") "<Topology TopologyType='3DSMesh' NumberOfElements='", &
807  nx3, nx2, nx1, "'/>"
808  write (file_id, "(a)") "<Geometry GeometryType='X_Y_Z'>"
809 
810  hdf5_file_format = .false.
811  if (present(file_format)) then
812  if (file_format == "HDF5") hdf5_file_format = .true.
813  end if
814  if (hdf5_file_format) then
815 
816  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
817  "' NumberType='Float' Precision='8' Format='HDF'>"
818  write (file_id, "(a)") array_name//".h5:/x1_values"
819  write (file_id, "(a)") "</DataItem>"
820  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
821  "' NumberType='Float' Precision='8' Format='HDF'>"
822  write (file_id, "(a)") array_name//".h5:/x2_values"
823  write (file_id, "(a)") "</DataItem>"
824  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
825  "' NumberType='Float' Precision='8' Format='HDF'>"
826  write (file_id, "(a)") array_name//".h5:/x3_values"
827  write (file_id, "(a)") "</DataItem>"
828 
829 #ifndef NOHDF5
830  call sll_s_hdf5_ser_file_create(array_name//".h5", hfile_id, error)
831  call sll_o_hdf5_ser_write_array(hfile_id, eta1, "/x1_values", error)
832  call sll_o_hdf5_ser_write_array(hfile_id, eta2, "/x2_values", error)
833  call sll_o_hdf5_ser_write_array(hfile_id, eta3, "/x3_values", error)
834 #endif
835 
836  else
837 
838  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
839  "' NumberType='Float' Precision='4' Format='XML'>"
840  call sll_o_ascii_write_array(file_id, eta1, error)
841  write (file_id, "(a)") "</DataItem>"
842  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
843  "' NumberType='Float' Precision='4' Format='XML'>"
844  call sll_o_ascii_write_array(file_id, eta2, error)
845  write (file_id, "(a)") "</DataItem>"
846  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
847  "' NumberType='Float' Precision='4' Format='XML'>"
848  call sll_o_ascii_write_array(file_id, eta3, error)
849  write (file_id, "(a)") "</DataItem>"
850 
851  end if
852 
853  write (file_id, "(a)") "</Geometry>"
854  if (present(iplot)) then
855  write (file_id, "(a)") &
856  "<Attribute Name='"//array_name//cplot//"' AttributeType='Scalar' Center='Node'>"
857  else
858  write (file_id, "(a)") &
859  "<Attribute Name='"//array_name//"' AttributeType='Scalar' Center='Node'>"
860  end if
861 
862 #ifndef NOHDF5
863  if (hdf5_file_format) then
864 
865  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
866  "' NumberType='Float' Precision='8' Format='HDF'>"
867  if (present(iplot)) then
868  write (file_id, "(a)") array_name//cplot//".h5:/node_values"
869  else
870  write (file_id, "(a)") array_name//".h5:/node_values"
871  end if
872 
873  call sll_o_hdf5_ser_write_array(hfile_id, array, "/node_values", error)
874  call sll_s_hdf5_ser_file_close(hfile_id, error)
875 
876  else
877 #endif
878 
879  write (file_id, "(a,3i5,a)") "<DataItem Dimensions='", nx3, nx2, nx1, &
880  "' NumberType='Float' Precision='4' Format='XML'>"
881  call sll_o_ascii_write_array(file_id, array, error)
882 
883 #ifndef NOHDF5
884  end if
885 #endif
886 
887  write (file_id, "(a)") "</DataItem>"
888  write (file_id, "(a)") "</Attribute>"
889  call sll_s_xml_file_close(file_id, error)
890 
891  end subroutine sll_s_xdmf_curv3d_nodes
892 
894  subroutine sll_s_xdmf_close(file_id, error)
895  sll_int32, intent(in) :: file_id
896  sll_int32, intent(out) :: error
897 
898  write (file_id, "(a)") "</Grid>"
899  write (file_id, "(a)") "</Domain>"
900  write (file_id, "(a)") "</Xdmf>"
901  close (file_id)
902  error = 0
903  end subroutine sll_s_xdmf_close
904 
916  subroutine sll_s_plot_f_cartesian(iplot, &
917  f, &
918  vec_x1, &
919  nnodes_x1, &
920  vec_x2, &
921  nnodes_x2, &
922  array_name, &
923  time)
924 
925  sll_int32, intent(in) :: iplot
926  sll_real64, dimension(:, :), intent(in) :: f
927  sll_real64, dimension(:), intent(in) :: vec_x1
928  sll_int32, intent(in) :: nnodes_x1
929  sll_real64, dimension(:), intent(in) :: vec_x2
930  sll_int32, intent(in) :: nnodes_x2
931  character(len=*), intent(in) :: array_name
932  sll_real64 :: time
933 
934  sll_int32 :: file_id
935  sll_int32 :: error
936  sll_real64, dimension(:, :), allocatable :: x1
937  sll_real64, dimension(:, :), allocatable :: x2
938  sll_int32 :: i
939  sll_int32 :: j
940  character(len=4) :: cplot
941 
942 #ifndef NOHDF5
943  type(sll_t_hdf5_ser_handle) :: hfile_id
944 #endif
945 
946  sll_assert_always(iplot > 0)
947  if (iplot == 1) then
948 
949  sll_allocate(x1(nnodes_x1, nnodes_x2), error)
950  sll_allocate(x2(nnodes_x1, nnodes_x2), error)
951  do j = 1, nnodes_x2
952  do i = 1, nnodes_x1
953  x1(i, j) = vec_x1(i) !x1_min+real(i-1,f32)*dx1
954  x2(i, j) = vec_x2(j) !x2_min+real(j-1,f32)*dx2
955  end do
956  end do
957 
958 #ifndef NOHDF5
959  call sll_s_hdf5_ser_file_create("cartesian_mesh-x1.h5", hfile_id, error)
960  call sll_o_hdf5_ser_write_array(hfile_id, x1, "/x1", error)
961  call sll_s_hdf5_ser_file_close(hfile_id, error)
962  call sll_s_hdf5_ser_file_create("cartesian_mesh-x2.h5", hfile_id, error)
963  call sll_o_hdf5_ser_write_array(hfile_id, x2, "/x2", error)
964  call sll_s_hdf5_ser_file_close(hfile_id, error)
965 #endif
966 
967  deallocate (x1)
968  deallocate (x2)
969 
970  end if
971 
972  call sll_s_int2string(iplot, cplot)
973  call sll_o_xdmf_open(trim(array_name)//cplot//".xmf", "cartesian_mesh", &
974  nnodes_x1, nnodes_x2, file_id, error)
975 
976  write (file_id, "(a,f8.3,a)") "<Time Value='", time, "'/>"
977 
978  call sll_o_xdmf_write_array(trim(array_name)//cplot, f, "values", &
979  error, file_id, "Node")
980  call sll_s_xdmf_close(file_id, error)
981 
982  end subroutine sll_s_plot_f_cartesian
983 
997 
998  subroutine sll_s_plot_f( &
999  iplot, &
1000  f, &
1001  nnodes_x1, &
1002  nnodes_x2, &
1003  array_name, &
1004  mesh_name, &
1005  time, &
1006  x1, &
1007  x2)
1008 
1009  sll_int32, intent(in) :: iplot
1010  sll_real64, intent(in) :: f(:, :)
1011  sll_int32, intent(in) :: nnodes_x1
1012  sll_int32, intent(in) :: nnodes_x2
1013  character(len=*), intent(in) :: array_name
1014  character(len=*), intent(in) :: mesh_name
1015  sll_real64, intent(in) :: time
1016  sll_real64, intent(in), optional :: x1(:, :)
1017  sll_real64, intent(in), optional :: x2(:, :)
1018 
1019  sll_int32 :: file_id
1020  sll_int32 :: ierr
1021  character(len=4) :: cplot
1022 
1023 #ifndef NOHDF5
1024  type(sll_t_hdf5_ser_handle) :: hfile_id
1025 
1026  if (present(x1) .and. present(x2)) then
1027 
1028  call sll_s_hdf5_ser_file_create(trim(mesh_name)//"-x1.h5", hfile_id, ierr)
1029  call sll_o_hdf5_ser_write_array(hfile_id, x1, "/x1", ierr)
1030  call sll_s_hdf5_ser_file_close(hfile_id, ierr)
1031  call sll_s_hdf5_ser_file_create(trim(mesh_name)//"-x2.h5", hfile_id, ierr)
1032  call sll_o_hdf5_ser_write_array(hfile_id, x2, "/x2", ierr)
1033  call sll_s_hdf5_ser_file_close(hfile_id, ierr)
1034 
1035  end if
1036 
1037 #endif
1038 
1039  call sll_s_int2string(iplot, cplot)
1040  call sll_o_xdmf_open( &
1041  trim(array_name)//cplot//".xmf", &
1042  trim(mesh_name), &
1043  nnodes_x1, &
1044  nnodes_x2, &
1045  file_id, &
1046  ierr)
1047 
1048  write (file_id, "(a,f8.3,a)") "<Time Value='", time, "'/>"
1049 
1050  call sll_o_xdmf_write_array( &
1051  trim(array_name)//cplot, &
1052  f, &
1053  "values", &
1054  ierr, &
1055  file_id, &
1056  "Node")
1057 
1058  call sll_s_xdmf_close(file_id, ierr)
1059 
1060  end subroutine sll_s_plot_f
1061 
1062 end module sll_m_xdmf
Write nD array in ascii file.
Write a nD array in a binary file.
Write nD array of double precision floats or integers into HDF5 file.
Create a xmf file.
Definition: sll_m_xdmf.F90:76
Write the field in xdmf format.
Definition: sll_m_xdmf.F90:82
write a data attribute in the xml file
write grid description in the xml file
Module that contains routines to write data in ASCII format file.
Implements the functions to write binary file to store heavy data.
subroutine, public sll_s_binary_file_close(file_id, error)
Open binary file.
subroutine, public sll_s_binary_file_create(filename, file_id, error)
Create binary file.
Implements the functions to write hdf5 file to store heavy data.
subroutine, public sll_s_hdf5_ser_file_create(filename, handle, error)
Create new HDF5 file.
subroutine, public sll_s_hdf5_ser_file_close(handle, error)
Close existing HDF5 file.
Some common numerical utilities.
subroutine, public sll_s_int2string(istep, cstep)
Convert an integer < 9999 to a 4 characters string.
Implements the functions to write xdmf file plotable by VisIt.
Definition: sll_m_xdmf.F90:24
subroutine sll_xdmf_open_2d(file_name, mesh_name, nnodes_x1, nnodes_x2, file_id, error)
Open a XDMF format file for a 2d plot.
Definition: sll_m_xdmf.F90:107
subroutine, public sll_s_plot_f(iplot, f, nnodes_x1, nnodes_x2, array_name, mesh_name, time, x1, x2)
Plot 2d distribution function for VisIt.
subroutine, public sll_s_plot_f_cartesian(iplot, f, vec_x1, nnodes_x1, vec_x2, nnodes_x2, array_name, time)
Plot 2d distribution function for VisIt.
Definition: sll_m_xdmf.F90:924
subroutine sll_xdmf_open_3d(file_name, mesh_name, nnodes_x1, nnodes_x2, nnodes_x3, file_id, error)
Open a XDMF format file for a 3d plot.
Definition: sll_m_xdmf.F90:129
subroutine, public sll_s_xdmf_rect3d_nodes(file_name, array, array_name, eta1, eta2, eta3, file_format, iplot)
Subroutine to write a 3D array in xdmf format. The field is describe on a cartesian mesh....
Definition: sll_m_xdmf.F90:539
subroutine, public sll_s_xdmf_close(file_id, error)
Close the XML file and finish to write last lines.
Definition: sll_m_xdmf.F90:895
subroutine sll_xdmf_array_2d(mesh_name, array, array_name, error, xmffile_id, center)
Write 2d array in binary or hdf5 file and the matching line in XDMF file.
Definition: sll_m_xdmf.F90:146
subroutine sll_xdmf_array_3d(mesh_name, array, array_name, error, xmffile_id, center)
Write 3d array in binary or hdf5 file and the matching line in XDMF file.
Definition: sll_m_xdmf.F90:197
subroutine, public sll_s_xdmf_curv2d_nodes(file_name, array, array_name, eta1, eta2, file_format, iplot)
Subroutine to write a 2D array in xdmf format. The field is describe on a cartesian mesh....
Definition: sll_m_xdmf.F90:635
subroutine, public sll_s_xdmf_corect2d_nodes(file_name, array, array_name, eta1_min, delta_eta1, eta2_min, delta_eta2, file_format, iplot, time)
Subroutine to write a 2D array in xdmf format The field is describe on a cartesian mesh Axis are perp...
Definition: sll_m_xdmf.F90:265
subroutine sll_xdmf_set_time(file_id, time)
Add the the good value of time in VisIt plot.
Definition: sll_m_xdmf.F90:91
subroutine, public sll_s_xdmf_rect2d_nodes(file_name, array, array_name, eta1, eta2, file_format, iplot, time)
Subroutine to write a 2D array in xdmf format. The field is describe on a cartesian mesh....
Definition: sll_m_xdmf.F90:447
subroutine, public sll_s_xdmf_corect3d_nodes(file_name, array, array_name, eta1_min, delta_eta1, eta2_min, delta_eta2, eta3_min, delta_eta3, file_format, iplot)
Subroutine to write a 3D array in xdmf format The field is describe on a cartesian mesh Axis are perp...
Definition: sll_m_xdmf.F90:356
subroutine, public sll_s_xdmf_curv3d_nodes(file_name, array, array_name, eta1, eta2, eta3, file_format, iplot)
Subroutine to write a 3D array in xdmf format. The field is describe on a cartesian mesh....
Definition: sll_m_xdmf.F90:769
Implements the functions to write xml file to store light data.
subroutine, public sll_s_xml_file_close(file_id, error)
Close the XML file and finish to write last lines. You give the file unit number.
subroutine, public sll_s_xml_file_create(filename, file_id, error)
Create the XML file and begin to write first lines. You get the file unit number.
    Report Typos and Errors