Revision 322f2b874d2ed979ad50cfd281b6b4b0cbd82918 authored by Rajesh Kommu on 25 September 2014, 21:09:29 UTC, committed by Rajesh Kommu on 25 September 2014, 21:09:29 UTC
1 parent ff03ffe
Checkpoints.c
/*
*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
*
*<LicenseText>
*
* CitcomS by Louis Moresi, Shijie Zhong, Lijie Han, Eh Tan,
* Clint Conrad, Michael Gurnis, and Eun-seo Choi.
* Copyright (C) 1994-2005, California Institute of Technology.
*
* This program is free software; you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation; either version 2 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program; if not, write to the Free Software
* Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
*
*</LicenseText>
*
*~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
*/
#include <sys/file.h>
#include <unistd.h>
#include "global_defs.h"
#include "composition_related.h"
/* Private function prototypes */
static void backup_file(const char *output_file);
static void write_sentinel(FILE *fp);
static void read_sentinel(FILE *fp, int me);
static void general_checkpoint(struct All_variables *E, FILE *fp);
static void tracer_checkpoint(struct All_variables *E, FILE *fp);
static void composition_checkpoint(struct All_variables *E, FILE *fp);
static void energy_checkpoint(struct All_variables *E, FILE *fp);
static void momentum_checkpoint(struct All_variables *E, FILE *fp);
static void read_general_checkpoint(struct All_variables *E, FILE *fp);
static void read_tracer_checkpoint(struct All_variables *E, FILE *fp);
static void read_composition_checkpoint(struct All_variables *E, FILE *fp);
static void read_energy_checkpoint(struct All_variables *E, FILE *fp);
static void read_momentum_checkpoint(struct All_variables *E, FILE *fp);
void myerror(struct All_variables *, char *);
void output_checkpoint(struct All_variables *E)
{
char output_file[255];
FILE *fp1;
sprintf(output_file, "%s.chkpt.%d.%d", E->control.data_file,
E->parallel.me, E->monitor.solution_cycles);
/* Disable the backup since the filename is unique. */
/* backup_file(output_file); */
fp1 = fopen(output_file, "wb");
/* checkpoint for general information */
/* this must be the first to be checkpointed */
general_checkpoint(E, fp1);
/* checkpoint for energy equation */
energy_checkpoint(E, fp1);
/* checkpoint for momentum equation */
momentum_checkpoint(E, fp1);
/* checkpoint for tracer/composition */
if(E->control.tracer) {
tracer_checkpoint(E, fp1);
if(E->composition.on)
composition_checkpoint(E, fp1);
}
fclose(fp1);
return;
}
void read_checkpoint(struct All_variables *E)
{
void initialize_material(struct All_variables *E);
void initial_viscosity(struct All_variables *E);
char output_file[255];
FILE *fp;
/* open the checkpoint file */
snprintf(output_file, 254, "%s.chkpt.%d.%d", E->control.old_P_file,
E->parallel.me, E->monitor.solution_cycles_init);
fp = fopen(output_file, "rb");
if(fp == NULL) {
fprintf(stderr, "Cannot open file: %s\n", output_file);
exit(-1);
}
if(E->parallel.me == 0)
fprintf(stderr,"read_checkpoint: restarting from %s\n",output_file);
/* check mesh information in the checkpoint file */
read_general_checkpoint(E, fp);
/* init E->mat */
initialize_material(E);
/* read energy information in the checkpoint file */
read_energy_checkpoint(E, fp);
/* read momentum information in the checkpoint file */
read_momentum_checkpoint(E, fp);
/* read tracer/composition information in the checkpoint file */
if(E->control.tracer) {
if(E->trace.ic_method_for_flavors == 99){
if(E->parallel.me == 0)
fprintf(stderr,"ic_method_for_flavors = 99 will override checkpoint restart\n");
}else{
read_tracer_checkpoint(E, fp);
if(E->composition.on)
read_composition_checkpoint(E, fp);
}
}
fclose(fp);
/* finally, init viscosity */
initial_viscosity(E);
return;
}
static void backup_file(const char *output_file)
{
char bak_file[255];
int ierr;
/* check the existence of output_file */
if(access(output_file, F_OK) == 0) {
/* if exist, renamed it to back up */
sprintf(bak_file, "%s.bak", output_file);
ierr = rename(output_file, bak_file);
if(ierr != 0) {
fprintf(stderr, "Warning, cannot backup checkpoint files\n");
}
}
return;
}
static void write_sentinel(FILE *fp)
{
int a[4] = {0, 0, 0, 0};
fwrite(a, sizeof(int), 4, fp);
}
static void read_sentinel(FILE *fp, int me)
{
int i, a[4];
int nonzero = 0;
fread(a, sizeof(int), 4, fp);
/* check whether a[i] are all zero */
for(i=0; i<4; i++)
nonzero |= a[i];
if(nonzero) {
fprintf(stderr, "Error in reading checkpoint file: wrong sentinel, "
"me=%d\n", me);
exit(-1);
}
return;
}
static void general_checkpoint(struct All_variables *E, FILE *fp)
{
/* write mesh information */
fwrite(&(E->lmesh.nox), sizeof(int), 1, fp);
fwrite(&(E->lmesh.noy), sizeof(int), 1, fp);
fwrite(&(E->lmesh.noz), sizeof(int), 1, fp);
fwrite(&(E->parallel.nprocx), sizeof(int), 1, fp);
fwrite(&(E->parallel.nprocy), sizeof(int), 1, fp);
fwrite(&(E->parallel.nprocz), sizeof(int), 1, fp);
fwrite(&(E->sphere.caps_per_proc), sizeof(int), 1, fp);
/* write timing information */
fwrite(&(E->monitor.solution_cycles), sizeof(int), 1, fp);
fwrite(&(E->monitor.elapsed_time), sizeof(float), 1, fp);
fwrite(&(E->advection.timestep), sizeof(float), 1, fp);
fwrite(&(E->control.start_age), sizeof(float), 1, fp);
return;
}
static void read_general_checkpoint(struct All_variables *E, FILE *fp)
{
int tmp[7];
double dtmp;
/* read mesh information */
fread(tmp, sizeof(int), 7, fp);
if((tmp[0] != E->lmesh.nox) ||
(tmp[1] != E->lmesh.noy) ||
(tmp[2] != E->lmesh.noz) ||
(tmp[3] != E->parallel.nprocx) ||
(tmp[4] != E->parallel.nprocy) ||
(tmp[5] != E->parallel.nprocz) ||
(tmp[6] != E->sphere.caps_per_proc)) {
fprintf(stderr, "Error in reading checkpoint file: mesh parameters mismatch, me=%d\n",
E->parallel.me);
fprintf(stderr, "%d %d %d %d %d %d %d\n",
tmp[0], tmp[1], tmp[2], tmp[3],
tmp[4], tmp[5], tmp[6]);
exit(-1);
}
/* read timing information */
tmp[0] = fread(&(E->monitor.solution_cycles), sizeof(int), 1, fp);
tmp[0]+= fread(&(E->monitor.elapsed_time), sizeof(float), 1, fp);
tmp[0]+= fread(&(E->advection.timestep), sizeof(float), 1, fp);
tmp[0]+= fread(&(E->control.start_age), sizeof(float), 1, fp);
if(tmp[0] != 4)
myerror(E,"read_general_checkpoint: header error");
E->advection.timesteps = E->monitor.solution_cycles;
return;
}
static void tracer_checkpoint(struct All_variables *E, FILE *fp)
{
int m, i;
write_sentinel(fp);
fwrite(&(E->trace.number_of_basic_quantities), sizeof(int), 1, fp);
fwrite(&(E->trace.number_of_extra_quantities), sizeof(int), 1, fp);
fwrite(&(E->trace.nflavors), sizeof(int), 1, fp);
fwrite(&(E->trace.ilast_tracer_count), sizeof(int), 1, fp);
fwrite(&(E->trace.ntracers), sizeof(int), 1, fp);
/* the 0-th element of basicq/extraq/ielement is not init'd
* and won't be used when read it. */
for(i=0; i<6; i++) {
fwrite(E->trace.basicq[i], sizeof(double),
E->trace.ntracers+1, fp);
}
for(i=0; i<E->trace.number_of_extra_quantities; i++) {
fwrite(E->trace.extraq[i], sizeof(double),
E->trace.ntracers+1, fp);
}
fwrite(E->trace.ielement, sizeof(int),
E->trace.ntracers+1, fp);
}
static void read_tracer_checkpoint(struct All_variables *E, FILE *fp)
{
void count_tracers_of_flavors(struct All_variables *E);
void allocate_tracer_arrays();
int m, i, itmp;
read_sentinel(fp, E->parallel.me);
fread(&itmp, sizeof(int), 1, fp);
if (itmp != E->trace.number_of_basic_quantities) {
fprintf(stderr, "Error in reading checkpoint file: tracer basicq, me=%d\n",
E->parallel.me);
fprintf(stderr, "%d\n", itmp);
exit(-1);
}
fread(&itmp, sizeof(int), 1, fp);
if (itmp != E->trace.number_of_extra_quantities) {
fprintf(stderr, "Error in reading checkpoint file: tracer extraq, me=%d\n",
E->parallel.me);
fprintf(stderr, "%d\n", itmp);
exit(-1);
}
fread(&itmp, sizeof(int), 1, fp);
if (itmp != E->trace.nflavors) {
fprintf(stderr, "Error in reading checkpoint file: tracer nflavors, me=%d\n",
E->parallel.me);
fprintf(stderr, "%d\n", itmp);
exit(-1);
}
fread(&itmp, sizeof(int), 1, fp);
E->trace.ilast_tracer_count = itmp;
/* # of tracers, allocate memory */
fread(&itmp, sizeof(int), 1, fp);
allocate_tracer_arrays(E, itmp);
E->trace.ntracers = itmp;
/* read tracer data */
for(i=0; i<6; i++) {
fread(E->trace.basicq[i], sizeof(double),
E->trace.ntracers+1, fp);
}
for(i=0; i<E->trace.number_of_extra_quantities; i++) {
fread(E->trace.extraq[i], sizeof(double),
E->trace.ntracers+1, fp);
}
fread(E->trace.ielement, sizeof(int),
E->trace.ntracers+1, fp);
/* init E->trace.ntracer_flavor */
count_tracers_of_flavors(E);
}
static void composition_checkpoint(struct All_variables *E, FILE *fp)
{
int i, m;
write_sentinel(fp);
fwrite(&(E->composition.ncomp), sizeof(int), 1, fp);
fwrite(E->composition.bulk_composition, sizeof(double),
E->composition.ncomp, fp);
fwrite(E->composition.initial_bulk_composition, sizeof(double),
E->composition.ncomp, fp);
/* the 0-th element of comp_el is not init'd
* and won't be used when read it. */
for(i=0; i<E->composition.ncomp; i++)
fwrite(E->composition.comp_el[i], sizeof(double),
E->lmesh.nel+1, fp);
}
static void read_composition_checkpoint(struct All_variables *E, FILE *fp)
{
double tmp;
int m, i, itmp;
read_sentinel(fp, E->parallel.me);
fread(&itmp, sizeof(int), 1, fp);
if (itmp != E->composition.ncomp) {
fprintf(stderr, "Error in reading checkpoint file: ncomp, me=%d\n",
E->parallel.me);
fprintf(stderr, "%d\n", itmp);
exit(-1);
}
fread(E->composition.bulk_composition, sizeof(double),
E->composition.ncomp, fp);
fread(E->composition.initial_bulk_composition, sizeof(double),
E->composition.ncomp, fp);
for(i=0; i<E->composition.ncomp; i++)
fread(E->composition.comp_el[i], sizeof(double),
E->lmesh.nel+1, fp);
/* init E->composition.comp_node */
map_composition_to_nodes(E);
/* preventing uninitialized access */
E->trace.istat_iempty = 0;
for (i=0; i<E->composition.ncomp; i++) {
E->composition.error_fraction[i] = E->composition.bulk_composition[i]
/ E->composition.initial_bulk_composition[i] - 1.0;
}
}
static void energy_checkpoint(struct All_variables *E, FILE *fp)
{
int m;
write_sentinel(fp);
fwrite(E->T, sizeof(double), E->lmesh.nno+1, fp);
fwrite(E->Tdot, sizeof(double), E->lmesh.nno+1, fp);
}
static void read_energy_checkpoint(struct All_variables *E, FILE *fp)
{
int m;
read_sentinel(fp, E->parallel.me);
/* the 0-th element of T/Tdot is not init'd
* and won't be used when read it. */
if(fread(E->T, sizeof(double), E->lmesh.nno+1, fp)!= E->lmesh.nno+1)
myerror(E,"read_energy_checkpoint: error at T");
if(fread(E->Tdot, sizeof(double), E->lmesh.nno+1, fp)!=E->lmesh.nno+1)
myerror(E,"read_energy_checkpoint: error at Tdot");
}
static void momentum_checkpoint(struct All_variables *E, FILE *fp)
{
int m;
float junk[2];
junk[0] = junk[1] = 0;
write_sentinel(fp);
/* for backward compatibility */
fwrite(junk, sizeof(float), 2, fp);
/* the 0-th element of P/NP/EVI/VI is not init'd
* and won't be used when read it. */
/* Pressure at equation points */
fwrite(E->P, sizeof(double), E->lmesh.npno, fp);
/* velocity at equation points */
fwrite(E->U, sizeof(double), E->lmesh.neq, fp);
}
static void read_momentum_checkpoint(struct All_variables *E, FILE *fp)
{
void v_from_vector();
void p_to_nodes();
double global_v_norm2(), global_p_norm2();
int m;
int lev = E->mesh.levmax;
float junk[2];
read_sentinel(fp, E->parallel.me);
/* for backward compatibility */
if(fread(junk, sizeof(float), 2, fp)!=2)
myerror(E,"read_momentum_checkpoint: error at vdotv");
/* Pressure at equation points */
if(fread(E->P, sizeof(double), E->lmesh.npno, fp) != E->lmesh.npno+1)
myerror(E,"read_momentum_checkpoint: error at P");
/* velocity at equation points */
if(fread(E->U, sizeof(double), E->lmesh.neq, fp) != E->lmesh.neq)
myerror(E,"read_momentum_checkpoint: error at U");
E->monitor.vdotv = global_v_norm2(E, E->U);
E->monitor.pdotp = global_p_norm2(E, E->P);
/* update velocity array */
v_from_vector(E);
/* init E->NP */
p_to_nodes(E, E->P, E->NP, lev);
}
Computing file changes ...