Gnuplot equation of state fitting loop

Question

I am trying to make a script to fit 50 files of free energy vs volume using the Birch Murnaghan equation of state. I'm fitting for 3 parameters, emin, v0 and k0. I can't seem to get a loop to work which would loop over each file fitting it and then making 50 new files of the three fitting parameters.

   birch(x) = emin + 1.5*v0*k0/1602 * ( 0.75*(1+2*xsi)*(v0/x)**(4.0/3.0) \
-xsi/2*(v0/x)**2 - 3.0/2.0*(1+xsi)*(v0/x)**(2.0/3.0) + 0.5*(xsi+1.5) )
emin=-8
k0=2000
v0=10
xsi=0

do for [n=100:5000:100]{
   file = sprintf('free.%.f',n)
   fit birch(x) file via emin,v0,k0
   A = A.sprintf(%f,emin)
   B = B.sprintf(%f,v0)
   C = C.sprintf(%f,k0)
}

''' The data is saved in files going free.100 free.200 and so on till free.5000, does anybody know how I could do this I only really want to make a script so I don't have to manually fit all the data

user_7295599 · Answer

(While I was coding @maij already posted an answer. I'll nevertheless post my example)

You don't specify how the output files should look like. In order to write something to a file in gnuplot you can use set print , check help set print. Check the following example:

Code:

### store fitting parameters into file
reset session

myFile(n)  = sprintf("free.%d",n)
myParamFile(n) = sprintf("free.%d.param",n)
set key left

# create some test data in files
do for [i=100:600:100] {
    set print myFile(i)
        m = int(rand(0)*10)+0.1
        c = int(rand(0)*10)
        do for [j=1:20] {
            print sprintf("%g %g",j, m*j + 5*m*rand(0) + c)
        }
    set print
}

set fit brief nolog
f(x) = a*x + b

set multiplot layout 3,2
    do for [i=100:600:100] {
        fit f(x) myFile(i) via a,b
        set print myParamFile(i)
            print sprintf("a = %g",a)
            print sprintf("b = %g",b)
        set print
        plot myFile(i) w p pt 7 lc "blue" ti myFile(i), \
             f(x) w l lc "red" 
    }
unset multiplot
### end of code

Result:

free.100.param

a = 3.29395
b = 8.13496

free.200.param

a = 1.99086
b = 7.94427

free.300.param

a = 2.11392
b = 12.3286

free.400.param

a = 3.97004
b = 13.3108

free.500.param

a = 7.60895
b = 15.8168

free.600.param

a = 0.0982883
b = 8.2272

user_7010554 · Answer

I would use the set print command to specify the output file for the print command.

do for [n=100:5000:100]{
   input_file = sprintf('free.%.f',n)
   fit birch(x) input_file via emin,v0,k0

   # set output file for print command
   output_file = sprintf('free.%.f.result', n)
   set print output_file

   # write output, choose format as you need it
   print(input_file)
   output_line = sprintf('emin = %f', emin)
   print(output_line)
   print(v0)
   print(k0)
}

This will create 50 output files free.100.result, ... in the current directory. The output of the following print command goes to these files.

Update

For printing values of a function, you can use two approaches:


do for [n=100:5000:100]{

   ...
   fit ...
   ...

   # with 'set table' command into a separate file
   table_file = sprintf('free.%.f.table', n)
   set table table_file
   set xrange [6:7]
   set samples 11
   plot birch(x)
   unset table

   # print manually into normal print output file
   # please note, step size must be an integer
   do for [i = 60:70:1] { 
      x = i/10.0
      output_line = sprintf("%f  %f", x, birch(x))
      print(output_line)
   } 
}

Gnuplot equation of state fitting loop

2 Answers