From: hackbard Date: Mon, 18 Feb 2008 18:42:42 +0000 (+0100) Subject: added diffusion code X-Git-Url: https://hackdaworld.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=30fafe0cb366808315649bd07bfaa7133ca0ecd2;p=physik%2Fposic.git added diffusion code --- diff --git a/Makefile b/Makefile index b25487e..1b875fe 100644 --- a/Makefile +++ b/Makefile @@ -4,7 +4,7 @@ CFLAGS = -Wall #CFLAGS += -Wextra -pedantic CFLAGS += -O3 CFLAGS += -g -#CFLAGS += -ffloat-store +CFLAGS += -ffloat-store CFLAGS += -DALBE #CFLAGS += -DTERSOFF_ORIG @@ -22,7 +22,7 @@ DEPS = moldyn.o random/random.o list/list.o DEPS += potentials/lennard_jones.o potentials/harmonic_oscillator.o DEPS += potentials/tersoff.o potentials/albe.o -all: posic sic fluctuation_calc postproc pair_correlation_calc +all: posic sic fluctuation_calc postproc pair_correlation_calc diffusion_calc posic: $(DEPS) @@ -32,6 +32,8 @@ postproc: $(DEPS) pair_correlation_calc: $(DEPS) +diffusion_calc: $(DEPS) + .PHONY:clean clean: rm -vf sic postproc fluctuation_calc pair_correlation_calc *.o */*.o diff --git a/diff_calc_script b/diff_calc_script new file mode 100755 index 0000000..272c7bd --- /dev/null +++ b/diff_calc_script @@ -0,0 +1,74 @@ +#!/bin/sh + +# +# calculate diffusion coefficient within a certain time period +# frank.zirkelbach@physik.uni-augsburg.de +# + +do_it() { + echo "$3 `./diffusion_calc $1 | grep coefficients | \ + awk '{ print $3 " " $4 " " $5 }'`" >> $2 +} + +if [ -d $1 ]; then + rm -f $1/diff_coeff.txt + echo "processing $1 ..." + for file in $1/*.save; do + time="`basename $file | awk -F- '{ print $2 }' | \ + sed 's/\.save//'`" + if [ ! -z "$3" ]; then + [ "$time" -gt "$4" -o "$time" -lt "$3" ] && continue + fi + do_it $file $1/diff_coeff.txt $time + done +else + echo "not a valid directory -> $1" + exit +fi + +# gnuplot + +ddir=$1 +dfile=$ddir/diff_coeff.scr + +cat > $dfile <<-EOF +set autoscale +unset log +unset label +set xtic auto +set ytic auto +set title 'Diffusion coefficients' +set xlabel 't [fs]' +set ylabel 'D [cm cm / s]' +set terminal postscript eps enhanced color solid lw 1 'Helvetica' 14 +set output '$ddir/diff_coeff.eps' +EOF + +echo -en "plot " >> $dfile + +komma=0 + +if [ ! -z `echo $2 | grep a` ]; then + [ "$komma" = "1" ] && + echo -en ", " >> $dfile + echo -en "\"$1/diff_coeff.txt\" u 1:2 w l t \"a\"" >> $dfile + komma=1 +fi + +if [ ! -z `echo $2 | grep b` ]; then + [ "$komma" = "1" ] && + echo -en ", " >> $dfile + echo -en "\"$1/diff_coeff.txt\" u 1:3 w l t \"b\"" >> $dfile + komma=1 +fi + +if [ ! -z `echo $2 | grep t` ]; then + [ "$komma" = "1" ] && + echo -en ", " >> $dfile + echo -en "\"$1/diff_coeff.txt\" u 1:4 w l t \"t\"" >> $dfile +fi + +echo -en "\n" >> $dfile + +gnuplot $dfile + diff --git a/diffusion_calc.c b/diffusion_calc.c new file mode 100644 index 0000000..576064a --- /dev/null +++ b/diffusion_calc.c @@ -0,0 +1,58 @@ +/* + * calculation of diffusion coefficient + * + * author: frank.zirkelbach@physik.uni-augsburg.de + * + */ + +#define _GNU_SOURCE +#include +#include +#include +#include +#include +#include +#include + +#include "moldyn.h" + +int usage(char *prog) { + + printf("\nusage:\n"); + printf(" %s \n\n",prog); + + return -1; +} + +int main(int argc,char **argv) { + + t_moldyn moldyn; + int ret; + double dc[3]; + + if(argc!=2) { + usage(argv[0]); + return -1; + } + + memset(&moldyn,0,sizeof(t_moldyn)); + + printf("[diffusion calc] reading save file ...\n"); + ret=moldyn_read_save_file(&moldyn,argv[1]); + if(ret) { + printf("[diffusion calc] exit!\n"); + return ret; + } + + calculate_diffusion_coefficient(&moldyn,dc); + + dc[0]*=(100*100*SECOND/(METER*METER)); + dc[1]*=(100*100*SECOND/(METER*METER)); + dc[2]*=(100*100*SECOND/(METER*METER)); + + printf("diffusion coefficients: %.10f %.10f %.10f [cm^2/s]\n", + dc[0],dc[1],dc[2]); + + return 0; +} +