[cig-commits] r18392 - mc/3D/CitcomCU/trunk/src
becker at geodynamics.org
becker at geodynamics.org
Wed May 18 12:22:22 PDT 2011
Author: becker
Date: 2011-05-18 12:22:22 -0700 (Wed, 18 May 2011)
New Revision: 18392
Added:
mc/3D/CitcomCU/trunk/src/io.c
mc/3D/CitcomCU/trunk/src/io.h
mc/3D/CitcomCU/trunk/src/output_ascii.c
mc/3D/CitcomCU/trunk/src/output_ascii.h
mc/3D/CitcomCU/trunk/src/output_vtk.c
mc/3D/CitcomCU/trunk/src/output_vtk.h
Modified:
mc/3D/CitcomCU/trunk/src/Makefile
Log:
Added contributed io output_ascii and output_vtk files to SVN, but
made sure that Makefile doesn't reference them yet (because functionality
isn't fully supported yet).
Modified: mc/3D/CitcomCU/trunk/src/Makefile
===================================================================
--- mc/3D/CitcomCU/trunk/src/Makefile 2011-05-17 22:20:11 UTC (rev 18391)
+++ mc/3D/CitcomCU/trunk/src/Makefile 2011-05-18 19:22:22 UTC (rev 18392)
@@ -146,11 +146,8 @@
Global_operations.c\
Stokes_flow_Incomp.c\
Instructions.c\
- io.c \
Nodal_mesh.c\
Output.c\
- output_ascii.c \
- output_vtk.c \
Pan_problem_misc_functions.c\
Parallel_related.c\
Parsing.c\
Added: mc/3D/CitcomCU/trunk/src/io.c
===================================================================
--- mc/3D/CitcomCU/trunk/src/io.c (rev 0)
+++ mc/3D/CitcomCU/trunk/src/io.c 2011-05-18 19:22:22 UTC (rev 18392)
@@ -0,0 +1,220 @@
+#include <assert.h>
+#include <errno.h>
+#include <math.h>
+#include <mpi.h>
+
+#include <stdbool.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+
+#include <sys/time.h>
+#include <sys/types.h>
+#include <sys/stat.h>
+#include <time.h>
+
+#include <unistd.h>
+
+#include "io.h"
+#include "output_ascii.h"
+#include "output_vtk.h"
+
+
+#define MAXPATHLENGTH 256
+#define MAXFILENAMELENGTH 256 /* TODO: sync it with global_defs.h */
+
+
+static bool io_setup_filename_root( struct All_variables * );
+
+/**
+ * create a directory
+ * @param const char *, directory name
+ * @return true or false
+ */
+bool io_directory_create( const char * dir_name )
+{
+ int i = 0;
+ int length = 0;
+ int status = 0;
+
+ char * tmp = NULL;
+
+ assert( dir_name );
+ length = strlen( dir_name );
+ tmp = (char *) malloc( (length+1) * sizeof( char ) );
+ assert( tmp );
+ memset( tmp, '\0', (length+1) );
+
+
+ for( i = 0; i != length; ++i )
+ {
+ tmp[i] = dir_name[i];
+
+ if( dir_name[i] == '/' )
+ if( !io_directory_exists( tmp ) )
+ status = mkdir( tmp, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH );
+
+ if( status != 0 )
+ {
+ io_error( );
+ return false;
+ }
+ }
+
+ free( tmp );
+ /*
+ * if dirName does not contain trailing / above
+ * loop does not create final directory
+ */
+ if( !io_directory_exists( dir_name ) )
+ {
+ status = mkdir( dir_name, S_IRWXU | S_IRWXG | S_IROTH | S_IXOTH );
+ }
+
+ if( status == 0 )
+ {
+ return true;
+ }
+
+ io_error( );
+
+ return false;
+}
+
+
+/**
+ * check if a directory exists
+ * @param const char *, directory name
+ * @return true or false
+ */
+bool io_directory_exists( const char * dir_name )
+{
+ struct stat dir_info;
+
+ assert( dir_name );
+ if( stat( dir_name, &dir_info ) < 0 )
+ {
+ return false;
+ }
+ else
+ {
+ return ( dir_info.st_mode & S_IFDIR ); /* S_IFDIR requires -D_GNU_SOURCE */
+ }
+}
+
+/**
+ *
+ */
+void io_error( )
+{
+}
+
+/**
+ * inferface for the output of results
+ * @param struct All_variables * E
+ * @return success or failure
+ */
+bool io_results( struct All_variables * E )
+{
+ static bool first_time = true;
+ bool result = false;
+
+ if( first_time )
+ {
+ int proceed = 1;
+ if( E -> parallel.me == 0 )
+ {
+ if( !io_setup_path( E ) )
+ {
+ fprintf( stdout, "failed to setup output directory %s\n", E -> control.output_path );
+ fflush( stdout );
+
+ proceed = 0;
+ }
+ }
+ MPI_Bcast( &proceed, 1, MPI_INT, 0, MPI_COMM_WORLD );
+
+ if( !proceed )
+ {
+ MPI_Finalize( );
+ exit( EXIT_FAILURE );
+ }
+
+ io_setup_filename_root( E );
+ first_time = false;
+ }
+
+ if( E -> control.output_ascii )
+ {
+ result = output_ascii( E );
+ assert( result == true );
+ }
+
+ if( E -> control.output_vtk )
+ {
+ result = output_vtk( E );
+ assert( result == true );
+ }
+
+ return result;
+}
+
+/**
+ * setup the output path for results
+ * @param struct All_variables * E
+ * @retrun success or failure
+ */
+bool io_setup_path( struct All_variables * E )
+{
+ char output_path[ MAXPATHLENGTH ];
+ int length = 0;
+
+ if( E -> control.output_path[0] != '/' )
+ {
+ sprintf( output_path, "%s", "./" );
+ }
+
+ length = strlen( E -> control.output_path );
+ if( length > MAXPATHLENGTH - 3 )
+ {
+ fprintf( stderr, "path length (%d), too long\n", length );
+ exit( EXIT_FAILURE );
+ }
+
+ strncat( output_path, E -> control.output_path, length );
+
+ if( !io_directory_exists( output_path ) )
+ {
+ if( !io_directory_create( output_path ) )
+ {
+ return false;
+ }
+ }
+
+ sprintf( E -> control.output_path, "%s", output_path );
+
+ return true;
+}
+
+static bool io_setup_filename_root( struct All_variables * E )
+{
+ char time_string[13];
+
+ struct timeval tv;
+ time_t curtime;
+
+ /* 15 = "_" + "201011312359\0" + "_" */
+ assert( MAXFILENAMELENGTH > 15 ); /* sanity check, in case it is unintentionally modified */
+ assert( strlen( E -> control.experiment_name ) < MAXFILENAMELENGTH - 15 );
+
+ strcat( E -> control.experiment_name, "_" );
+
+ gettimeofday( &tv, NULL );
+ curtime = tv.tv_sec;
+ strftime( time_string, 13, "%Y%m%d%H%M", localtime( &curtime ) );
+ strcat( E -> control.experiment_name, time_string );
+
+ strcat( E -> control.experiment_name, "_" );
+
+ return true;
+}
Added: mc/3D/CitcomCU/trunk/src/io.h
===================================================================
--- mc/3D/CitcomCU/trunk/src/io.h (rev 0)
+++ mc/3D/CitcomCU/trunk/src/io.h 2011-05-18 19:22:22 UTC (rev 18392)
@@ -0,0 +1,42 @@
+#ifndef IO_HEADER_H
+#define IO_HEADER_H
+
+#include <stdbool.h>
+
+#include "global_defs.h"
+
+/**
+ * create a directory
+ * @param const char *, directory name
+ * @return true or false
+ */
+bool io_directory_create( const char * );
+
+/**
+ * check if a directory exists
+ * @param const char *, directory name
+ * @return true or false
+ */
+bool io_directory_exists( const char * );
+
+/**
+ *
+ */
+void io_error( );
+
+/**
+ * inferface for the output of results
+ * @param struct All_variables * E
+ * @return success or failure
+ */
+bool io_results( struct All_variables * );
+
+/**
+ * setup the output path for results
+ * @param struct All_variables * E
+ * @retrun success or failure
+ */
+bool io_setup_path( struct All_variables * );
+
+
+#endif
Added: mc/3D/CitcomCU/trunk/src/output_ascii.c
===================================================================
--- mc/3D/CitcomCU/trunk/src/output_ascii.c (rev 0)
+++ mc/3D/CitcomCU/trunk/src/output_ascii.c 2011-05-18 19:22:22 UTC (rev 18392)
@@ -0,0 +1,396 @@
+#include <assert.h>
+#include <stdio.h>
+#include <string.h>
+
+#include "output_ascii.h"
+
+#define MAX_FILENAME_LENGTH 256
+
+
+static bool output_ascii_ave( struct All_variables * E );
+static bool output_ascii_coords( struct All_variables * E );
+static bool output_ascii_temp( struct All_variables * E );
+static bool output_ascii_th_b( struct All_variables * E );
+
+static bool output_ascii_th_t( struct All_variables * E );
+static bool output_ascii_traces( struct All_variables * E );
+static bool output_ascii_velocity( struct All_variables * E );
+static bool output_print_coords( FILE * fd, int number_of_coords, float * x0, float * x1, float * x2 );
+
+
+/**
+ * interface for the output of results in ascii format
+ * @param E pointer to struct All_variables
+ * @return success or failure
+ */
+bool output_ascii( struct All_variables * E )
+{
+ bool written = true;
+ static bool first_time = true;
+
+ if( first_time )
+ {
+ first_time = false;
+ written = output_ascii_coords( E );
+ assert( written );
+ }
+
+ if( E -> parallel.me < E -> parallel.nprocz )
+ {
+ written = output_ascii_ave( E );
+ assert( written );
+ }
+
+ if( ( E -> monitor.solution_cycles % E -> control.record_every ) == 0 )
+ {
+ written = output_ascii_temp( E );
+ assert( written );
+
+ written = output_ascii_velocity( E );
+ assert( written );
+
+
+ if( E -> parallel.me_loc[3] == E -> parallel.nprocz - 1 )
+ {
+ written = output_ascii_th_t( E );
+ assert( written );
+ }
+
+ if( E -> parallel.me_loc[3] == 0 )
+ {
+ written = output_ascii_th_b( E );
+ assert( written );
+ }
+ }
+
+ if( E -> control.composition && E -> monitor.solution_cycles % ( 10 * E -> control.record_every ) == 0 )
+ {
+ written = output_ascii_traces( E );
+ assert( E );
+ }
+
+ return written;
+}
+
+
+static bool output_ascii_ave( struct All_variables * E )
+{
+ /* TODO: ave == averages? of what? */
+ bool written = true;
+
+ int i;
+ char output_filename[ MAX_FILENAME_LENGTH ];
+ FILE *fd = NULL;
+
+ memset( output_filename, '\0', MAX_FILENAME_LENGTH );
+ sprintf( output_filename, "%s/%s%05d_%08d.ave",
+ E -> control.output_path,
+ E -> control.experiment_name,
+ E -> parallel.me,
+ E -> monitor.solution_cycles );
+ assert( strlen( output_filename ) < MAX_FILENAME_LENGTH );
+
+ if( ( fd = fopen( output_filename, "w" ) ) != NULL )
+ {
+ fprintf( fd, "%6d %6d %.5e %.5e %.5e %.4e %.4e %.5e %.5e\n",
+ E -> lmesh.noz,
+ E -> advection.timesteps, E -> monitor.elapsed_time,
+ E -> slice.Nut, E -> slice.Nub,
+ E -> data.T_adi0, E -> data.T_adi1,
+ E -> monitor.Sigma_interior, E -> monitor.Sigma_max );
+
+ for( i = 1; i <= E -> lmesh.noz; i++ )
+ {
+ fprintf( fd, "%.4e %.4e %.4e %.4e %.4e %.4e %.4e %.4e\n",
+ E -> Have.T[i],
+ E -> Have.vrms[i],
+ E -> Have.Vi[i],
+ E -> Have.Rho[i],
+ E -> Have.F[i],
+ E -> Have.f[i],
+ E -> Have.C[i],
+ E -> Have.Tadi[i] );
+ }
+ }
+ else
+ {
+ written = false;
+ fprintf( stderr, "failed to open file %s\n", output_filename );
+ fflush( stderr );
+ }
+
+ fclose( fd );
+ return written;
+}
+
+static bool output_ascii_coords( struct All_variables * E )
+{
+ bool written = true;
+
+ char output_filename[ MAX_FILENAME_LENGTH ];
+ FILE *fd = NULL;
+
+ memset( output_filename, 0, MAX_FILENAME_LENGTH );
+ sprintf( output_filename, "%s/%s%05d.coordinates",
+ E -> control.output_path,
+ E -> control.experiment_name,
+ E -> parallel.me );
+ assert( strlen( output_filename ) < MAX_FILENAME_LENGTH );
+
+ if( ( fd = fopen( output_filename, "w" ) ) != NULL )
+ {
+ fprintf( fd, "%6d\n", E -> lmesh.nno );
+ if( E -> control.CART3D )
+ {
+ written = output_print_coords( fd, E -> lmesh.nno, E -> X[1], E -> X[2], E -> X[3] );
+ }
+ else
+ {
+ written = output_print_coords( fd, E -> lmesh.nno, E -> SX[1], E -> SX[2], E -> SX[3] );
+ }
+ }
+ else
+ {
+ written = false;
+ fprintf( stderr, "failed to open file %s\n", output_filename );
+ fflush( stderr );
+ }
+
+ fclose( fd );
+ return written;
+}
+
+static bool output_print_coords( FILE * fd, int number_of_coords, float * x_0, float * x_1, float * x_2 )
+{
+ int i = 1;
+ int result = 0;
+
+ for( i = 1; i <= number_of_coords; ++i )
+ {
+ result = fprintf( fd, "%.5e %.5e %.5e\n", x_0[i], x_1[i], x_2[i] );
+ }
+
+ return ( result > 0 ? true : false );
+}
+
+static bool output_ascii_temp( struct All_variables * E )
+{
+ bool written = true;
+
+ int i;
+ char output_filename[ MAX_FILENAME_LENGTH ];
+ FILE *fd = NULL;
+
+ memset( output_filename, 0, MAX_FILENAME_LENGTH );
+ sprintf( output_filename, "%s/%s%05d_%08d.temperature",
+ E -> control.output_path,
+ E -> control.experiment_name,
+ E -> parallel.me,
+ E -> monitor.solution_cycles );
+ assert( strlen( output_filename ) < MAX_FILENAME_LENGTH );
+
+ if( ( fd = fopen( output_filename, "w" ) ) != NULL )
+ {
+ fprintf( fd, "%6d %6d %.5e\n", E -> lmesh.nno, E -> advection.timesteps, E -> monitor.elapsed_time );
+
+ if( E -> monitor.solution_cycles % ( 20 * E -> control.record_every ) == 0 )
+ {
+ if( E -> control.composition == 0 )
+ {
+ for( i = 1; i <= E -> lmesh.nno; i++ )
+ {
+ /* note in original, the second arg to this fprintf was E -> V[3][i], not sure why */
+ fprintf( fd, "%.5e %.4e %.4e\n", E -> T[i], E -> heatflux[i], E -> heatflux_adv[i] );
+ }
+ }
+ else if( E -> control.composition )
+ {
+ for( i = 1; i <= E -> lmesh.nno; i++ )
+ {
+ /* note in original, the second arg to this fprintf was E -> V[3][i], not sure why */
+ fprintf( fd, "%.5e %.4e %.4e %.4e\n", E -> T[i], E -> heatflux[i], E -> heatflux_adv[i], E -> C[i] );
+ }
+ }
+ }
+ else
+ {
+ if( E -> control.composition == 0 )
+ {
+ for( i = 1; i <= E -> lmesh.nno; i++ )
+ {
+ fprintf( fd, "%.5e %.4e %.4e\n", E -> T[i], E -> V[3][i], E -> heatflux_adv[i] );
+ }
+ }
+ else if( E -> control.composition )
+ {
+ for( i = 1; i <= E -> lmesh.nno; i++ )
+ {
+ fprintf( fd, "%.5e %.4e %.4e %.4e\n", E -> T[i], E -> V[3][i], E -> heatflux_adv[i], E -> C[i] );
+ }
+ }
+ }
+ }
+ else
+ {
+ written = false;
+ fprintf( stderr, "failed to open file %s\n", output_filename );
+ fflush( stderr );
+ }
+
+ fclose( fd );
+ return written;
+}
+
+static bool output_ascii_th_b( struct All_variables * E )
+{
+ /* TODO: find out what th_b is */
+ bool written = true;
+
+ int i;
+ char output_filename[ MAX_FILENAME_LENGTH ];
+ FILE *fd = NULL;
+
+ memset( output_filename, 0, MAX_FILENAME_LENGTH );
+ sprintf( output_filename, "%s/%s%05d_%08d.th_b",
+ E -> control.output_path,
+ E -> control.experiment_name,
+ E -> parallel.me,
+ E -> monitor.solution_cycles );
+ assert( strlen( output_filename ) < MAX_FILENAME_LENGTH );
+
+ if( ( fd = fopen( output_filename, "w" ) ) != NULL )
+ {
+ fprintf( fd, "%6d %6d %.5e %.5e\n", E -> lmesh.nsf, E -> advection.timesteps, E -> monitor.elapsed_time, E -> slice.Nub );
+
+ for( i = 1; i <= E -> lmesh.nsf; i++ )
+ {
+ fprintf( fd, "%.5e %.5e %.5e %.5e\n", E -> slice.tpgb[i], E -> slice.bhflux[i], E -> Fas410_b[i], E -> Fas670_b[i] );
+ }
+ }
+ else
+ {
+ written = false;
+ fprintf( stderr, "failed to open file %s\n", output_filename );
+ fflush( stderr );
+ }
+
+ fclose( fd );
+ return written;
+}
+
+static bool output_ascii_th_t( struct All_variables * E )
+{
+ /* TODO: find out what th_t is */
+ bool written = true;
+
+ int i;
+ char output_filename[ MAX_FILENAME_LENGTH ];
+ FILE *fd = NULL;
+
+ memset( output_filename, 0, MAX_FILENAME_LENGTH );
+ sprintf( output_filename, "%s/%s%05d_%08d.th_t",
+ E -> control.output_path,
+ E -> control.experiment_name,
+ E -> parallel.me,
+ E -> monitor.solution_cycles );
+ assert( strlen( output_filename ) < MAX_FILENAME_LENGTH );
+
+ if( ( fd = fopen( output_filename, "w" ) ) != NULL )
+ {
+ fprintf( fd, "%6d %6d %.5e %.5e\n", E -> lmesh.nsf, E -> advection.timesteps, E -> monitor.elapsed_time, E -> slice.Nut );
+
+ for( i = 1; i <= E -> lmesh.nsf; i++ )
+ {
+ fprintf( fd, "%.5e %.5e %.5e %.5e\n", E -> slice.tpg[i], E -> slice.shflux[i], E -> Fas410_b[i], E -> Fas670_b[i] );
+ }
+ }
+ else
+ {
+ written = false;
+ fprintf( stderr, "failed to open file %s\n", output_filename );
+ fflush( stderr );
+ }
+
+ fclose( fd );
+ return written;
+}
+
+
+static bool output_ascii_traces( struct All_variables * E )
+{
+ bool written = true;
+
+ int i;
+ char output_filename[ MAX_FILENAME_LENGTH ];
+ FILE *fd = NULL;
+
+ memset( output_filename, 0, MAX_FILENAME_LENGTH );
+ sprintf( output_filename, "%s/%s%05d_%08d.traces",
+ E -> control.output_path,
+ E -> control.experiment_name,
+ E -> parallel.me,
+ E -> monitor.solution_cycles );
+ assert( strlen( output_filename ) < MAX_FILENAME_LENGTH );
+
+ if( ( fd = fopen( output_filename, "w" ) ) != NULL )
+ {
+ fprintf( fd, "%6d %6d %.5e\n", E -> advection.markers, E -> advection.timesteps, E -> monitor.elapsed_time );
+
+ for( i = 1; i <= E -> advection.markers; i++ )
+ {
+ fprintf( fd, "%g %g %g %d %d\n", E -> XMC[1][i], E -> XMC[2][i], E -> XMC[3][i], E -> CElement[i], E -> C12[i] );
+ }
+
+ for( i = 1; i <= E -> lmesh.nel; i++ )
+ {
+ fprintf( fd, "%g\n", E -> CE[i] );
+ }
+
+ }
+ else
+ {
+ written = false;
+ fprintf( stderr, "failed to open file %s\n", output_filename );
+ fflush( stderr );
+ }
+
+ fclose( fd );
+ return written;
+}
+
+static bool output_ascii_velocity( struct All_variables * E )
+{
+ bool written = true;
+
+ int i;
+ char output_filename[ MAX_FILENAME_LENGTH ];
+ FILE *fd = NULL;
+
+ memset( output_filename, 0, MAX_FILENAME_LENGTH );
+ sprintf( output_filename, "%s/%s%05d_%08d.velocity",
+ E -> control.output_path,
+ E -> control.experiment_name,
+ E -> parallel.me,
+ E -> monitor.solution_cycles );
+ assert( strlen( output_filename ) < MAX_FILENAME_LENGTH );
+
+ if( ( fd = fopen( output_filename, "w" ) ) != NULL )
+ {
+
+ fprintf( fd, "%6d %6d %.5e\n", E -> lmesh.nno, E -> advection.timesteps, E -> monitor.elapsed_time );
+
+ for( i = 1; i <= E -> lmesh.nno; i++ )
+ {
+ fprintf( fd, "%.6e %.6e %.6e\n", E -> V[1][i], E -> V[2][i], E -> V[3][i] );
+ }
+ }
+ else
+ {
+ written = false;
+ fprintf( stderr, "failed to open file %s\n", output_filename );
+ fflush( stderr );
+ }
+ fclose( fd );
+
+ return written;
+}
Added: mc/3D/CitcomCU/trunk/src/output_ascii.h
===================================================================
--- mc/3D/CitcomCU/trunk/src/output_ascii.h (rev 0)
+++ mc/3D/CitcomCU/trunk/src/output_ascii.h 2011-05-18 19:22:22 UTC (rev 18392)
@@ -0,0 +1,15 @@
+#ifndef OUTPUT_ASCII_H
+#define OUTPUT_ASCII_H
+
+#include <stdbool.h>
+
+#include "global_defs.h"
+
+/**
+ * interface for the output of results in ascii format
+ * @param E pointer to struct All_variables
+ * @return success or failure
+ */
+bool output_ascii( struct All_variables * E );
+
+#endif
Added: mc/3D/CitcomCU/trunk/src/output_vtk.c
===================================================================
--- mc/3D/CitcomCU/trunk/src/output_vtk.c (rev 0)
+++ mc/3D/CitcomCU/trunk/src/output_vtk.c 2011-05-18 19:22:22 UTC (rev 18392)
@@ -0,0 +1,506 @@
+#include <assert.h>
+#include <mpi.h>
+#include <stdbool.h>
+#include <string.h>
+
+#include "output_vtk.h"
+
+
+#define LONG_VTK_LINE 256
+
+
+static void output_vtk_calculate_extents( struct All_variables * E, int ** extent );
+static void output_vtk_print_vts_coords( struct All_variables * E, FILE * data_file, float * w1, float * w2, float * w3 );
+
+static inline void output_vtk_sort_vts_coords_w1( struct All_variables * E, float * restrict x1, float * restrict xcoords );
+static inline void output_vtk_sort_vts_coords_w2( struct All_variables * E, float * restrict x2, float * restrict ycoords );
+static inline void output_vtk_sort_vts_coords_w3( struct All_variables * E, float * restrict x3, float * restrict zcoords );
+
+static bool output_vtk_write_pvts( struct All_variables * E );
+static void output_vtk_write_vts_coords( struct All_variables * E, FILE * data_file, float * w1, float * w2, float * w3 );
+static void output_vtk_write_vts_data( struct All_variables * E, FILE * data_file );
+
+static bool output_vtk_write_vts_output( struct All_variables *E );
+static void output_vtk_write_vts_epilog( FILE *data_file );
+static void output_vtk_write_vts_prolog( struct All_variables *E, FILE *data_file );
+
+__attribute__((unused)) static void output_vtk_write_vts_scalar_d( struct All_variables * E, FILE * data_file, char * name, double * w );
+static void output_vtk_write_vts_scalar_f( struct All_variables * E, FILE * data_file, char * name, float * w );
+__attribute__((unused)) static void output_vtk_write_vts_scalar_i( struct All_variables * E, FILE * data_file, char * name, int * w );
+
+__attribute__((unused)) static void output_vtk_write_vts_vector_d( struct All_variables * E, FILE * data_file, char * name, double * w1, double * w2, double * w3 );
+static void output_vtk_write_vts_vector_f( struct All_variables * E, FILE * data_file, char * name, float * w1, float * w2, float * w3 );
+__attribute__((unused)) static void output_vtk_write_vts_vector_i( struct All_variables * E, FILE * data_file, char * name, int * w1, int * w2, int * w3 );
+
+
+bool output_vtk( struct All_variables * E )
+{
+ bool written = false;
+ written = output_vtk_write_vts_output( E ); assert( written );
+ written = output_vtk_write_pvts( E ); assert( written );
+
+ return written;
+}
+
+
+/* static functions below */
+
+
+static void output_vtk_calculate_extents( struct All_variables * E, int ** extent )
+{
+ int current_proc;
+ int elx, ely, elz;
+ int w1, w2, w3;
+
+ elx = E -> mesh.elx / E -> parallel.nprocx;
+ ely = E -> mesh.ely / E -> parallel.nprocy;
+ elz = E -> mesh.elz / E -> parallel.nprocz;
+
+ current_proc = 0;
+ for( w3 = 0; w3 < E -> parallel.nprocz; ++w3 )
+ {
+ for( w2 = 0; w2 < E -> parallel.nprocy; ++w2 )
+ {
+ for( w1 = 0; w1 < E -> parallel.nprocx; ++w1 )
+ {
+ extent[ current_proc ][0] = elx * ( w1 );
+ extent[ current_proc ][1] = elx * ( w1 + 1 );
+ extent[ current_proc ][2] = ely * ( w2 );
+ extent[ current_proc ][3] = ely * ( w2 + 1 );
+ extent[ current_proc ][4] = elz * ( w3 );
+ extent[ current_proc ][5] = elz * ( w3 + 1 );
+
+ ++current_proc;
+ }
+ }
+ }
+}
+
+
+static void output_vtk_print_vts_coords( struct All_variables *E, FILE *data_file, float * w1, float * w2, float * w3 )
+{
+ int i, j, k;
+
+ fprintf( data_file, "%s\n", "<Points>" );
+ fprintf( data_file, "%s\n", "<DataArray type=\"Float32\" NumberOfComponents=\"3\" format=\"ascii\">" );
+
+ for( k = 0; k != E -> lmesh.noz; ++k )
+ {
+ for( j = 0; j != E -> lmesh.noy; ++j )
+ {
+ for( i = 0; i != E -> lmesh.nox; ++i )
+ {
+ fprintf( data_file, "%.5e %.5e %.5e\n", w1[i], w2[j], w3[k] );
+ }
+ }
+ }
+
+ fprintf( data_file, "%s\n", "</DataArray>" );
+ fprintf( data_file, "%s\n", "</Points>" );
+}
+
+
+static inline void output_vtk_sort_vts_coords_w1( struct All_variables * E, float * restrict w1, float * restrict xcoords )
+{
+ int i, j;
+
+ /* x-coords */
+ for( i = 1, j = 0; i <= E -> lmesh.noz * E -> lmesh.nox; ++i )
+ {
+ if( i % E -> lmesh.noz == 0 )
+ {
+ xcoords[j] = w1[i];
+ ++j;
+ }
+ }
+}
+
+
+static inline void output_vtk_sort_vts_coords_w2( struct All_variables * E, float * restrict w2, float * restrict ycoords )
+{
+ int i, j;
+
+ /* y-coords */
+ for( i = 1, j = 0; i <= E -> lmesh.nno; ++i )
+ {
+ if( i % ( E -> lmesh.noz * E -> lmesh.nox ) == 0 )
+ {
+ ycoords[j] = w2[i];
+ ++j;
+ }
+ }
+}
+
+
+static inline void output_vtk_sort_vts_coords_w3( struct All_variables * E, float * restrict w3, float * restrict zcoords )
+{
+ int i, j;
+
+ /* z-coords */
+ for( i = 1, j = 0; i <= E -> lmesh.noz; ++i, ++j )
+ {
+ zcoords[j] = w3[i];
+ }
+}
+
+
+static bool output_vtk_write_pvts( struct All_variables * E )
+{
+ int i = 0;
+ int ret_i = 0;
+ bool written = false;
+
+ char output_file[ 512 ];
+ FILE * pvts_file = NULL;
+ int ** extent = NULL;
+
+ if( E -> parallel.me == 0 )
+ {
+ extent = ( int ** ) malloc( ( E -> parallel.nproc ) * sizeof( int * ) ); assert( extent != NULL );
+ for( i = 0; i < E -> parallel.nproc; ++i )
+ {
+ extent[i] = NULL;
+ extent[i] = ( int * ) malloc( 6 * sizeof( int ) ); assert( extent[i] != NULL );
+ }
+ output_vtk_calculate_extents( E, extent );
+
+ sprintf( output_file, "%s/%s%08d.pvts",
+ E -> control.output_path,
+ E -> control.experiment_name,
+ E -> monitor.solution_cycles );
+
+ if( ( pvts_file = fopen( output_file, "w" ) ) != NULL )
+ {
+ fprintf( pvts_file, "<?xml version=\"1.0\"?>\n" );
+ fprintf( pvts_file, "<VTKFile type=\"PStructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">\n" );
+ fprintf( pvts_file, "<PStructuredGrid GhostLevel=\"0\" WholeExtent=\"0 %d 0 %d 0 %d\">\n",
+ E -> mesh.nox - 1,
+ E -> mesh.noy - 1,
+ E -> mesh.noz - 1 );
+
+ fprintf( pvts_file, " <PCellData>\n" );
+ fprintf( pvts_file, " </PCellData>\n" );
+
+ fprintf( pvts_file, " <PPoints>\n" );
+ fprintf( pvts_file, " <PDataArray type=\"Float32\" NumberOfComponents=\"3\"/>\n" );
+ fprintf( pvts_file, " </PPoints>\n" );
+
+ fprintf( pvts_file, " <PPointData Scalars=\"Temperature\" Vectors=\"Velocity\">\n" );
+ fprintf( pvts_file, " <PDataArray type=\"Float32\" Name=\"Temperature\" format=\"ascii\" NumberOfComponents=\"1\"/>\n" );
+ fprintf( pvts_file, " <PDataArray type=\"Float32\" Name=\"Heatflux\" format=\"ascii\" NumberOfComponents=\"1\"/>\n" );
+ fprintf( pvts_file, " <PDataArray type=\"Float32\" Name=\"Heatflux_adv\" format=\"ascii\" NumberOfComponents=\"1\"/>\n" );
+ fprintf( pvts_file, " <PDataArray type=\"Float32\" Name=\"Velocity\" format=\"ascii\" NumberOfComponents=\"3\"/>\n" );
+ fprintf( pvts_file, " </PPointData>\n" );
+
+ for( i = 0; i < E -> parallel.nproc; ++i )
+ {
+ fprintf( pvts_file, " <Piece Extent=\" %d %d %d %d %d %d\" Source=\"%s%05d_%08d.vts\"/>\n",
+ extent[i][0], extent[i][1],
+ extent[i][2], extent[i][3],
+ extent[i][4], extent[i][5],
+ E -> control.experiment_name,
+ i,
+ E -> monitor.solution_cycles );
+ }
+ fprintf( pvts_file, " </PStructuredGrid>\n" );
+ fprintf( pvts_file, "</VTKFile>\n" );
+
+ fclose( pvts_file );
+
+ for( i = 0; i < E -> parallel.nproc; ++i )
+ {
+ free( extent[i] );
+ extent[i] = NULL;
+ }
+ free( extent );
+
+ written = true;
+ }
+
+ ret_i = written ? 1 : 0;
+ }
+
+ MPI_Bcast( &ret_i, 1, MPI_INT, 0, MPI_COMM_WORLD );
+ written = ( ret_i == 1 ? true : false );
+
+ return written;
+}
+
+
+static void output_vtk_write_vts_coords( struct All_variables * E, FILE * data_file, float * w1, float * w2, float * w3 )
+{
+ /* v1, v2, v3 will contain the sorted (for vtk/vts) cartesian or spherical coordinates */
+ float * v1 = NULL;
+ float * v2 = NULL;
+ float * v3 = NULL;
+
+ v1 = malloc( sizeof( float ) * ( E -> lmesh.nox ) ); assert( v1 != NULL );
+ v2 = malloc( sizeof( float ) * ( E -> lmesh.noy ) ); assert( v2 != NULL );
+ v3 = malloc( sizeof( float ) * ( E -> lmesh.noz ) ); assert( v3 != NULL );
+
+ output_vtk_sort_vts_coords_w1( E, w1, v1 );
+ output_vtk_sort_vts_coords_w2( E, w2, v2 );
+ output_vtk_sort_vts_coords_w3( E, w3, v3 );
+
+ output_vtk_print_vts_coords( E, data_file, v1, v2, v3 );
+
+ free( v3 );
+ free( v2 );
+ free( v1 );
+}
+
+
+static void output_vtk_write_vts_data( struct All_variables * E, FILE * data_file )
+{
+ fprintf( data_file, "%s\n", "<PointData Scalars = \"Temperature\" Vectors = \"Velocity\">" );
+ output_vtk_write_vts_scalar_f( E, data_file, "Temperature", E -> T );
+ output_vtk_write_vts_scalar_f( E, data_file, "Heatflux", E -> heatflux );
+ output_vtk_write_vts_scalar_f( E, data_file, "Heatflux_adv", E -> heatflux_adv );
+ output_vtk_write_vts_vector_f( E, data_file, "Velocity", E -> V[1], E -> V[2], E -> V[3] );
+ fprintf( data_file, "%s\n", "</PointData>" );
+}
+
+
+static bool output_vtk_write_vts_output( struct All_variables *E )
+{
+ bool written = true;
+ char output_file[ 512 ];
+ FILE *data_file = NULL;
+
+ sprintf( output_file, "%s/%s%05d_%08d.vts",
+ E -> control.output_path,
+ E -> control.experiment_name,
+ E -> parallel.me,
+ E -> monitor.solution_cycles );
+
+ if( ( data_file = fopen( output_file, "w" ) ) != NULL )
+ {
+ output_vtk_write_vts_prolog( E, data_file );
+ if( E -> control.CART3D )
+ {
+ output_vtk_write_vts_coords( E, data_file, E -> X[1], E -> X[2], E -> X[3] );
+ }
+ else
+ {
+ output_vtk_write_vts_coords( E, data_file, E -> SX[1], E -> SX[2], E -> SX[3] );
+ }
+
+ output_vtk_write_vts_data( E, data_file );
+
+ output_vtk_write_vts_epilog( data_file );
+
+ fclose( data_file );
+ }
+ else
+ {
+ written = false;
+ }
+
+ return written;
+}
+
+
+static void output_vtk_write_vts_epilog( FILE *data_file )
+{
+ fprintf( data_file, "%s\n", "</Piece>" );
+ fprintf( data_file, "%s\n", "</StructuredGrid> " );
+ fprintf( data_file, "%s\n", "</VTKFile>" );
+}
+
+
+static void output_vtk_write_vts_prolog( struct All_variables *E, FILE *data_file )
+{
+ fprintf( data_file, "%s\n", "<?xml version=\"1.0\"?>" );
+ fprintf( data_file, "%s\n", "<VTKFile type=\"StructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">" );
+
+ fprintf( data_file, "%s%d %d %d %d %d %d%s\n",
+ "<StructuredGrid WholeExtent=\"",
+ ( E -> mesh.nox - 1 ) * 0,
+ ( E -> mesh.nox - 1 ) * 1,
+ ( E -> mesh.noy - 1 ) * 0,
+ ( E -> mesh.noy - 1 ) * 1,
+ ( E -> mesh.noz - 1 ) * 0,
+ ( E -> mesh.noz - 1 ) * 1,
+ "\">" );
+
+ fprintf( data_file, "%s%d %d %d %d %d %d%s\n",
+ "<Piece Extent=\"",
+ ( E -> lmesh.nox - 1 ) * ( E -> parallel.me_loc[1] ),
+ ( E -> lmesh.nox - 1 ) * ( E -> parallel.me_loc[1] + 1 ),
+ ( E -> lmesh.noy - 1 ) * ( E -> parallel.me_loc[2] ),
+ ( E -> lmesh.noy - 1 ) * ( E -> parallel.me_loc[2] + 1 ),
+ ( E -> lmesh.noz - 1 ) * ( E -> parallel.me_loc[3] ),
+ ( E -> lmesh.noz - 1 ) * ( E -> parallel.me_loc[3] + 1 ),
+ "\">" );
+
+ fprintf( data_file, "%s\n", "<CellData></CellData>" );
+}
+
+
+
+__attribute__((unused)) static void output_vtk_write_vts_scalar_d( struct All_variables * E, FILE * data_file, char * name, double * w )
+{
+ int i, j, k, l;
+ char line[ LONG_VTK_LINE ];
+
+ sprintf( line, "%s", "<DataArray type=\"Float32\" Name=\"" );
+ strcat( line, name );
+ strcat( line, "\" format=\"ascii\" NumberOfComponents=\"1\">" );
+ fprintf( data_file, "%s\n", line );
+
+ /* data is stored in zxy direction, paraview wants it in xyz */
+ for( k = 0; k < E -> lmesh.noz; ++k )
+ {
+ for( j = 0; j < E -> lmesh.noy; ++j )
+ {
+ for( i = 0; i < E -> lmesh.nox; ++i )
+ {
+ l = k + j * E -> lmesh.nox * E -> lmesh.noz + i * E -> lmesh.noz;
+ l = ( l % ( E -> lmesh.nno ) ) + 1;
+ fprintf( data_file, "%lf\n", w[l] );
+ }
+ }
+ }
+
+ fprintf( data_file, "%s\n", "</DataArray>" );
+}
+
+
+static void output_vtk_write_vts_scalar_f( struct All_variables * E, FILE * data_file, char * name, float * w )
+{
+ int i, j, k, l;
+ char line[ LONG_VTK_LINE ];
+
+ sprintf( line, "%s", "<DataArray type=\"Float32\" Name=\"" );
+ strcat( line, name );
+ strcat( line, "\" format=\"ascii\" NumberOfComponents=\"1\">" );
+ fprintf( data_file, "%s\n", line );
+
+ /* data is stored in zxy direction, paraview wants it in xyz */
+ for( k = 0; k < E -> lmesh.noz; ++k )
+ {
+ for( j = 0; j < E -> lmesh.noy; ++j )
+ {
+ for( i = 0; i < E -> lmesh.nox; ++i )
+ {
+ l = k + j * E -> lmesh.nox * E -> lmesh.noz + i * E -> lmesh.noz;
+ l = ( l % ( E -> lmesh.nno ) ) + 1;
+ fprintf( data_file, "%f\n", w[l] );
+ }
+ }
+ }
+
+ fprintf( data_file, "%s\n", "</DataArray>" );
+}
+
+
+__attribute__((unused)) static void output_vtk_write_vts_scalar_i( struct All_variables * E, FILE * data_file, char * name, int * w )
+{
+ int i, j, k, l;
+ char line[ LONG_VTK_LINE ];
+
+ sprintf( line, "%s", "<DataArray type=\"Float32\" Name=\"" );
+ strcat( line, name );
+ strcat( line, "\" format=\"ascii\" NumberOfComponents=\"1\">" );
+ fprintf( data_file, "%s\n", line );
+
+ /* data is stored in zxy direction, paraview wants it in xyz */
+ for( k = 0; k < E -> lmesh.noz; ++k )
+ {
+ for( j = 0; j < E -> lmesh.noy; ++j )
+ {
+ for( i = 0; i < E -> lmesh.nox; ++i )
+ {
+ l = k + j * E -> lmesh.nox * E -> lmesh.noz + i * E -> lmesh.noz;
+ l = ( l % ( E -> lmesh.nno ) ) + 1;
+ fprintf( data_file, "%d\n", w[l] );
+ }
+ }
+ }
+
+ fprintf( data_file, "%s\n", "</DataArray>" );
+}
+
+
+__attribute__((unused)) static void output_vtk_write_vts_vector_d( struct All_variables * E, FILE * data_file, char * name, double * w1, double * w2, double * w3 )
+{
+ int i, j, k, l;
+ char line[ LONG_VTK_LINE ];
+
+ sprintf( line, "%s", "<DataArray type=\"Float32\" Name=\"" );
+ strcat( line, name );
+ strcat( line, "\" format=\"ascii\" NumberOfComponents=\"3\">" );
+ fprintf( data_file, "%s\n", line );
+
+ /* data is stored in zxy direction, paraview wants it in xyz */
+ for( k = 0; k < E -> lmesh.noz; ++k )
+ {
+ for( j = 0; j < E -> lmesh.noy; ++j )
+ {
+ for( i = 0; i < E -> lmesh.nox; ++i )
+ {
+ l = k + j * E -> lmesh.nox * E -> lmesh.noz + i * E -> lmesh.noz;
+ l = ( l % ( E -> lmesh.nno ) ) + 1;
+ fprintf( data_file, "%.6e %.6e %.6e\n", w1[l], w2[l], w3[l] );
+ }
+ }
+ }
+
+ fprintf( data_file, "%s\n", "</DataArray>" );
+}
+
+
+__attribute__((unused)) static void output_vtk_write_vts_vector_f( struct All_variables * E, FILE * data_file, char * name, float * w1, float * w2, float * w3 )
+{
+ int i, j, k, l;
+ char line[ LONG_VTK_LINE ];
+
+ sprintf( line, "%s", "<DataArray type=\"Float32\" Name=\"" );
+ strcat( line, name );
+ strcat( line, "\" format=\"ascii\" NumberOfComponents=\"3\">" );
+ fprintf( data_file, "%s\n", line );
+
+ /* data is stored in zxy direction, paraview wants it in xyz */
+ for( k = 0; k < E -> lmesh.noz; ++k )
+ {
+ for( j = 0; j < E -> lmesh.noy; ++j )
+ {
+ for( i = 0; i < E -> lmesh.nox; ++i )
+ {
+ l = k + j * E -> lmesh.nox * E -> lmesh.noz + i * E -> lmesh.noz;
+ l = ( l % ( E -> lmesh.nno ) ) + 1;
+ fprintf( data_file, "%.6e %.6e %.6e\n", w1[l], w2[l], w3[l] );
+ }
+ }
+ }
+
+ fprintf( data_file, "%s\n", "</DataArray>" );
+}
+
+
+__attribute__((unused)) static void output_vtk_write_vts_vector_i( struct All_variables * E, FILE * data_file, char * name, int * w1, int * w2, int * w3 )
+{
+ int i, j, k, l;
+ char line[ LONG_VTK_LINE ];
+
+ sprintf( line, "%s", "<DataArray type=\"Float32\" Name=\"" );
+ strcat( line, name );
+ strcat( line, "\" format=\"ascii\" NumberOfComponents=\"3\">" );
+ fprintf( data_file, "%s\n", line );
+
+ /* data is stored in zxy direction, paraview wants it in xyz */
+ for( k = 0; k < E -> lmesh.noz; ++k )
+ {
+ for( j = 0; j < E -> lmesh.noy; ++j )
+ {
+ for( i = 0; i < E -> lmesh.nox; ++i )
+ {
+ l = k + j * E -> lmesh.nox * E -> lmesh.noz + i * E -> lmesh.noz;
+ l = ( l % ( E -> lmesh.nno ) ) + 1;
+ fprintf( data_file, "%d %d %d\n", w1[l], w2[l], w3[l] );
+ }
+ }
+ }
+
+ fprintf( data_file, "%s\n", "</DataArray>" );
+}
Added: mc/3D/CitcomCU/trunk/src/output_vtk.h
===================================================================
--- mc/3D/CitcomCU/trunk/src/output_vtk.h (rev 0)
+++ mc/3D/CitcomCU/trunk/src/output_vtk.h 2011-05-18 19:22:22 UTC (rev 18392)
@@ -0,0 +1,10 @@
+#ifndef OUTPUT_VTR_H
+#define OUTPUT_VTR_H
+
+#include <stdbool.h>
+
+#include "global_defs.h"
+
+bool output_vtk( struct All_variables * );
+
+#endif
More information about the CIG-COMMITS
mailing list