[cig-commits] r7846 - mc/3D/CitcomS/trunk/lib
becker at geodynamics.org
becker at geodynamics.org
Sun Aug 19 16:54:47 PDT 2007
Author: becker
Date: 2007-08-19 16:54:46 -0700 (Sun, 19 Aug 2007)
New Revision: 7846
Modified:
mc/3D/CitcomS/trunk/lib/Instructions.c
mc/3D/CitcomS/trunk/lib/Output.c
mc/3D/CitcomS/trunk/lib/Output_gzdir.c
mc/3D/CitcomS/trunk/lib/Output_h5.c
mc/3D/CitcomS/trunk/lib/Parsing.c
Log:
- merged changes with Eh's additions
- added coor=2 option (for full citcom) to assign finer radial node
spacing to top and lower layers of shell. The
coor_refine=0.1,0.15,0.1,0.2 parameter specifies the radius fraction
of the bottom layer [0], the fraction of the nodes in this layer
[1], the top layer fraction [2], and the top layer node fraction
[3]. I.e. the defaults will put 15% of all nz nodes into the 10%
lower layer, 20% in the top 10% upper layer, and the rest in
between.
- renamed gzipped output option ascii-gz
- built in restart facilities for temperature and tracers when using
gzdir I/O
- added a composition viscosity function, CDEPV, based on two tracer
flavors
- for this to work, I had to move viscosity_input() *behind*
tic_input() and tracer_input() in instructions
- added tracer_enriched option. If tracer = on and tracer_enriched =
on, will reader Q0_enriched and vary the element heat production
between Q0 for C = 0 and Q0_enriched for C = 1. I.e. this only works
if C varies between 0 and 1.
- added an experimental option to write to a single VTK file (note
caveats!) if ascii-gz is activated
gzdir_vtkio = 2, will try to write VTK straight (experimental)
the VTK output is of the "legacy", serial type, and
requires that all processors see the same
filesystem. This will likely lead to a bottleneck
for larg CPU computations as each processor has to
wait til the previous is done.
Modified: mc/3D/CitcomS/trunk/lib/Instructions.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Instructions.c 2007-08-18 03:14:06 UTC (rev 7845)
+++ mc/3D/CitcomS/trunk/lib/Instructions.c 2007-08-19 23:54:46 UTC (rev 7846)
@@ -74,9 +74,8 @@
void set_sphere_harmonics (struct All_variables*);
void set_starting_age(struct All_variables*);
void tracer_initial_settings(struct All_variables*);
+void get_vtk_filename(char *,int,struct All_variables *,int);
-
-
void initial_mesh_solver_setup(struct All_variables *E)
{
@@ -1285,7 +1284,7 @@
void output_finalize(struct All_variables *E)
{
- char message[255];
+ char message[255],files[255];
if (E->fp)
fclose(E->fp);
@@ -1296,18 +1295,20 @@
fclose(E->fp_out);
/*
- remove VTK geo file in case we used that for IO (we'll only do
- this for one processor, since this IO option requires shared
- filesystems anyway
+ remove VTK geo file in case we used that for IO
*/
- if((E->parallel.me == 0) && (E->output.gzdir.vtk_io == 2)&&
+ if((E->output.gzdir.vtk_io != 1) &&
(strcmp(E->output.format, "ascii-gz") == 0)){
- /* delete the vtk geo pre-file */
- snprintf(message,255,"rm -f %s/vtk_geo",
- E->control.data_dir);
- system(message);
- /* close the log */
- fclose(E->output.gzdir.vtk_fp);
+ if((E->output.gzdir.vtk_io == 3)||(E->parallel.me == 0)){
+ /* delete the geo files */
+ get_vtk_filename(files,1,E,0);
+ sprintf(message,"rm %s",files);system(message);
+ fprintf(stderr,"%s",message);
+ if(E->parallel.me == 0){
+ /* close the log */
+ fclose(E->output.gzdir.vtk_fp);
+ }
+ }
}
}
Modified: mc/3D/CitcomS/trunk/lib/Output.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Output.c 2007-08-18 03:14:06 UTC (rev 7845)
+++ mc/3D/CitcomS/trunk/lib/Output.c 2007-08-19 23:54:46 UTC (rev 7846)
@@ -70,7 +70,8 @@
if(strcmp(E->output.format, "ascii-gz") == 0){
/*
vtk_io = 1: write files for post-processing into VTK
- vtk_io = 2: write legacy VTK file straight
+ vtk_io = 2: write serial legacy VTK file
+ vtk_io = 3: write paralle legacy VTK file
*/
input_int("gzdir_vtkio",&(E->output.gzdir.vtk_io),"0",m);
@@ -510,7 +511,8 @@
for(i=1;i<=E->lmesh.nno;i++) {
for(k=0;k<E->composition.ncomp;k++) {
- fprintf(fp1,"%.6e ",E->composition.comp_el[j][k][i]);
+ fprintf(fp1,"%.6e ",
+ E->composition.comp_el[j][k][i]);
}
fprintf(fp1,"\n");
}
Modified: mc/3D/CitcomS/trunk/lib/Output_gzdir.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Output_gzdir.c 2007-08-18 03:14:06 UTC (rev 7845)
+++ mc/3D/CitcomS/trunk/lib/Output_gzdir.c 2007-08-19 23:54:46 UTC (rev 7846)
@@ -28,19 +28,23 @@
/* Routine to process the output of the finite element cycles
and to turn them into a coherent suite files */
-/*
+/*
this version uses gzipped, ascii output to subdirectories for the
ascii-gz option
if, additionally, gzdir.vtk_io = 1, will write different format files
- for post-processing into VTK
+ for later post-processing into VTK
- gzdir.vtk_io = 2, will try to write VTK straight (experimental)
+ gzdir.vtk_io = 2, will try to write legacy serial VTK (experimental)
+ gzdir.vtk_io = 3, will try to write to legacy parallel VTK (experimental)
+
+
+
the VTK output is the "legacy" type, requires that
all processors see the same filesystem, and will
- likely lead to a bottleneck for larg CPU
+ likely lead to a bottleneck for large CPU
computations as each processor has to wait til the
previous is done.
@@ -75,8 +79,8 @@
/* */
FILE* output_open_mode(char *, char *);
-void get_vtk_filename(char *,int,struct All_variables *);
+void get_vtk_filename(char *,int,struct All_variables *,int);
gzFile *gzdir_output_open(char *,char *);
void gzdir_output(struct All_variables *, int );
@@ -119,20 +123,21 @@
char output_dir[255];
if (cycles == 0) {
/* initial I/O */
-
+
gzdir_output_coord(E);
/*gzdir_output_mat(E);*/
}
- /*
- make a new directory for all the other output
-
+ /*
+ make a new directory for all the other output
+
(all procs need to do that, because we might be using a local tmp
dir)
-
+
*/
/* make a directory */
snprintf(output_dir,255,"%s/%d",E->control.data_dir,cycles);
+
mkdatadir(output_dir);
@@ -160,7 +165,7 @@
gzdir_output_horiz_avg(E, cycles);
if(E->control.tracer){
- if(E->output.tracer ||
+ if(E->output.tracer ||
(cycles == E->advection.max_timesteps))
gzdir_output_tracer(E, cycles);
}
@@ -192,7 +197,12 @@
return fp1;
}
+/*
+initialization output of geometries, only called once
+
+
+ */
void gzdir_output_coord(struct All_variables *E)
{
int i, j, offset,ix[9],out;
@@ -200,31 +210,28 @@
float x[3];
gzFile *gz1;
FILE *fp1;
-
- /* for dealing with several processors */
MPI_Status mpi_stat;
int mpi_rc, mpi_inmsg, mpi_success_message = 1;
- if(E->output.gzdir.vtk_io == 2){
+ if((E->output.gzdir.vtk_io == 2)||(E->output.gzdir.vtk_io == 3)){
+ /*
+ direct VTK file output
+ */
+ if(E->output.gzdir.vtk_io == 2) /* serial */
+ parallel_process_sync(E);
+ /*
+
+ start geometry pre-file, to which data will get appended later
- parallel_process_sync(E);
- /*
- new VTK output, write directly to a VTK file
- this will get deleted from
-
- output_finalize(struct All_variables *E)
-
-
*/
E->output.gzdir.vtk_ocount = -1;
- snprintf(output_file,255,"%s/vtk_geo",E->control.data_dir);
+ get_vtk_filename(output_file,1,E,0); /* geometry file */
if(E->parallel.me == 0){
-
- /* log times */
- snprintf(message,255,
- "%s/vtk_time.log",E->control.data_dir);
+ /* start log file */
+ snprintf(message,255,"%s/vtk_time.log",E->control.data_dir);
E->output.gzdir.vtk_fp = output_open_mode(message,"w");
-
-
+ }
+ if((E->parallel.me == 0) || (E->output.gzdir.vtk_io == 3)){
+ /* either in first CPU or parallel output */
/* start geo file */
fp1 = output_open_mode(output_file,"w");
myfprintf(fp1,"# vtk DataFile Version 2.0\n");
@@ -235,12 +242,16 @@
myfprintf(fp1,"BINARY\n");
#endif
myfprintf(fp1,"DATASET UNSTRUCTURED_GRID\n");
- sprintf(message,"POINTS %i float\n", /* total number of nodes */
- E->lmesh.nno * E->parallel.nproc *
- E->sphere.caps_per_proc);
+ if(E->output.gzdir.vtk_io == 2) /* serial */
+ sprintf(message,"POINTS %i float\n", /* total number of nodes */
+ E->lmesh.nno * E->parallel.nproc *
+ E->sphere.caps_per_proc);
+ else /* parallel */
+ sprintf(message,"POINTS %i float\n",
+ E->lmesh.nno * E->sphere.caps_per_proc);
myfprintf(fp1,message);
- }else{
- /* if not first, wait for previous */
+ }else{ /* serial output */
+ /* if not first CPU, wait for previous before appending */
mpi_rc = MPI_Recv(&mpi_inmsg, 1, MPI_INT, (E->parallel.me-1), 0, MPI_COMM_WORLD, &mpi_stat);
/* open for append */
fp1 = output_open_mode(output_file,"a");
@@ -250,26 +261,34 @@
for(j=1;j <= E->sphere.caps_per_proc;j++) {
for(i=1;i <= E->lmesh.nno;i++) {
x[0]=E->x[j][1][i];x[1]=E->x[j][2][i];x[2]=E->x[j][3][i];
- if(be_write_float_to_file(x,3,fp1) != 3)BE_WERROR;
+ if(be_write_float_to_file(x,3,fp1) != 3)
+ BE_WERROR;
out++;
}
}
- fclose(fp1);fflush(fp1); /* close file and flush buffer */
- if(E->parallel.me < E->parallel.nproc-1){/* send to next if not last*/
- mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 0, MPI_COMM_WORLD);
+ if(E->output.gzdir.vtk_io == 2){ /* serial output, close and have
+ next one write */
+ fclose(fp1);fflush(fp1); /* close file and flush buffer */
+ if(E->parallel.me < E->parallel.nproc-1){/* send to next if not last*/
+ mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 0, MPI_COMM_WORLD);
+ }
+ /*
+ node numbers for all the elements
+ */
+ parallel_process_sync(E);
}
- /*
- node numbers for all the elements
- */
- parallel_process_sync(E);
- if(E->parallel.me == 0){
- fp1 = output_open_mode(output_file,"a");
- j = E->parallel.nproc * E->lmesh.nel *
- E->sphere.caps_per_proc; /* total number of elements */
+ if((E->output.gzdir.vtk_io == 3) || (E->parallel.me == 0)){ /* in first CPU, or parallel output */
+ if(E->output.gzdir.vtk_io == 2){ /* need to reopen, serial */
+ fp1 = output_open_mode(output_file,"a");
+ j = E->parallel.nproc * E->lmesh.nel *
+ E->sphere.caps_per_proc; /* total number of elements */
+ }else{ /* parallel */
+ j = E->lmesh.nel * E->sphere.caps_per_proc;
+ }
sprintf(message,"CELLS %i %i\n", /* number of elements
total number of int entries
-
- */
+
+ */
j,j*(enodes[E->mesh.nsd]+1));
myfprintf(fp1,message);
}else{
@@ -277,16 +296,20 @@
mpi_rc = MPI_Recv(&mpi_inmsg, 1, MPI_INT, (E->parallel.me-1), 0, MPI_COMM_WORLD, &mpi_stat);
fp1 = output_open_mode(output_file,"a");
}
- /*
- write CELL element nodes
+ /*
+ write CELL element nodes
*/
if(enodes[E->mesh.nsd] != 8)
myerror(E,"vtk error, only eight node hexes supported");
- offset = E->lmesh.nno * E->parallel.me - 1;
+ if(E->output.gzdir.vtk_io == 2){ /* serial, global node numbers */
+ offset = E->lmesh.nno * E->parallel.me - 1;
+ }else{ /* parallel, only use local node numbers? */
+ offset = -1;
+ }
ix[0] = enodes[E->mesh.nsd];
for(j=1;j <= E->sphere.caps_per_proc;j++) {
for(i=1;i <= E->lmesh.nel;i++) {
- /*
+ /*
need to add offset according to the processor for global
node numbers
*/
@@ -294,42 +317,46 @@
ix[3]= E->ien[j][i].node[3]+offset;ix[4] = E->ien[j][i].node[4]+offset;
ix[5]= E->ien[j][i].node[5]+offset;ix[6] = E->ien[j][i].node[6]+offset;
ix[7]= E->ien[j][i].node[7]+offset;ix[8] = E->ien[j][i].node[8]+offset;
- if(be_write_int_to_file(ix,9,fp1)!=9)BE_WERROR;
+ if(be_write_int_to_file(ix,9,fp1)!=9)
+ BE_WERROR;
}
}
- fclose(fp1);fflush(fp1); /* close file and flush buffer */
- if(E->parallel.me < E->parallel.nproc-1)
- mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 0, MPI_COMM_WORLD);
- /*
- add the cell type flags
- */
- parallel_process_sync(E);
-
- if(E->parallel.me == 0){
- fp1 = output_open_mode(output_file,"a");
- j=E->parallel.nproc*E->lmesh.nel*E->sphere.caps_per_proc;
+ if(E->output.gzdir.vtk_io == 2){ /* serial IO */
+ fclose(fp1);fflush(fp1); /* close file and flush buffer */
+ if(E->parallel.me < E->parallel.nproc-1)
+ mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 0, MPI_COMM_WORLD);
+ parallel_process_sync(E);
+ }
+ if((E->output.gzdir.vtk_io==3) || (E->parallel.me == 0) ){
+ if(E->output.gzdir.vtk_io == 2){ /* serial */
+ fp1 = output_open_mode(output_file,"a");
+ j=E->parallel.nproc*E->lmesh.nel*E->sphere.caps_per_proc;
+ }else{ /* parallel */
+ j = E->lmesh.nel*E->sphere.caps_per_proc;
+ }
sprintf(message,"CELL_TYPES %i\n",j); /* number of elements*/
myfprintf(fp1,message);
ix[0] = 12;
for(i=0;i<j;i++)
if(be_write_int_to_file(ix,1,fp1)!=1)BE_WERROR;
- fclose(fp1);fflush(fp1); /* close file and flush buffer */
- report(E,"vtk_io: vtk geometry done\n");
+ fclose(fp1);fflush(fp1); /* all procs close file and flush buffer */
+ if(E->parallel.me == 0)
+ fprintf(stderr,"vtk_io: vtk geometry done for %s\n",output_file);
}
- /* done straight VTK */
+ /* done straight VTK output, geometry part */
}else{
- /*
-
+ /*
+
either zipped regular, or old VTK type for post-processing
-
+
*/
- /*
+ /*
don't use data file name
*/
snprintf(output_file,255,"%s/coord.%d.gz",
E->control.data_dir,E->parallel.me);
gz1 = gzdir_output_open(output_file,"w");
-
+
/* nodal coordinates */
for(j=1;j<=E->sphere.caps_per_proc;j++) {
gzprintf(gz1,"%3d %7d\n",j,E->lmesh.nno);
@@ -337,16 +364,16 @@
gzprintf(gz1,"%.6e %.6e %.6e\n",
E->sx[j][1][i],E->sx[j][2][i],E->sx[j][3][i]);
}
-
+
gzclose(gz1);
- if(E->output.gzdir.vtk_io){
- /*
-
+ if(E->output.gzdir.vtk_io == 1){
+ /*
+
output of Cartesian coordinates and element connectivitiy for
vtk visualization
-
+
*/
- /*
+ /*
nodal coordinates in Cartesian
*/
snprintf(output_file,255,"%s/vtk_ecor.%d.gz",
@@ -359,7 +386,7 @@
}
}
gzclose(gz1);
- /*
+ /*
connectivity for all elements
*/
offset = E->lmesh.nno * E->parallel.me - 1;
@@ -373,7 +400,7 @@
gzprintf(stderr,"gzdir: Output: error, only eight node hexes supported");
parallel_process_termination();
}
- /*
+ /*
need to add offset according to the processor for global
node numbers
*/
@@ -385,23 +412,24 @@
}
}
gzclose(gz1);
- } /* end vtkio = 1 */
+ } /* end vtkio = 1 (pre VTK) */
}
return;
}
-/*
+/*
-this needs to be called first before any of the other stuff if VTK
-straight output is chosen
+this needs to be called after the geometry files have been
+established, and before any of the other stuff if VTK straight output
+is chosen
*/
void gzdir_output_velo_temp(struct All_variables *E, int cycles)
{
int i, j, k,os;
- char output_file[255],output_file2[255],message[255];
+ char output_file[255],output_file2[255],message[255],geo_file[255];
float cvec[3];
gzFile *gzout;
FILE *fp1;
@@ -414,7 +442,7 @@
if(E->output.gzdir.vtk_io){ /* all VTK modes need basis vectors */
os = E->lmesh.nno*9;
-
+
if((!E->output.gzdir.vtk_base_init) ||(!E->output.gzdir.vtk_base_save)){
if(!E->output.gzdir.vtk_base_init)
E->output.gzdir.vtk_base = (float *)safe_malloc(sizeof(float)*os*E->sphere.caps_per_proc);
@@ -427,36 +455,38 @@
E->output.gzdir.vtk_base_init = 1;
}
}
+ if((E->output.gzdir.vtk_io == 2) || (E->output.gzdir.vtk_io == 3)){
+ /*
- if(E->output.gzdir.vtk_io == 2){
- parallel_process_sync(E);
+ direct VTK
- /*
-
- straight VTK output, always need to start here!
-
*/
- E->output.gzdir.vtk_ocount++;
- get_vtk_filename(output_file,E->output.gzdir.vtk_ocount,E);
- /*
+ if(E->output.gzdir.vtk_io == 2)
+ parallel_process_sync(E); /* serial needs sync */
+ E->output.gzdir.vtk_ocount++; /* regular output file name */
+ get_vtk_filename(geo_file,1,E,cycles);
+ get_vtk_filename(output_file,0,E,cycles);
+ /*
+
start with temperature
*/
- if(E->parallel.me == 0){
+ if((E->parallel.me == 0) || (E->output.gzdir.vtk_io == 3)){
/* copy geo file over to start out vtk file */
- snprintf(output_file2,255,"cp %s/vtk_geo %s",
- E->control.data_dir,output_file);
+ snprintf(output_file2,255,"cp %s %s",geo_file,output_file);
system(output_file2);
/* should we do something to check if this has worked? */
-
- /* write a time log */
- fprintf(E->output.gzdir.vtk_fp,"%12i %12i %12.6e %s\n",
- E->output.gzdir.vtk_ocount,cycles,E->monitor.elapsed_time,output_file);
-
- /* */
+ if(E->parallel.me == 0){
+ /* write a time log */
+ fprintf(E->output.gzdir.vtk_fp,"%12i %12i %12.6e %s\n",
+ E->output.gzdir.vtk_ocount,cycles,E->monitor.elapsed_time,output_file);
+ }
fp1 = output_open_mode(output_file,"a");
- sprintf(message,"POINT_DATA %i\n",E->lmesh.nno*E->parallel.nproc*E->sphere.caps_per_proc);
+ if(E->output.gzdir.vtk_io == 2) /* serial */
+ sprintf(message,"POINT_DATA %i\n",E->lmesh.nno*E->parallel.nproc*E->sphere.caps_per_proc);
+ else /* parallel */
+ sprintf(message,"POINT_DATA %i\n",E->lmesh.nno*E->sphere.caps_per_proc);
myfprintf(fp1,message);
myfprintf(fp1,"SCALARS temperature float 1\n");
myfprintf(fp1,"LOOKUP_TABLE default\n");
@@ -468,21 +498,26 @@
}
for(j=1; j<= E->sphere.caps_per_proc;j++) /* print the temperatures */
for(i=1;i<=E->lmesh.nno;i++){
- cvec[0] = E->T[j][i];
- if(be_write_float_to_file(cvec,1,fp1)!=1)BE_WERROR;
+ cvec[0] = E->T[j][i];
+ if(be_write_float_to_file(cvec,1,fp1)!=1)
+ BE_WERROR;
}
- fclose(fp1);fflush(fp1); /* close file and flush buffer */
- if(E->parallel.me < E->parallel.nproc-1){
- mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 7, MPI_COMM_WORLD);
- }else{
- mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, 0, 6, MPI_COMM_WORLD); /* tell m=0 to go ahead */
+ if(E->output.gzdir.vtk_io == 2){
+ fclose(fp1);fflush(fp1); /* close file and flush buffer */
+ if(E->parallel.me < E->parallel.nproc-1){
+ mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 7, MPI_COMM_WORLD);
+ }else{
+ mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, 0, 6, MPI_COMM_WORLD); /* tell m=0 to go ahead */
+ }
}
- /*
+ /*
velocities second
*/
- if(E->parallel.me == 0){
- mpi_rc = MPI_Recv(&mpi_inmsg, 1, MPI_INT, E->parallel.nproc-1 , 6, MPI_COMM_WORLD, &mpi_stat);
- fp1 = output_open_mode(output_file,"a"); /* append velocities */
+ if((E->output.gzdir.vtk_io == 3) || (E->parallel.me == 0)){
+ if(E->output.gzdir.vtk_io == 2){
+ mpi_rc = MPI_Recv(&mpi_inmsg, 1, MPI_INT, E->parallel.nproc-1 , 6, MPI_COMM_WORLD, &mpi_stat);
+ fp1 = output_open_mode(output_file,"a"); /* append velocities */
+ }
sprintf(message,"VECTORS velocity float\n");myfprintf(fp1,message);
}else{
mpi_rc = MPI_Recv(&mpi_inmsg, 1, MPI_INT, (E->parallel.me-1), 5, MPI_COMM_WORLD, &mpi_stat);
@@ -496,30 +531,35 @@
}
}
fclose(fp1);fflush(fp1); /* close file and flush buffer */
- if(E->parallel.me < E->parallel.nproc-1){
- mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 5, MPI_COMM_WORLD);
+ if(E->output.gzdir.vtk_io == 2){
+ if(E->parallel.me < E->parallel.nproc-1){
+ mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 5, MPI_COMM_WORLD);
+ }else{
+ fprintf(stderr,"vtk_io: geo, temp, & vel writtend to %s\n",output_file);
+ }
}else{
- fprintf(stderr,"vtk_io: geo, temp, & vel writtend to %s\n",output_file);
+ if(E->parallel.me == 0)
+ fprintf(stderr,"vtk_io: geo, temp, & vel written to %s\n",output_file);
}
/* new VTK velo and temp done */
}else{
- /*
+ /*
modified zipped output
*/
/*
-
-
+
+
temperatures are printed along with velocities for old type of
output
-
+
if VTK is selected, will generate a separate temperature file
-
+
*/
- if(E->output.gzdir.vtk_io) {
- /*
- for VTK, only print temperature
+ if(E->output.gzdir.vtk_io == 1) {
+ /*
+ for VTK, only print temperature
*/
snprintf(output_file2,255,"%s/%d/t.%d.%d",
E->control.data_dir,
@@ -530,7 +570,7 @@
E->parallel.me,cycles);
}
snprintf(output_file,255,"%s.gz",output_file2); /* add the .gz */
-
+
gzout = gzdir_output_open(output_file,"w");
gzprintf(gzout,"%d %d %.5e\n",
cycles,E->lmesh.nno,E->monitor.elapsed_time);
@@ -538,21 +578,21 @@
gzprintf(gzout,"%3d %7d\n",j,E->lmesh.nno);
if(E->output.gzdir.vtk_io){
/* VTK */
- for(i=1;i<=E->lmesh.nno;i++)
- gzprintf(gzout,"%.6e\n",E->T[j][i]);
- } else {
+ for(i=1;i<=E->lmesh.nno;i++)
+ gzprintf(gzout,"%.6e\n",E->T[j][i]);
+ } else {
/* old */
- for(i=1;i<=E->lmesh.nno;i++)
+ for(i=1;i<=E->lmesh.nno;i++)
gzprintf(gzout,"%.6e %.6e %.6e %.6e\n",
E->sphere.cap[j].V[1][i],
E->sphere.cap[j].V[2][i],
- E->sphere.cap[j].V[3][i],E->T[j][i]);
+ E->sphere.cap[j].V[3][i],E->T[j][i]);
}
}
gzclose(gzout);
if(E->output.gzdir.vtk_io){
- /*
- write Cartesian velocities to file
+ /*
+ write Cartesian velocities to file
*/
snprintf(output_file,255,"%s/%d/vtk_v.%d.%d.gz",
E->control.data_dir,cycles,E->parallel.me,cycles);
@@ -570,7 +610,7 @@
}
}
gzclose(gzout);
-
+
}
} /* end gzipped and old VTK out */
if(E->output.gzdir.vtk_io){ /* all VTK modes */
@@ -581,8 +621,8 @@
return;
}
-/*
- viscosity
+/*
+ viscosity
*/
void gzdir_output_visc(struct All_variables *E, int cycles)
{
@@ -597,7 +637,8 @@
int mpi_rc;
int mpi_inmsg, mpi_success_message = 1;
- if(E->output.gzdir.vtk_io != 2){ /* old output */
+
+ if(E->output.gzdir.vtk_io < 2){
snprintf(output_file,255,
"%s/%d/visc.%d.%d.gz", E->control.data_dir,
cycles,E->parallel.me, cycles);
@@ -607,14 +648,14 @@
for(i=1;i<=E->lmesh.nno;i++)
gzprintf(gz1,"%.4e\n",E->VI[lev][j][i]);
}
-
+
gzclose(gz1);
}else{
- parallel_process_sync(E);
-
- /* new legacy VTK */
- get_vtk_filename(output_file,E->output.gzdir.vtk_ocount,E);
- if(E->parallel.me == 0){
+ if(E->output.gzdir.vtk_io == 2)
+ parallel_process_sync(E);
+ /* new legacy VTK */
+ get_vtk_filename(output_file,0,E,cycles);
+ if((E->parallel.me == 0) || (E->output.gzdir.vtk_io == 3)){
fp1 = output_open_mode(output_file,"a");
myfprintf(fp1,"SCALARS log10(visc) float 1\n");
myfprintf(fp1,"LOOKUP_TABLE default\n");
@@ -624,16 +665,17 @@
/* open for append */
fp1 = output_open_mode(output_file,"a");
}
- for(j=1; j<= E->sphere.caps_per_proc;j++)
+ for(j=1; j<= E->sphere.caps_per_proc;j++)
for(i=1;i<=E->lmesh.nno;i++){
- ftmp = log10(E->VI[lev][j][i]);
+ ftmp = log10(E->VI[lev][j][i]);
if(fabs(ftmp) < 5e-7)ftmp = 0.0;
if(be_write_float_to_file(&ftmp,1,fp1)!=1)BE_WERROR;
}
fclose(fp1);fflush(fp1); /* close file and flush buffer */
- if(E->parallel.me < E->parallel.nproc-1){
- mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 0, MPI_COMM_WORLD);
- }
+ if(E->output.gzdir.vtk_io == 2)
+ if(E->parallel.me < E->parallel.nproc-1){
+ mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 0, MPI_COMM_WORLD);
+ }
}
return;
}
@@ -824,7 +866,7 @@
int mpi_rc;
int mpi_inmsg, mpi_success_message = 1;
- if(E->output.gzdir.vtk_io != 2){ /* old */
+ if(E->output.gzdir.vtk_io < 2){ /* old */
snprintf(output_file,255,"%s/%d/pressure.%d.%d.gz", E->control.data_dir,cycles,
E->parallel.me, cycles);
gz1 = gzdir_output_open(output_file,"w");
@@ -836,9 +878,10 @@
}
gzclose(gz1);
}else{/* new legacy VTK */
- parallel_process_sync(E);
- get_vtk_filename(output_file,E->output.gzdir.vtk_ocount,E);
- if(E->parallel.me == 0){
+ if(E->output.gzdir.vtk_io == 2)
+ parallel_process_sync(E);
+ get_vtk_filename(output_file,0,E,cycles);
+ if((E->parallel.me == 0) || (E->output.gzdir.vtk_io == 3)){
fp1 = output_open_mode(output_file,"a");
myfprintf(fp1,"SCALARS pressure float 1\n");
myfprintf(fp1,"LOOKUP_TABLE default\n");
@@ -846,15 +889,16 @@
mpi_rc = MPI_Recv(&mpi_inmsg, 1, MPI_INT, (E->parallel.me-1), 0, MPI_COMM_WORLD, &mpi_stat);
fp1 = output_open_mode(output_file,"a");
}
- for(j=1; j<= E->sphere.caps_per_proc;j++)
+ for(j=1; j<= E->sphere.caps_per_proc;j++)
for(i=1;i<=E->lmesh.nno;i++){
- ftmp = E->NP[j][i];
+ ftmp = E->NP[j][i];
if(be_write_float_to_file(&ftmp,1,fp1)!=1)BE_WERROR;
}
fclose(fp1);fflush(fp1); /* close file and flush buffer */
- if(E->parallel.me < E->parallel.nproc-1){
- mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 0, MPI_COMM_WORLD);
- }
+ if(E->output.gzdir.vtk_io == 2)
+ if(E->parallel.me < E->parallel.nproc-1){
+ mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 0, MPI_COMM_WORLD);
+ }
}
return;
}
@@ -867,7 +911,7 @@
char output_file[255];
gzFile *fp1;
- snprintf(output_file,255,"%s/%d/tracer.%d.%d.gz",
+ snprintf(output_file,255,"%s/%d/tracer.%d.%d.gz",
E->control.data_dir,cycles,
E->parallel.me, cycles);
fp1 = gzdir_output_open(output_file,"w");
@@ -902,7 +946,7 @@
void gzdir_output_comp_nd(struct All_variables *E, int cycles)
{
int i, j, k;
- char output_file[255];
+ char output_file[255],message[255];
gzFile *gz1;
FILE *fp1;
float ftmp;
@@ -910,53 +954,52 @@
MPI_Status mpi_stat;
int mpi_rc;
int mpi_inmsg, mpi_success_message = 1;
-
- if(E->output.gzdir.vtk_io != 2){
- snprintf(output_file,255,"%s/%d/comp_nd.%d.%d.gz",
+
+ if(E->output.gzdir.vtk_io < 2){
+ snprintf(output_file,255,"%s/%d/comp_nd.%d.%d.gz",
E->control.data_dir,cycles,
E->parallel.me, cycles);
gz1 = gzdir_output_open(output_file,"w");
for(j=1;j<=E->sphere.caps_per_proc;j++) {
- gzprintf(gz1,"%3d %7d %.5e %d\n",
+ gzprintf(gz1,"%3d %7d %.5e %.5e %.5e\n",
j, E->lmesh.nel,
- E->monitor.elapsed_time, E->composition.ncomp);
- for(i=0;i<E->composition.ncomp;i++) {
- gzprintf(fp1," %.5e %.5e ",
- E->composition.initial_bulk_composition[i],
- E->composition.bulk_composition[i]);
- }
- gzprintf(fp1,"\n");
-
+ E->monitor.elapsed_time,
+ E->composition.initial_bulk_composition,
+ E->composition.bulk_composition);
for(i=1;i<=E->lmesh.nno;i++) {
- for(k=0;k<E->composition.ncomp;k++) {
+ for(k=0;k<E->composition.ncomp;k++)
gzprintf(gz1,"%.6e ",E->composition.comp_node[j][k][i]);
- }
- gzprintf(fp1,"\n");
+ gzprintf(gz1,"\n");
}
}
gzclose(gz1);
}else{/* new legacy VTK */
- parallel_process_sync(E);
- get_vtk_filename(output_file,E->output.gzdir.vtk_ocount,E);
- if(E->parallel.me == 0){
+ if(E->output.gzdir.vtk_io == 2)
+ parallel_process_sync(E);
+ get_vtk_filename(output_file,0,E,cycles);
+ if((E->output.gzdir.vtk_io==3) || (E->parallel.me == 0)){
fp1 = output_open_mode(output_file,"a");
- myfprintf(fp1,"SCALARS composition float %d\n",E->composition.ncomp);
+ if(E->composition.ncomp > 4)
+ myerror(E,"vtk out error: ncomp out of bounds (needs to be < 4)");
+ sprintf(message,"SCALARS composition float %d\n",E->composition.ncomp);
+ myfprintf(fp1,message);
myfprintf(fp1,"LOOKUP_TABLE default\n");
- }else{
+ }else{ /* serial wait */
mpi_rc = MPI_Recv(&mpi_inmsg, 1, MPI_INT, (E->parallel.me-1), 0, MPI_COMM_WORLD, &mpi_stat);
fp1 = output_open_mode(output_file,"a");
}
- for(j=1; j<= E->sphere.caps_per_proc;j++)
- for(k=0;k<E->composition.ncomp;k++) {
+ for(j=1; j<= E->sphere.caps_per_proc;j++)
for(i=1;i<=E->lmesh.nno;i++){
- ftmp = E->composition.comp_node[j][k][i];
- if(be_write_float_to_file(&ftmp,1,fp1)!=1)BE_WERROR;
+ for(k=0;k<E->composition.ncomp;k++){
+ ftmp = E->composition.comp_node[j][k][i];
+ if(be_write_float_to_file(&ftmp,1,fp1)!=1)BE_WERROR;
+ }
}
+ fclose(fp1);fflush(fp1); /* close file and flush buffer */
+ if(E->output.gzdir.vtk_io == 2) /* serial */
+ if(E->parallel.me < E->parallel.nproc-1){
+ mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 0, MPI_COMM_WORLD);
}
- fclose(fp1);fflush(fp1); /* close file and flush buffer */
- if(E->parallel.me < E->parallel.nproc-1){
- mpi_rc = MPI_Send(&mpi_success_message, 1, MPI_INT, (E->parallel.me+1), 0, MPI_COMM_WORLD);
- }
}
return;
}
@@ -964,7 +1007,7 @@
void gzdir_output_comp_el(struct All_variables *E, int cycles)
{
- int i, j;
+ int i, j, k;
char output_file[255];
gzFile *fp1;
@@ -973,31 +1016,24 @@
fp1 = gzdir_output_open(output_file,"w");
for(j=1;j<=E->sphere.caps_per_proc;j++) {
- gzprintf(fp1,"%3d %7d %.5e %d\n",
- j, E->lmesh.nel,
- E->monitor.elapsed_time, E->composition.ncomp);
+ gzprintf(fp1,"%3d %7d %.5e %.5e %.5e\n",
+ j, E->lmesh.nel,
+ E->monitor.elapsed_time,
+ E->composition.initial_bulk_composition,
+ E->composition.bulk_composition);
- for(i=0;i<E->composition.ncomp;i++) {
- gzprintf(fp1,"%.5e %.5e ",
- E->composition.initial_bulk_composition[i],
- E->composition.bulk_composition[i]);
+ for(i=1;i<=E->lmesh.nel;i++) {
+ for(k=0;k<E->composition.ncomp;k++)
+ gzprintf(fp1,"%.6e ",E->composition.comp_el[j][k][i]);
+ gzprintf(fp1,"\n");
}
- fprintf(fp1,"\n");
-
- for(i=1;i<=E->lmesh.nno;i++) {
- for(k=0;k<E->composition.ncomp;k++) {
- gzprintf(fp1,"%.6e ",E->composition.comp_node[j][k][i]);
- }
- gzprintf(fp1,"\n");
- }
-
}
gzclose(fp1);
return;
}
-/*
+/*
restart facility for zipped/VTK style , will init temperature
@@ -1013,34 +1049,43 @@
float v1, v2, v3, g;
ii = E->monitor.solution_cycles_init;
-
- if(E->output.gzdir.vtk_io) { /* VTK I/O */
- snprintf(output_file,255,"%s/%d/t.%d.%d",
- E->control.data_dir_old,
- ii,E->parallel.me,ii);
- }else{ /* old output */
+ switch(E->output.gzdir.vtk_io){
+ case 2:
+ case 3:
+ myerror(E,"sorry, restart with vtk_io 2 or 3 not implemented yet");
+ break;
+ case 1:
+ /* VTK I/O */
+ snprintf(output_file,255,"%s/%d/t.%d.%d",
+ E->control.data_dir_old,
+ ii,E->parallel.me,ii);
+ break;
+ default:
snprintf(output_file,255,"%s/%d/velo.%d.%d",
E->control.data_dir_old,ii,
E->parallel.me,ii);
+ break;
}
/* open file */
rezip = open_file_zipped(output_file,&fp,E);
if (E->parallel.me==0)
fprintf(E->fp,"restart_tic_from_gzdir_file: using %s for restarted temperature\n",
output_file);
-
+
if(fscanf(fp,"%i %i %f",&ll,&mm,&restart_elapsed_time) != 3)
myerror(E,"restart vtkl read error 0");
-
- if(E->output.gzdir.vtk_io) { /* VTK */
+
+ switch(E->output.gzdir.vtk_io) {
+ case 1: /* VTK */
for(m=1;m <= E->sphere.caps_per_proc;m++) {
if(fscanf(fp,"%i %i",&ll,&mm) != 2)
myerror(E,"restart vtkl read error 1");
- for(i=1;i<=E->lmesh.nno;i++)
+ for(i=1;i<=E->lmesh.nno;i++)
if(fscanf(fp,"%f",&(E->T[m][i]))!=1)
myerror(E,"restart vtkl read error 2");
}
- }else{ /* old style velo */
+ break;
+ default: /* old style velo */
for(m=1;m <= E->sphere.caps_per_proc;m++) {
fscanf(fp,"%i %i",&ll,&mm);
for(i=1;i<=E->lmesh.nno;i++) {
@@ -1052,18 +1097,19 @@
E->T[m][i] = g;
}
}
+ break;
}
fclose (fp);
if(rezip) /* rezip */
gzip_file(output_file);
-
+
temperatures_conform_bcs(E);
return;
}
-/*
-
+/*
+
tries to open 'name'. if name exists, out will be pointer to file and
return 0. if name doesn't exist, will check for name.gz. if this
exists, will unzip and open, and return 1
@@ -1077,14 +1123,14 @@
char mstring[1000];
*in = fopen(name,"r");
if (*in == NULL) {
- /*
+ /*
unzipped file not found
*/
snprintf(mstring,1000,"%s.gz",name);
*in= fopen(mstring,"r");
- if(*in != NULL){
- /*
- zipped version was found
+ if(*in != NULL){
+ /*
+ zipped version was found
*/
fclose(*in);
snprintf(mstring,1000,"gunzip -f %s.gz",name); /* brutal */
@@ -1095,7 +1141,7 @@
myerror("open_file_zipped: unzipping error",E);
return 1;
}else{
- /*
+ /*
no file, either zipped or unzipped
*/
snprintf(mstring,1000,"no files %s and %s.gz were found, exiting",
@@ -1104,8 +1150,8 @@
return 0;
}
}else{
- /*
- file was found unzipped
+ /*
+ file was found unzipped
*/
return 0;
}
@@ -1116,21 +1162,39 @@
{
char command_string[300];
snprintf(command_string,300,"gzip -f %s",output_file); /* brutal */
- system(command_string);
+ system(command_string);
}
-void get_vtk_filename(char *output_file,int cycles,struct All_variables *E)
+void get_vtk_filename(char *output_file,
+ int geo,struct All_variables *E,
+ int cycles)
{
- snprintf(output_file,255,"%s/d.%08i.vtk",E->control.data_dir,cycles);
+ if(E->output.gzdir.vtk_io == 2){ /* serial */
+ if(geo) /* geometry */
+ sprintf(output_file,"%s/vtk_geo",
+ E->control.data_dir);
+ else /* data part */
+ sprintf(output_file,"%s/d.%08i.vtk",
+ E->control.data_dir, E->output.gzdir.vtk_ocount);
+ }else{ /* parallel */
+ if(geo) /* geometry */
+ sprintf(output_file,"%s/vtk_geo.%i",
+ E->control.data_dir,E->parallel.me);
+ else /* data part */
+ sprintf(output_file,"%s/%d/d.%08i.%i.vtk",
+ E->control.data_dir,cycles,
+ E->output.gzdir.vtk_ocount,
+ E->parallel.me);
+ }
}
-/*
+/*
big endian I/O (needed for vtk)
@@ -1138,7 +1202,7 @@
*/
-/*
+/*
write the x[n] array to file, making sure it is written big endian
@@ -1150,13 +1214,13 @@
size_t bsize;
float ftmp;
#ifdef ASCII_DEBUG
- for(i=0;i<n;i++)
- fprintf(out,"%11g ",x[i]);
+ for(i=0;i<n;i++)
+ fprintf(out,"%11g ",x[i]);
fprintf(out,"\n");
nout = n;
#else
- /*
- do we need to flip?
+ /*
+ do we need to flip?
*/
if(be_is_little_endian()){
nout = 0;
@@ -1178,13 +1242,13 @@
size_t bsize;
int itmp;
#ifdef ASCII_DEBUG
- for(i=0;i<n;i++)
- fprintf(out,"%11i ",x[i]);
+ for(i=0;i<n;i++)
+ fprintf(out,"%11i ",x[i]);
fprintf(out,"\n");
nout = n;
#else
- /*
- do we need to flip?
+ /*
+ do we need to flip?
*/
if(be_is_little_endian()){
nout = 0;
@@ -1217,14 +1281,14 @@
return *(const unsigned char *)&a;
}
-/*
+/*
flip endian-ness
*/
-/*
+/*
flip endianness of x
Modified: mc/3D/CitcomS/trunk/lib/Output_h5.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Output_h5.c 2007-08-18 03:14:06 UTC (rev 7845)
+++ mc/3D/CitcomS/trunk/lib/Output_h5.c 2007-08-19 23:54:46 UTC (rev 7846)
@@ -1463,7 +1463,7 @@
status = set_attribute_float(input, "gruneisen",
(E->control.inv_gruneisen == 0)?
1.0/E->control.inv_gruneisen :
- E->control.inv_gruneisen));
+ E->control.inv_gruneisen);
status = set_attribute_float(input, "Q0", E->control.Q0);
status = set_attribute_int(input, "stokes_flow_only", E->control.stokes);
Modified: mc/3D/CitcomS/trunk/lib/Parsing.c
===================================================================
--- mc/3D/CitcomS/trunk/lib/Parsing.c 2007-08-18 03:14:06 UTC (rev 7845)
+++ mc/3D/CitcomS/trunk/lib/Parsing.c 2007-08-19 23:54:46 UTC (rev 7846)
@@ -114,9 +114,8 @@
exit(11);
}
- unique_copy_file(E,filename,"copy");
+
-
/* now the parameter file is open, read into memory */
while( fgets(line,MAXLINE,fp) != NULL )
@@ -168,6 +167,12 @@
DESCRIBE=j;
BEGINNER=k;
+ /* make this an optional behavior */
+ input_boolean("copy_input_file",&k,"on",m);
+ if(k)
+ unique_copy_file(E,filename,"copy");
+
+
}
void shutdown_parser(E)
More information about the cig-commits
mailing list