From d8c1238ec3ded1886970fa6ecacc1e5d823394b5 Mon Sep 17 00:00:00 2001 From: "W. Trevor King" Date: Tue, 12 Feb 2013 16:54:47 -0500 Subject: [PATCH] error_analysis.c: Implement error analysis Currently it is hard-coded to read column five, but adding a command line switch to change columns wouldn't be too difficult. --- error_analysis.c | 57 +++++++++++++++++++++++++++++++++++++++++++----- 1 file changed, 52 insertions(+), 5 deletions(-) diff --git a/error_analysis.c b/error_analysis.c index c242304..924af91 100644 --- a/error_analysis.c +++ b/error_analysis.c @@ -1,7 +1,7 @@ /* -A simple "hello world" example in C. +Calculate the maximum relative energy error in a trajectory. -Copyright (C) 2012 W. Trevor King +Copyright (C) 2013 W. Trevor King 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 @@ -17,9 +17,56 @@ You should have received a copy of the GNU General Public License along with this program. If not, see . */ -#include +#include /* for *printf(), getc(), fscanf(), EOF, stdin, stderr */ +#include /* for fabs() */ -int main() { - printf("hello, world\n"); +/* Read a double column from the current line in a tab-delimited + * stream. + */ +int read_column(FILE *stream, int column, size_t line, double *value) +{ + int c, col = 0; + + while ((c = getc(stream)) != EOF) { + if (c == '\t') { /* found a new column */ + col++; + if (col == column) { /* read our value */ + c = fscanf(stream, "%lf", value); + if (c != 1) { + fprintf(stderr, "could not parse column %d in line %ld\n", + column, line); + return -1; + } + while ((c = getc(stream)) != EOF) + if (c == '\n') + return 0; /* advanced to next line */ + } + } else if (c == '\n') { + fprintf(stderr, "only %d columns in line %ld\n", column, line); + return -1; + } + } + return -1; /* reached the end of the file */ +} + +int main(int argc, char **argv) { + int line, column = 5; + double val_i, val, error, max_error = 0; + + for (line=0; read_column(stdin, column, line, &val) == 0; line++) { + if (line == 0) { + val_i = val; + if (val_i == 0) { + fprintf(stderr, + "cannot calculate relative error when initial value is zero\n"); + return -1; + } + } else { + error = fabs((val - val_i)/(val_i)); + if (error > max_error) + max_error = error; + } + } + printf("%g\n", max_error); return 0; } -- 2.26.2